Search code examples
pythonfloating-pointbigintegerexpbessel-functions

Handling large numbers and precision of Ricean Fading PDF python


I'm trying to calculate Ricean Fading PDF using following equation. RIcean Fading PDF. where 'y' is normalized envelope and 'gamma' is SNR

if the K value is large, then

math.exp(-((1.+_gamma)*pow(_y,2.) + _gamma))

exp results in big floating point (e.q. 1.01e-5088). in python it will shows '0.0' as value

mpmath.besseli(0,2. * _y * np.sqrt(_gamma * (1. + _gamma)))

the value of Bessel Function shows big int value (e.q. 7.78e+5092). in python it will shows '**inf**' value

How can i store big integer and floating point value in python and calculate the pdf?

def rice_pdf(self, _y, _gamma):
   return 2. * _y * (1. + _gamma) * math.exp(-((1.+_gamma)*pow(_y,2.) + _gamma)) * special.i0(2. * _y * np.sqrt(_gamma * (1. + _gamma))) 

Thanks.


Solution

  • If you have a way to compute the logarithm of the bessel function, you can avoid the multiplication of very large and very small numbers by transforming into a summation with subsequent exonentiation, which should solve the numerical problems (exploit the fact that exp(a) * exp(b) == exp(a + b)).

    def rice_pdf(_y, _gamma):
        a = np.log(2. * _y * (1. + _gamma))
        b = -((1.+_gamma)*pow(_y,2.) + _gamma)
        c = lni(2. * _y * np.sqrt(_gamma * (1. + _gamma))) 
        return np.exp(a + b + c)
    

    This function assumes there exist an implementation of lni that computes log(i0(z)). However, I am aware of no existing implementation of such a function. You can work around this by using mpmath for intermediate results:

    def lni(z):
        i0 = mpmath.besseli(0, z)  # may become a big number
        logi0 = mpmath.log(i0)  # logarithm brings big number into sensible range
        return float(logi0)  # convert back to normal floating point