Что-то похожее кажется уже спрашивал на форумах, по моему до конкретики не дошло...
Имеется алгоритм обработки сигнала, исползуются различные фильтры (банк фильтров - один длинный FIR и много коротких простых IIRов второго порядка). Система рассчитана и просимулирована в Матлабе, теперь нужно имплементировать на 16-bit, fixed point DSP. Матлаб ессно вычислает в 64 бита floating point, т.е. коеффициенты фильтров имеют до 15 десятичных знаков после запятой что дает практически абсолютную точность - нет проблем. Однако, при переводе на 16 бит - имеем 4 десятичные цифры на все-провсе, т.е. проблема квантизации, плюс fixed point.
Вопросы:

1. Для подгонки под данного типа DSP, беру рассчетные (в Матлабе) коеффициенты, умножаю их на 10^3 после чего делаю round на результат. Таким образом получаю целые числа размером до 4 десятичных цифр которые и намереваюсь использовать как коеффициенты фильтров. После обработки делю результаты на 10^3. Так ли в реальной практике "сражаются" с переводом рассчетных параметров из Матлаба на 16-битные, fixed point процессоры ?

2. В процессе рутин фильтрации возможны overflow промежуточных и конечного результатов. Как боремся с этим ? Проверяем ли на overflow condition каждые промежуточные результаты ? Если да, и в случае обнаружения overflow - просто обрезаем результат на максимально-допустимое число (для данной разрядности данных) или предварительно делаем scaling ? Кроме того, насколько понял может потребоваться scaling входных и выходных сигналов во избежание overflow ? Каким образом это делаем ?

Спасибо.

По первому вопросу - Простое умножение коэффициентов и последующее округление или усечение результата
до целого может исковеркать АЧХ до неузнаваемости.
Умножать коэффициенты принято не на 10^..., а на 2^.... чтобы получить целочисленные значения
Загляните вот сюда -> http://www.analog.spb.ru/pub_dsp.htm
один из путей решения примерно такой ->
%--------------------------------------------------------------------------
function out=Create_IIR()
[bn,an] = butter(2,0.125,'low');
Hq=qfilt('df2',{bn,an},...
'CoefficientFormat', quantizer('fixed', 'round', 'saturate', [16 15]),...
'InputFormat', quantizer('fixed', 'round', 'saturate', [16 15]),...
'OutputFormat', quantizer('fixed', 'round', 'saturate', [16 15]),...
'MultiplicandFormat', quantizer('fixed', 'round', 'saturate', [16 15]),...
'ProductFormat', quantizer('fixed', 'round', 'saturate', [32 30]),...
'SumFormat', quantizer('fixed', 'round', 'saturate', [32 30]));
Hq=normalize(Hq)
figure;
freqz(Hq.ReferenceCoefficients{1},Hq.ReferenceCoefficients{2},1024,8000);
figure;
freqz(Hq.QuantizedCoefficients{1},Hq.QuantizedCoefficients{2},1024,8000);
%--------------------------------------------------------------------------
fp=fopen('d:\iir.h','w');
fprintf(fp,'const short IIR_Coefs_DF2[]=\n{\n');
%...
% выводим массив kоэффициентов в том виде как нам надо
%...
fprintf(fp,'//- Numerator coefs \n');
for i=length(Hq.QuantizedCoefficients{1}):-1:1
temp=num2hex(Hq.CoefficientFormat,Hq.QuantizedCoefficients{1}(i));
fprintf(fp,' 0x');
fprintf(fp,temp);
fprintf(fp,', // B%x- Numerator coefs \n',i-1);
end
fprintf(fp,'//- Denumerator coefs \n');
for i=length(Hq.QuantizedCoefficients{2}):-1:1
temp=num2hex(Hq.CoefficientFormat,Hq.QuantizedCoefficients{2}(i));
fprintf(fp,' 0x');
fprintf(fp,temp);
fprintf(fp,', // B%x- Numerator coefs \n',i-1);
end
fprintf(fp,'};\n\n');
fclose(fp)
out=Hq;
%----------------------------------------------------------------------
Ключевые слова для поиска в хелпе матлаба - qfilt,quantizer,normalize
%------------------------------------------------------------------------
По Второму - Действительно перегрузку аккумуляторов надо предотвращать
предварительно масштабируя сигнал, и используя насыщенную арифметику. Я Бы порекомендовал предварительно
написать функцию фильтра в матлабе и помоделировать, заодно увидите, где и
насколько возможна перегрузка во внутренних состояниях. Таким образом
сможете определить как надо масштабировать по входу и по выходу чтобы
достичь нужного результата.

Большое спасибо, как всегда ваша помощь весьма существенна.
Сейчас пытаюсь разбираться в вашем примере и релевантным функциям (qfilt, quantizer, normalize)...
У меня - МАТЛАБ 7: qfilt там уже obsolete, посему хелпа по нему не дает, но рекоммендует вместо него использовать dfilt структуру (они несколько отличаются полями).
Сейчас пытаюсь разобратся на вашем примере что и как произходит и зачем. Думаю будут еще вопросы...

Пока насчет второго: насчет маштабирования сигнала - проблема: маштабирование для выхода по входному сигналу (output scaling) ведь нелинейная операция - вносит искажения в сигнал, значит путь прямого маштабирования входного сигнала - видимо не самый лучший вариант. С другой стороны, моя "DSP библия" рекоммендует соотв. маштабирование коеффициентов фильтра с тем что-бы в результате конволюции избежать overflow как результат входного сигнала и заодно использовать полный диапазон разрядности. Там описываются несколько способов/подходов к этому - в time domain и 2-3 подхода во frequency domain. В результате рассчитывается scaling factors на коеффициенты фильтра.

Кроме того, насчет возможного overflow внутренних (промежуточных) результатов - есть следующие утверждения:
1. Произведения сигнала на константы - если оба маштабированы ниже 1 - их произведение всегда меньше 1 - значит нет опасности overflow.
2. Сумма сигналов служащая промежуточной суммой для последующих сложений. Есть потенциал на overflow, НО, если вычисления ведуться в two's complement арифметике (что обычно и есть на практике сегодня) то есть правило гласящее что ежели конечная сумма не дает overflow, то даже ежели и будут overflow среди промежуточных сумм ее формирующих, эти overflow cancel out друг друга. Вывод: при two's complement арифметике, overflow промежуточных сумм можно игнорировать.
3. Результат суммы 2х сигналов формуриющий конечную сумму которая требуется в фильтре в дальнейшем или вне фильтра. Тут может быть overflow с серьезными последствиями. Такой резыльтат должен быть маштабирован способом идентичным маштабированию для выхода (т.е. маштабированием impulse response фильтра): т.е. рассчитать transfer function от входа до проблемной суммы и маштабировать ее одним из способов упомянутых выше. Кроме того, как вы и упоминали, применять saturation по overflow. Вроде-бы у многих DSP процессоров есть спец. saturation mode в котором он сам определяет состояние overflow после каждой операции после сложения и делает saturation.

Подозреваю что некоторые операции которые вы применили в вашем примере (например normalize) имеют отношение к вышесказанному. Пытаюсь разобраться в деталях что и как...

Кстати - вы правы насчет фактора множителя коеффициентов 2^ а не 10^, 2^ ведет просто к бинарным сдвигам что ессно дает целочисленный результат...

Я почти так же делаю. Только коэффициенты сначала нормализую в смысле максимального модуля значения. То есть, умножаю их на такое число, чтобы максимальный или минимальный к-т были по модулю равны 32767 или -32768. После вычисления выходного отсчёта он домножается на обратное число. Потеря точности при таком подходе получается минимальной. Этот способ годится для КИХ-фильтров; для БИХ - не знаю, подумать надо, как отмасштабировать...
Впрочем, иногда достаточно просто "сдвинуть" к-ты (домножить на степень 2), а потом вых. отсчёт, в противоположную сторону. Это гораздо проще...

У DSP, как правило, есть "запасные" разряды аккумулятора, именно для того, чтобы избегать переполнения. А также режим "насыщения" (saturation), который при переполнении сохраняет хотя бы знак числа.

Спасибо Станислав. В принципе ежели подгонка коеффициентов фильтра делается в Матлабе (т.е. не в самом DSP) для переноса в реализацию в DSP, но тогда наверно не имеет значения умножать на 10^3 (для 16-бит ессно) или на 2^10, я ошибаюсь ?
С другой стороны, в плане нормализации до маскимального модуля коеффициентов - тут у меня есть некоторые НО: в принципе - говориться об маштабировании в целях избежания overflow выхода. Это можно делать либо маштабируя (ограничивая) входной сигнал что неприемлимо ибо вносит искажения (ограничение сигнала по амплитуде ведь нелинейная операция - вносит искажения), либо соотв. образом маштабируя коеффицеинты фильтра (т.е. его impulse response с которым и деламе конволюцию). Маштабирование фильтра не искажает сигнал, при том дает scaling factor такой что одновременно избегает overflow на выходе и более полное использование диапазона разрядности.
К чему я это ? Я имею ввиду что мне кажется в свете изложенного что нормализация коеффициентов (т.д. подбирая фактор) целесообразно из вышеприведенного расчета, но не просто с целью разтянуть сами коеффициенты до максимума разрядности (без учета максимума входного сигнала).
Насколько понял, обычно стараются нормализовать например входной сигнал вгоняя его в range [-1 1], что дает нам всегда знание максимально-возможного входного сигнала |Xin|max = 1, что в свою очеред позволяет определить scaling factor коеффициентов фильтра.

Насчет запасной разрядности DSP - тут вы правы. У того что у меня (TI C5402) - аккумулятор имет 40 бит (при том что сам процессор - 16-битный), старшие 8 бит - quard bits как раз для предотварщения overflow накопления. Т.е. что произошло его переполнение, нужно что было 256 последовательных сумм каждая из которых давала-бы overflow, что наверно редкость в реальности. И что тоже верно, я читал что у большинства DSP процессоров есть режим saturation, пока не в курсе как он включается у C5402 который у меня... (может он автоматически включен как default ?).
Наверно такая разрядность аккумулятора впридачу с режимом saturation позволяет избежать головной боли поиска и устранения overflow промежуточных результатов, но вероятно что не заменяет output scaling коеффициентами фильтра...хотя может я и не прав...

Насчет IIRов (БИХ это IIR ?), 4x и более порядковые IIRы обычно вроде как реализуются парными каскадами (4х порядковый IIR раскладывается на 2 последовательных 2х порядковых каскада), каждый такой каскад реализуется в coupled form что дает равномерное распределение полюсов по площади единичного круга). И тогда, делается маштабирование коеффициентов каждого их 2х порядковых фильтров, т.е. получаем scaling factor для каждого каскада. Затем вычисляются доп. факторы которые умножают выходы каждого каскада с тем чтобы получить общие gains всех каскадов равным тому что было в оригинале (до маштабирование коеффициентов). Т.е. отличия в подходе к FIR и к IIR не большие, но есть...

На самом деле - все это из теории, но вроде-бы по идее очень близкой к практике и ориенторованной на реализацию...

Для КИХ-фильтров есть простой, но довольно эффективный способ, называемый блочная плав. запятая. Имеющий смысл при более-менее приличных длинах.
БИХ-фильтры (IIR) 2-го порядка. Там все проблемы с рекурсивной частью. В ней единственный к-т может быть больше 1. Второй всегда меньше 1. Соответственно можно либо принять формат для к-тов 2.14(2 разряда, один из них по совместительству знаковый на целую часть, 14 на дробную), что достигается умножением ваших плавающих к-тов на 16384 (с дальнейшим округлением), либо принять формат для к-тов 1.15, он же q15, а к-т, который может быть больше 1, представить как сумму 2-х, в 2-раза меньших. Так точнее будет. С к-тами при нерекурсивной части проблем не будет. Там совершенно спокойно можно все отнормировать и к формату 2.14, и к 1.15, а компенсацию нормировки сделать по выходу. Дополнительно по БИХ. Промежуточные 32-х разрядные результаты и соответсвующая организация вычислений обычно идут во благо , если конечно есть ресурсы.
На переполнение прверок не делают. Устанавливается режим работы с насыщением (saturation). Он должен быть в Вашем DSP.

Спасибо.
Что такое блочная запятая ? Интересно расширить свой кругозор
Насчет IIRов - тут у меня в принципе меньше волнений - они у меня все простые - второго порядка, не более того, т.е. не нужно разбивать на каскады. Кроме того 40-битный (и даже 32х) аккумулятор с режимом насыщения думаю будет достаточно адекватной мерой для них...хотя это нужно еще понять более точно...
Тут в универских лаб. работах преподаватель просто домножал входной сигнал (который изначально floating 32 бита) на 2^10 (т.е. 10 эффективных бит целого) после чего - round, а коэффициенты домножал на 2^9 (т.е. беря 9 эффективных битов целого). В результате одного умножения получаем
до 19 бит, последующее сложение дает до 20 бит, аккумулятор имеет 32 реальных бита (+ 8 guard bits), что дает 2^12 операций сложения с переполнением для достижения overflow аккумулятора что наверно с случае не черезвычайно длинных фильтров статистически достаточно (т..е вероятность такого переполнения - очень мала). С другой стороны, для достижения приемлимой точности frequency response IIRов, на лабе получили необходимость минимум 5-6 десятичных знаков коеффициента (т.е. превышает 16 бит), если я это правильно понял.

Ну конечно, валяйте, расширяйте кругозор. Святое дело.
Что такое просто плавающая запятая, Вы в курсе? Знаете, что каждое число представляется в виде: мантисса*2^порядок? Т.е. для каждого числа есть и мантисса и порядок?
Довольно часто возникают ситуации, и КИХ фильтр (не всякий) одна из таких, когда для целой группы, блока, чисел можно взять один порядок, мантиссы остаются индивидуальные и в пределах порядка денормированные. В пределах этой группы все вычисления можно производить только с мантиссами. И вычисления эти будут в чистом виде вычислениями с фиксированной запятой. А по окончании вычислений подкорректировать порядок.
КИХ-фильтры, классического, что ли, типа (ФНЧ, ФВЧ, ПФ, РФ, Гильберт, дифференциаторы), как правило имеют большую по величине огибающую импульсной характеристики в своем центре и малую, стремящуюсю к 0, по краям.
Этим можно воспользоваться и разбить ИХ на несколько блоков. У каждого блока свой порядок.
В пределах каждого блока Вы обрабатываете сигнал используя арифметику с фиксированной запятой. А потом сшиваете несколько чисел (блочных результатов) используя принципы плавающей запятой.
Окончательный результат округляете, усекаете в любой нужный Вам формат. Все.
Число блоков, по опыту, 3-7. Потери, следовательно, мизерные по сравнению с чистой фмксир. запятой .
А бих я Вам на каскады разбивать не предлагал, он у Вас и так 2 порядка.
Основное, что я Вам говорил про бих 2-го порядка, у них единственный к-т в рекурсивной части может быть больше 1, но при этом он всегда меньше 2-х. Следовательно без проблем к-ты можно умножать на 16384=2^14.
Что, как Вы наверное должны догадываться, немного лучше умножения на 512, в смысле точности.
А чуть-чуть извернувшись легко можно умножать на 32768=2^15. И это я Вам тоже писал.
По поводу переполнений. Уже писал Вам про включение режима saturation.
Второй способ - можно заранее определить, грозят Вам переполнения, или не грозят.
Для КИХ фильтра это вообще очень просто определить. Достаточно просуммировать модули его к-тов. И сравнить с 1 (во избежании недоразумений - считается, что все к-ты по модулю меньше равны 1). Если сумма модулей к-тов КИХ фильтра меньше 1 - у Вас никогда не будет переполнения. Больше-равна 1 - есть вероятность этого.
Простейший способ борьбы, кроме включения насыщения, разбиение к-тов на несколько блоков, с гарантией отсутствия переполнения в каждом из блоков. Блочная плавающая запятая, кстати, как правило, решает заодно и эту задачу.

OK спасибо. Плавающая точка конечно-же знаю что это такое, как и fixed. Я просто плохо ориентируюсь в русскоязычных терминах, очень часто путаю в этом плане божий дар с яичницей ибо вся моя теория по которой учился, литература и 90% общения на профи темы - англоязычна как минимум в плане терминологии.
Кстати не совсем понял что такое огибающая импульсной характевристики... , как это в оригинале ? (в англоязычной терминологии) Вы горотите об impulse response shape ? Или это что-то касающееся frequency response ?
В целом насчет блочной точки - более-менее понятно.

Кстати, в лабе о которой упоминал, преподаватель просто ограничился более малой разрядностью (10 и 9 бит) так что-бы с гарантией избежать overflows, не утруждая себя целесообразностью использования по возможности более широкого динамического диапазона разрядности (т.е. вместо фактора 2^14 или 2^15 с последующим scalingом коеффициентов, он просто дал фактор 2^10/2^9, чем возможно ухудшил точность, но гарантированно избегая overflows). Это то что я понял на свой хлопский розум.

Кстати сейчас проверил согласно вашему совету сумму модулей коеффициентов моего FIRa = 2.76...
Значит при желании успользования как можно более широкого диапазона разрядности коеффициентов - нужно учетсь возможность overflows...

Елки-палки. Импульсная характеристика это импульсная характеристика, это импульс респонз, импульсный отклик, это отклик на дельта-импульс, последний иногда называют единичный импульс. Частотная - это частотная, это преобразование Фурье от импульсной, это фрекуенси респонз, частотный отклик. Коэффициенты КИХ фильтра и есть его импульсная характеристика. Огибающая это envelope.


Мотивы Вашего преподавателя мне мало интересны.
БИХ фильтры намного чуствительнее к к-там, чем КИХ. И даже фильтры 2 порядка могут потребовать разрядности к-тов, большей, чем 16. Не говоря уж о 9-10. Поэтому применение 2^14 или 2^15 во 2 порядке просто святое дело, это максимум, выжимаемый на 16 разрядах. И этот максимум практически ничего не стоит. Проблема переполнений - надуманная, решается чрезвычайно просто.
Куча проблем возникает именно из-за грубого квантования к-тов рекурсивных фильтров.
И еще одна куча из-за грубого представления промежуточных результатов.
Следствием этих грубостей очень часто являются дальнейшие вопросы из серии - у меня в подполье стук, помогите.


Да, должны учесть. При 16-ти разрядныхданных и к-тах неокругленный результат у Вас может быть 34 разряда.