Friday 25 August 2017

Moving Average Filter 3db


E eu uso a segunda resposta no meu algoritmo para calcular a frequência de corte 3dB do meu filtro, o que funciona de forma excelente, já que o tamanho do filtro geralmente está acima de 300. Verifiquei com a resposta de passo. Mas eu gostaria de ter uma fonte ou derivação para esta fórmula. Eu tentei à mão com a série taylor parando após o segundo e terceiro terços. Eu chego perto, mas não exatamente com a fórmula e o mapple me dá um resultado válido, mas extremamente complicado. Espero que vocês possam ajudar. E você não precisa aproximar qualquer soma disso com uma integral. Mas você precisa aproximar sin2 () com um número finito de termos da série Maclaurin. O que você precisa é uma solução exata para isso: 2 sin2 (omega0 N 2) N2 sin2 (omega0 2). E a resposta que tenho é, da melhor forma possível, a aproximação mais próxima, fazendo o menor pressuposto. Ndash robert bristow-johnson 13 de janeiro às 5:46 Considere uma média móvel de fase zero de comprimento N: filtros de comprimento par em operação em seqüências discretas com índices de tempo inteiro não podem ser fase zero. Nós contornamos isso, permitindo que os índices de tempo de saída sempre tenham uma fração de fração, no caso de N mesmo. Como um exemplo do mundo real, se a entrada fosse amostrada a cada meia-noite, a média móvel de fase zero de comprimento igual seria calculada para cada meio-dia. Esta indexação incomum fornece convenientemente a mesma forma de resposta de freqüência de fase zero FN (omega) para N ímpar e N mesmo: Infelizmente, a resposta de freqüência não possui uma solução simbólica para o omegac de freqüência de corte de -3 dB, de modo que: Strictly speaking sqrt is Cerca de -3.01 dB, mas acho que é isso que as pessoas querem dizer quando dizem -3 dB, porque, de outra forma, é apenas um número arbitrário. Uma resposta de freqüência aproximada hat N (omega) usa uma integral em vez de uma soma: os lobos principais das respostas de freqüência verdadeira (soma) e aproximada (integral) convergem em grande N: podemos provar a convergência através da introdução de funções GN (chi ) FN (omega) e chapéu N (chi) hat N (omega) com o argumento normalizado de forma que omega frac, trazendo o primeiro zero de ambas as funções para chi 1: GN (chi) é conhecido como N-periódico banda limitada Trem de impulso. Seu limite em geral N e o chapéu de função N (chi) são a função de texto. Infelizmente, a frequência de corte de -3 dB não tem solução simbólica no chapéu de aproximação N (omega). Para diferentes N, a aproximação só difere da aproximação N 1 por um mapeamento omega rightarrow omega N, por isso basta resolver a freqüência de corte -3 dB hatomegac (N) nummente para N 1: dando a freqüência de corte aproximada para N arbitrário : Isso parece ser outra, aproximação mais simples do que Massimos. Para o seu N gt 300, não deve haver problema em usá-lo. Massimos e as constantes das respostas estão relacionadas por: Eu olhei um pouco mais e descobriu que Massimo se aproxima de FN (omega) com o chapéu M (omega), escolhendo M de tal forma que os (limites do) segundo derivado da resposta de freqüência e aproximação combinam Em omega 0: isto melhora a aproximação em omega pequeno que inclui o ponto de corte de -3dB, especialmente na pequena aproximação de N: Massimos sempre superestima a freqüência de corte (veja a comparação de erro), deixando espaço para melhorá-lo alterando a constante 1. O O erro é o maior para N 2. Se o seu erro for limitado igual ao erro (atualmente o segundo maior) no N 3, obtemos uma aproximação ainda melhor, mas tão barata: Este e outros ajustes da constante, como a constante de Matts 0,863031932778066. Trabalhe surpreendentemente bem para o N grande (veja a comparação do erro). Para o grande N, o erro cai em um fator de 1000 para cada aumento de N por um fator 10. A explicação para essas coisas é que a verdadeira freqüência de corte como função de N tem uma série Laurent: e a aproximação e sua série Laurent São: de modo que: a1 a 2.78311475650302030063992 a3 aprox-frac Se a aproximação do N - term for feita com exatidão, o erro de aproximação deve diminuir em um fator de 105 para um aumento de N grande por um fator de 10. Os coeficientes ak Da série de Laurent f (x) soma fração de uma função f (x) como xrightarrowinfty pode ser encontrada iterativamente por: Quando não temos f (x) em forma simbólica, mas podemos resolvê-lo numericamente a qualquer precisão para x grande , Podemos fazer o equivalente ao procedimento acima numericamente. O seguinte script Python que usa SymPy e mpmath calculará um número dado (aqui 10) dos primeiros coeficientes ak na precisão desejada para a série Laurent da freqüência de corte real: no meu computador, o programa funciona por aproximadamente 7 minutos. Ele imprime o seguinte, mostrando que a série Laurent consiste apenas em poderes negativos estranhos: esses números, mostrados a 24 casas decimais, não são uma aproximação no sentido de que a série Laurent é única, não há nenhuma outra série Laurent que seja igual a omegac ( N). Usando apenas a1 e a3, pode-se construir uma simples aproximação truncada de duas séries de Laurent: e por c-frac a aproximação: ambos apresentam 1 decadência de erro N5 na grande N, ver comparação de erros. Colunas h) e i) respectivamente. Uma série de Laurent mais longa e truncada com mais termos do resultado dos scripts decai ainda mais rápido, 1 N para a aproximação de 5 períodos na coluna j) na comparação de erros. Acima de mim, Olli. Mas por algum motivo, acho que a resposta é muito mais simples. Normalmente eu gosto de projetar filtros simétricos de FIR acausal, porque eles são zero fase, mas geralmente eu me limito a um número ímpar de torneiras não-zero. Para fazer isso de forma mais geral, eu posso ficar com a média móvel causal da FIR. Digamos que o número de torneiras é N. aplicando mathcal - transform (e a fórmula de soma geométrica): substituindo z leftarrow e para obter o DTFT: normalmente chamamos a coisa que multiplica X (z) a função de transferência e a coisa que multiplica X (E), a resposta de frequência o fator e significa uma fase linear, atraso constante de amostras de fração. Não altera o ganho. O fator de fracasso é o fator de ganho. A freqüência de -3 dB, omegac, (normalmente significamos a freqüência de -3,0103 dB porque isso corresponde à metade da freqüência de potência) é tal que 2 sin2 (omegac N 2) N2 sin2 (omegac 2). Então, dado o número de torneiras N, você tem que resolver omegac. Isso pode não ser tão fácil de fazer para uma forma fechada, mas você pode cavar sua calculadora e conectar e chug até obter uma resposta que tenha precisão suficiente. Ou você pode fazer MATLAB para fazê-lo. Pode-se ter uma aproximação decente para omegac para N grande utilizando uma identidade trigonométrica (uma que costumo usar quando estou mexendo com a transformação bilinear) e os três primeiros termos da série Maclaurin para cos (). Se você conecta essa aproximação para sin2 () na equação anterior e resolva. (Ignorando muitos passos porque eu sou muito preguiçoso para sair daqui.) Olli, quão bem isso se compara aos seus resultados fazendo isso melhor com outro termo para a aproximação de sin2 (), é realizável, requer apenas uma solução quadrática para Omega02. A aproximação para usar (manter os primeiros quatro termos da expansão cos ()) é: tente essa aproximação e resolva omegac2. A resposta mais consistente que recebo é com a opção que parece e com a opção - parece que está muito mais perto da aproximação de primeira ordem que fiz acima. Então eu acho que eu tomaria a opção -. Então, embora eu não possa dizer de forma analítica por que a opção deve ser rejeitada, acho que minha resposta mais precisa seria qual o limite, para o N grande, mostrado acima. Alguém mais tem uma maneira melhor de olhar para uma boa solução de forma fechada aproximada para esta última semana antes de se retirar. A aproximação sin2 (theta) aproximadamente theta2 esquerda (1 - frac theta2 frac theta4 right) realmente deve ser boa para todos os 0 the theta le frac, para que isso aconteça e para tornar o comportamento continuar sendo real bom em theta ll 1, nós Deve fudge o último coeficiente, frac, para ser fraco, de modo que a aproximação é boa para sin2left (frac direita). Isso não aumenta a complexidade, mas pode fazer uma melhor resposta. Ndash robert bristow-johnson 13 de janeiro às 6:27 frac)) é na verdade um trem de impulso de banda limitada, de modo que aproximar-se com uma função de texto, como na minha resposta é exato (dentro da precisão de 2.78311475650302030063992) no limite de grande N onde Sua fração omega0 dá cerca de 0,88 vezes o limite verdadeiro e seu omega0 sqrt direito) dá cerca de 1,035 vezes o verdadeiro corte. Eu acho que se você quiser fazer uma melhor aproximação, você deve incluir essa longa constante. Ndash Olli Niemitalo 13 de janeiro às 8:46 Robert, você precisa usar o sinal 39-39 em sua fórmula de equação quadrática, porque isso dá a solução onde a série Taylor ainda se aproxima da função original. A outra solução é válida apenas para o polinômio de Taylor, mas não para a função original, pois, para esse valor maior, o polinômio de Taylor não se aproxima mais da função original. Então, para uma expansão de Taylor em torno de x00, você normalmente precisa escolher a menor solução (em magnitude), porque essa é a que a aproximação funciona melhor. Ndash Matt L. 13 de janeiro às 14:23 Vamos comparar os erros numéricos reais para diferentes aproximações da freqüência de corte. O erro dado na tabela é calculado subtraindo o omegac de freqüência de corte verdadeiro (com resolução numérica) de -3 dB da aproximação. Notas: A aproximação e) não permite N2. Alguns dos erros estão listados como 0, mas significa apenas que sua magnitude é inferior a aproximadamente 1E-17. Essa e outras imprecisões possíveis são do uso de aritmética de ponto flutuante de dupla precisão no cálculo da aproximação e do erro. Sinta-se livre para editar, adicione outra aproximação. OK, isso é divertido. Eu vou adicionar meus próprios pensamentos e aproximações, o primeiro se revela idêntico ao dado por Massimo nesta resposta. E aquele derivado por Olli neste tópico. Continuo incluído aqui porque a sua derivação é diferente. Então, eu mostro uma melhor aproximação, que tem um erro relativo máximo de 0,002 para N2 (para o qual, é claro, temos uma solução analítica para a freqüência de corte exata: omegacpi 2) e para o qual o erro relativo é menor que 1,2 Cdot 10 para Nge 10. É bem conhecido, e foi mostrado por Olli e Robert em suas respostas, que a função de amplitude de valor real de um filtro de média móvel N de comprimento é dada pela freqüência de corte de 3 dB Omegac satisfaz o tempo todo Como eu sei, não existe uma solução analítica razoavelmente simples para a Eq. (2). A chave para as aproximações apresentadas aqui é - não surpreendentemente - uma aproximação de Taylor. A diferença entre a série de Taylor usada na resposta de Roberts é que não me aproximo separadamente das funções seno (ou dos seus valores quadrados como nas respostas de Roberts), mas eu me aproximo diretamente da função de amplitude completa dada em (1). O pecado aproximado (Nomega 2) (ou seu valor quadrado) resultará em erros maiores do que quando a função completa é aproximada, porque o argumento Nomega 2 nunca se aproxima de zero, mesmo para valores grandes de N. Aproximando apenas o denominador sin (omega 2) (Ou o seu valor quadrado) está OK, porque o seu argumento omegaomegac aproxima-se de zero para o N. grande. De qualquer forma, não usarei nenhuma das duas aproximações, mas vou usar a série Taylor de HN (omega). Para uma notação mais simples, eu uso xomega 2 e FN (x) HN (omega). A série Taylor de FN (x) em torno de x00 é dada por Para grandes valores de N, essa aproximação é legítima porque a frequência de corte omegac tende a valores pequenos. Para a primeira aproximação, eu uso apenas os dois primeiros termos em (3): a resolução (4) fornece uma primeira solução aproximada: o problema com esta solução é que é tendencioso, o que significa que seu erro não converge para zero para N. grande No entanto, verifica-se que por uma simples escala de (5), este tendencioso pode ser removido. Para uma polarização zero, exigimos onde usei a notação omega (N) para enfatizar sua dependência de N. Resolvendo (6) com a expressão geral nos leva à equação que deve ser resolvida numericamente para a solução agora famosa. A aproximação ( 7) com um dado por (9) é idêntico à fórmula de Massimos (você precisa dividir por 2pi para obter sua constante de magia), e também é o mesmo que o derivado por Olli de uma maneira diferente neste tópico. Nós vemos que uma aproximação de Taylor nos deu a forma correta da equação, mas a constante teve que ser determinada por um processo limite para obter uma fórmula com viés zero. Para a maioria dos propósitos práticos, esta fórmula é suficientemente precisa com um erro relativo máximo de 6.9cdot 10 para Nge 10. Usar todos os termos na aproximação (3) nos dará uma aproximação ainda melhor. O processo é exatamente o mesmo que antes: definimos a aproximação de Taylor de FN (x) igual a 1 sqrt e resolvemos para xc (existem apenas poderes de x, então precisamos apenas resolver uma equação quadrática). Isso nos dá a seguinte fórmula: Note que das quatro soluções da equação quartic, precisamos escolher o menor dos dois positivos, porque esse é o valor em que a série Taylor se aproxima de FN (x). A outra solução positiva é um artefato em uma faixa em que a série de Taylor diverge de FN (x). A aproximação (10) tem o mesmo pequeno problema que a primeira versão da aproximação anterior dada por (5) na medida em que possui um pequeno viés. Este viés pode ser removido exatamente da mesma forma que antes considerando o limite (6), desta vez com omega (N). Minha aproximação final com base em (10), mas com preconceito zero é dada por onde b também pode ser obtido resolvendo uma equação semelhante a (8). Pode realmente ser escrito em termos de um dado por (9): bfrac sqrt -1) 0,997314251642175tag. Eu computei os valores exatos de omegac numericamente para N na faixa de 2,100, então eu poderia calcular o erro relativo que permite a comparação de diferentes Aproximações omega. Eu apenas discuto as aproximações com preconceito zero: omega dado por (7) com um dado por (9) e omega dado por (11) (e (10)), com b dado por (12). A figura abaixo mostra os erros relativos conforme definido por (13) como uma função de N. A curva vermelha é o erro relativo de aproximação (7), e a curva verde é o erro de aproximação (11). Ambas as aproximações têm polarização zero (convergem para os valores exatos para N grande), mas a curva verde converge significativamente mais rapidamente. As fórmulas de polarização zero mostradas acima são aproximações decentes às freqüências de corte reais, mas a melhor (fórmulas (10,11,12)) é muito estranha. Olli teve a ótima idéia para ajustar o denominador constante na fórmula simples (7). Enquanto usarmos o valor ideal de um determinado por (9), podemos mudar a constante do denominador sem perder a propriedade do desvio zero. Então, obtemos uma nova fórmula com uma constante c para ser otimizada. Se entendi corretamente, Olli baseou sua otimização de c no valor de erro para N2. No entanto, acho que o valor N2 não é muito relevante porque, para N2, omegac pode ser calculado de forma analítica: omegac (2) pi 2. Portanto, não precisamos necessariamente otimizar a fórmula (14) para o caso N2 se ocorrer às custas Da aproximação em valores maiores de N. Eu otimizei a constante c em (14) da seguinte maneira. Se omegac (N) são as frequências de corte exatas para um determinado conjunto de comprimentos de filtro N, temos um sistema de equações sobredeterminado: onde podemos escolher qualquer conjunto razoável de valores para N. Reorganizando (15) dá outro conjunto de equações , Desta vez linear no desconhecido c: A solução ideal de mínimos quadrados de (16) é onde L é o número de valores diferentes para N usado na soma. Se você usar todos os valores inteiros de N na faixa de 2,100, você obtém o que está próximo do valor de Ollis, mas que dá uma aproximação ainda melhor para todos os Nge 3. Eu adicionei os valores de erro a esta tabela. Coluna f). Em sua resposta, Robert estava perguntando por que ele deve descartar a segunda (maior) solução positiva para omegac ao usar uma série Taylor de quarta ordem para sin2 (x). A figura abaixo mostra o motivo. A função de amplitude quadrada original é mostrada em azul (para N10). A linha 3dB está em vermelho. A função verde é a aproximação de Taylor, que cruza a linha vermelha duas vezes. Estas são as duas soluções positivas para omegac. Uma vez que a função é igual, também obtemos as mesmas duas soluções com sinais negativos, o que o torna quatro, como deveria ser o caso de um polinômio de quarta ordem. No entanto, é óbvio que a maior das duas soluções positivas é um artefato devido à divergência da aproximação de Taylor para argumentos maiores. Portanto, é apenas a solução menor que faz sentido, a outra não. Eu forneço outra resposta porque essa abordagem é completamente diferente no sentido de não tentar aproximar a função de amplitude dos filtros para calcular uma aproximação da freqüência de corte, mas uso uma abordagem de ajuste de dados pura, dado as freqüências de corte exatas , Que foram computados numericamente (e também são fornecidos para um conjunto de comprimentos de filtro na coluna mais à esquerda desta tabela). Com o ajuste de dados, muitas vezes o problema mais difícil é encontrar uma parametrização apropriada da função de aproximação. Uma vez que das outras respostas neste tópico sabemos que com as constantes adequadamente escolhidas a e c é uma aproximação surpreendentemente boa para uma ampla gama de valores de N, e como o Wolfram Alpha nos diz que a expansão da série Laurent de (1) na Ninfty tem Apenas termos com poderes ímpares de 1 N, parece razoável parametrizar a freqüência de corte por uma série Laurent com poderes ímpares de 1 N: podemos calcular o valor exato de a1 em (2) a partir da exigência de que a estimativa c (N) tem polarização zero, ou seja, converge para a verdadeira freqüência de corte para o N. grande. Isto é explicado na minha outra resposta. Seu valor é As outras constantes em (2) podem ser calculadas pelo ajuste de mínimos quadrados de (2) aos dados, que são as freqüências de corte exatas. O ajuste de mínimos quadrados pode ser calculado pelo seguinte script Matlab Octave (assumindo que a variável wc é um vetor com frequências de corte exatas pré-calculadas para o conjunto desejado de comprimentos de filtro): os coeficientes resultantes começam a3amp1.201014809636180 a5amp0. 768735238011194 a7amp0.514237034990353 a9amp0.548681885931852end com a1 dado por (3). Esta aproximação vem extremamente próxima dos valores exatos de omegac. O erro de aproximação pode ser encontrado nesta tabela (coluna g). Preciso projetar um filtro médio móvel que tenha uma freqüência de corte de 7,8 Hz. Eu usei filtros de média móvel antes, mas, na medida em que estou ciente, o único parâmetro que pode ser alimentado é o número de pontos a serem calculados. Como isso se relaciona com uma freqüência de corte O inverso de 7,8 Hz é de 130 ms, e estou trabalhando com dados que são amostrados a 1000 Hz. Isso implica que eu deveria estar usando uma média móvel de tamanho de janela de filtro de 130 amostras, ou há algo mais que estou faltando aqui? 18 de julho 13 às 9:52 O filtro de média móvel é o filtro usado no domínio do tempo para remover O som adicionado e também para fins de alisamento, mas se você usar o mesmo filtro de média móvel no domínio da freqüência para a separação de freqüência, o desempenho será pior. Então, nesse caso, use filtros de domínio de freqüência ndash user19373 3 de fevereiro às 5:53 O filtro de média móvel (às vezes conhecido coloquialmente como um filtro de caixa) tem uma resposta de impulso retangular: Ou, afirmado de forma diferente: lembrando que uma resposta de freqüência de sistemas discretos é Igual à transformação de Fourier de tempo discreto de sua resposta de impulso, podemos calculá-lo da seguinte maneira: o que mais interessou para o seu caso é a resposta de magnitude do filtro, H (omega). Usando algumas manipulações simples, podemos obter isso de forma mais fácil de entender: isso pode não parecer mais fácil de entender. No entanto, devido à identidade do Eulers. Lembre-se que: Portanto, podemos escrever o acima como: Como eu disse anteriormente, o que você realmente está preocupado é a magnitude da resposta de freqüência. Então, podemos tomar a magnitude do acima para simplificá-lo ainda mais: Nota: Podemos soltar os termos exponenciais porque eles não afetam a magnitude do resultado e 1 para todos os valores de ômega. Como xy xy para dois números complexos finitos x e y, podemos concluir que a presença dos termos exponenciais não afeta a resposta de magnitude global (em vez disso, eles afetam a resposta de fase dos sistemas). A função resultante dentro dos suportes de magnitude é uma forma de um kernel Dirichlet. Às vezes, é chamado de função periódica sinc, porque se parece com a função sinc em um pouco de aparência, mas é periodicamente. De qualquer forma, uma vez que a definição de frequência de corte é um pouco sub-especificada (ponto -3 dB -6 dB ponto primeiro sidelobe nulo), você pode usar a equação acima para resolver o que você precisa. Especificamente, você pode fazer o seguinte: Ajuste H (omega) para o valor correspondente à resposta do filtro que você deseja na freqüência de corte. Defina omega 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 omega 2pi frac, onde fs é a taxa de amostragem. Encontre o valor de N que lhe dá o melhor acordo entre os lados esquerdo e direito da equação. Esse deve ser o comprimento da sua média móvel. Se N é o comprimento da média móvel, então uma frequência de corte aproximada F (válida para N gt 2) na frequência normalizada Ff fs é: O inverso disso é Esta fórmula é assintoticamente correta para N grande e tem cerca de 2 Erro para N2 e menos de 0,5 para N4. P. S. Depois de dois anos, aqui, finalmente, qual era a abordagem seguida. O resultado baseou-se na aproximação do espectro de amplitude MA em torno de f0 como uma parábola (série de 2ª ordem) de acordo com MA (Omega) aprox. 1 (frac - frac) Omega2 que pode ser feita mais exatamente perto do cruzamento zero de MA (Omega) Frac, multiplicando Omega por um coeficiente que obteve MA (Omega) aproximadamente 10.907523 (frac - frac) Omega2 A solução de MA (Omega) - frac 0 dá os resultados acima, onde 2pi F Omega. Todos os itens acima dizem respeito à frequência de corte -3dB, o assunto desta publicação. Às vezes, é interessante obter um perfil de atenuação na banda de parada que é comparável ao de um filtro de passagem baixa IIR de 1ª ordem (LPF de um único pólo) com uma freqüência de corte de -3dB dada (como um LPF também é chamado de integrador vazado, Tendo um pólo não exatamente na DC, mas perto disso). De fato, tanto o MA quanto o LPR de 1ª ordem IIR têm declive de década de 20 dB na faixa de parada (um precisa de um N maior que o usado na figura, N32, para ver isso), mas enquanto o MA tem nulos espectrales no Fk N E um 1 f evelope, o filtro IIR possui apenas um perfil de 1 f. Se alguém quiser obter um filtro MA com capacidades semelhantes de filtragem de ruído como este filtro IIR e corresponder às freqüências de corte 3dB para serem iguais, ao comparar os dois espectros, ele perceberia que a ondulação da faixa de parada do filtro MA termina 3dB abaixo do 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) como o filtro IIR, as fórmulas podem ser modificadas da seguinte maneira: encontrei o script Mathematica onde eu calculava o corte para vários filtros, incluindo o MA. O resultado foi baseado em aproximar o espectro de MA em torno de f0 como uma parábola de acordo com MA (Omega) Sin (OmegaN 2) Sin (Omega 2) Omega 2piF MA (F) aproximadamente N1 6F2 (N-N3) pi2. E derivando o cruzamento com 1 sqrt a partir daí. Ndash Massimo 17 de janeiro às 2:08

No comments:

Post a Comment