Um filtro FIR simples

Há três métodos básicos para construir-se um filtro digital: FIR, IIR e DFT. Neste artigo demonstrarei o filtro FIR (resposta impulsional finita).

Primeiro, precisamos visitar um conceito da matemática, a convolução, que pode ser entendida como uma espécie de "multiplicação gráfica", "borrão", ou média móvel.

Considere a seguinte equação:

média móvel do
índice da Bolsa no dia n = 0.5 * índice(dia n) +
                           0.3 * índice(dia n-1) +
                           0.2 * índice(dia n-2)

Esta é uma média móvel do índice da Bolsa de Valores, que leva em conta os três últimos dias, com diferentes pesos para cada dia: 0.5, 0.3 e 0.2.

Figura 1: Gráfico com o índice da Bolsa de Valores (azul) e a média móvel (vermelha)

A média móvel é claramente amortecida em relação ao índice no qual ela é baseada. A linha vermelha ainda segue o mercado, porém "amacia" ou "borra" picos e vales.

Dado o movimento do índice de uma Bolsa de Valores, eu posso divisar uma "linha de tendência", que é mais fácil de enxergar se usarmos uma média móvel em vez do índice cru que é muito volátil. Muitos métodos de análise técnica lançam mão de médias móveis.

A linha vermelha pode ser matematicamente entendida como a convolução de duas seqüências: a primeira seqüência é a linha azul, cujo comprimento pode ser muito grande. A segunda seqüência são os pesos da média móvel, cujo comprimento é exatamente três: [0.5, 0.3, 0.2].

Calcular uma média móvel não é barato, em termos computacionais. Se há três pesos, precisamos fazer três multiplicações por amostra. Uma convolução de 60 pesos em cima de uma seqüência de 5000 amostras exige 300 mil multiplicações.

Médias móveis são filtros passa-baixas naturais. Filtro passa-baixas é justamente o que vamos construir neste artigo, para um sinal de áudio.

A convolução estabelece uma relação poderosa entre domínio do tempo e domínio da freqüência. O que estamos fazendo aqui é relacionado ao trabalho que desenvolvemos no artigo sobre filtros DFT.

Converter um sinal no domínio do tempo para o domínio da freqüência tem o poder de "converter" uma convolução em uma simples multiplicação, e vice-versa. Isto significa que, se podemos definir um filtro no domínio da freqüência em termos de multiplicações, podemos definir o mesmo filtro em termos de uma convolução no domínio do tempo.

Olhando no código do filtro DFT, temos uma seção onde definimos uma "máscara" de freqüência. Esta máscara é multiplicada pelo espectro, com o objetivo de filtrar uma banda:

mask = []
for f in range(0, FFT_LENGTH / 2 + 1):
    rampdown = 1.0
    if f > LOWPASS:
        rampdown = 0.0
    elif f < HIGHPASS:
        rampdown = 0.0
    mask.append(rampdown)

Então, em princípio, se convertermos esta máscara para o domínio do tempo com IFFT, armazenar o resultado num WAV, e fazer a convolução desse WAV com nosso áudio, obteremos o mesmo tipo de filtro (tirante errors de arredondamento e quantização).

O arquivo WAV gerado pela máscara não é audível. Ele é chamado de "resposta impulsional", e é simplesmente uma lista de pesos, tal qual a lista (0.5, 0.3, 0.2) da média móvel da Bolsa do início do artigo. A diferença é que são muito mais pesos, já que filtragem de áudio é trabalho mais delicado.

Escolhemos gravar esses pesos num arquivo WAV apenas para demonstrar que, uma vez de posse deles, nunca mais precisaremos de FFT para fazer a filtragem. A resposta impulsional de um filtro passa-baixas de 3000Hz é algo semelhante a isto:

Figura 2: Resposta impulsional de um filtro FIR

Tudo que temos de fazer agora é convolver um sinal de áudio com esta resposta impulsional, e a filtragem estará feita. Aqui está o código:

#!/usr/bin/env python

import wave, struct, numpy
SAMPLE_RATE = 44100

original = wave.open("NubiaCantaDalva.wav", "r")
filter = wave.open("fir_filter.wav", "r")
filtered = wave.open("NubiaFilteredFIR.wav", "w")
filtered.setnchannels(1)
filtered.setsampwidth(2)
filtered.setframerate(SAMPLE_RATE)

n = original.getnframes()
original = struct.unpack('%dh' % n, original.readframes(n))
original = [s / 2.0**15 for s in original]

n = filter.getnframes()
filter = struct.unpack('%di' % n, filter.readframes(n))
filter = [s / 2.0**31 for s in filter]

'''
result = []
# original content is prepended with zeros to simplify convol. logic
original = [0 for i in range(0, len(filter)-1)] + original

for pos in range(len(filter)-1, len(original)):
    # convolution of original * filter
    # i.e. a weighted average of len(filter) past samples
    sample = 0
    for cpos in range(0, len(filter)):
        sample += original[pos - cpos] * filter[cpos]
    result.append(sample)
'''

result = numpy.convolve(original, filter)

result = [ sample * 2.0**15 for sample in result ]
filtered.writeframes(struct.pack('%dh' % len(result), *result))

Agora, você está curioso para ouvir os resultados. Aqui está o áudio original e a versão filtrada:

(link to audio)

Talvez você tenha notado que há uma seção do código que está comentada (desativada). Ela implementa a convolução de forma "manual", usando dois loops e aritmética básica. Isto funciona perfeitamente bem, e deixei comentado no código para deixar claro que convolução é fundamentalmente simples. O único senão é a lentidão. Usar o método numpy.convolve() é muitíssimo mais rápido. De fato, este filtro FIR é mais rápido que a versão FFT.

É claro, o filtro FIR é rápido, mas precisamos ter os pesos, ou seja, a resposta impulsional à mão. O código a seguir gera a resposta impulsional e só precisa ser executado uma vez (ou quando você quiser mudar a freqüência de corte do filtro):

#!/usr/bin/env python

LOWPASS = 3000 # Hz
SAMPLE_RATE = 44100 # Hz

import wave, struct, math
from numpy import fft
FFT_LENGTH = 512

NYQUIST_RATE = SAMPLE_RATE / 2.0
LOWPASS /= (NYQUIST_RATE / (FFT_LENGTH / 2.0))

# Builds filter mask. Note that this sharp-cut filter is BAD
mask = []
negatives = []
l = FFT_LENGTH / 2
for f in range(0, l+1):
    rampdown = 1.0
    if f > LOWPASS:
        rampdown = 0.0
    mask.append(rampdown)
    if f > 0 and f < l:
        negatives.append(rampdown)

negatives.reverse()
mask = mask + negatives

fir = wave.open("fir_filter.wav", "w")
fir.setnchannels(1)
fir.setsampwidth(4)
fir.setframerate(SAMPLE_RATE)

# Convert filter from frequency domain to time domain
impulse_response = fft.ifft(mask).real.tolist()

# swap left and right sides
left = impulse_response[0:FFT_LENGTH / 2]
right = impulse_response[FFT_LENGTH / 2:]
impulse_response = right + left

# write in a normal WAV file
impulse_response = [ sample * 2**31 for sample in impulse_response ]
fir.writeframes(struct.pack('%di' % len(impulse_response),
                                 *impulse_response))

Apesar da simplicidade, um filtro FIR nem sempre é mais rápido que um filtro FFT. Nos meus testes, o "empate" ocorre em torno de 64 pesos. Mais pesos que isso, e o filtro FFT passa a ser mais barato. Esse ponto pode variar de acordo com o tipo de processador ou microcontrolador que executa a filtragem.

O filtro FIR é menos poderoso que FFT, no sentido de que ele não pode fazer todos os truques que o FFT pode. (Exemplo: aqueles efeitos sonoros "auto tuner" que fazem a voz ficar mais aguda ou mais grave sem alterar a velocidade de dicção.) FIR está para FFT assim como um abridor de latas está para um canivete suíço.