Janela gaussiana média móvel

Variações na média móvel O filtro de média móvel é mais ou menos perfeito para suavizar dados na presença de ruído, se as informações úteis em seus dados estiverem completamente no domínio do tempo. Nesse caso, você não se preocupa com seu desempenho bastante fraco no domínio da freqüência. A Figura 1 mostra as respostas de impulso, passo e freqüência do filtro básico de média móvel (com três amostras extras em ambos os lados que não fazem parte das respostas de impulso e passo, para maior clareza). Às vezes, no entanto, você precisa trabalhar com dados para os quais ambos os domínios são importantes. Para esses casos, existem versões ponderadas da média móvel que são mais ou menos equivalentes no domínio do tempo, mas que têm um desempenho muito melhor no domínio da freqüência. Média móvel repetida A primeira coisa que você pode fazer para melhorar a resposta de freqüência da média móvel é aplicá-la várias vezes. Após duas repetições, isso equivale a uma ponderação triangular dos coeficientes (Figura 2). Uma vez que a aplicação do mesmo filtro duplica duas vezes seu efeito, o primeiro lóbulo lateral da resposta de freqüência é apenas metade tão alto quanto o da Figura 1. O motivo da forma triangular é que a média móvel é uma convolução com um pulso retangular. A aplicação de duas vezes provoca uma convolução desse pulso retangular, resultando em uma janela triangular para o filtro combinado. Observe que Ive tirou o mesmo comprimento de filtro na Figura 2 como na Figura 1, deslocando assim o primeiro zero da resposta de freqüência. Uma verdadeira convolução do filtro retangular original resultaria em um filtro mais longo e teria mantido os zeros exatamente no mesmo lugar, é claro. Se o filtro de média móvel for repetido várias vezes, seus coeficientes convergem para uma janela gaussiana (Figura 3) devido ao teorema do limite central. Claro, um Gaussiano real se estende infinitamente em ambas as direções, então não há outra opção do que cortá-la em algum ponto (ou talvez multiplicá-la com uma segunda janela). Além disso, o desvio padrão do gaussiano deve ser escolhido. Para esta ilustração (e para a implementação do Designer de filtros), adotei as configurações padrão do MATLAB. Na prática, você pode querer simplesmente aplicar repetidamente a média móvel em vez de aplicar uma janela gaussiana. Quando implementado de forma recursiva, a média móvel é muito eficiente. Enquanto a janela gaussiana deve ser implementada através da convolução. Blackman Window Outra possibilidade é escolher uma das funções de janela clássicas que são usadas para filtros windowed-sinc, e usá-lo como kernel de filtro (veja a excelente página da Wikipedia em funções de janela). Por exemplo, eu escolhi a janela Blackman (Figura 4). Isso melhora a atenuação da faixa de parada ainda mais, enquanto ainda mostra uma resposta de domínio do tempo suave sem qualquer toque ou superação. Em conclusão, se você precisa alisar dados, mas precisa de um melhor desempenho de freqüência do que a média móvel básica tem para oferecer, várias alternativas estão disponíveis. Ferramenta de design de filtro Este artigo é complementado com uma ferramenta de design de filtro. Experimente as diferentes funções da janela e o comprimento do filtro e veja o efeito na resposta de freqüência. Experimente agora O Guia de cientistas e engenheiros para processamento de sinal digital Por Steven W. Smith, Ph. D. Capítulo 15: Filtros médios móveis Parentes do filtro de média móvel Em um mundo perfeito, os designers de filtros só precisam lidar com informações codificadas no domínio do tempo ou no domínio da freqüência, mas nunca uma mistura dos dois no mesmo sinal. Infelizmente, existem algumas aplicações em que ambos os domínios são simultaneamente importantes. Por exemplo, sinais de televisão se enquadram nesta categoria desagradável. As informações de vídeo são codificadas no domínio do tempo, ou seja, a forma da forma de onda corresponde aos padrões de brilho na imagem. No entanto, durante a transmissão, o sinal de vídeo é tratado de acordo com sua composição de freqüência, como sua largura de banda total, como as ondas de suporte para cor de amplificador de som são adicionadas, restauração de amplificação de eliminação do componente de CC, etc. Como outro exemplo, interferência eletromagnética É melhor entendido no domínio de freqüência, mesmo que a informação de sinais seja codificada no domínio do tempo. Por exemplo, o monitor de temperatura em uma experiência científica pode estar contaminado com 60 hertz das linhas de energia, 30 kHz de uma fonte de alimentação de comutação ou 1320 kHz de uma estação de rádio AM local. Parentes do filtro de média móvel têm melhor desempenho de domínio de freqüência e podem ser úteis nestas aplicações de domínio misto. Os filtros médios móveis de passagem múltipla envolvem passar o sinal de entrada através de um filtro médio móvel duas ou mais vezes. A Figura 15-3a mostra o kernel geral do filtro resultante de uma, duas e quatro passagens. Duas passagens equivalem a usar um kernel de filtro triangular (um kernel de filtro retangular convolvido com ele próprio). Após quatro ou mais passagens, o kernel de filtro equivalente parece um Gaussiano (lembre-se do Teorema do Limite Central). Conforme mostrado em (b), as passagens múltiplas produzem uma resposta de passo em forma de S, em comparação com a linha reta da única passagem. As respostas de freqüência em (c) e (d) são dadas pela Eq. 15-2 multiplicado por si mesmo por cada passagem. Ou seja, cada vez que a convolução do domínio resulta em uma multiplicação dos espectros de freqüência. A Figura 15-4 mostra a resposta de freqüência de dois outros familiares do filtro de média móvel. Quando um Gaussiano puro é usado como um kernel de filtro, a resposta de freqüência também é gaussiana, conforme discutido no Capítulo 11. O gaussiano é importante porque é a resposta de impulso de muitos sistemas naturais e manmade. Por exemplo, um breve pulso de luz entrando em uma longa linha de transmissão de fibra óptica sairá como um pulso gaussiano, devido aos diferentes caminhos captados pelos fótons dentro da fibra. O kernel de filtro gaussiano também é usado extensivamente no processamento de imagens porque possui propriedades únicas que permitem rápidas ondulações bidimensionais (ver Capítulo 24). A segunda resposta de freqüência na Fig. 15-4 corresponde ao uso de uma janela Blackman como kernel de filtro. (A janela do termo não tem significado aqui é simplesmente parte do nome aceito desta curva). A forma exata da janela Blackman é dada no Capítulo 16 (Eq. 16-2, Fig. 16-2) no entanto, parece muito com um gaussiano. Como estes parentes do filtro de média móvel melhor do que o filtro de média móvel em si. Três maneiras: primeiro e mais importante, esses filtros têm melhor atenuação de parada do que o filtro de média móvel. Em segundo lugar, os grãos de filtro se afilam a uma amplitude menor perto das extremidades. Lembre-se de que cada ponto no sinal de saída é uma soma ponderada de um grupo de amostras da entrada. Se o kernel do filtro diminui, as amostras no sinal de entrada que estão mais distantes recebem menos peso do que as próximas. Em terceiro lugar, as respostas passo a passo são curvas suaves, em vez da linha direta abrupta da média móvel. Estes últimos dois geralmente são de benefício limitado, embora você possa encontrar aplicativos onde eles são vantagens genuínas. O filtro de média móvel e seus parentes são quase iguais em reduzir o ruído aleatório enquanto mantém uma resposta passo a passo. A ambigüidade reside na forma como o tempo de subida da resposta passo é medido. Se o tempo de subida for medido de 0 a 100 da etapa, o filtro médio móvel é o melhor que você pode fazer, como mostrado anteriormente. Em comparação, medir o tempo de subida de 10 a 90 torna a janela Blackman melhor do que o filtro de média móvel. O argumento é que isso é apenas dificuldades teóricas consideram esses filtros iguais neste parâmetro. A maior diferença nesses filtros é a velocidade de execução. Usando um algoritmo recursivo (descrito em seguida), o filtro de média móvel funcionará como um raio em seu computador. Na verdade, é o filtro digital mais rápido disponível. Várias passagens da média móvel serão correspondentemente mais lentas, mas ainda muito rápidas. Em comparação, os filtros gaussianos e negros são incrivelmente lentos, porque devem usar convolução. Acho um fator de dez vezes o número de pontos no kernel do filtro (com base na multiplicação sendo cerca de 10 vezes mais lento do que a adição). Por exemplo, espere que um gaussiano de 100 pontos seja 1000 vezes mais lento do que uma média móvel usando recursão. Esta versão HTML é fornecida por conveniência, mas não é o melhor formato para o livro. Em particular, alguns dos símbolos não são processados ​​corretamente. Você pode preferir ler a versão em PDF. CapítuloXA08XA0XA0Filtração e convolução Neste capítulo, apresento uma das idéias mais importantes e úteis relacionadas ao processamento do sinal: o teorema da convolução. Mas antes que possamos entender o teorema da convolução, temos que entender a convolução. IX2019ll começa com um exemplo simples, alisamento e weX2019ll vão de lá. O código para este capítulo está em chap08.ipynb. Que está no repositório deste livro (veja SectionXA0 0.2). Você também pode vê-lo em tinyurlthinkdsp08. 8.1XA0XA0Smoothing Figura 8.1: Preço de fechamento diário do estoque do Facebook e uma média móvel de 30 dias. O suavização é uma operação que tenta remover variações de curto prazo de um sinal para revelar tendências a longo prazo. Por exemplo, se você traçar mudanças diárias no preço de uma ação, pareceria barulhento que um operador de suavização talvez tornasse mais fácil ver se o preço geralmente aumentava ou diminuía ao longo do tempo. Um algoritmo de suavização comum é uma média móvel, que calcula a média dos valores n anteriores, para algum valor de n. Por exemplo, o FigureXA0 8.1 mostra o preço de fechamento diário do Facebook de 17 de maio de 2017 a 8 de dezembro de 2017. A linha cinza é o dado bruto, a linha mais escura mostra a média móvel de 30 dias. O alisamento remove as mudanças mais extremas e torna mais fácil ver as tendências a longo prazo. As operações de suavização também se aplicam aos sinais de som. Como exemplo, IX2019ll começa com uma onda quadrada a 440 Hz. Como vimos na Seção XA0 2.2. Os harmônicos de uma onda quadrada caem devagar, portanto, ele contém muitos componentes de alta freqüência. O primeiro IX2019ll constrói o sinal e duas ondas: a onda é uma fatia de 1 segundo do segmento de sinal é uma fatia mais curta IX2019ll para plotar. Para calcular a média móvel deste sinal, IX2019ll usa uma janela semelhante àquelas na seção XA0 3.7. Anteriormente usamos uma janela Hamming para evitar vazamentos espectrales causados ​​pela descontinuidade no início e no final de um sinal. Mais geralmente, podemos usar o Windows para calcular a soma ponderada de amostras em uma onda. Por exemplo, para calcular uma média móvel, IX2019ll crie uma janela com 11 elementos e normalize-a para que os elementos se somem a 1. Agora eu posso calcular a média dos primeiros 11 elementos, multiplicando a janela pela matriz de ondas: preenchido é um Versão da janela com zeros adicionados ao final por isso é o mesmo comprimento que segment. ys. Adicionar zeros como este é chamado de preenchimento. Prod é o produto da janela e da matriz de ondas. A soma dos produtos elementwise é a média dos primeiros 11 elementos da matriz. Como esses elementos são todos -1, sua média é -1. Figura 8.2: Um sinal quadrado a 400 Hz (cinza) e uma média móvel de 11 elementos. Para calcular o próximo elemento da média móvel, rolaremos a janela, que desloca os da direita e envolve um dos zeros do fim ao início. Quando multiplicamos a janela enrolada e a matriz de ondas, obtemos a média dos próximos 11 elementos da matriz de ondas, começando pela segunda. O resultado é -1 novamente. Podemos calcular o resto dos elementos da mesma maneira. A seguinte função envolve o código que vimos até agora em um loop e armazena os resultados em uma matriz. Suavizado é a matriz que irá conter os resultados preenchidos é uma matriz que contém a janela e zeros suficientes para ter o comprimento N e rolado é uma cópia do acolchoado que é deslocado para a direita por um elemento cada vez que o loop. No interior do loop, multiplicamos os dados por rolado para selecionar 11 elementos e adicioná-los. FiguraXA0 8.2 mostra o resultado para uma onda quadrada. A linha cinza é o sinal original, a linha mais escura é o sinal suavizado. O sinal suavizado começa a acelerar quando a borda principal da janela atinge a primeira transição, e se estabiliza quando a janela atravessa a transição. Como resultado, as transições são menos abruptas, e os cantos são menos afiados. Se você ouvir o sinal suavizado, ele parece menos buzzy e ligeiramente abafado. 8.2XA0XA0Convolução A operação que acabamos de executar X2017 aplicando uma função de janela para cada segmento sobreposto de uma onda X2017 é chamada de convolução. Convolução é uma operação tão comum que o NumPy fornece uma implementação mais simples e rápida do que a minha versão: np. convolve calcula a convolução da matriz de ondas e a janela. A bandeira de modo válida indica que ele deve calcular apenas valores quando a janela e a matriz de ondas se sobrepõem completamente, então ele pára quando a borda direita da janela atinge o final da matriz de ondas. Além disso, o resultado é o mesmo que no FigureXA0 8.2. Na verdade, há outra diferença. O loop na seção anterior calcula a correlação cruzada. (F XA0X22C6XA0 g) n XA0XA0 O suavização torna as transições em um sinal quadrado menos abrupto e faz o som um pouco abafado. LetX2019s vêem o efeito que esta operação tem no espectro. O primeiro IX2019ll traça o espectro da onda original: então a onda suavizada: a mesma bandeira do modo indica que o resultado deve ter o mesmo comprimento que a entrada. Neste exemplo, ele incluirá alguns valores que X201Cargue aroundX201D, mas thatX2019s está certo por enquanto. FiguraXA0 8.3 mostra o resultado. A freqüência fundamental é quase inalterada, os primeiros harmônicos são atenuados e os harmônicos mais altos são quase eliminados. Portanto, o alisamento tem o efeito de um filtro passa-baixa, que vimos na Seção XA0 1.5 e SectionXA0 4.4. Para ver o quanto cada componente foi atenuado, podemos calcular a proporção dos dois espectros: a proporção é a proporção da amplitude antes e depois do suavização. Quando os amplificadores são pequenos, essa relação pode ser grande e barulhenta, então, por simplicidade, eu ajuste a razão para 0, exceto onde os harmônicos são. Figura 8.4: Ratio de espectros para a onda quadrada, antes e depois do alisamento. FiguraXA0 8.4 mostra o resultado. Como esperado, a proporção é alta para baixas freqüências e cai em uma freqüência de corte perto de 4000 Hz. Mas há outra característica que não esperávamos: acima do ponto de corte, a proporção salta em torno de 0 a 0,2. WhatX2019 com esse 8.4XA0XA0. Teorema da convolução Figura 8.5: Relação dos espectros para a onda quadrada, antes e depois do suavização, juntamente com a DFT da janela de suavização. A resposta é o teorema da convolução. Datado matematicamente: DFT (f XA0X2217XA0 g) XA0XA0 DFT (f) XA0XB7XA0 DFT (g) XA0 onde f é uma matriz de ondas e g é uma janela. Em palavras, o Teorema da Convolução diz que se convolver f e g. E depois calculamos a DFT, obtemos a mesma resposta que a computação da DFT de f e g. E, em seguida, multiplicando os resultados por elementos. Quando aplicamos uma operação como convolução a uma onda de onda, dizemos que estamos trabalhando no domínio do tempo. Porque a onda é uma função do tempo. Quando aplicamos uma operação como multiplicação para a DFT, estamos trabalhando no domínio da freqüência. Porque a DFT é uma função de freqüência. Usando esses termos, podemos indicar o teorema da convolução de forma mais concisa: a convolução no domínio do tempo corresponde à multiplicação no domínio da freqüência. E isso explica FigureXA0 8.4. Porque quando convolvemos uma onda e uma janela, multiplicamos o espectro da onda com o espectro da janela. Para ver como isso funciona, podemos calcular o DFT da janela: acolchoado contém a janela de suavização, preenchida com zeros para ter o mesmo comprimento que a onda dftwindow contém o DFT de preenchido. FiguraXA0 8.5 mostra o resultado, juntamente com as proporções que calculamos na seção anterior. As proporções são exatamente as amplitudes em dftwindow. Matemmaticamente: abs (DFT (f XA0X2217XA0 g)) XA0XA0 abs (DFT (f)) XA0XA0 abs (DFT (g)) XA0 Neste contexto, o DFT de uma janela é chamado de filtro. Para qualquer janela de convolução no domínio do tempo, existe um filtro correspondente no domínio da frequência. E para qualquer filtro que possa ser expresso pela multiplicação de elementos no domínio da freqüência, há uma janela correspondente. 8.5XA0XA0Gaussian filter Na Seção XA0 8.2 Apresentamos definições de correlação cruzada e convolução, e vimos que são quase iguais, exceto que em convolução a janela é invertida. Agora que temos um algoritmo eficiente para convolução, também podemos usá-lo para calcular correlações cruzadas e autocorrelações. Usando os dados da seção anterior, podemos calcular os preços das ações de autocorrelação no Facebook: Com modeX2019sameX2019. O resultado tem o mesmo comprimento que o próximo. Correspondendo a atrasos de X2212 N 2 para N 2X22121. A linha cinza no Figura 8.8 mostra o resultado. Exceto em lag0. Não há picos, então não existe um comportamento periódico aparente nesse sinal. No entanto, a função de autocorrelação cai devagar, sugerindo que esse sinal se assemelha ao ruído rosa, como vimos na Seção XA0 5.3. Para calcular a autocorrelação usando a convolução, temos que ligar o sinal para dobrar o comprimento. Este truque é necessário porque a FFT baseia-se no pressuposto de que o sinal é periódico, que se envolve do final para o início. Com dados de séries temporais como este, essa suposição é inválida. Adicionando zeros, e depois cortando os resultados, remove os valores falsos. Além disso, lembre-se de que a convolução inverte a direção da janela. Para cancelar esse efeito, invertimos a direção da janela antes de chamar fftconvolve. Usando np. flipud. Que virou uma matriz NumPy. O resultado é uma visão da matriz, não uma cópia, então esta operação é rápida. O resultado de fftconvolve tem o comprimento 2 N. Daquelas, o primeiro e o último N 2 são válidos, o resto é o resultado de zero-padding. Para selecionar o elemento válido, rolo os resultados e selecione o primeiro N. correspondente aos atrasos de X2212 N 2 para N 2X22121. Conforme mostrado em FiguraXA0 8.8, os resultados de fftautocorr e np. correlate são idênticos (com cerca de 9 dígitos de precisão). Observe que as correlações em FigureXA0 8.8 são números grandes que podemos normalizá-los (entre -1 e 1) como mostrado na Seção XA0 5.6. A estratégia que usamos aqui para auto-correlação também funciona para correlação cruzada. Mais uma vez, você precisa preparar os sinais, folheando um e preenchendo ambos, e então você deve cortar as partes inválidas do resultado. Esse preenchimento e aparar é um incômodo, mas isso é o que as plataformas como o NumPy fornecem funções para fazê-lo por você. 8.8XA0XA0Exercícios As soluções para esses exercícios estão em chap08soln. ipynb. ExerciseXA01 XA0XA0 O notebook para este capítulo é chap08.ipynb. Leia-o e execute o código. Ele contém um widget interativo que permite experimentar os parâmetros da janela gaussiana para ver o efeito que eles têm na freqüência de corte. O que dá errado quando você aumenta a largura do Gaussian, std, sem aumentar o número de elementos na janela, M ExerciseXA02 XA0XA0 Neste capítulo, afirmei que a transformada de Fourier de uma curva gaussiana é também uma curva Gaussiana. Para Transformações de Fourier Discretas, essa relação é aproximadamente verdadeira. Experimente alguns exemplos. O que acontece com a transformada de Fourier à medida que você varia. Exercício XA03 XA0XA0 Se você fez os exercícios no Capítulo XA0 3, você viu o efeito da janela Hamming, e algumas das outras janelas fornecidas pelo NumPy, em vazamento espectral. Podemos obter alguma visão sobre o efeito dessas janelas, observando as DFTs. Além da janela Gaussiana usada neste capítulo, crie uma janela Hamming com o mesmo tamanho. Zero-pad as janelas e traçar seus DFTs. Qual janela atua como um melhor filtro de passagem baixa Você pode achar útil traçar os DFTs em uma escala de log. Experimente com alguns diferentes janelas e alguns tamanhos diferentes. Você está usando um dos nossos livros em uma classe, Wed, como para saber sobre isso. Por favor, considere preencher esta pequena pesquisa.