Qual é a frequência de corte de um filtro de média móvel?

Preciso projetar um filtro de média móvel que tenha uma frequência de corte de 7,8 Hz. Eu já usei filtros de média móvel antes, mas, até onde sei, o único parâmetro que pode ser alimentado é o número de pontos a serem calculados em média ... Como isso pode estar relacionado a uma frequência de corte?

O inverso de 7,8 Hz é ~ 130 ms e estou trabalhando com dados que são amostrados em 1000 Hz. Isso implica que eu deveria estar usando um tamanho médio de janela de filtro em movimento de 130 amostras, ou há algo mais que estou perdendo aqui?

O filtro de média móvel (às vezes conhecido coloquialmente como filtro de vagão ) tem uma resposta retangular ao impulso:

$$h[n] = \frac{1}{N}\sum_{k=0}^{N-1} \delta[n-k]$$

Ou, dito diferentemente:

$$h[n] = \begin{cases} \frac{1}{N}, && 0 \le n < N \\ 0, && \text{otherwise} \end{cases}$$

Lembrando que a resposta de frequência de um sistema de tempo discreto é igual à transformação de Fourier em tempo discreto de sua resposta de impulso, podemos calculá-lo da seguinte forma:

$$\begin{align} H(\omega) &= \sum_{n=-\infty}^{\infty} x[n] e^{-j\omega n} \\ &= \frac{1}{N}\sum_{n=0}^{N-1} e^{-j\omega n} \end{align}$$

Para simplificar isso, podemos usar a fórmula conhecida para a soma dos primeiros termos de uma série geométrica$N$ :

$$\sum_{n=0}^{N-1} e^{-j\omega n} = \frac{1-e^{-j \omega N}}{1 - e^{-j\omega}}$$

O que mais interessa para o seu caso é a resposta de magnitude do filtro,. Usando algumas manipulações simples, podemos obtê-lo de uma forma mais fácil de compreender:$|H(\omega)|$

$$\begin{align} H(\omega) &= \frac{1}{N}\sum_{n=0}^{N-1} e^{-j\omega n} \\ &= \frac{1}{N} \frac{1-e^{-j \omega N}}{1 - e^{-j\omega}} \\ &= \frac{1}{N} \frac{e^{-j \omega N/2}}{e^{-j \omega/2}} \frac{e^{j\omega N/2} - e^{-j\omega N/2}}{e^{j\omega /2} - e^{-j\omega /2}} \end{align}$$

Isso pode não parecer mais fácil de entender. No entanto, devido à identidade de Euler , lembre-se de que:

$$\sin(\omega) = \frac{e^{j\omega} - e^{-j\omega}}{j2}$$

Portanto, podemos escrever o acima como:

$$\begin{align} H(\omega) &= \frac{1}{N} \frac{e^{-j \omega N/2}}{e^{-j \omega/2}} \frac{j2 \sin\left(\frac{\omega N}{2}\right)}{j2 \sin\left(\frac{\omega}{2}\right)} \\ &= \frac{1}{N} \frac{e^{-j \omega N/2}}{e^{-j \omega/2}} \frac{\sin\left(\frac{\omega N}{2}\right)}{\sin\left(\frac{\omega}{2}\right)} \end{align}$$

Como afirmei antes, você está realmente preocupado com a magnitude da resposta em frequência. Portanto, podemos aproveitar a magnitude do exposto para simplificá-lo ainda mais:

$$|H(\omega)| = \frac{1}{N} \left|\frac{\sin\left(\frac{\omega N}{2}\right)}{\sin\left(\frac{\omega}{2}\right)}\right|$$

Nota: podemos descartar os termos exponenciais porque eles não afetam a magnitude do resultado; para todos os valores de . Desdepara quaisquer dois números complexos finitos e , podemos concluir que a presença dos termos exponenciais não afeta a resposta de magnitude geral (em vez disso, eles afetam a resposta de fase do sistema).$|e^{j\omega}| = 1$$\omega$$|xy| = |x||y|$$x$$y$

A função resultante dentro dos parênteses de magnitude é uma forma de um núcleo Dirichlet . Às vezes, é chamada de função sinc periódica , porque se assemelha à função sinc um pouco na aparência, mas é periódica.

De qualquer forma, como a definição de frequência de corte é um pouco subespecificada (-3 dB ponto? -6 dB ponto? Primeiro sidelobe nulo?), Você pode usar a equação acima para resolver o que precisar. Especificamente, você pode fazer o seguinte:

1. Definirpara o valor correspondente à resposta do filtro que você deseja na frequência de corte.$|H(\omega)|$

2. Defina igual à frequência de corte. Para mapear uma frequência de tempo contínuo para o domínio de tempo discreto, lembre-se de que , em que é a sua taxa de amostragem.$\omega$$\omega = 2\pi \frac{f}{f_s}$$f_s$

3. Encontre o valor de que oferece a melhor concordância entre os lados esquerdo e direito da equação. Essa deve ser a duração da sua média móvel.$N$

Se é o comprimento da média móvel, uma frequência de corte aproximada (válida para ) na frequência normalizada é:$N$$F_{co}$$N >= 2$$F=f/fs$

$F_{co} = \frac {0.442947} {\sqrt{N^2-1}}$

O inverso disso é

$N = \frac {\sqrt{0.196202 + F_{co}^2}}{F_{co}}$

Essa fórmula é assintoticamente correta para N grande e possui cerca de 2% de erro para N = 2 e menos de 0,5% para N> = 4.

PS: Depois de dois anos, aqui finalmente, qual foi a abordagem seguida. O resultado foi baseado na aproximação do espectro de amplitude MA em torno de como uma parábola (Série de 2ª ordem) de acordo com$f=0$

$MA(\Omega) = \frac {Sin(\Omega∗N/2)}{Sin(\Omega/2)}$

$MA(\Omega) \approx 1+(\frac{1}{24}-\frac{N^2}{24})\Omega^2$

que pode ser mais exato perto do cruzamento zero de multiplicando por um coeficiente$MA(\Omega)-\frac{\sqrt2}{2}$$\Omega$

$\alpha=0.95264$

obtendo $MA(\Omega) \approx 1+0.907523(\frac{1}{24}-\frac{N^2}{24})\Omega^2$

A solução de fornece os resultados acima, onde .$MA(\Omega)-\frac{\sqrt2}{2}=0$$2\pi F_{co}=\Omega_{co}$

Tudo acima se refere à frequência de corte de -3dB, o assunto deste post.

Às vezes, embora seja interessante obter um perfil de atenuação na banda de parada que seja comparável ao de um filtro passa-baixa IIR de 1ª ordem (LPF de polo único) com uma determinada frequência de corte de -3dB (esse LPF também é chamado de integrador com vazamento, ter um poste não exatamente em DC, mas perto dele).

De fato, o MA e o IIR LPF de 1ª ordem têm uma inclinação de -20dB / década na banda de parada (é necessário um N maior que o usado na figura, N = 32, para ver isso), mas enquanto o MA possui nulos espectrais em e um envelope , o filtro IIR possui apenas um perfil .$F=k/N$$1/f$$1/f$

$H_{IIR} = \frac{1-Exp(-\Omega_{co})}{1-Exp(-\Omega_{co})*Exp(j \Omega)}$

Se alguém quiser obter um filtro MA com recursos de filtragem de ruído semelhantes aos do filtro IIR e combinar as frequências de corte 3dB iguais, comparando os dois espectros, ele perceberá que a ondulação da banda de parada do filtro MA acaba ~ 3dB abaixo do filtro IIR.

Para obter a mesma ondulação de banda de parada (ou seja, a mesma atenuação de potência de ruído) que o filtro IIR, as fórmulas podem ser modificadas da seguinte maneira:

$F_{co,IIR} = \frac {0.32} {\sqrt{N^2-1}}$

$N = \frac {\sqrt{0.1024 + F_{co,IIR}^2}}{F_{co,IIR}}$