Alguém está familiarizado com o algoritmo de Gustafson para minimizar os transientes na filtragem reversa para frente [1]? Estou tentando implementá-lo e meu primeiro palpite foi verificar o filtfilt.m do Matlab, já que eles estão fazendo referência ao artigo. Na função Matlab também é resolvido um sistema de equações lineares para encontrar condições iniciais zi que minimizem os transientes de inicialização, mas a relação entre referência e código não é óbvia para mim. As únicas linhas de código relacionadas à minimização são (nfilt é o comprimento dos vetores do coeficiente):
zi = ( eye(nfilt-1) - [-a(2:nfilt), [eye(nfilt-2); zeros(1,nfilt-2)]] ) \...
( b(2:nfilt) - b(1)*a(2:nfilt) );
Alguém pode me indicar a direção certa de como essas linhas se relacionam com o algoritmo descrito no artigo de Gustafson?
[1] Gustafsson, F. "Determinando os estados iniciais na filtragem para frente e para trás". Transações IEEE® no processamento de sinais. Vol. 44, abril de 1996, pp. 988–992.
filtfilt()
não consigo ver o porquê. eu não li o artigo de Gustafson (eu não sou o IEEE e não posso obtê-lo de graça, qualquer pessoa que tenha uma cópia pode enviar e-mail com um arquivo .pdf). ao usar o conceito defiltfilt
, pode-se fazer isso com um arquivo inteiro de amostras (para mim seria um arquivo de áudio ou som, como um .wav) primeiro filtrar o som para a frente, filtrando o som com o preenchimento zero no final, desde que você espera que a resposta de impulso do filtro direto seja. isso aumenta o arquivo, mas a saída fica praticamente zero. depois, execute o arquivo resultante pelo filtro para trás.filtfilt
em tempo real dividindo a entrada em blocos de amostras, preenchendo com zero cada bloco, filtrando os blocos para trás, mas mantendo as "caudas" voltando-as para a direção para frente e sobreposição de adição. Powell-Chau não fez isso, mas acho que essa é uma boa aplicação dos filtros IIR truncados, para que você saiba quando a saída do bloco decadente termina.Respostas:
Para quem estiver interessado, encontrei por acaso um documento descrevendo o método implementado no filtfilt.m do matlab. Um link para o artigo está anexado. Pelo menos para meu entendimento, o filtfilt.m do matlab não implementa o algoritmo Gustafson.
Sadovsky, P .; Bartusek, K: Otimização da resposta transitória de um filtro digital, Radioengineering vol. 9, n. 2, 2000
fonte
scipy
documentação delfilter_zi
, que é o padrãoscipy.signal.filtfilt
para determinar as condições iniciais, como você pode ver na fonte . Nesse caso, o preenchimento 'ímpar' é usado por padrão, mas pode usar o método de Gustafsson como uma opção (veja a definição_filtfilt_gust
na fonte).A
zi = (...)\(...)
linha na pergunta do OP determina o estado inicial do filtro. Eu acredito que exatamente a mesma abordagem é usada pelo Pythonscipy
. De acordo com os documentos scipy :(minha ênfase)
Observe que
filtfilt
calcula o estado inicialzi
, dimensiona-o pelo primeiro valor de amostra e passa para ofilter
que realmente o aplica ( docs ).Um exemplo básico de como esse estado inicial
zi
pode ser aplicado em um filtro, usando a representação do espaço de estado, é fornecido aqui .fonte