Курс лекций дисциплины вариативной части математического и естественнонаучного цикла (855805), страница 19
Текст из файла (страница 19)
4.30подавленные участки спектраX2(f)АЧХ ФНЧf-0.5fs- fcfc0.5(f2- f1)Рис.4.30 – Спектр сдвинутого и отфильтрованного сигнала0.5fs154 выполнить децимацию результата на M;x3 (n ) =↓ M [x2 (n )] – рис. 4.31копии спектраX3(f)f-0.5fs-0.5 fд0.5 fдfд= fs /MРис.4.31 – Спектр сигнала в выделенном частотномдиапазоне0.5fs расчет ДПФ с высоким разрешением, т.к. fд<fs.Для оценки эффективности метода по сравнению с прямым вычислениемДПФ положим:− частота дискретизации исходного сигнала 512 кГц;− частоты f1=120 кГц и f2=136 кГц (полоса полезного участка спектра 16 кГц);− требуемое разрешение по частоте при выполнении ДПФ 1 кГц;− частоты прочих спектральных компонент относительно полезного участкаспектра находятся на удалении не менее 40 кГц.Тогда прямое вычисление ДПФ с заданным разрешением требуетвыборкисигналадлинойN=512отсчетовисоответственноN log 2 N = 4608 комплексных умножений.В соответствии с рассмотренным методом необходимо сдвинуть сигнал– 512 комплексных умножений, выполнить фильтрацию сигнала – Pfкомплексных умножений, выполнить децимацию сигнала на M < f s / f д , гдеfд>16 кГц.
Полоса пропускания ФНЧ должна быть 16 кГц, положим ширинупереходной полосы фильтра менее 16 кГц. Тогда получим: fд=32 кГц и M=16.Для реализации ДПФ с разрешением 1 кГц потребуется выборка N`=32отсчетов, т.е. 160 комплексных умножений. Если ФНЧ реализовать в видеКИХ фильтра с линейной фазой порядка K, то фильтрация потребуетPf=0.5(K+1)N` комплексных умножений. Для реализации ФНЧ с полосой 16кГц и переходной полосой 16 кГц и частоте дискретизации сигнала 512 кГцнеобходим фильтр порядка K≥27.155Итогоалгоритмтребует512+160+0.5·28·32=1120комплексныхумножений. Поэтому выигрыш (в смысле вычислительной эффективности поколичеству умножений) по сравнению с прямым вычислением ДПФтребуемого разрешения составляет более 75.6%, т.е.
решение данной задачи спомощью рассмотренного метода (структура представлена на рис. 4.32)эффективнее более чем в 4.1 раза.x(t)АЦПfs=512кГцx(n)x1(n)×БПФФНЧ↓16x2(n)x3(n)fд=32кГцX(k)e − j 2 f c nt sРис.4.32 – Структурная схема устройства спектрального анализаучастка спектра сигнала1564.5.3 Алгоритм локализации спектральных пиковРасчет ДПФ сигнала с неизвестной частотой (даже с применениемоконных функций или с передискретизацией сигнала) фактически всегдасодержит эффект размытия. Это означает, что при попытке точногоопределения частоты гармоники или ее амплитуды неизбежно возникаетпроблема точности измерения (поскольку максимум, который ассоциируется сдельта-функцией, отображающей в спектре гармоническую компонентусигнала, оказывается смещенным по частоте и искаженным по амплитуде) –рис. 4.33 (синусоидальный сигнал с амплитудой 1 В, частотой 100 Гц, частотадискретизации Fs=2 кГц, длина ДПФ N=128 отсчетов: разрешение по частоте=Fs= 15.625 Гц ).Nамплитуда93.75 Гц109.375 Гцf, ГцРис.
4.33 – Проблема оценки частоты и амплитуды синусоидальной компонентыОчевидно, частота синусоидальной компоненты лежит между спектральнымибинами с частотами 93.75 Гц и 109.375 Гц (или номера частотных отсчетов 6 и7 соответственно, нумерация с нуля).
Конечно, можно применить методы поуменьшению размытия: увеличить длину ДПФ (уменьшив разрешение по частоте) – не гарантируетрешения проблемы и требует значительных вычислительных затрат; применить метод дополнения нулями – не решает задачи, т.к. частотыспектральных пиков при этом смещаются, а амплитуды искажаются;157 применить подходящее окно анализа – размытие уменьшается, но недостаточно – рис.
4.34.1амплитуда0.80.60.40.20f, Гц050100150200250300Рис. 4.34 – Фрагмент спектра сигнала (рис. 4.33) при применении окнаРешение (без увеличения выборки сигнала, применения окно и т.п.)может быть найдено с помощью аппарата аппроксимации функции.
При этомнеобходимо быть уверенными, что в окрестности аппроксимируемогоспектрального пика (в данном примере 100 Гц) нет других частотныхсоставляющих сигнала.Положим, что индекс (номер отсчета), соответствующий максимумуспектрального пика (т.е. гармонической составляющей сигнала с частотойFp =Fsm p ) имеет вид: m p = m + Re{},Nгде m – индекс (номер отсчета), соответствующий максимуму в рассчитанномДПФ (содержащим размытие спектра) – является целым числом;δ – характеризует смещение спектрального пика (гармоники сигнала) оттекущего максимума.Определив положение максимума – индекс m по рассчитанному ДПФ(результат подвержен размытию спектра), можно выделить три комплексныхотсчета: X(m-1), X(m), X(m+1), т.е. значения спектра в точке максимума изначения в соседних точках.
По этим трем равноотстоящим точкам можнопостроить полином второй степени. Затем по найденным коэффициентамуравнения определить положение его единственного экстремума (максимума).158Положим коэффициенты a, b, c – коэффициенты полинома второйстепени am 2 + bm + c = X (m ) . Тогда справедливы уравнения:am 2 + bm + c = X (m )2a(m + 1) + b(m + 1) + c = X (m + 1)2a(m − 1) + b(m − 1) + c = X (m − 1)Решая систему уравнений, определим коэффициенты a, b и вычислимположение экстремума: = mopt = −bX (m + 1) − X (m − 1)=.2a X (m + 1) + X (m − 1) − 2 X (m )Тогда, для примера рис.
4.33, получаем следующие результаты:m p = m + Re{}=6+0.401825Fp =Fsm p =100.0285 ГцN1594.6 Эффективная реализация КИХ фильтровКакужеговорилось,применениеКИХфильтровимеетрядпреимуществ, однако при решении практических задач порядок фильтраполучается большим (десятки, сотни). Это требует решения вопросовэффективной реализации свертки. Одним из простых способов является метод,в основе которого лежит свойство ДПФ: спектр свертки сигналов – естьпроизведение спектров исходных сигналов. Тогда выполнив ДПФ от сигнала иимпульсной характеристики фильтра (т.к.
у КИХ фильтра она не изменяется,ее ДПФ вычисляют однократно и хранят в памяти устройства обработкисигнала), выполняют перемножение полученных спектров и затем вычисляютобратное ДПФ, которое будет представлять собой свертку импульснойхарактеристики фильтра и сигнала.Положим, дана импульсная характеристика фильтра h(k), порядокфильтра K. Требуется отфильтровать сигнал x(n), длина выборки N.x(n)x(n)*h(k)y(n)X(k)x(n)БПФ×ОБПФy(n)H(k)Рис.4.35 – Прямая (слева) и эффективная (справа) реализация КИХ фильтрацииПрямое решение задачи требует (K+1)N операций умножения.Реализация посредством алгоритма БПФ требует вычисления прямогоДПФ ( N log 2 N операций умножения), перемножения спектров (4N операцийумножения), вычисления обратного ДПФ ( N log 2 N операций умножения).Итого: 2 N log 2 N + 4 N (отметим, величина не зависит от порядка фильтра).Рассчитаем коэффициент эффективности как отношение числа операцийумножения при прямом решении задачи к числу операций умножения прирешении задачи посредством применения алгоритма БПФ: =(K + 1)2 log 2 N + 4.Зададимся вопросом, при каком порядке фильтра K имеет место улучшение,т.е.
η>1. Это выполняется если K > 2 log 2 N + 3 . Например: при N=1024, K>23.1604.7 Вычисление свертки секционированиемБылопоказано,чтосложностьарифметическихоперацийприреализации КИХ фильтров с помощью ДПФ не зависит от порядка фильтра, вто время как сложность реализации с непосредственным вычислением сверткипропорциональна порядку фильтра. Если порядок фильтра мал, можноожидать, что более эффективной будет реализация методом прямой свертки,однако по мере увеличения порядка фильтра реализация с ДПФ станет болееэффективной. Для фильтров очень высоких порядков выигрыш можетдостигать десятков и сотен раз. С другой стороны, реализация с ДПФ требуетзначительных объемов памяти и соответственно растет время обработки блока(требуетсянакопитьбуферотсчетовзаданнойдлины).Методысекционирования свертки являются компромиссным решением.
Суть ихзаключается в том, что операция свертки выполняется над секциями илиблоками данных с использованием БПФ для ускорения процесса вычислений.Ограничение размера секций уменьшает объем требуемой памяти, аиспользование алгоритмов БПФ сохраняет вычислительную эффективностьпроцедуры,компенсируя«усилия»наорганизациюпроцессаидополнительную вычислительную нагрузку.4.7.1 Метод перекрытия с суммированиемРазделим сигнал x(n) (длина выборки N может быть и бесконечной) насекции длиной по M отсчетов. Мысленно выделим секцию с индексамиmM≤n≤(m+1)M-1 (m – произвольный номер секции) – рис. 4.36.Отсчеты сигнала x(n)0 … M-1 …секция0секция1x`m(i)mM……секцияmmM+1секцияm+1mM+2……N-1секцияN/M-1mM+i…mM+M-1Рис.4.36 – Разбиение сигнала на непересекающиеся секции длиной M отсчетов161Расчет линейной свертки импульсной характеристики фильтра длиной Kотсчетов и секцией сигнала длиной M отсчетов дает M+K-1 отсчетов(рис.4.37):mM ≤ n ≤ (m + 1)M − 1 x(n),,xm′ (i ) = иначе0,K −1y m′ (m1 ) = ∑ h(k ) xm′ (m1 − k ) ,0≤i≤M-1.n = mM + m1 ,k =0m1 ∈ [0, M − 1] – смещение в секции,m ∈ [0, N / M − 1] – номер секции.Отсчеты сигнала x(n)0 … M-1 … 2M-1секция1секция2………секцияmсекцияm+1…N-1секцияN/Mh(k) K отсчетовВыходные отсчеты y`(n)K-1отсчетотклик 0отклик 1…M+K-1отсчетотклик mотклик m+1сложить… откликN/M-1Рис.4.37 – Метод перекрытия с суммированиемВ силу линейности операции свертки выходной сигнал есть сумма откликовфильтра h(k) на каждую секцию:y (n ) =N / M −1∑ ym′ (n mod mM ) =m =0N / M −1 K −1∑ ∑ h(k ) xm′ ((nm =0k =0mod mM ) − k ) .Для достижения вычислительной эффективности свертку в каждой секциивыполняют с помощью эффективных алгоритмов, например: с помощью БПФ.А пересекающиеся результаты вычисления сверток суммируются по границаминдексов исходных секций.1624.7.2 Метод перекрытия с накоплениемПри подаче возмущения на вход системы (фильтра) с конечнойимпульсной характеристикой h(n) длиной K отсчетов ее отклик имеет длинуK-1 отсчет.