[Speex-dev] A technical question about the speex preprocessor.

[Jean-Marc Valin]

It's been a while since I did the maths. M(-.5;1;-x) optimises something
else, though it's not too far. I think it may be [M(-.25;1;-x)]^2 (or is
it the sqrt?) that is supposed to be there. In any case, if there's a
mismatch between the doc and the code, the code is the one likely to be
correct.

	Jean-Marc

John Ridges a écrit :
> I got the approximation from a Google book:
> 
> http://books.google.com/books?id=2CAqsF-RebgC&pg=PA385
> 
> Page 392, formula (10.33)
> 
> Using this formula, you're right, hypergeom_gain() would *not* converge
> to 1 for large x, but would instead be gamma(1.25)/sqrt(sqrt(x)) which
> would approach zero. Now if the formula for the hypergeometric gain were
> instead gamma(1.5) * M(-.5;1;-x) / sqrt(x) that *would* approach 1, but
> that's just me noodling around with the formula to get something that
> approaches 1. Since I don't know how the hypergoemetric gain was derived
> (or even really what it means) I don't know if that's useful or not. Can
> you tell me what the source was for the original table values?
> 
> John Ridges
> 
> 
> Jean-Marc Valin wrote:
>> Something looks odd without your values (or the doc) because
>> hypergeom_gain()
>> should really approach 1 as x goes to infinity. But in the end, an
>> approximation is probably OK because denoising is anything but an
>> exact science
>> :-)
>>
>>    Jean-Marc
>>
>> Quoting John Ridges <jridges at masque.com>:
>>
>>  
>>> By my reckoning the confluent hypergoemetric functions should have the
>>> following values:
>>>
>>> M(-.25;1;-.5) = 1.11433
>>> M(-.25;1;-1) = 1.21088
>>> M(-.25;1;-1.5) = 1.29385
>>> M(-.25;1;-2) = 1.36627
>>> M(-.25;1;-2.5) = 1.43038
>>> M(-.25;1;-3) = 1.48784
>>> M(-.25;1;-3.5) = 1.53988
>>> M(-.25;1;-4) = 1.58747
>>> M(-.25;1;-4.5) = 1.63134
>>> M(-.25;1;-5) = 1.67206
>>> M(-.25;1;-5.5) = 1.71009
>>> M(-.25;1;-6) = 1.74579
>>> M(-.25;1;-6.5) = 1.77947
>>> M(-.25;1;-7) = 1.81136
>>> M(-.25;1;-7.5) = 1.84167
>>> M(-.25;1;-8) = 1.87056
>>> M(-.25;1;-8.5) = 1.89818
>>> M(-.25;1;-9) = 1.92466
>>> M(-.25;1;-9.5) = 1.95009
>>> M(-.25;1;-10) = 1.97456
>>>
>>> Which would give table values of:
>>>
>>>    static const float table[21] = {
>>>       0.82157f, 0.91549f, 0.99482f, 1.06298f, 1.12248f, 1.17515f,
>>> 1.22235f,
>>>       1.26511f, 1.30421f, 1.34025f, 1.37371f, 1.40495f, 1.43428f,
>>> 1.46195f,
>>>       1.48815f, 1.51305f, 1.53679f, 1.55948f, 1.58123f, 1.60212f,
>>> 1.62223f};
>>>
>>> There is also a formula which asymptotically approaches M(a;b;-x) for
>>> high values of x that is:
>>>
>>> (gamma(b)/gamma(b-a))*(x^-a)
>>>
>>> Which for M(-.25;1;-x) is sqrt(sqrt(x))*1.10326
>>>
>>> at x=10 this gives a value of 1.96191 (vs. what I think is the true
>>> value of 1.97456).
>>>
>>> The reason I've gotten into all this is I'm trying to vectorize with SSE
>>> intrinsics some of the slower loops in the preprocessor, and the
>>> hypergeom_gain function is the only thing stopping me from removing all
>>> the branches in the loops. I don't know how critical the accuracy of the
>>> function is to the performance of the preprocessor, but if that
>>> aforementioned approximation was good enough for all the values of x, it
>>> would really speed the loops up.
>>>
>>> Let me know what you think. I hope I'm helping out here (and not just
>>> confusing things).
>>>
>>> John Ridges
>>>
>>>
>>> Jean-Marc Valin wrote:
>>>    
>>>> OK, so the problem is that the table you see does not match the
>>>> definition?
>>>> y = gamma(1.25)^2 * M(-.25;1;-x) / sqrt(x)
>>>> Note that the table data has an interval of .5 for the x axis.
>>>>
>>>> How far are your results from the data in the table?
>>>>
>>>> Cheers,
>>>>
>>>>    Jean-Marc
>>>>
>>>> Quoting John Ridges <jridges at masque.com>:
>>>>
>>>>
>>>>      
>>>>> Thanks for the confirmation  Jean-Marc. I kind of suspected from the
>>>>> comments that it was the confluent hypergoemetric function, which I
>>>>> was
>>>>> trying to evaluate using Kummer's equation, namely:
>>>>>
>>>>> M(a;b;x) is the sum from n=0 to infinity of (a)n*x^n / (b)n*n!
>>>>> where (a)n = a(a+1)(a+2) ... (a+n-1)
>>>>>
>>>>> But when I use Kummer's equation, I don't get the values in the
>>>>> "hypergeom_gain" table. Did you use a different solution to the
>>>>> confluent hypergoemetric function when you created the table?
>>>>>
>>>>> John Ridges
>>>>>
>>>>>
>>>>> Jean-Marc Valin wrote:
>>>>>
>>>>>        
>>>>>> Hi John,
>>>>>>
>>>>>> M(;;) is the confluent hypergeometric function.
>>>>>>
>>>>>>     Jean-Marc
>>>>>>
>>>>>> John Ridges a écrit :
>>>>>>
>>>>>>
>>>>>>          
>>>>>>> Hi,
>>>>>>>
>>>>>>> I've been trying to re-create the table in the function
>>>>>>> "hypergeom_gain"
>>>>>>> in preprocess.c, and I just simply can't get the same values. I
>>>>>>> get the
>>>>>>> same value for the first element, so I know I'm computing
>>>>>>> gamma(1.25)^2
>>>>>>> correctly, but I can't get the same numbers for M(-.25;1;-x),
>>>>>>> which I
>>>>>>> assume is Kummer's function. Is it possible that the comment is
>>>>>>> out of
>>>>>>> date and the values of Kummer's function used to make the table were
>>>>>>> different? Any help would be greatly appreciated.
>>>>>>>
>>>>>>> John Ridges
>>>>>>>
>>>>>>>
>>>>>>> _______________________________________________
>>>>>>> Speex-dev mailing list
>>>>>>> Speex-dev at xiph.org
>>>>>>> http://lists.xiph.org/mailman/listinfo/speex-dev
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>>             
>>>>>>           
>>>>>         
>>>>
>>>>
>>>>
>>>>       
>>>
>>>     
>>
>>
>>
>>
>>   
> 
> 
>