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

[Jean-Marc Valin]

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
> >>>>
> >>>>
> >>>>
> >>>>
> >>>
> >>>
> >>
> >>
> >
> >
> >
> >
> >
>
>
>