I have come across the Levy's theorem which finds the difference in distributions as follows:
$$\frac{F(x+h) - F(x-h)}{2h} = \frac{1}{2\pi}\int_{-\infty}^{\infty} \frac{\sin(ht)}{ht}\exp(-itx)\psi_X(t)\mathrm{d}t,$$ where $F(\cdot)$ denotes the cumulative distribution function of random variable $X$ and $\psi_X(\cdot)$ denotes the characteristic function of $X$.
I am unable to understand how this theorem allows us to calculate the distribution given the characteristic function, as it actually gives the difference of distribution functions.
Secondly, if $X$ is a positive random variable, how Levy's theorem is more suitable than the Gil-Pelaez theorem which is as follows?
$$F_X(x) = \frac{1}{2} - \frac{1}{\pi}\int_{0}^{\infty}\frac{\text{Im}(\exp(-itx)\psi_X(t))}{t}\mathrm{d}t,$$ where $\text{Im}(\cdot)$ denotes the imaginary part.
As per this paper, Levy's theorem allows us to compute the distribution function numerically, which is much more efficient than the calculation of the Gil-Pelaez inversion formula. According to the paper, the distribution function can be approximated as
$$F(x) \sim \frac{\lambda x}{\pi} + \frac{2}{\pi}\sum\limits_{\nu = 1}^{N}\frac{\sin \lambda \nu x}{\nu}C(\nu/N)Re(\psi(\lambda\nu)),$$ where one needs to choose $\lambda$ and $N$ carefully as per the guidelines in the paper to get accurate results. $C(\nu/N)$ takes care of the approximation error.