Почему результаты свертки имеют разную длину при выполнении во временной области по сравнению с частотной областью? - PullRequest
10 голосов
/ 03 октября 2010

Я не эксперт по DSP, но я понимаю, что есть два способа применения дискретного фильтра во временной области к дискретному сигналу во временной области. Первый состоит в том, чтобы свести их во временную область, а второй - взять БПФ обоих, умножить оба комплексных спектра и взять ОБПФ результата. Одним из ключевых отличий в этих методах является то, что второй подход подвержен круговой свертке.

В качестве примера, если фильтр и формы волны имеют длину N точек, первый подход (то есть свертка) дает результат длиной N + N-1 точек, где первая половина этого отклика представляет собой заполнение фильтра и вторая половина - фильтр опустошения. Чтобы получить устойчивый отклик, фильтр должен иметь меньше точек, чем сигнал, подлежащий фильтрации.

Продолжая этот пример со вторым подходом и предполагая, что все дискретные данные формы сигнала во временной области являются действительными (не сложными), БПФ фильтра и формы сигнала оба создают БПФ длиной N точек. Умножение обоих спектров IFFT на результат дает результат во временной области также длиной N точек. Здесь ответ, где фильтр заполняется и очищается, перекрывает друг друга во временной области, и нет никакого ответа устойчивого состояния. Это эффект круговой свертки. Чтобы избежать этого, обычно размер фильтра будет меньше размера формы волны, и оба будут дополнены нулями, чтобы дать возможность пространственной свертке увеличиваться во времени после ОБПФ произведения двух спектров.

У меня такой вопрос, я часто вижу в литературе работы известных экспертов / компаний, у которых есть дискретная (реальная) форма сигнала во временной области (N точек), они БПФ, умножают ее на некоторый фильтр (также N баллы), и ОБЪЯВЛЯЕТ результат для последующей обработки. Мое наивное мышление заключается в том, что этот результат не должен содержать ответа в стационарном состоянии и, следовательно, должен содержать артефакты от заполнения / опорожнения фильтра, которые могут привести к ошибкам при интерпретации полученных данных, но я должен что-то упустить. При каких обстоятельствах это может быть правильным подходом?

Любое понимание будет с благодарностью

Ответы [ 4 ]

8 голосов
/ 03 октября 2010

Основная проблема заключается не в заполнении нулями в сравнении с предполагаемой периодичностью, а в том, что анализ Фурье разбивает сигнал на синусоидальные волны, которые на самом базовом уровне предполагаются бесконечными по степени. Оба подхода верны в том смысле, что IFFT, использующий полное БПФ, вернет точную форму входного сигнала, а оба подхода неверны в том смысле, что использование меньшего, чем полный спектр, может привести к эффектам накрая (которые обычно простираются на несколько длин волн).Единственное различие заключается в деталях того, что, как вы предполагаете, заполняет остальную часть бесконечности, а не в том, делаете ли вы предположение.

Назад к вашему первому абзацу: Обычно в DSP самая большая проблема, с которой я сталкиваюсьс БПФ заключается в том, что они не являются причинно-следственными, и по этой причине я часто предпочитаю оставаться во временной области, используя, например, фильтры FIR и IIR.

Обновление:

В постановке вопроса ОП правильно указывает на некоторые проблемы, которые могут возникнуть при использовании БПФ для фильтрации сигналов, например, краевые эффекты, которые могут быть особенно проблематичными при выполнении свертки, сравнимой по длине (во временной области) к выбранной форме волны.Однако важно отметить, что не вся фильтрация выполняется с использованием БПФ, и в статье, цитируемой ФП, они не используют фильтры БПФ, и проблемы, которые могут возникнуть при реализации фильтра БПФ, не возникают при использовании их подхода..

Рассмотрим, например, фильтр, который реализует простое среднее значение по 128 точкам выборки с использованием двух разных реализаций.

БПФ : При подходе БПФ / свертка можно получить выборку, скажем, из 256 точек, и свернуть ее с wfm, который является постоянным для первой половины и стремится к нулю ввторая половинаВопрос здесь (даже после того, как эта система провела несколько циклов), что определяет значение первой точки результата?БПФ предполагает, что wfm является круговым (то есть бесконечно периодическим), так что либо: первая точка результата определяется последними 127 (то есть будущими) выборками wfm (пропуская по середине wfm), или на 127 нулей, если вы нули.Ни один из них не является правильным.

FIR : Другой подход заключается в реализации среднего с помощью фильтра FIR.Например, здесь можно использовать среднее значение в 128-регистровой очереди FIFO.То есть по мере поступления каждой точки выборки: 1) поставить ее в очередь, 2) удалить из очереди самый старый элемент, 3) усреднить все 128 элементов, оставшихся в очереди;и это ваш результат для этого образца.Этот подход работает непрерывно, обрабатывая одну точку за раз и возвращая отфильтрованный результат после каждой выборки, и не имеет никаких проблем, возникающих из БПФ, поскольку он применяется к конечным фрагментным выборкам.Каждый результат представляет собой только среднее значение текущей выборки и 127 выборок, которые были до него.

В статье, на которую ссылается ОП, используется подход, гораздо более похожий на КИХ-фильтр, чем на БПФ-фильтр (обратите внимание, чтоФильтр в статье более сложный, и вся статья в основном представляет собой анализ этого фильтра.) См., например, эту бесплатную книгу , в которой описывается, как анализировать и применять различные фильтры, и обратите внимание также, чтоПодход Лапласа к анализу КИХ и БИХ-фильтров весьма схож с тем, что можно найти в цитируемой статье.

7 голосов
/ 19 ноября 2010

Вот пример свертки без дополнения нулями для DFT (круговая свертка) против линейной свертки.Это свертка последовательности длиной M = 32 с последовательностью длины L = 128 (с использованием Numpy / Matplotlib):

f = rand(32); g = rand(128)
h1 = convolve(f, g)
h2 = real(ifft(fft(f, 128)*fft(g)))
plot(h1); plot(h2,'r')
grid()

alt text Первые точки M-1 отличаются, и они короткиеM-1 балл, так как он не был заполнен нулями.Эти различия являются проблемой, если вы делаете свертку блоков, но для преодоления этой проблемы используются такие методы, как наложение и сохранение или наложение и добавление.В противном случае, если вы просто вычисляете операцию одноразовой фильтрации, действительный результат начнется с индекса M-1 и закончится индексом L-1 с длиной L-M + 1.

Что касаетсяВ цитируемой статье я посмотрел на их код MATLAB в приложении A. Я думаю, что они допустили ошибку, применив передаточную функцию Hfinal к отрицательным частотам без предварительного сопряжения.В противном случае вы можете увидеть на их графиках, что дрожание тактового сигнала является периодическим сигналом, поэтому использование круговой свертки хорошо для анализа в устойчивом состоянии.

Редактировать: Что касается сопряжения передаточной функции, ФАПЧ имеют действительную импульсную характеристику, а каждый действительный сигнал имеет сопряженный симметричный спектр.В коде вы можете видеть, что они просто используют Hfinal [Ni], чтобы получить отрицательные частоты, не принимая конъюгат.Я построил их передаточную функцию от -50 МГц до 50 МГц:

N = 150000                    # number of samples. Need >50k to get a good spectrum. 
res = 100e6/N                 # resolution of single freq point  
f = res * arange(-N/2, N/2)   # set the frequency sweep [-50MHz,50MHz), N points
s = 2j*pi*f                   # set the xfer function to complex radians 

f1 = 22e6       # define 3dB corner frequency for H1 
zeta1 = 0.54    # define peaking for H1 
f2 = 7e6        # define 3dB corner frequency for H2 
zeta2 = 0.54    # define peaking for H2    
f3 = 1.0e6      # define 3dB corner frequency for H3 

# w1 = natural frequency   
w1 = 2*pi*f1/((1 + 2*zeta1**2 + ((1 + 2*zeta1**2)**2 + 1)**0.5)**0.5)  
# H1 transfer function 
H1 = ((2*zeta1*w1*s + w1**2)/(s**2 + 2*zeta1*w1*s + w1**2))            

# w2 = natural frequency 
w2 = 2*pi*f2/((1 + 2*zeta2**2 + ((1 + 2*zeta2**2)**2 + 1)**0.5)**0.5)  
# H2 transfer function  
H2 = ((2*zeta2*w2*s + w2**2)/(s**2 + 2*zeta2*w2*s + w2**2))            

w3 = 2*pi*f3        # w3 = 3dB point for a single pole high pass function. 
H3 = s/(s+w3)       # the H3 xfer function is a high pass

Ht = 2*(H1-H2)*H3   # Final transfer based on the difference functions

subplot(311); plot(f, abs(Ht)); ylabel("abs")
subplot(312); plot(f, real(Ht)); ylabel("real")
subplot(313); plot(f, imag(Ht)); ylabel("imag")

alt text

Как видите, реальный компонент имеет четную симметрию, а мнимый компонент имеет нечетноесимметрии.В своем коде они только вычислили положительные частоты для графика журнала (достаточно разумно).Однако для расчета обратного преобразования они использовали значения для положительных частот для отрицательных частот путем индексации Hfinal [Ni], но забыли спрягать его.

1 голос
/ 05 октября 2010

Я могу пролить некоторый свет на причину, по которой «применение окон» применяется до применения БПФ.

Как уже указывалось, БПФ предполагает, что у нас бесконечный сигнал.Когда мы берем выборку за конечное время T, это математически эквивалентно умножению сигнала на прямоугольную функцию.

Умножение во временной области становится сверткой в ​​частотной области.Частотная характеристика прямоугольника - это функция синхронизации, то есть sin (x) / x.Символом «x» в числителе является кикер, поскольку он гаснет O (1 / N).

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

Однако, если у вас есть синус, который находится между 2 точками, вы увидите функцию синхронизации, отобранную для точки частоты.Это выглядит как увеличенная версия функции синхронизации, и «побочные» сигналы, вызванные сверткой, затухают с 1 / N или 6 дБ / октава.Если у вас есть сигнал на 60 дБ выше минимального уровня шума, вы не увидите шум на 1000 частот слева и справа от основного сигнала, он будет затоплен «юбками» функции синхронизации.

Если выиспользуя другое временное окно, вы получаете другую частотную характеристику, например, косинус исчезает с 1 / x ^ 2, есть специальные окна для различных измерений.Окно Ханнинга часто используется как окно общего назначения.

Дело в том, что прямоугольное окно, используемое без применения какой-либо «оконной функции», создает гораздо худшие артефакты, чем хорошо выбранное окно.то есть, «искажая» выборки времени, мы получаем намного лучшую картину в частотной области, которая ближе напоминает «реальность», или, скорее, «реальность», которую мы ожидаем и хотим видеть.

1 голос
/ 03 октября 2010

Несмотря на наличие артефактов, связанных с предположением, что прямоугольное окно данных является периодическим при ширине апертуры БПФ, что является одной из интерпретаций того, что делает круговая свертка без достаточного заполнения нулями, различия могут быть или не быть достаточно большими, чтобы затопитьанализ данных в вопросе.

...