what is the application of lavinson and durbin,, i mean where do u
imolement
>
Reply by Rune Allnor●January 21, 20112011-01-21
On Jan 20, 9:55=A0pm, John McDermick <johnthedsp...@gmail.com> wrote:
> Hi,
>
> I keep on getting a divide-by-zero exception when I execute the
> following code:
>
> Can somebody please tell me where my fixed-point math is "broken" ???
Don't know about the fixed-point math, but the *algorithm*
contains a step where a divide-by-0 type of error might
occur.
You need to identify that step and test if the result is 0.
If it is, you are dealing with a perfect AR(p) model with
no noise, and can terminate the function.
Rune
Reply by Randy Yates●January 20, 20112011-01-20
On 01/20/2011 04:09 PM, Randy Yates wrote:
> On 01/20/2011 03:55 PM, John McDermick wrote:
>> Hi,
>>
>>
>> I keep on getting a divide-by-zero exception when I execute the
>> following code:
>>
>> Can somebody please tell me where my fixed-point math is "broken" ???
>> Thank you!
>>
>> void levinsonDurbin(SINT32 *a, SINT32 *k, SINT32 *r, SINT32 n)
>> {
>> SINT32 a_temp[10+1], alpha, epsilon;
>> int i, j;
>> SINT32 tmpval;
>>
>> k[0] = 0;
>> a[0] = 32768;
>> a_temp[0] = 32768;
>> alpha = r[0];
>>
>> for(i = 1; i<= n; i++)
>> {
>> epsilon = 0;
>> for(j = 0; j< i; j++)
>> {
>> epsilon += ((a[j]*r[i-j])>> 0);
>> }
>>
>> epsilon = epsilon>> 15;
>>
>> a[i] = ((-epsilon)<< 10);
>> a[i] = a[i] / alpha;
>
> Probably here ^^^^^
>
>> a[i] = (a[i]<< 5);
>> k[i] = a[i];
>> tmpval = 32768 - ((k[i]*k[i])>> 15);
>> alpha = ((alpha*tmpval)>> 15);
>
> Very possibly, alpha is going to zero here ^^^^
>
>>
>>
>> /* update a[ ] array into temporary array */
>> for(j=1; j<i; j++)
>> {
>> a_temp[j] = a[j] + ((k[i]*a[i-j])>> 15);
>> }
>>
>> /* update a[ ] array */
>> for(j=1; j<i; j++)
>> {
>> a[j] = a_temp[j];
>> }
>> }
>> }
>
>
John,
If you were simulating this under cygwin or linux, gdb would be your
friend. In fact, the running gdb under emacs would be your bestest
buddy!
--
Randy Yates % "My Shangri-la has gone away, fading like
Digital Signal Labs % the Beatles on 'Hey Jude'"
yates@digitalsignallabs.com %
http://www.digitalsignallabs.com % 'Shangri-La', *A New World Record*, ELO
Reply by Randy Yates●January 20, 20112011-01-20
On 01/20/2011 03:55 PM, John McDermick wrote:
> Hi,
>
>
> I keep on getting a divide-by-zero exception when I execute the
> following code:
>
> Can somebody please tell me where my fixed-point math is "broken" ???
> Thank you!
>
> void levinsonDurbin(SINT32 *a, SINT32 *k, SINT32 *r, SINT32 n)
> {
> SINT32 a_temp[10+1], alpha, epsilon;
> int i, j;
> SINT32 tmpval;
>
> k[0] = 0;
> a[0] = 32768;
> a_temp[0] = 32768;
> alpha = r[0];
>
> for(i = 1; i<= n; i++)
> {
> epsilon = 0;
> for(j = 0; j< i; j++)
> {
> epsilon += ((a[j]*r[i-j])>> 0);
> }
>
> epsilon = epsilon>> 15;
>
> a[i] = ((-epsilon)<< 10);
> a[i] = a[i] / alpha;
--
Randy Yates % "My Shangri-la has gone away, fading like
Digital Signal Labs % the Beatles on 'Hey Jude'"
yates@digitalsignallabs.com %
http://www.digitalsignallabs.com % 'Shangri-La', *A New World Record*, ELO
Reply by John McDermick●January 20, 20112011-01-20
Hi,
I keep on getting a divide-by-zero exception when I execute the
following code:
Can somebody please tell me where my fixed-point math is "broken" ???
Thank you!
void levinsonDurbin(SINT32 *a, SINT32 *k, SINT32 *r, SINT32 n)
{
SINT32 a_temp[10+1], alpha, epsilon;
int i, j;
SINT32 tmpval;
k[0] = 0;
a[0] = 32768;
a_temp[0] = 32768;
alpha = r[0];
for(i = 1; i <= n; i++)
{
epsilon = 0;
for(j = 0; j < i; j++)
{
epsilon += ((a[j]*r[i-j]) >> 0);
}
epsilon = epsilon >> 15;
a[i] = ((-epsilon) << 10);
a[i] = a[i] / alpha;
a[i] = (a[i] << 5);
k[i] = a[i];
tmpval = 32768 - ((k[i]*k[i]) >> 15);
alpha = ((alpha*tmpval) >> 15);
/* update a[ ] array into temporary array */
for(j=1; j<i; j++)
{
a_temp[j] = a[j] + ((k[i]*a[i-j]) >> 15);
}
/* update a[ ] array */
for(j=1; j<i; j++)
{
a[j] = a_temp[j];
}
}
}