Tuesday 28 November 2017

Moving Average Filter Design In Matlab


Criado em quarta-feira, 08 de Outubro de 2008 20:04 Última atualização em Quinta, 14 Março 2013 01:29 Escrito por Batuhan Osmanoglu Hits: 38931 Média móvel Em Matlab Muitas vezes eu me encontro na necessidade de calcular a média dos dados que tenho para reduzir o ruído um pouco pouco. Eu escrevi funções de casal para fazer exatamente o que eu quero, mas matlabs construído em função de filtro funciona muito bem também. Aqui Ill escrever sobre 1D e 2D média dos dados. 1D filtro pode ser realizado usando a função de filtro. A função de filtro requer pelo menos três parâmetros de entrada: o coeficiente de numerador para o filtro (b), o coeficiente do denominador para o filtro (a) e os dados (X), é claro. Um filtro de média em execução pode ser definido simplesmente por: Para dados 2D, podemos usar a função Matlabs filter2. Para obter mais informações sobre como o filtro funciona, você pode digitar: Aqui está uma implementação rápida e suja de um filtro de média móvel 16 por 16. Primeiro precisamos definir o filtro. Uma vez que tudo o que queremos é a contribuição igual de todos os vizinhos, podemos apenas usar a função uns. Nós dividimos tudo com 256 (1616) desde que nós não queremos mudar o nível geral (amplitude) do sinal. Para aplicar o filtro podemos simplesmente dizer o seguinte Abaixo estão os resultados para a fase de um interferograma SAR. Neste caso, Range está no eixo Y e Azimuth é mapeado no eixo X. O filtro tinha 4 pixels de largura em Alcance e 16 pixels de largura em Azimute. Login SearchUsing MATLAB, como posso encontrar a média móvel de 3 dias de uma coluna específica de uma matriz e anexar a média móvel para a matriz que eu estou tentando calcular a média móvel de 3 dias de baixo para cima da matriz. Eu forneci o meu código: Dada a seguinte matriz a e máscara: Eu tentei implementar o comando conv, mas estou recebendo um erro. Aqui está o comando conv que eu tenho tentado usar na segunda coluna da matriz a: A saída que desejo é dada na seguinte matriz: Se você tiver alguma sugestão, eu gostaria muito. Obrigado Para a coluna 2 da matriz a, estou computando a média móvel de 3 dias da seguinte maneira e colocando o resultado na coluna 4 da matriz a (I renomeado matriz a como 39desiredOutput39 apenas para ilustração). A média de 3 dias de 17, 14, 11 é 14 a média de 3 dias de 14, 11, 8 é 11 a média de 3 dias de 11, 8, 5 é 8 ea média de 3 dias de 8, 5, 2 é 5. Não há valor nas 2 linhas inferiores para a 4ª coluna porque a computação para a média móvel de 3 dias começa na parte inferior. A saída 39valid39 não será mostrada até pelo menos 17, 14 e 11. Espero que isso faz sentido ndash Aaron Jun 12 13 em 1:28 1 Resposta Em geral, seria útil se você mostrar o erro. Neste caso você está fazendo duas coisas erradas: Primeiro, sua convolução precisa ser dividida por três (ou o comprimento da média móvel) Segundo, observe o tamanho de c. Você não pode apenas caber c em um. A maneira típica de obter uma média móvel seria usar o mesmo: mas isso não se parece com o que você quer. Em vez disso, você é forçado a usar um par de linhas: filtros FIR, filtros IIR e equação de diferença de coeficiente constante linear Filtros de média móvel causal (FIR) Nós discutimos sistemas em que cada amostra da saída é uma soma ponderada de As) amostras da entrada. Vamos tomar um sistema de soma ponderada causal, onde causal significa que uma dada amostra de saída depende apenas da amostra de entrada atual e outros insumos mais cedo na seqüência. Nem os sistemas lineares em geral, nem os sistemas finitos de resposta ao impulso em particular, precisam ser causais. No entanto, a causalidade é conveniente para um tipo de análise que iria explorar em breve. Se simbolizamos as entradas como valores de um vetor x. E as saídas como valores correspondentes de um vetor y. Então tal sistema pode ser escrito como onde os valores de b são quotweights aplicados às amostras de entrada atuais e anteriores para obter a amostra de saída atual. Podemos pensar na expressão como uma equação, com o sinal de igual signo igual a, ou como uma instrução processual, com o sinal de igual significação atribuição. Vamos escrever a expressão para cada amostra de saída como um loop MATLAB de instruções de atribuição, onde x é um vetor N-comprimento de amostras de entrada, e b é um vetor M-comprimento de pesos. A fim de lidar com o caso especial no início, vamos incorporar x em um vetor mais longo xhat cujas primeiras M-1 amostras são zero. Vamos escrever a soma ponderada para cada y (n) como um produto interno, e faremos algumas manipulações das entradas (como inverter b) para este fim. Esse tipo de sistema é muitas vezes chamado de filtro de média móvel, por razões óbvias. De nossas discussões anteriores, deve ser óbvio que tal sistema é linear e invariante ao deslocamento. Claro, seria muito mais rápido usar a convolução de função MATLAB conv () em vez do nosso mafilt (). Em vez de considerar as primeiras M-1 amostras da entrada de ser zero, poderíamos considerá-los a ser o mesmo que as últimas M-1 amostras. Isso é o mesmo que tratar a entrada como periódica. Bem, use cmafilt () como o nome da função, uma pequena modificação da função mafilt () anterior. Na determinação da resposta de impulso de um sistema, não há geralmente nenhuma diferença entre estes dois, desde que todas as amostras não-iniciais da entrada são zero: Uma vez que um sistema deste tipo é linear e shift-invariante, sabemos que seu efeito em qualquer Sinusoid será apenas a escala e deslocá-lo. Aqui é importante que usemos a versão circular A versão circularmente convoluta é deslocada e escalada um pouco, enquanto a versão com convolução ordinária é distorcida no início. Vamos ver o que a escala exata e deslocamento é usando um fft: Tanto a entrada ea saída têm amplitude apenas nas freqüências 1 e -1, que é como deveria ser, uma vez que a entrada era uma sinusoid eo sistema era linear. Os valores de saída são maiores em uma proporção de 10,6251 / 8 1,3281. Este é o ganho do sistema. E quanto à fase Nós só precisamos olhar onde a amplitude é diferente de zero: A entrada tem uma fase de pi / 2, como nós pedimos. A fase de saída é deslocada por um 1,0594 adicional (com sinal oposto para a freqüência negativa), ou cerca de 1/6 de um ciclo à direita, como podemos ver no gráfico. Agora vamos tentar uma sinusoid com a mesma freqüência (1), mas em vez de amplitude 1 e fase pi / 2, vamos tentar amplitude 1,5 e fase 0. Sabemos que apenas a freqüência 1 e -1 terá amplitude não-zero, então vamos Basta olhar para eles: Novamente a relação de amplitude (15.9377 / 12.0000) é 1.3281 - e quanto à fase é novamente deslocado por 1.0594 Se esses exemplos são típicos, podemos prever o efeito do nosso sistema (resposta ao impulso .1.2 .3 .4 .5) em qualquer sinusoide com frequência 1 - a amplitude será aumentada em um fator de 1,3281 e a fase (freqüência positiva) será deslocada em 1,0594. Poderíamos continuar a calcular o efeito desse sistema sobre sinusóides de outras freqüências pelos mesmos métodos. Mas há uma maneira muito mais simples, e uma que estabelece o ponto geral. Dado que a circunvolução (circular) no domínio do tempo significa a multiplicação no domínio da frequência, daí decorre que, por outras palavras, a DFT da resposta de impulso é a razão da DFT da saída para a DFT da entrada. Nesta relação os coeficientes de DFT são números complexos. Desde abs (c1 / c2) abs (c1) / abs (c2) para todos os números complexos c1, c2, esta equação nos diz que o espectro de amplitude da resposta de impulso será sempre a relação do espectro de amplitude da saída para que Da entrada. No caso do espectro de fase, ângulo (c1 / c2) ângulo (c1) - ângulo (c2) para todos os c1, c2 (com a ressalva de que as fases que diferem por n2pi são considerados iguais). Portanto, o espectro de fase da resposta ao impulso será sempre a diferença entre os espectros de fase da saída e da entrada (com quaisquer correções de 2pi são necessárias para manter o resultado entre - pi e pi). Podemos ver os efeitos de fase mais claramente se desempacotar a representação da fase, isto é, se adicionarmos vários múltiplos de 2pi conforme necessário para minimizar os saltos que são produzidos pela natureza periódica da função ângulo (). Embora a amplitude e a fase sejam normalmente utilizadas para apresentação gráfica e mesmo tabular, uma vez que são uma maneira intuitiva de pensar sobre os efeitos de um sistema sobre os vários componentes de frequência de sua entrada, os coeficientes de Fourier complexos são mais úteis algébricamente, A expressão simples da relação A abordagem geral que acabamos de ver funcionará com filtros arbitrários do tipo esboçado, em que cada amostra de saída é uma soma ponderada de algum conjunto de amostras de entrada. Como mencionado anteriormente, estes são freqüentemente chamados filtros de resposta de impulso finito, porque a resposta ao impulso é de tamanho finito, ou às vezes filtros de média móvel. Podemos determinar as características de resposta de freqüência de tal filtro da FFT de sua resposta de impulso e também podemos projetar novos filtros com características desejadas por IFFT a partir de uma especificação da resposta de freqüência. Filtros Autoregressivos (IIR) Não haveria nenhum ponto em ter nomes para filtros FIR a menos que houvesse algum outro tipo de distinção, de modo que aqueles que estudaram pragmática não ficarão surpresos ao saber que existe de fato outro tipo principal Do filtro tempo-invariante linear. Estes filtros são às vezes chamados recursivos porque o valor das saídas anteriores (assim como entradas anteriores) importa, embora os algoritmos sejam geralmente escritos usando construções iterativas. Eles também são chamados filtros Infinite Impulse Response (IIR), porque em geral sua resposta a um impulso continua para sempre. Eles também são às vezes chamados de filtros auto-regressivos, porque os coeficientes podem ser considerados como o resultado de fazer uma regressão linear para expressar valores de sinal em função de valores de sinal anteriores. A relação dos filtros FIR e IIR pode ser vista claramente numa equação de diferença de coeficiente constante linear, isto é, estabelecendo uma soma ponderada de saídas igual a uma soma ponderada de entradas. Isto é como a equação que damos anteriormente para o filtro causal FIR, exceto que, além da soma ponderada de insumos, também temos uma soma ponderada de saídas. Se quisermos pensar nisso como um procedimento para gerar amostras de saída, precisamos reorganizar a equação para obter uma expressão para a amostra de saída corrente y (n), Adotando a convenção de que a (1) 1 (por exemplo, escalando outros como E bs), podemos nos livrar do termo 1 / a (1): y (n) b (1) x (n) b (2) x (n-1). B (Nb1) x (n-nb) - a (2) y (n-1) -. - a (Na1) y (n-na) Se todos os a (n) diferentes de a (1) são zero, isso reduz a nosso velho amigo o filtro FIR causal. Este é o caso geral de um filtro (causal) LTI, e é implementado pelo filtro de função MATLAB. Vejamos o caso em que os coeficientes b diferentes de b (1) são zero (em vez do caso FIR, onde a (n) são zero): Neste caso, a amostra de saída corrente y (n) é calculada como um (N-1), y (n-2), etc. Para ter uma idéia do que acontece com esses filtros, vamos começar com o caso em que: Ou seja, a amostra de saída atual é a soma da amostra de entrada atual e metade da amostra de saída anterior. Bem, tome um impulso de entrada através de alguns passos de tempo, um de cada vez. Deve ficar claro neste ponto que podemos facilmente escrever uma expressão para o n-ésimo valor de amostra de saída: é apenas (se MATLAB contado a partir de 0, isso seria simplesmente .5n). Como o que estamos calculando é a resposta ao impulso do sistema, demonstrámos por exemplo que a resposta ao impulso pode de fato ter infinitas amostras diferentes de zero. Para implementar esse filtro trivial de primeira ordem no MATLAB, poderíamos usar o filtro. A chamada será assim: eo resultado é: Este negócio é realmente ainda linear Podemos olhar para isto empiricamente: Para uma abordagem mais geral, considere o valor de uma amostra de saída y (n). Por substituição sucessiva poderíamos escrever isto como Isto é exatamente como nosso velho amigo a forma convolução-soma de um filtro FIR, com a resposta de impulso fornecida pela expressão .5k. E o comprimento da resposta ao impulso é infinito. Assim, os mesmos argumentos que usamos para mostrar que os filtros FIR eram lineares agora se aplicam aqui. Até agora isso pode parecer um monte de barulho por não muito. O que é toda esta linha de investigação bom para Bem responder esta questão em etapas, começando com um exemplo. Não é uma grande surpresa que possamos calcular uma amostra exponencial por multiplicação recursiva. Vamos olhar para um filtro recursivo que faz algo menos óbvio. Este tempo bem torná-lo um filtro de segunda ordem, de modo que a chamada para filtro será da forma Permite definir o segundo coeficiente de saída a2 para -2cos (2pi / 40), eo terceiro coeficiente de saída a3 para 1, e olhar para A resposta ao impulso. Não muito útil como um filtro, na verdade, mas ele gera uma onda senoidal amostrada (de um impulso) com três multiplicações por amostra. Para entender como e por que ele faz isso, e como os filtros recursivos podem ser projetados e analisados ​​em O caso mais geral, precisamos dar um passo atrás e dar uma olhada em algumas outras propriedades de números complexos, no caminho para a compreensão da transformada z. filtro de filtro de documentação y (b, a, x) filtra os dados de entrada x usando um racional Função de transferência definida pelos coeficientes numerador e denominador b e a. Se a (1) não for igual a 1. então o filtro normaliza os coeficientes do filtro por a (1). Portanto, a (1) deve ser diferente de zero. Se x é um vetor, então o filtro retorna os dados filtrados como um vetor do mesmo tamanho que x. Se x é uma matriz, o filtro age ao longo da primeira dimensão e retorna os dados filtrados para cada coluna. Se x é uma matriz multidimensional, então o filtro age ao longo da primeira dimensão da matriz cujo tamanho não é igual a 1. y filter (b, a, x, zi) usa as condições iniciais zi para os atrasos do filtro. O comprimento de zi deve ser igual a max (comprimento (a), comprimento (b)) - 1. Y (b, a, x, zi, dim) atua ao longo da dimensão dim. Por exemplo, se x é uma matriz, então filtro (b, a, x, zi, 2) retorna os dados filtrados para cada linha. Y, zf filter () também retorna as condições finais zf dos atrasos do filtro, usando qualquer uma das sintaxes anteriores. Função de Transferência Racional A descrição entrada-saída da operação de filtro em um vetor no domínio de transformação Z é uma função de transferência racional. Uma função de transferência racional é da forma, Y (z) b (1) b (2) z x2212 1. B (n b 1) z x 2212 n b 1 a (2) z x 2212 1. A (n a 1) z x 2212 n a X (z). Que manipula os filtros FIR e IIR 1. n a é a ordem do filtro de realimentação e n b é a ordem do filtro feedforward. Você também pode expressar a função de transferência racional como a seguinte equação de diferença, a (1) y (n) b (1) x (n) b (2) x (n x2212 1). B (n b 1) x (n x 2212 n b) x 2212 a (2) y (n x 2212 1) x 2212. X2212 a (n a 1) y (n x2212 n a). Além disso, você pode representar a função de transferência racional usando sua implementação de transposição direta da forma II, como no diagrama a seguir. Devido à normalização, assumir a (1) 1. A operação do filtro na amostra m é dada pelas equações de diferença no domínio do tempo y (m) b (1) x (m) z 1 (m x 2212 1) z 1 (m) b (2) x (m) z 2 m x2212 1) x2212 um (2) Y (m) x00A0x00A0 x22EE x00A0x00A0x00A0x00A0x00A0 x00A0x00A0x00A0x00A0x00A0 x22EE x00A0x00A0x00A0x00A0x00A0x00A0x00A0x00A0x00A0x00A0x00A0x00A0 x22EE Zn x2212 2 (m) b (n x2212 1) x (m) x2212 Zn 1 (m x2212 1) x2212 um (N x2212 1 ) Y (m) zn x2212 1 (m) b (n) x (m) x2212 a (n) y (m). Dicas Se você tiver o Processamento de Sinal Toolboxx2122, você pode projetar um filtro, d. Usando designfilt. Em seguida, você pode usar o filtro Y (d, X) para filtrar seus dados. Selecione seu PaísUsando MATLAB, como posso encontrar a média móvel de 3 dias de uma coluna específica de uma matriz e anexar a média móvel àquela matriz Estou tentando calcular a média móvel de 3 dias de baixo para cima da matriz. Eu forneci o meu código: Dada a seguinte matriz a e máscara: Eu tentei implementar o comando conv, mas estou recebendo um erro. Aqui está o comando conv que eu tenho tentado usar na segunda coluna da matriz a: A saída que desejo é dada na seguinte matriz: Se você tiver alguma sugestão, eu gostaria muito. Obrigado Para a coluna 2 da matriz a, estou computando a média móvel de 3 dias da seguinte maneira e colocando o resultado na coluna 4 da matriz a (I renomeado matriz a como 39desiredOutput39 apenas para ilustração). A média de 3 dias de 17, 14, 11 é 14 a média de 3 dias de 14, 11, 8 é 11 a média de 3 dias de 11, 8, 5 é 8 ea média de 3 dias de 8, 5, 2 é 5. Não há valor nas 2 linhas inferiores para a 4ª coluna porque a computação para a média móvel de 3 dias começa na parte inferior. A saída 39valid39 não será mostrada até pelo menos 17, 14 e 11. Espero que isso faz sentido ndash Aaron Jun 12 13 em 1:28 1 Resposta Em geral, seria útil se você mostrar o erro. Neste caso você está fazendo duas coisas erradas: Primeiro, sua convolução precisa ser dividida por três (ou o comprimento da média móvel) Segundo, observe o tamanho de c. Você não pode apenas caber c em um. A maneira típica de obter uma média móvel seria usar o mesmo: mas isso não se parece com o que você quer. Em vez disso, você é forçado a usar um par de linhas:

No comments:

Post a Comment