Filtri IIR digitali e cascata di biquad

Cosa fa lo strumento

Progetta un filtro IIR digitale a partire da una specifica (tipo, famiglia, ordine, frequenza di campionamento e di taglio) e restituisce i coefficienti biquad, la risposta in frequenza (ampiezza e fase), la mappa poli/zeri e l’export per CMSIS-DSP. Il calcolo gira nel browser su un’implementazione propria, senza dipendenze runtime: prototipo analogico → trasformazione LP→LP / LP→HP → bilineare → SOS, con l’algoritmo di Orfanidis per le funzioni ellittiche di Jacobi dell’ellittico.

Un filtro IIR di ordine N si realizza come cascata di sezioni del secondo ordine (SOS, second-order sections): è la forma numericamente robusta, perché una singola sezione di ordine alto accumula errori di quantizzazione sui coefficienti. Ogni sezione è un biquad:

$$H_k(z)=\frac{b_0+b_1{z}^{-1}+b_2{z}^{-2}}{1+a_1{z}^{-1}+a_2{z}^{-2}},H(z)=\prod _{k=1}^{N/2}{H}_k(z)$$

Dal continuo al discreto: trasformata bilineare

Il prototipo analogico (Butterworth o Bessel) viene mappato nel piano z con la trasformata bilineare. Per far coincidere esattamente la frequenza di taglio si applica il pre-warping:

$$s\leftarrow \frac{2}{T}\frac{1-z^{-1}}{1+z^{-1}},{\omega }_a=\frac{2}{T}\tan \left(\frac{{\omega }_{d}T}{2}\right)$$

dove T = 1/Fs. La bilineare conserva la stabilità e comprime l’asse delle frequenze: tutto l’asse analogico 0…∞ finisce nell’intervallo 0…Fs/2 (Nyquist).

Stabilità: si legge dai poli

I poli sono le radici del denominatore di ciascun biquad; gli zeri quelle del numeratore. La bilineare manda il semipiano sinistro analogico (stabile) nell’interno del cerchio unitario. Quindi:

$$\text{stabile}⟺|p_i|<1\forall i$$

Lo strumento calcola i poli risolvendo z² + a₁z + a₂ = 0 per ogni sezione e segnala se anche uno solo esce dal cerchio unitario.

Le cinque approssimazioni: quale scegliere

Le famiglie sono i diversi modi di approssimare la risposta ideale “a muro”. Ognuna è un compromesso tra ripidità della transizione, ripple (ondulazione) e linearità di fase (fedeltà nel tempo). Lo strumento le copre tutte.

Butterworth — modulo massimamente piatto

Nessun ripple in banda passante, attenuazione asintotica di 20·N dB/decade, roll-off il più dolce del gruppo.

$$|H_{\text{BW}}(j\omega ){|}^2=\frac{1}{1+{\left(\omega /{\omega }_c\right)}^{2N}}$$

Quando: banda passante pulita, ordine non critico. Default ragionevole quando non hai requisiti stringenti. Casi d’uso: condizionamento di sensori, anti-aliasing non aggressivo, livellamento generico.

Chebyshev I — ripple in banda passante

Ammette un ripple Ap in banda passante (parametro ε) in cambio di una transizione più ripida di Butterworth a parità di ordine. Stopband monotona.

$$|H_{\text{C1}}(j\omega ){|}^2=\frac{1}{1+{\epsilon }^2{T}_N^2(\omega /{\omega }_c)},\epsilon =\sqrt{10^{A_p/10}-1}$$

dove Tₙ è il polinomio di Chebyshev. Quando: serve più selettività e si tollera un po’ di ondulazione in banda. Casi d’uso: filtri di canale, separazione di bande dove l’ampiezza in banda non deve essere perfettamente piatta.

Chebyshev II (inverso) — ripple in stopband

Banda passante piatta (come Butterworth), ripple in stopband (parametro As). La transizione è ripida e l’attenuazione equiripple garantisce un minimo di reiezione costante.

$$|H_{\text{C2}}(j\omega ){|}^2=\frac{1}{1+\frac{1}{{\epsilon }^2{T}_N^2({\omega }_s/\omega )}}$$

Quando: vuoi banda passante pulita ma accetti ondulazioni nell’attenuata, con un livello minimo di reiezione garantito. Casi d’uso: anti-aliasing dove conta la planarità della banda utile e basta un floor di stopband.

Ellittico / Cauer — il più ripido

Ripple in entrambe le bande. A parità di ordine ha la transizione più ripida in assoluto: nessuna altra approssimazione classica lo batte in selettività per dato N. Il prezzo è la fase più non lineare (ritardo di gruppo molto variabile vicino al taglio) e maggiore overshoot/ringing nel tempo.

$$|H_{\text{ell}}(j\omega ){|}^2=\frac{1}{1+{\epsilon }^2{R}_N^2(\xi ,\omega /{\omega }_c)}$$

dove Rₙ è la funzione razionale ellittica (costruita su funzioni ellittiche jacobiane) e ξ il fattore di selettività che lega ripple in banda passante e attenuata. Casi limite (inquadrano la sua generalità):

• ripple stopband → 0 degenera in Chebyshev I
• ripple passband → 0 degenera in Chebyshev II
• entrambi → 0 degenera in Butterworth

Quindi l’ellittico è la forma generale; le altre tre sono suoi casi particolari. Quando: serve la massima ripidità con l’ordine (e quindi il costo computazionale) più basso, e puoi tollerare ripple e fase non lineare. Casi d’uso: anti-aliasing aggressivi, reiezione stretta di toni/bande, applicazioni dove conta solo il modulo e non la forma d’onda.

Il ritardo di gruppo, misurato

Dire che il Bessel ha «ritardo di gruppo quasi costante» è vero ma non dice quanto. Calcolato per differenze finite sulla fase, in campioni, per i quattro filtri a ordine 4 con Fc = 100 Hz e fs = 1000 Hz:

L’ultima colonna è l’escursione percentuale del ritardo fra 10 e 100 Hz, e mette in fila le famiglie in un ordine preciso: Bessel 12 %, Butterworth 61 %, Chebyshev I 214 %, ellittico 404 %. Un fronte che attraversa l’ellittico esce con le sue componenti spettrali sfasate di un fattore quattro nel tempo di transito: è esattamente l’origine del ringing. Se il segnale è un impulso, un fronte, una forma d’onda da misurare nel tempo, quella colonna conta più della ripidezza.

Si noti anche che il ritardo cala a 150 Hz per tutte: fuori banda il filtro non ritarda, attenua, e il ritardo di gruppo perde significato pratico.

Quanto resta della linearità di fase in digitale

La linearità del Bessel è una proprietà del prototipo analogico, e la bilineare la deforma. Quanto? Dipende da dove sta Fc rispetto a fs:

Sotto fc/fs ≈ 0,1 il Bessel digitale conserva il comportamento analogico: l’escursione del 18–21 % è quella intrinseca del prototipo a ordine 4, non un danno della trasformazione. Attorno a 0,2 la bilineare inizia a mordere, e a fc/fs = 0,4 l’escursione arriva al 750 %: il filtro ha perso qualunque somiglianza con un Bessel. La regola pratica è quindi netta — se si sceglie il Bessel per la fase, si sovracampiona: la fase lineare digitale si compra con un fs almeno dieci volte la banda utile, oppure si passa a un FIR, dove è esatta per costruzione.

Bessel / Thomson — fase massimamente lineare

Ritardo di gruppo quasi costante (fase lineare): preserva la forma d’onda nel tempo, niente overshoot sul gradino. Roll-off molto dolce. Nota: la linearità di fase del Bessel è una proprietà analogica; la trasformazione bilineare la distorce, quindi in digitale è meno netta (per fase lineare esatta in digitale si usano i FIR).

Quando: conta la fedeltà nel dominio del tempo (impulsi, transitori, forme d’onda) più dell’attenuazione. Casi d’uso: catene audio/strumentazione dove l’overshoot è inaccettabile, ricostruzione di segnali a fronte ripido.

Le cinque famiglie non sono confrontabili a frequenza fissa

Prima di ogni confronto va sciolto un equivoco che manda fuori strada quasi tutti. Il campo Fc dello strumento non significa la stessa cosa per le cinque famiglie: ognuna normalizza il proprio prototipo su un bordo diverso, e all’ampiezza a Fc questo si vede subito.

Quindi la domanda «quanto attenua a 2·Fc?» non ha una risposta confrontabile: si stanno misurando cinque filtri con l’origine in cinque posti diversi. Chi confronta le famiglie a frequenza fissa conclude, per esempio, che il Chebyshev II è fenomenale a ordine 2 — quando in realtà sta solo leggendo il valore che gli è stato imposto per definizione.

Le due metriche che invece si confrontano

La prima è la larghezza di transizione: a che multiplo di Fc si raggiungono i 40 dB di attenuazione.

Letta così, la gerarchia è quella che ci si aspetta: l’ellittico vince (1,02·Fc a ordine 8, praticamente un muro), poi Chebyshev I, poi Butterworth e Bessel appaiati. Il Chebyshev II segna 1,00 a ogni ordine perché la sua Fc è il bordo di stopband: la riga non dice che è il più ripido, dice che la domanda è mal posta per quella famiglia.

La seconda è il pavimento di stopband, cioè la peggiore attenuazione in tutta la banda oscura:

Qui si vede la differenza strutturale fra le due metà della tabella. Butterworth, Bessel e Chebyshev I non hanno una specifica di stopband, quindi il pavimento scende con l’ordine — il Chebyshev I arriva a −88 dB a ordine 8. Chebyshev II ed ellittico hanno una stopband equiripple fissata dal parametro Rs: aumentare l’ordine non approfondisce l’attenuazione, la lascia inchiodata a −40 dB e stringe la transizione. Sono due filosofie di progetto diverse, e sceglierle a caso è il modo più comune di sprecare ordine.

Un caso di frontiera vale la nota: l’ellittico a ordine 2 non arriva a 40 dB (segna −13,8) perché con due poli non esiste una soluzione che soddisfi contemporaneamente 1 dB di ripple e 40 dB di stopband. Non è un difetto del calcolo: è l’equazione del grado che dice che quella specifica chiede almeno l’ordine 4.

Export CMSIS-DSP

I coefficienti sono esportati per arm_biquad_cascade_df2T_f32 di Arm: 5 valori per stadio nell’ordine {b0, b1, b2, a1, a2}, con a1 e a2 negati rispetto alla funzione di trasferimento. CMSIS implementa infatti:

$$y[n]=b_{0}x[n]+b_{1}x[n-1]+b_{2}x[n-2]-a_{1}y[n-1]-a_{2}y[n-2]$$

Il segno meno è già incorporato nei coefficienti generati: si copia l’array float32_t direttamente nel firmware.

Validazione e limiti

La risposta in ampiezza di tutte e cinque le famiglie (LP e HP) è stata confrontata con scipy.signal (butter, cheby1, cheby2, ellip, bessel → sosfreqz), con tutti i poli dentro il cerchio unitario (test automatici inclusi nel repository).

I riferimenti sono trascritti a due decimali, quindi un confronto con questa tabella non può dire più di ±0,005 dB: i millesimi di differenza sono il residuo della troncatura, non un errore del progetto.

Dove esiste una forma chiusa il riferimento si calcola invece di trascriverlo, e la verifica diventa esatta. Per il Butterworth bilineare con prewarping

$$|H(f)|=\frac{1}{\sqrt{1+{\left(\frac{\tan (\pi f/f_s)}{\tan (\pi f_c/f_s)}\right)}^{2N}}}$$

e a f = 2·f_c il rapporto delle tangenti vale esattamente √5, quindi il valore atteso è −10·log₁₀(1 + 5^N) = −27,965743 dB per N = 4 — la cifra che compare in tabella. La cascata di biquad generata coincide con la forma chiusa a precisione macchina su tutta la banda e per tutti gli ordini 2–8: lo verifica il test Butterworth contro la forma chiusa, con tolleranza 10⁻⁹ dB.

Un avvertimento sul metodo, se volete rifare la verifica: il dB va valutato esattamente alla frequenza, non letto dal punto più vicino della griglia della risposta. A 512 punti quella lettura sbaglia di 0,09 dB a 200 Hz, venti volte la differenza che si sta cercando.

L’ellittico usa l’algoritmo di Orfanidis (trasformazione di Landen per le funzioni ellittiche di Jacobi).

Limiti v0.2: passa-basso e passa-alto, ordini pari 2–8. Bandpass/bandstop in arrivo. Cutoff Fc = bordo della banda passante (−3 dB per Butterworth; bordo di ripple per Chebyshev I ed Ellittico; bordo di stopband per Chebyshev II; normalizzazione del ritardo per Bessel). Strumento dimostrativo: non sostituisce gli strumenti di progetto e verifica di una toolchain DSP di produzione.