Да, я тоже считаю, что происходит самовозбуждение из-за недостаточной точности вычислений, хотя, чтобы утверждать точно, нужно моделировать, а на это нет времени. Для уменьшения эффекта нужно отодвинуть полюса фильтра внутрь единичной окружности (в Z-плоскости). Для этого достаточно модифицировать к-ты фильтра: A1 умножить на число 0,9<b<0.98, а А2 - на b^2. Сначала выберите b близким к 0,9. Если устойчивость не будет достигнута - причина ошибки в другом. Возможно, эффекты округления промежуточных данных и коэффициентов фильтра. Так что рекомендую прислушаться к совету - использовать матлаб.

ЗЫ. Из-за постоянки такая беда тоже может быть. Попробуйте перед укладыванием числа в буфер, вычесть из него постоянку.

Сообщение отредактировал Stanislav - Jan 25 2006, 16:13

Полагаю, что х - это входной сигнал, а у - выходной.

Из написанного уравнения следует: y[t]=(x[t]+a1*x[t-1]+a2*x[t-2]-b2*y[t-2])/b0

поскольку сумма коэффициентов при х не равна 0, это значит, что данный фильтр как минимум пропускает постоянную составляющую.

Есть большое подозрение, что запись исходного уравнения (для приведенных выше коэффициентов) неверна. Оно должно выглядеть так:

y[t] +a1*y[t-1] + a2*y[t-2] = b0*x[t] +b2*x[t-2]

а конечная рассчетная формула фильтра:

y[t] = b0*x[t]+b2*x[t-2]-a1*y[t-1]-a2*y[t-2]

В схеме как мне кажется - та же ошибка со знаками - задержанный на два такта сигнал должен не суммироваться, а вычитаться. Уравнение
2. Uad = 4 * b0 * abs(D0-D2)

Еще раз о коэффициентах фильтра.

Если рассчитанные коэффициенты:
A1=1.584584
A2=-0.987981
B0=-B2=0.006009

То для рассчетных формул вида:


1. D0=x[t]+a1*D1[t]+a2*D2[t]
2. y[t]=b0*(D0-D2)

где D1[t]=D0[t-1], D2[t]=D1[t-1], x[t] - входной сигнал, y[t] - выходной.

Коэффициенты в целочисленных процессорах коэффициенты обычно сдвигают так, чтобы максимальный из них занимал все значащие разряды.

В нашем случае, для 16-разрядной арифметики сдвиг на 14 разрядов,

a1=25961 (6569h)
a2=-16187 (С0С5h)
b0=98 (62h)
конечный результат надо сдвинуть на 14 разрядов вправо

Надо еще внимательно посмотреть на разрядность D0,D1,D2, чтобы не было переполнений. В любом cлучае, результат умножения a1*D1 и a2*D2 - как минимум 32-разрядный, для дальнейших вычислений младшие 16 бит отбрасываются.

Возможно, для оптимизации разрядности промежуточных результатов, сдвигать надо не на 14 бит, а на 13 или 12 (соответственно рассчитанные целочисленные коэффициенты сдвинуть на 1 или 2 бита вправо).

Сообщение отредактировал Johny - Jan 26 2006, 08:18

Почему-то не смог отредактировать. У результата умножения a1*D1 и a2*D2 не отбрасываются младшие 16 бит, он сдвигается вправо на 14 бит.

Кстати, надо следить за разрядностью входного сигнала. Если хотите, чтобы сумматоры были 16-разрядными, входной сигнал навскидку должен быть не более 14 разрядов. Иначе надо переходить на более высокоразрядную арифметику (например на 32-разрядную)

Как я понял, при расчете в таком виде

реально вы получаете фильтр с коэффициентами b0=x, b1=0, b2=0, а не b0=x, b1=0, b2=-x. В результате АЧХ, конечно, изменяется, но описываемых вами чудес не наблюдается.

Посмотрел еще в FDATool как сказываются эффекты квантования на характеристики фильтра на 430 Гц. При формате коэффициентов 4.12 все нормально при условии использования режима переполнения wrap. Если же использованить saturate фильтр становится неустойчивым. А если верить FDATool переполнения случаются. Возможно, причина глюков в этом.

Если есть сомнения в программной реализации, попробуйте проверить ее в Simulink. Сделайте ваш фильтр как S-function и расчет делайте точно также как вы делаете в AVR'е, т.е. с теми же типами переменных, форматами и прочее.

Спасибо за совет. Но удивительное в том, что первоначально использовалась арифметика 8.8. В этом случае из-за эффектов квантования B0=-B2=0.006009 на самом деле равнялось 0. Но устройство, в основе которог лежит эта программа, работает и серийно выпускается и насколько я знаю неплохо справляется со своими функциями. Странно что при небольшом увеличении разрядности, фильтры перестали работать..

Stanislav: Постоянку уже вычитаю. Все картинки - из модифицированной программы.

Johny: похоже разностное уравнение я немного переврал, а вот за "2. Uad = 4 * b0 * abs(D0-D2)" я вам благодарен. Мне тоже не нравится исходное уравнение, попробую изменить

Вопрос: может ли отсутствие ФНЧ так сильно оказывать влияние, при условии, что сигнал амплитуды -21дб идет из звуковой карты компьютера?



На самом деле все скорее вот так:

x[kT]+a1*x[(k-1)T]+a2*x[(k-2)T]=bo*y[kT]+b1*y[(k-1)T]+b2*y[(k-2)T]
уравнение для классического биквадратного звена

в моем случае b1=0
поэтому уравнение запишется как:

x[kT] + a1*x[(k-1)T] + a2*x[(k-2)T] = bo*y[kT] + b2*y[(k-2)T], где x - входной сигнал, y - выходной.

Наоборот, при такой записи уравнения

x-ВЫХОД, y-ВХОД

Кроме того, если коэффициенты b0=b1=b2=0, то получается, что выход не зависит от входа, на выходе фильтра всегда будет 0. Тогда не понятно как вообще ваше устройство работало.

так в том-то и дело что этих b0 и b2 в уравнениях нет
есть только левая часть, а правая почему-то вынесена за сумматор. Выход левой части
1. y=x[t]+a1*y[t-1]+a2*y[t-2].
2. берется абсолютное значение y.
И человек, который это придумал и написал предыдущую версию программы (мой шеф надо заметить), утверждает что с точки зрения математики это преобразование верно. И он доказал это работающим устройством.

Если это не так, мне нужно по крайней мере аргументированно доказать, раз уж не получается воспользоваться его мат. аппаратом.

И кстати, если убрать правую часть уравнения, получится цифровой резонатор

Есть тут b0, только b0=1, а вот b2, действительно, равно 0. Просто при такой записи выполнено масштабирование на b0.
Математика вообщем-то правильная, только вот, как я уже говорил, у меня подозрение, что при описанных вами уравнениях b2=0 и фильтр получается измененный.

Попробуйте реализовать обычную биквадратную форму без этих хитростей, которые непонятно как работают.