Movimento média filtro matlab conv no Brasil


29 Setembro, 2017 Média móvel por convolução O que é média móvel e para que é bom Como a média móvel é feita usando a convolução Média móvel é uma operação simples usada geralmente para suprimir o ruído de um sinal: ajustamos o valor de cada ponto para a Média dos valores em sua vizinhança. Por uma fórmula: Aqui x é a entrada ey é o sinal de saída, enquanto o tamanho da janela é w, suposto ser ímpar. A fórmula acima descreve uma operação simétrica: as amostras são tomadas de ambos os lados do ponto real. Abaixo está um exemplo da vida real. O ponto em que a janela é colocada realmente é vermelho. Valores fora de x são supostos ser zeros: Para brincar e ver os efeitos da média móvel, dê uma olhada nesta demonstração interativa. Como fazê-lo por convolução Como você pode ter reconhecido, o cálculo da média móvel simples é semelhante à convolução: em ambos os casos, uma janela é deslizada ao longo do sinal e os elementos na janela são resumidos. Então, dar-lhe uma tentativa de fazer a mesma coisa usando convolução. Use os seguintes parâmetros: A saída desejada é: Como primeira aproximação, vamos tentar o que obtemos ao converter o sinal x pelo k kernel seguinte: A saída é exatamente três vezes maior do que o esperado. Também pode ser visto que os valores de saída são o resumo dos três elementos na janela. É porque durante a convolução a janela é deslizada ao longo, todos os elementos nele são multiplicados por um e, em seguida, resumido: yk 1 cdot x 1 cdot x 1 cdot x Para obter os valores desejados de y. A saída deve ser dividida por 3: Por uma fórmula incluindo a divisão: Mas não seria ótimo para fazer a divisão durante convolução Aqui vem a idéia, reorganizando a equação: Então vamos usar o k kernel seguinte: Desta forma, vamos Obter a saída desejada: Em geral: se queremos fazer a média móvel por convolução tendo um tamanho de janela de w. Usaremos o k kernel seguinte: Uma função simples que faz a média móvel é: Um exemplo de uso é: Preciso calcular uma média móvel em uma série de dados, dentro de um loop for. Eu tenho que obter a média móvel em N9 dias. O array Im computing in é 4 séries de 365 valores (M), que são valores médios de outro conjunto de dados. Eu quero traçar os valores médios dos meus dados com a média móvel em um gráfico. Eu pesquisei um pouco sobre as médias móveis eo comando conv e encontrei algo que eu tentei implementar no meu código. Então, basicamente, eu computo o meu médio e plotá-lo com uma média móvel (errada). Eu escolhi o valor de wts fora do site mathworks, de modo que está incorreto. (Fonte: mathworks. nlhelpeconmoving-average-trend-estimation. html) Meu problema, porém, é que eu não entendo o que este wts é. Alguém poderia explicar Se tem algo a ver com os pesos dos valores: que é inválido neste caso. Todos os valores são ponderados da mesma forma. E se eu estou fazendo isso inteiramente errado, eu poderia obter alguma ajuda com ele Meus mais sinceros agradecimentos. September 23 14 at 19:05 Usando conv é uma excelente maneira de implementar uma média móvel. No código que você está usando, wts é o quanto você está pesando cada valor (como você adivinhou). A soma desse vetor deve ser sempre igual a um. Se você deseja pesar cada valor uniformemente e fazer um filtro de tamanho N em movimento, então você gostaria de fazer Usando o argumento válido em conv resultará em ter menos valores em Ms do que você tem em M. Use o mesmo se você não se importa os efeitos de Zero preenchimento. Se você tiver a caixa de ferramentas de processamento de sinal, você pode usar o cconv se quiser experimentar uma média móvel circular. Algo como Você deve ler a documentação conv e cconv para obter mais informações se você já não. Você pode usar o filtro para encontrar uma média em execução sem usar um loop for. Este exemplo encontra a média em execução de um vetor de 16 elementos, usando um tamanho de janela de 5. 2) suave como parte da Caixa de Ferramentas de Ajuste de Curva (que está disponível na maioria dos casos) yy suave (y) suaviza os dados no vetor de coluna Y usando um filtro de média móvel. Os resultados são retornados no vetor de coluna yy. O intervalo padrão para a média móvel é 5.Time-domain filtragem em Matlab Este tutorial é sobre o comando Matlab ldquoconvrdquo. Para descobrir o que ldquoconvrdquo faz, você pode tentar ler a ajuda online do Matlab. E confira a explicação. (Você pode achar mais fácil copiar as linhas de comando em itálico para a janela de comandos do Matlab usando copiar e colar.) No entanto, a menos que você saiba muito sobre matemática de engenharia, provavelmente encontrará o ldquohelprdquo Matlab muito mais ldquoconvolutedrdquo do que ldquohelpfulrdquo. No caso improvável de que você não tenha esse problema, continue com o resto desses exercícios (caso contrário, veja se uma das pessoas que você gosta nessa classe parece perdida e poderia se beneficiar de suas idéias e seu charme, ou simplesmente ir para A praia por um pouco). Agora, ldquoconvolvingrdquo é algo que fazemos para ldquosignalsrdquo, e em sinais Matlab são vetores. Itrsquos muitas vezes uma boa idéia para jogar com muito pouco artificial sinal vetores em primeiro lugar para entender o que os comandos de processamento de sinal, antes de lançá-los em grandes sinais ldquorealrdquo. Em Matlab, vamos fazer um brinquedo ldquodatardquo vetor d com os seguintes valores: então tente o seguinte O que você começa O que você começa Convolvendo um vetor com um ldquoscalarrdquo (ou seja, um número simples ou vetor com apenas um elemento), ldquoscalesrdquo o vetor. (Agora você sabe por que os chamam de escalares). Agora, para algo um pouco mais misterioso. Tente cada um destes, por sua vez: Whatrsquos acontecendo aqui Pense em um vetor que é todos os zeros, exceto um único como um impulserdquo ldquounit. Agora, se ldquoconvolverdquo um sinal com um impulso unitário, então não mudamos o sinal (exceto, talvez, por ldquozero paddingrdquo o fim dele), mas podemos atrasá-lo (efetivamente colocando zeros na frente dele). Então, você pode adivinhar o que aconteceria se você tentar este Pense sobre isso primeiro, em seguida, experimentá-lo para ver se você estava correto. O que você deve ter encontrado é que você pode combinar o dimensionamento e as propriedades de atraso do comando conv, você deve ter sinal duas vezes maior que d, atrasado por amostras (ou ldquotapsrdquo). Antes que possamos entender o que tudo isso tem a ver com a filtragem, precisamos de mais uma coisa. Nada de novo aqui, apenas atribuímos nosso ldquoimpulsesrdquo de variável atraso para variáveis ​​Mas agora, se eu fizer E o que, portanto, de O que você deve ver é que ab tem dois impulsos, e que a convolução, portanto, contém duas cópias do sinal d, Um com atraso zero, e um atrasado por seis torneiras. Se fizermos o atraso menor e menor, podemos fazer as cópias de reboque do sinal ldquocolliderdquo, ou sobrepor. Tente cada um destes, por sua vez: Você deve ver que, como as cópias atrasadas do sinal ldquocolliderdquo, eles acrescentam. Em resumo, convolução combina escala, atraso e adição em um pacote acessível. Agora, o que tem tudo isso a ver com a filtragem Letrsquos considerar um exemplo relativamente simples: o filtro de média móvel. Imagine que eu quero obter um sinal filtrado f por filtragem de média móvel de dois pontos no vetor d d (1) d (2) d (3). . O que isto significa, quero que o primeiro ponto de f seja a média dos dois primeiros pontos de d, isto é, f (1) d (1) 2 d (2) 2. Similarmente, eu quero f (2) d (2) 2d (3) 2 e assim por diante até que f (n) d (n) 2d (n1) 2. Se você pensar nisso, você pode notar que f é simplesmente a soma de duas cópias de d, com a segunda cópia adiada por uma torneira, e as duas cópias escaladas por uma metade, e em vez de fazer isso um elemento de cada vez, nós Deixe conv fazê-lo para nós dar-lhe-á a média móvel de dois pontos filtrada (suavizada) versão de d. Agora, supostamente, este alisamento pode ajudar a reduzir ldquonoiserdquo. Vamos tentar isso. Primeiro, letrsquos fazer um sinal sinusoidal agradável: Agora nós contaminar o sinal com algum ruído gaussian feio Alisamento com conv recuperar o sinal limpo Letrsquos experimentá-lo: Bem, olhando para o enredo é claro que o filtro de média móvel de dois pontos melhorou matérias Apenas um pouquinho. Como sobre um filtro ldquomore aggressiverdquo Tente a média móvel de cinco pontos (By the way, para filtros de média móvel de n pontos mais longos do que n5, escrevê-los para fora torna-se explicitamente tedioso, e você vai achar mais fácil escrever uns (1, n) n Ao invés do explícito. Isso removeu um pouco do ruído (não significa que tudo isso, mas não é direto para fazer muito melhor, é por isso que é bom tentar ter o menor ruído Como possível para começar), mas se você olhar cuidadosamente você verá que o sinal filtrado é realmente tempo deslocado em relação ao sinal original (e é um pouco mais). Os engenheiros elétricos diriam que nosso filtro de média de cinco movimentos teve um retardo de duas torneiras e, tendo em vista os exercícios anteriores com impulsos retardados, você pode ter um pouco de intuição a respeito de onde esses atrasos vêm. Estes atrasos podem ser um incômodo. Imagine que você queria medir os tempos quando os picos ocorrem. A filtragem tornaria os picos mais limpos, mas os atrasaria. Você poderia se livrar dos atrasos, jogando fora os primeiros n2 pontos do sinal de filtragem de média móvel de n pontos. Matlab tem realmente uma variante inteligente em conv chamado filtfilt que evita os atrasos executando a filtragem ldquotwicerdquo, uma vez ldquoforwardrdquo então ldquobackwardrdquo. A filtragem é, então, duas vezes mais ldquoeffectiverdquo muito embora. Tente isso: Não há atraso eo sinal filtrado é bastante suave, mas itrsquos realmente começam a ser reduzidos em amplitude. O filtro está começando a comer no sinal, bem como o ruído. Talvez cinco pontos de média móvel é muito se aplicada duas vezes. Não há uma única resposta correta para isso, e para obter uma adequada apreciação de quanto um determinado filtro afetará quais freqüências requerem uma compreensão do processamento de sinal que está além desses exercícios introdutórios simples. Você pode, no entanto, obter um pouco de uma compreensão intuitiva se você considerar o seguinte. Anteriormente você viu que a convolução fconv (s, a), onde a é n torneiras longas, significa que f (i) s (i) a (1) s (i1) a (2) s (i2) a (3) . S (em-1) a (n). Assim f (i) é o produto ponto ou o produto interno dos vetores a e s (i: em-1). Se você estudou álgebra linear, você pode lembrar que o produto ponto mede o ldquoprojectionrdquo de s (i: em-1) em a. Mesmo se você não estudou a álgebra linear, se você pensar nela, provavelmente verá que f (i) é efetivamente uma medida de quão bem o segmento de sinal s (i: em-1) está correlacionado com o filtro a. Assim, o i-ésimo ponto de tempo do sinal filtrado f (i) indica-nos o quão bem o segmento de sinal original s (i: em-1) corresponde à forma de onda do filtro a. Nos exemplos que vimos até agora, nossos filtros de média móvel eram todos os segmentos de linha reta. Claramente, flutuações rápidas em um sinal correspondem a uma linha reta mais mal do que flutuações lentas, então os segmentos de linha reta que usamos na filtragem média de corrida são exemplos de filtros passa-baixos. Eles são um pouco como flat ldquoironsrdquo, usado para passar pequenas rugas no sinal. E se você pensa sobre eles dessa forma, deve tornar-se intuitivamente claro que: um filtro de média móvel de n pontos irá suprimir fortemente flutuações cujas frequências têm um comprimento de onda menor que n, mas terá muito pouco efeito em freqüências com comprimentos de onda muito maiores que n . A convolução pode ser usada não apenas para filtragem de baixa passagem. Se fizermos um filtro que contém muitas flutuações rápidas (altas freqüências), então ele irá amplificar as flutuações rápidas e suprimir as baixas. Vamos tentar isso. Faça um n-ponto ldquohigh passar filterrdquo: O filtro a agora oscila o mais rápido possível na taxa de amostragem atual. O que ele vai fazer: o que você deve ver é que a amplitude da sinusoide é muito reduzida, mas as flutuações rápidas ldquosort ofrdquo permanecem (embora eles não são afetados) Como eu cheguei a1 -1 0 13 Um palpite cego, que é Não algo que você sempre fez se você quisesse filtrar sinais para a pesquisa científica (isso era apenas um ldquoproof de principlerdquo). Para fazer filtros ldquorealrdquo que você pode usar para sinais high-pass, band-pass ou low-pass ldquorealrdquo, o Matlab oferece um comando acessível chamado ldquofir2rdquo. Letrsquos ilustram seu uso com alguns dados reais do córtex auditivo. (Você pode precisar baixar o arquivo auditory. mat abaixo e copiá-lo para o seu diretório de trabalho Matlab) Apenas observando o sinal, vemos grandes flutuações em uma escala de cerca de 120 de segundo, provavelmente devido a entradas sinápticas, bem como Uma pequena quantidade de ruído de alta freqüência. Digamos que queremos passar o filtro de baixa passagem para que todas as freqüências abaixo de 1 kHz não sejam afetadas, mas as freqüências mais altas são suprimidas. O comando fir2 nos dará esse filtro, mas precisa de dois vetores, A e F, que lhe dizem que Amplitudes queremos em que Freqüências. Infelizmente, não podemos especificar apenas F em kHz, mas temos de fornecê-la como uma fração da freqüência de ldquoNyquistrdquo, que é igual a metade da nossa taxa de amostragem. Neste exemplo, a taxa de amostragem é 12000, portanto, se quisermos manter nossas freqüências até 2kHz plana, então queremos A1 para todos os 0ltFlt20006000 e A0 para 20006000ltFlt1. Por razões teóricas, não é uma boa idéia para tentar fazer o filtro demasiado ldquosteeprdquo, então em vez de dizer F0 20006000 20016000 1 A1 1 0 0 seria melhor aconselhados a usar agora letrsquos fazer um 20 ponto ldquofinite impulser responserdquo filtro para estes (Time, vfiltered, 39r39) hold off Utilize o comando zoom na janela de figura para plotar. Experimente com diferentes frequências de corte. Como você pode usar o comando fir2 para fazer um filtro passa-alto Como você pode fazer um filtro passa bandaFIR, filtros IIR ea equação de diferença de coeficiente constante linear Filtros de média móvel causal (FIR) Discutimos sistemas em que cada amostra Da saída é uma soma ponderada de (algumas das) 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 amostras M-1 da entrada como sendo zero, poderíamos considerá-las iguais às últimas amostras M-1. 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 numa razão de 10,62518 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 pi2, como nós pedimos. A fase de saída é deslocada por um adicional 1.0594 (com sinal oposto para a freqüência negativa), ou cerca de 16 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 pi2, 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 apenas olhar Para eles: Novamente a razão de amplitude (15.937712.0000) é 1.3281 - e quanto à fase é novamente deslocada por 1.0594 Se esses exemplos são típicos, podemos prever o efeito do nosso sistema (resposta de impulso .1 .2 .3 .4 .5) em qualquer sinusoide com freqüê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 (c1c2) 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 entre o espectro de amplitude da saída para a da entrada . No caso do espectro de fase, ângulo (c1c2) ângulo (c1) - ângulo (c2) para todos os c1, c2 (com a condição de que as fases diferentes 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 desempacotarmos 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 forma intuitiva de pensar 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 muitas vezes 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 a partir 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 de 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 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 como uma 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 atual y (n), Adotando a convenção de que a (1) 1 (por exemplo, escalando outros como E bs), podemos nos livrar do termo 1a (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: Isto é, a amostra de saída atual é a soma da amostra de entrada corrente 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). Uma vez que o que estamos calculando é a resposta ao impulso do sistema, temos demonstrado 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 o nosso velho amigo a forma convolução-soma de um filtro FIR, com a resposta ao 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. Desta vez, torná-lo um filtro de segunda ordem, de modo que a chamada para filtrar será da forma Vamos definir o segundo coeficiente de saída a2 para -2cos (2pi40) eo terceiro coeficiente de saída a3 para 1, e olhar para o impulso resposta. Não é muito útil como um filtro, na verdade, mas gera uma onda senoidal amostrada (de um impulso) com três multiplicações por amostra. Para entender como e por que faz isso, e como 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.

Comments