Курс лекций дисциплины вариативной части математического и естественнонаучного цикла (855805), страница 21
Текст из файла (страница 21)
На рис. 4.45 представлены кубическиеполиномы Лагранжа (N=4). Полином, выделенный красным цветом, равенединице в точке t=-2 и равен нулю в других точках. Полином, выделенныйсиним цветом, равен единице в точке t=-1 и равен нулю в других точках и т.д.L3(t)tРис. 4.45 – Полиномы Лагранжа (N=4)Для получения эффективной вычислительной процедуры распишемполиномы Лагранжа для случая N=4 и подставим необходимые моментывремени: t0=-2; t1=-1; t2=0; t3=1:(t − t1 )(t − t 2 )(t − t3 ) =(t + 1)t (t − 1)1= − (t 3 − t )(t0 − t1 )(t0 − t 2 )(t0 − t3 ) (− 2 + 1)(− 2)(− 2 − 1) 6(t − t0 )(t − t 2 )(t − t3 ) = (t + 2)(t )(t − 1) = 0.5t 3 + 0.5t 2 − tL13 (t ) =(t1 − t0 )(t1 − t 2 )(t1 − t3 ) (− 1 + 2)(− 1)(− 1 − 1)(t − t0 )(t − t1 )(t − t3 ) = (t + 2)(t + 1)(t − 1) = −0.5t 3 − t 2 + 0.5t + 1L32 (t ) =(t 2 − t0 )(t 2 − t1 )(t 2 − t3 )(2)(1)(− 1)(t − t0 )(t − t1 )(t − t2 ) = (t + 2)(t + 1)(t ) = 1 t 3 + 0.5t 2 + 1 tL33 (t ) =(t3 − t0 )(t3 − t1 )(t3 − t2 ) (1 + 2)(1 + 1)(1) 63L30 (t ) =171Произвольный отсчет интерполируемого сигнала является взвешенной суммой3кубических полиномов: y (t ) = ∑ x(n )L3n (t ) и при подстановке полученныхn =0выражений для полиномов Лагранжа получим:1 3(t − t )x(0) + 1 t 3 + 1 t 2 − t x(1) + − 1 t 3 − t 2 + 1 t + 1 x(2) +6222 211 1+ t 3 + t 2 + t x(3)23 6y (t ) = −Раскрыв скобки и сгруппировав слагаемые относительно степеней t, получимискомыйполином:y (t ) = a3t 3 + a2 t 2 + a1t + a0 ,осуществляющийрасчетзначения сигнала в любой момент времени с погрешностью определяемойпорядком интерполяции, т.е.
решающей задачу передискретизации илизадержки сигнала на момент времени некратный периоду дискретизации:1111a3 = − x(0 ) + x(1) − x(2 ) + x(3)622611a2 = x(1) − x(2 ) + x(3)22111a1 = x(0 ) − x(1) + x(2 ) + x(3)623a0 = x(2 )Каждый коэффициент полинома является взвешенной суммой отсчетовсигнала, как и процесс расчета значения полинома можно представить в видеКИХ фильтров – рис. 4.46.
Полученный фильтр известен как фильтр Фарроу(аналогично рассмотренному примеру можно получить фильтр любогопорядка). В представленной структурной схеме всего 9 умножений на расчеткоэффициентов полинома (умножения на ноль, на плюс/минус единицу несчитаем) и 3 умножения на расчет значения сигнала (по сравнению с 67умножениями это достижение). Однако нетрудно видеть, что отсчетыисходного сигнала умножаются на одни и те же коэффициенты – значит, естьвозможность группировки слагаемых перед умножением, что позволитснизить их количество.172x(n)-1/6x(0)0×z-1×z-10.5x(1)×z-1+×z-1+x(3)×z-10.5+×a3a2××z-1××+0+×a1++1+1/3+×z-1×z-10+0.5+0.5z-1×z-1-1×-1+×z-11/60××z-1-0.5x(2)1/6+a0+×+ty(t)Рис.
4.46 – Структурная схема фильтра Фарроу (N=4)Заметим, что коэффициенты:111111a3 = − x(0 ) + x(1) − x(2 ) + x(3) = [x(3) − x(0 )] + [x(1) − x(2 )]622662a1 + a3 =1[x(3) − x(1)]2a1 + a 2 + a3 = x(3) − x(2 )→a1 =1[x(3) − x(1)] − a32→a 2 = x(3) − x(2 ) − a1 − a3a0 = x(2 )Таким образом, модифицированный фильтр Фарроу содержит только 3умножения для вычисления коэффициентов полинома еще 3 умножения длявычисления значения отсчета. Структурная схема модифицированногофильтра представлена на рис.
4.47.173x(0)x(n)Σx(1)z-1+ΣΣ-+Σ--+Σ0.5 0.5+a0+-Σ-a3x(3)z-1+1/6+x(2)z-1- Σa1+a2×+×+×+ty(t)Рис. 4.47 – Структурная схема модифицированного фильтра Фарроу (N=4)В заключение, необходимо рассмотреть какие значения может приниматьпеременная t. При синтезе фильтра рассматривались четыре момента времени:t0=-2; t1=-1; t2=0; t3=1, поэтому значение переменной t должно лежать внутриинтервала [-2; 1]. Однако учитывая конечную степень интерполирующегополинома, максимальная точность интерполяции (минимальная погрешность)обеспечивается в интервале [-1; 0], т.е. в середине отрезка.Продолжая рассмотрение примера обработки сигнала с частотойдискретизации fs=44.1 кГц, изменим его частоту дискретизации до fsn=48 кГц.Поскольку требуется немедленное воспроизведение гипотетического контента,то данная задача относится к задачам реального времени, т.е. выполнениеопераций по обработке сигнала выполняется потоковым методом, при этомположим без потери общности, что сигнал монофонический.Обозначим:- шкала времени исходного сигнала: t x = nt s ;- шкала времени передискретизированного сигнала: t y = kt sn .174Введем для удобства обозначений коэффициент R =f s t sn== 0.91875 .
Тогдаf sn t sдля расчета отсчета выходного сигнала с произвольным индексом kпонадобятся четыре отсчета исходного сигнала с индексами: n-2, n-1, n, n+1(взяты для удобства нумерации в соответствии с моментами времени t0=-2;t1=-1; t2=0; t3=1). Рассчитать индекс n можно сопоставив шкалы временисигнала, при этом ближайшим отсчетом к требуемому моменту времени kt snбудет являться отсчет с индексом n = kR + 1 , где односторонние скобочкиx означают отбрасывание дробной части величины x.
Таким образом, точкеt0=-2 соответствует отсчет x(n-2), а точке t3=1 – отсчет x(n+1). В этом случаенетрудно пересчитать шкалу nts в интервал [-2; 1]: t = kR − n .Рассмотрим пример – запись звучания ноты «ля» (частота порядка 440Гц) с частотой дискретизации fs=44.1 кГц (показан на рис. 4.48) требуетсяпередискретизировать до частоты fsn=48 кГц.x(t)t, мсРис. 4.48 – Затухающая синусоида с частотой 440 Гц (модель ноты «ля»)В соответствии со спроектированным фильтром Фарроу для каждогоотсчета выходного сигнала y(k) рассчитаем:− индексы отсчетов исходного сигнала n = kR + 1 , т.е. из исходного сигналаделаем выборку из N=4 отсчетов: x(n-2), x(n-1), x(n), x(n+1);− переменную t = kR − n175− коэффициенты интерполятора a3 =a1 =1[x(3) − x(0)] + 1 [x(1) − x(2)], a0 = x(2) ,621[x(3) − x(1)] − a3 , a2 = x(3) − x(2) − a1 − a3 ,2− собственно интерполированные значения – исходный сигнал с другойчастотой дискретизации y (k ) = [(a3t + a2 )t + a1 ]t + a0Фрагмент результата работы фильтра Фарроу представлен на рис.
4.49(зеленые кружочки – исходные отсчеты, красные крестики – отсчетыпередискретизированного сигнала).y(k)t, мсРис. 4.49 – Результат передискретизации (фильтром Фарроу N=4)1764.9 Передискретизация как задача регрессии данныхС математической точки зрения задача передискретизации сигнала, какужеотмечалось,относитсякзадачаминтерполяциифункциипоравноотстоящим точкам (узлам интерполяции).
Однако реальные сигналычасто зашумлены, а ограниченный набор отсчетов сигнала не позволяет вполной мере говорить о статистических характеристиках шума. ПоэтомуинтерполированиесигналасмалымчувствительностикоэффициентовОСШприводитинтерполирующегокбольшоймногочленаотреализации шума (помех). Конечно, можно применить фильтрацию сигнала допередискретизации или после нее. Однако это требует дополнительныхвычислительных затрат.Альтернативой может быть «подбор» такой функции, которая научастках максимально близка именно к полезной составляющей сигнала – этозадачаоптимальнойаппроксимации.Этообеспечитпроцесспередискретизации и подавление шума одновременно. В этом контекстезадачу можно сформулировать в следующем виде: минимизация ошибкимежду сигналом x(t) и аппроксимирующей функцией d(t):pmin ∫ [x(t ) − d (t )] dt , где p – показатель нормы.Если p=1 – мера наименьших модулей, например, метод Лагранжа даетхорошиерезультатыприобработкеслучайныхсигналовсзакономраспределения близким к распределению Лапласа.Еслиp=2–квадратичнаямераобеспечиваетмаксимальноеправдоподобие функции приближения при нормальном законе распределенияслучайной составляющей зависимой переменной.При p=∞ – норма ошибки равна максимальному абсолютномуотклонению сигнала от аппроксимирующей функции.
Минимизация этойнормы соответствует минимаксной (равноволновой) аппроксимации (известнатакже как Чебышевская аппроксимация).177Аппроксимация данных с учетом их статистических параметровотносится к задачам регрессии, которые обычно возникают при обработкеэкспериментальныхданных,полученныхврезультатеизмеренийстатистических по своей природе сигналов. Как правило, при регрессионноманализе усреднение данных производится методом наименьших квадратов(p=2).Рассмотрим стационарный сигнал x(n), который содержит шум,распределенный, как правило, по нормальному закону.
Требуется определитьтакую функцию d(n), часто в виде полинома: d (n ) =минимизирует выражение:∑ [x(n ) − d (n )]2M −1∑ am n m ,котораяm =0, обеспечивающее минимальнуюnпогрешность приближения в среднеквадратическом смысле. Для определениякоэффициентов полинома минимизируемое выражение дифференцируютотносительно неизвестных величин am. Полученные частные производныеприравнивают нулю и решают систему уравнений, из которой находяткоэффициенты полинома (это соответствует отысканию экстремума функции).Виды регрессии обычно называются по типу аппроксимирующихфункций: полиномиальная, экспоненциальная, логарифмическая и т.п. Приприменении регрессии большое значение имеет адекватность модели (вид,порядок) полезной составляющей сигнала.4.9.1 Линейная регрессияАппроксимирующийN −1полином:d (n ) = an + b ,минимизируетсяфункционал: F = ∑ [x(n ) − an − b] , что приводит к системе уравнений:2n =0N −1 ∂F=−2∑ [x(n ) − an − b]n = 0 ∂an =0 ∂FN −1= −2 ∑ [x(n ) − an − b] = 0n =0 ∂bРешение системы уравнений:178a=N −16 2 N −1()nxn−x(n )∑∑N ( N + 1) N − 1 n=0n =0b=1NN −1∑ x(n ) −n =0N −1a24.9.2 Квадратичная регрессияАппроксимирующий полином: d (n ) = an 2 + bn + c ,[минимизируется]N −1функционал: F = ∑ x(n ) − an 2 − bn − c , что приводит к системе уравнений:n =02[[]][]N −1 ∂F=−2∑ x(n ) − an 2 − bn − c n 2 = 0 ∂an =0 ∂FN −1=−2∑ x(n ) − an 2 − bn − c n = 0 ∂bn =0N −1 ∂F = −2 x(n ) − an 2 − bn − c = 0∑ ∂cn =0Решение системы уравнений:N −1a = 30N −1N −1n =0n =06 ∑ n x(n ) − 6( N − 1)∑ nx(n ) − ( N − 1)( N − 2 )∑ x(n )2n =0(N − 2)(N − 1)N (N + 1)(N + 2)N −1N −1126b=∑ nx(n ) − N (N + 1) ∑ x(n ) − a(N − 1)N ( N − 1)( N + 1) n=0n =0c=1NN −1∑ x(n ) − an =0(N − 1)(2 N − 1) − b (N − 1)624.9.3 Экспоненциальная регрессияАппроксимирующий[полином:d (n ) = Ae n ,минимизируется]N −1функционал: F = ∑ x(n ) − Ae n , что приводит к системе уравнений:n =02[]N −1 ∂F=−2∑ x(n ) − Aen en = 0 ∂An =0 ∂FN −1= −2 ∑ x(n ) − Ae n Ae n n = 0 ∂n =0[]179Решение системы уравнений при x(n)>0:=6N −1N −1n =0n =02 ∑ n ln ( x(n )) − ( N − 1)∑ ln ( x(n ))1A = exp NN ( N − 1)( N + 1)N −1∑ ln(x(n )) − n =0(N − 1)24.9.4 Применение регрессии к зашумленному сигналуРассмотрим несколько сигналов (при n<0 сигнал отсутствует, N=128):0.5 + 2(1 − e n ), 0 ≤ n ≤ N − 1x1 (n ) = , = −0.0625 – сигнал имеет излом; (n− N )0.5+2e,N≤n≤2N−1x 2 (n ) = 0.4 sin (n ) + 0.6 cos(n ), =2, n ∈ [0, N − 1] – гладкий сигнал.NОценим возможности задержки этих сигналов на половину периодадискретизации при воздействии стационарного гауссова шума (n) спараметрами: M=0, =0.4.Таким образом, по имеющемуся дискретному сигналу x(n) необходимовычислить значение отсчета, расположенного посередине каждой парыисходного сигнала посредствомследующих методов (для возможностисравнения выберем значения апертур M=4 и M=8): регрессия линейная; регрессия квадратичная; регрессия экспоненциальная; фильтрация по Фарроу (апертура M=4).Сопоставим результаты обработки сигнала с помощью различныхрегрессий и величин апертуры, а также результатов кубической регрессии ифильтрации по Фарроу 3-го порядка: y (k ) = [(a3t + a2 )t + a1 ]t + a0 , t=0.5a3 =1[x(3) − x(0)] + 1 [x(1) − x(2)], a0 = x(2) ,62180a1 =1[x(3) − x(1)] − a3 , a2 = x(3) − x(2) − a1 − a32Вычислим СКО разности обработанного сигнала и эталона в каждомслучае.