I assume we know how to simulate a Poisson point process with constant intensity in an interval (e.g. by considering partial sums of i.i.d. exponential variables.)
That allows you to simulate a standard Poisson point process in a rectangle $[a,b] \times [0,d]$ by simulating an intensity $d$ Poisson process in $[a,b]$ and assigning the resulting points $\{x_i\}_{i=1}^T$ uniform heights $\{y_i\}_{i=1}^T$ in $[0,d]$ to obtain $\{(x_i,y_i)\}_{i=1}^T$.
To simulate (with a small error $<\epsilon$ in total variation) a Poisson point process with slowly decreasing integrable intensity $f(x)dx$ in $[0,\infty)$, first pick an increasing sequence $\{N_k\}_{k \ge 0}$ with $N_0=0$, e.g. $N_k=2^k-1$. Then pick $m$ large enough so that $\int_{N_m}^\infty f(x) \, dx <\epsilon$, so we can neglect points beyond $N_m$. Then for each $k=0,1,\ldots m-1$ independently simulate a standard Poisson process in a rectangle $[N_k,N_{k+1}] \times [0,f(N_k)]$. For each $k$ this yields points $(x_i(k),y_i(k))$
for $1 \le i \le T_k$. Reject those points for which $y_i(k)>f(x_i(k))$.
Note that with $f$ decreasing polynomially and $N_k$ increasing geometrically, the rejection probability is bounded away from 1. The $x$ coordinates of the surviving points will yield a sample from a Poisson process with intensity $f(x)\, dx $ in $[0,N_m]$.
