Elaborazione Dati e Segnali Biomedici II - PDF

Segnale Elettroencefalografico Elaborazione Dati e Segnali Biomedici II A.A. 2023-24  Pietro Aricò, PhD [email protected] 27/09/2023 Summary Richia...

Segnale Elettroencefalografico Elaborazione Dati e Segnali Biomedici II A.A. 2023-24  Pietro Aricò, PhD [email protected] 27/09/2023 Summary Richiami di fisiologia alla base del segnale EEG Analisi nel tempo Analisi in frequenza Catena di processamento di un segnale EEG Strumentazione per l’acquisizione del segnale EEG Live Demo  Sistema Nervoso Il sistema nervoso comprende il sistema nervoso centrale e il sistema nervoso periferico. Il sistema nervoso centrale è costituito dal cervello e dal midollo spinale. Il sistema nervoso periferico è costituito dal sistema nervoso somatico e da quello autonomo.  Sistema Nervoso Centrale Il Sistema Nervoso Centrale (SNC) raccoglie informazioni sull'ambiente attraverso le sensazioni, controlla il pensiero e il controllo motorio. Il cervello è organizzato in lobi, preposti a specifiche funzioni: ✓ Lobo Frontale ✓ Lobo Parietale ✓ Lobo Occipitale ✓ Lobo Temporale  Sistema Nervoso Centrale L'area del lobo frontale del cervello è coinvolta in diverse attività importanti, tra cui la funzione motoria, la risoluzione dei problemi, la memoria, il linguaggio, il giudizio, il controllo degli impulsi emotivi ed il comportamento sociale.  Sistema Nervoso Centrale Questa parte del cervello svolge diverse funzioni. In primo luogo, le funzioni cognitive di sensazione e percezione. Questi input sensoriali vengono poi integrati per formare un sistema di coordinate spaziali corrispondente all'ambiente.  Sistema Nervoso Centrale La capacità di elaborare le immagini visive si trova nel lobo occipitale. Questa parte del cervello gestisce la percezione del movimento, la discriminazione dei colori e l'elaborazione visiva e spaziale.  Sistema Nervoso Centrale I lobi temporali sono altamente associati alle capacità di memoria e sono coinvolti nell'organizzazione primaria degli input sensoriali (Read, 1981). Quest'area del cervello è coinvolta nella risposta emotiva, nella memoria e nel riconoscimento del linguaggio. La responsabilità di questi lobi comprende anche funzioni linguistiche come la denominazione e la comprensione verbale. Neuroimaging techniques  Neurone I neuroni comunicano tra loro There are 60,000,000,000 neurons in the human brain attraverso segnali elettrochimici. La connessione tra due Dendrite neuroni è chiamata sinapsi. Axon I neuroni presinaptici Terminal rilasciano neurotrasmettitori che creano potenziali Myelin Nucleus all'interno del neurone Axon Post-synaptic -- -- Neuron + ++ postsinaptico. - --- + - Pre-synaptic Neuron - ++++ Se il neurone postsinaptico ---- + diventa più positivo, si parla di Gerstner & Kistler, 2002 EPSP (Potenziale postsinaptico eccitatorio).  Neurone There are 60,000,000,000 neurons in the human brain I neuroni presinaptici rilasciano neurotrasmettitori che creano potenziali Dendrite all'interno del neurone Axon postsinaptico. Terminal Se il neurone postsinaptico Myelin diventa più negativo, si parla Nucleus di IPSP (potenziali postsinaptici Axon Post-synaptic Neuron inibitori). Pre-synaptic Neuron Gerstner & Kistler, 2002  Neurone È possibile misurare l'attivazione (EPSP e IPSP) delle cellule piramidali giganti. Le cellule piramidali giganti che sono allineate perpendicolarmente alla superficie! Adjamian (2014)  Segnale EEG L'EEG è una registrazione dell'attività elettrica del cervello, effettuata da elettrodi sulla superficie del cuoio capelluto (EEG di superficie) o da elettrodi ad ago inseriti nel cervello (EEG invasivo). L'EEG viene misurato utilizzando elettrodi che registrano la differenza di potenziale elettrico tra un elettrodo con un segnale neurale attivo e un elettrodo posizionato su una regione presumibilmente inattiva. Queste registrazioni sono i potenziali di campo risultanti di molti neuroni attivi.  Segnale EEG La larghezza di banda di questo segnale va da meno di 1 Hz a oltre 100 Hz. La Federazione Internazionale di Elettroencefalografia e Neurofisiologia Clinica ha adottato un sistema di misurazione EEG a 10-20 (Jasper, 1958). Tale sistema prevede il posizionamento fisico standardizzato degli elettrodi sullo scalpo. Gli elettrodi sono etichettati in base alle aree cerebrali adiacenti: F (frontale), C (centrale), T (temporale), P (posteriore) e O (occipitale), con numeri dispari sul lato sinistro e numeri pari sul lato destro.  Potenziali evocati I potenziali evocati sono componenti temporali dell'EEG che sorgono in risposta a uno stimolo, che può essere elettrico, visivo, uditivo, tattile, ecc. Tali segnali sono spesso valutati mediando un certo numero di prove per migliorare il rapporto segnale/rumore (prossime lezioni vedremo come)  Potenziali evento correlati (ERPs) Il potenziale evento-correlato (ERP) P300 è una deflessione positiva del segnale EEG che si manifesta nel processo decisionale. Visual Oddball Brain-Computer Interface «A BCI is a system that measures CNS activity and converts it into artificial output that replaces, restores, enhances, supplements, or improves natural CNS output and thereby changes the ongoing interactions between the CNS and its external on internal environment» Wolpaw & Wolpaw, 2012 P300 Speller: Brain-Computer Interface P300 Speller Active BCI (communication & control, restoration, rehab) provaFFT.m MATLAB Dominio della frequenza La frequenza di una sinusoide è l'inverso del periodo ed è pari al numero di volte che l'onda si ripete in un secondo (si misura in Hertz, Hz). Nella figura sono riportate una sinusoide con frequenza 1 Hz, una a 3 Hz e una a 1/3 Hz. Più l'oscillazione è veloce nel tempo e maggiore è la frequenza risultante. provaFFT.m MATLAB Dominio della frequenza Se sommo questi segnali ottengo una onda un po’ più complicata, perché nel tempo i segnali si sovrappongono. Nel dominio della frequenza invece questi tre segnali rimangono sempre separati.  Dominio della frequenza In parole povere, un grafico del dominio del tempo mostra come un segnale cambia nel tempo, mentre un grafico del dominio della frequenza mostra quanta parte del segnale si trova all'interno di ciascuna banda di frequenza in un intervallo di frequenze. La trasformata di Fourier converte la rappresentazione della funzione nel dominio del tempo, mostrata in rosso, nella rappresentazione della funzione nel dominio della frequenza, mostrata in blu. Le frequenze componenti, distribuite nello spettro di frequenza, sono rappresentate come picchi nel dominio della frequenza.  Bande di frequenza L'attività oscillatoria dell'EEG spontaneo viene tipicamente classificata in cinque diverse bande di frequenza: delta (0- 4 Hz), theta (4-7), alfa (8-12), beta (12-30) e gamma (30- 100 Hz).  Bande di frequenza L'attività delta è caratterizzata da un'elevata ampiezza e da una bassa frequenza. Nella ricerca sul sonno è solitamente associata al sonno a onde lente. Si suggerisce che le onde delta rappresentino l'inizio delle fasi di sonno profondo negli adulti sani.  Bande di frequenza I ritmi theta sono importanti per le funzioni di apprendimento e memoria (Sammer et al., 2007), concentrazione (Hagemann, 2008). Meccanismo di controllo dell'attenzione (Kubota et al., 2001; Smith et al., 2001) e spesso è dimostrato che aumentano in presenza di una maggiore richiesta di compiti cognitivi (carico mentale).  Bande di frequenza L'attività della banda alfa si trova a livello della corteccia visiva (lobo occipitale) durante i periodi di rilassamento o di ozio (occhi chiusi ma svegli). È caratterizzata da oscillazioni regolari e di elevata ampiezza, con un massimo sugli elettrodi parietali e occipitali. Viene inibita durante compiti attivi.  Bande di frequenza L'onda beta è predominante quando l'uomo è sveglio. Spazialmente, predomina nelle aree frontali e centrali del cervello. È stato descritto che l'elevata potenza della banda beta è associata all'aumento dell'eccitazione e dell'attività.  Bande di frequenza La banda gamma è l'attività più veloce dell'EEG e si ritiene che sia poco frequente durante gli stati di veglia (Dooley, 2009). Tallon-Baudry et al. (2005) hanno rivelato che le aree della corteccia occipitale laterale svolgono un ruolo importante nella codifica degli stimoli visivi e mostrano ampie oscillazioni gamma diversamente influenzate dalla modulazione attenzionale. Studi recenti rivelano che la gamma è collegata a molte altre funzioni cognitive come l'attenzione, l'apprendimento, la memoria (Jensen et al., 2007) e la percezione del linguaggio (Eulitz et al., 1996).  Passive Brain-Computer Interface Nell'ultimo decennio il concetto di BCI si è evoluto, in particolare non è più l'utente a modulare la propria attività cerebrale per eseguire un controllo sull'ambiente circostante, ma è l'ambiente stesso che è in grado di decodificare lo stato mentale ed emotivo dell'utente e di utilizzare queste informazioni per adattare il suo comportamento. Aricò, P. et al. 2017. Passive BCI in operational environments: Insights, recent advances, and future trends. IEEE Transactions on Biomedical Engineering Aricò, P. et al. 2017. Human Factors and Neurophysiological Metrics in Air Traffic Control: A Critical Review. IEEE Reviews in Biomedical Engineering Aricò, P. et al. 2018. Passive BCI beyond the lab: Current trends and future directions. Physiological Measurement Aricò, P. et al. 2020. Brain–computer interfaces: Toward a daily life employment. Brain Sciences  Human Factor assessment: Why Oltre 1,2 milioni di persone muoiono ogni anno sulle strade del mondo, mentre altri milioni riportano gravi lesioni e vivono con conseguenze negative a lungo termine sulla salute. L'errore umano è la causa principale del 57% degli incidenti stradali e una concausa in oltre il 90% di essi. Oltre il 70% degli incidenti aerei è dovuto a errori umani, la maggior parte dei quali causati dal sovraccarico dei piloti o da alterazioni dello stato mentale. Gli errori medici causano un'elevata mortalità, circa 100.000 persone all'anno. Inoltre, circa il 10% dei pazienti ricoverati in ospedale ha subito complicazioni durante il trattamento a causa di errori medici. Il fattore umano è è il fattore più importante ma meno controllabile negli ambienti operativi. (Boeing Report, 2011) (Feyer and Williamson, 2011) (WHO Report, 2015) Spare Cognitive Capacity  The concept Cognitive Neuroscience applied to operational environments VIGILANCE WORKL EEG OAD STRESS EYES GAZE … BIOSIGNAL ECG PROCESSING OPERATOR’S GSR PSYCHOPHYSIOLOGICAL STATE … BIOSIGNALS Neurophysiological metrics of Human Factor Components How to do for an effective research in the passive BCI field ML applied to Processing and biosignals features classification extraction Neurophysiological EEG signal quality characterization and related HW Typical Neurometrics processing chain from EEG Buffer [X} sec, shift [Y] sec (EEG) Artifacts Correction PSD Frequency bands … Eye Blink Removing (e.g. ICA, REBLINCA, MWF) selection method Threshold criteria; Individual Alpha Trend estimation; Frequency (IAF) Channels Sample-to-sample difference Theta, Alpha, Beta, GFP Gamma Samples PRO: CONS: PRO: CONS: ✓ Calibration data Needed ✓ Group analysis ✓ Single subject analysis ✓ No calibration data ✓ Lower temporal resolution ✓ High temporal resolution ✓ High Computational (minute) (few seconds) effort ✓ Low Computational effort ML/DL method Discriminant function (y) ML-based Neurometric Formulas Neurometric (MA [Z]) (MA [Z]) Classification result E.g. WORKLOAD: Frontal Hard Theta/Parietal Alpha ytest(t) = ∑iwi train ∙ fi test(t) + btrain Low MACHINE LEARNING APPROACH STANDARD PHYSIOLOGICAL APPROACH  Adaptive automation in ATM domain 12 ATCOs (23 ± 2 years). Aricò, P. et al. 2016. Adaptive automation triggered by EEG-based mental workload index: A passive brain-computer interface application in realistic air traffic control environment. Frontiers in Human Neuroscience  HMI evaluation and ergonomy Which is the most effective technology to use for the operator? Augmented Technology for Remote Tower Operators In the MOTO project, neurophysiological measures have been used to compare different technologies (i.e. Augmentation vs No Augmentation) in terms of induced workload. Aricò, P. et al. 2019. How Neurophysiological Measures Can be Used to Enhance the Evaluation of Remote Tower Solutions. Frontiers in Human Neuroscience Marucci, M., … and Aricò, P. 2021. The impact of multisensory integration and perceptual load in virtual reality settings on performance, workload and presence. Scientific Reports  Stress on Air Traffic Controllers 16 ATCOs (23.8 ± 1.3 years). Borghini, G., Di Flumeri G., Aricò, P. et al. 2020. A multimodal and signals fusion approach for assessing the impact of stressful events on Air Traffic Controllers. Scientific Reports Stress on Air Traffic Controllers  Stress on Air Traffic Controllers The neurophisiological results are coherent with literature: stress is not a rapidly reversible phenomenon, otherwise it has long-time effects. Current effect: F(3, 45)=6,6606, p=,00081 4,2 4,0 Efficiency Evaluation (ISA Score) 3,8 3,6 The introduction of stressful events induces the reduction of the efficiency. 3,4 3,2 3,0 2,8 2,6 2,4 2,2 Var1 Var2 NewVar1 NewVar2 ATM Scenario Acrobatic sky divers Emergency landing Strumentazione EEG  Strumentazione EEG Properties Type of EEG electrode Sponge with Gel Dry Solid Gel liquid electrolyte Signal quality +++ --- + ++ Preparation / --- +++ 0 0 Cleaning time Comfort / +++ --- + ++ Recording time Residuals in hair --- +++ +++ ++ Strumentazione EEG Strumentazione EEG https://mindtooth-eeg.com/  Strumentazione EEG 2023 2014 Strumentazione EEG Measure of Neurometrics of a team of surgeons during a real laparoscopic operation. Mindtooth EEG headset  Experimental task Increasing difficulty (6 levels) ✓ Workload index EEG Frontal Theta and Parietal Alpha features Gevins et al.“High-resolution EEG mapping of cortical activation related to working memory: effects of task difficulty, type of processing, and practice.,” Cereb. cortex, 1997 ✓ Stress index 1. Select the 2. Direct the airplane EEG Parietal Beta airplane with to the landing strip by features the mouse the arrow keys Seo and Lee. “Stress and EEG,” Converg. Hybrid Inf. Technol., 2010  Elaborazione Dati e Segnali Biomedici II A.A. 2023-24  Dr. Pietro Aricò, PhD [email protected]  Biosegnali Elaborazione Dati e Segnali Biomedici II A.A. 2022-23  Pietro Aricò, PhD [email protected] 02/10/2023  Biosegnali – tipico sistema di misura biomedico Conversione in segnale elettrico attraverso il trasduttore. È poi necessaria un'elaborazione del segnale analogico, che spesso comprende l'amplificazione e il filtraggio passa-basso (o passa-banda). Poiché la maggior parte dell'elaborazione del segnale è più facile da implementare con metodi digitali, il segnale analogico viene convertito in formato digitale (ADC). Il segnale viene poi memorizzato o bufferizzato nella memoria per facilitare la successiva elaborazione. In alternativa, nelle applicazioni in tempo reale, i dati in ingresso vengono elaborati man mano che arrivano, spesso con un buffering minimo, e possono non essere memorizzati in modo permanente. Al segnale digitalizzato possono poi essere applicati algoritmi di elaborazione digitale del segnale, classificazione (es. seizures) ed eventuale visualizzazione. Summary Sistemi di misura Trasduttori Amplificatori/Detettori Filtri analogici Conversione ADC Quantizzazione Campionamento  Biosegnali Gran parte dell'attività dell'ingegneria biomedica, sia clinica che di ricerca, comporta la misurazione, l'elaborazione, l'analisi, la visualizzazione e/o la generazione di biosegnali. I segnali sono variazioni di energia che trasportano informazioni. L’elaborazione del segnale attraverso opportuno processamento, consente di estrarre tali informazioni, talvolta completamente nascoste all’interno del segnale stesso. Biosegnali  Biosegnali La variabile che trasporta l'informazione (la fluttuazione energetica specifica) dipende dal tipo di energia coinvolta. I segnali biologici sono solitamente codificati in variazioni di energia elettrica, chimica o meccanica, anche se, occasionalmente, sono interessanti le variazioni di energia termica. Energia Variabili Misure comunemente usate Attività chimica e/o Ioni nel sangue, O2, CO2, pH, Chimica concentrazione concentrazioni ormonali… Posizione, Forza, Coppia, Movimento muscolare, Meccanica Pressione pressioni cardiovascolari… Elettrica Tensione, Corrente EEG, ECG, GSR… Temperatura corporea, Termica Temperatura termografia  Trasduttori Un trasduttore, nella definizione del termine, è un dispositivo che converte l'energia da una forma all'altra (es. un motore). Nelle applicazioni di elaborazione dei segnali, l'energia viene utilizzata per trasportare informazioni; lo scopo della conversione energetica è trasferire tali informazioni, non trasformare l'energia come nel caso di un motore. Nei sistemi di misura, la maggior parte dei i trasduttori converte l'energia non elettrica in un segnale elettronico. Caso particolare l’elettrodo, un trasduttore che converte l'energia elettrica ionica in energia elettrica elettronica. EEG ECG  Trasduttori L'energia che viene convertita dal trasduttore di ingresso può essere generata dal processo fisiologico stesso, può essere energia indirettamente correlata al processo fisiologico o può essere energia prodotta da una fonte esterna (es. raggi X vengono assorbiti dal tessuto, e la misurazione dell’assorbimento serve a costruire un’immagine). Molti processi fisiologici producono energia che può essere rilevata direttamente (es. la misurazione dell'attività elettrica nel cuore, nei muscoli o nel cervello fornisce esempi di misurazione diretta dell'energia fisiologica. Per queste misurazioni, l'energia è già elettrica e deve solo essere convertita da corrente ionica a corrente elettronica utilizzando un "elettrodo". A queste fonti viene solitamente attribuito il termine ExG, dove "x" rappresenta il processo fisiologico che produce energia elettrica, ECG-elettrocardiogramma, EEG-elettroencefalogramma, EMG-elettromiogramma, EOG-elettrooculogramma. Un'eccezione a questa terminologia è la risposta galvanica cutanea (GSR), l'attività elettrica generata dalla pelle.  Amplificatore/detettore La progettazione del primo stadio analogico di un sistema di misura biomedico dipende dal funzionamento di base del trasduttore. Se il trasduttore si basa su una variazione della proprietà elettrica, il primo stadio deve convertire tale variazione in una variazione di tensione. Se nel trasduttore è presente un solo elemento sensibile, il segnale del trasduttore si dice «single ended» (ad esempio, il termistore). Se il trasduttore produce un'uscita in corrente, viene utilizzato un amplificatore da corrente a tensione per produrre un'uscita in tensione (ad esempio, i rilevatori di luce). In alcune circostanze, può essere necessaria un'amplificazione aggiuntiva oltre al primo stadio.  Amplificatore biomedico È un amplificatore operazionale, definito per strumentazione, particolarmente adatto all'amplificazione di segnali deboli provenienti da trasduttori nel campo delle misure elettroniche e degli strumenti professionali. La sua struttura è derivata dall'amplificatore differenziale: rispetto a questo ha due amplificatori operazionali in più che migliorano (aumentandola) l'impedenza d'ingresso e permettono di variare l'amplificazione del segnale differenziale d'ingresso Vin variando solo un componente di guadagno (Rgain). V1 R2 R3 La relazione tra tensione di uscita e tensione differenziale di ingresso è data da: R1 Rgain Vout R1 V2 R2 R3  Processamento del segnale analogico Sebbene l'elaborazione del segnale più completa venga eseguita su dati digitalizzati utilizzando algoritmi implementati nel software, è solitamente necessaria un'elaborazione del segnale analogico. Una delle elaborazioni più frequenti del segnale analogico utilizza filtri analogici, limitando la gamma di frequenze o la larghezza di banda del segnale. È questo filtraggio che di solito stabilisce la larghezza di banda del sistema di misura complessivo. Poiché la larghezza di banda del segnale ha un impatto diretto sul processo di conversione di un segnale analogico in un segnale digitale equivalente, è spesso un elemento essenziale in qualsiasi sistema di misura biomedico.  Filtri analogici I filtri analogici sono dispositivi elettronici che eliminano frequenze selezionate. I filtri vengono solitamente definiti in base alla gamma di frequenze che non sopprimono. Pertanto, i filtri passa-basso consentono il passaggio delle basse frequenze con un'attenuazione minima, mentre le frequenze più alte vengono attenuate. Al contrario, i filtri passa-alto fanno passare le frequenze alte, ma attenuano quelle basse. I filtri passa-banda rifiutano le frequenze al di sopra e al di sotto di una regione di banda passante. Un'eccezione a questa terminologia è rappresentata dai filtri band-stop, che fanno passare le frequenze su entrambi i lati di una gamma di frequenze attenuate (es. Notch).  Filtri analogici Spettro Il guadagno del filtro è il rapporto tra l'ampiezza del segnale di uscita e l'ampiezza del segnale di ingresso in funzione della frequenza (dB). Ad esempio, Il filtro passa-basso ha un guadagno di 1,0 per le frequenze più basse. Ciò significa che l'uscita è uguale all'ingresso a quelle frequenze. Tuttavia, all'aumentare della frequenza, il guadagno diminuisce, indicando che anche il segnale di uscita diminuisce per un dato segnale di ingresso.  Filtri analogici Quando il rapporto uscita/ingresso è dato analiticamente come funzione della frequenza, viene definito funzione di trasferimento. La rappresentazione grafica della funzione di trasferimento viene definita diagramma di Bode. La larghezza di banda di un filtro è definita dalla gamma di frequenze che non vengono attenuate. Queste frequenze non attenuate vengono anche chiamate frequenze della banda passante.  Filtri analogici: banda passante I filtri ideali hanno banda passante perfettamente piatta e una pendenza di attenuazione infinita. I filtri reali possono essere abbastanza piatti nella regione della banda passante, ma attenuano con una pendenza più dolce  Filtri analogici: banda passante L’ampiezza di banda di un filtro ideale è facile da determinare. Per un filtro reale determinare la larghezza di banda è più complesso. Si deve identificare una frequenza che definisca il confine tra le porzioni attenuate e non. Questo limite è stato definito in modo un po' arbitrario come la frequenza in cui l'attenuazione è di 3 dB.  Filtri analogici: Ordine La pendenza della curva di attenuazione di un filtro è correlata alla complessità del filtro: i filtri più complessi hanno una pendenza maggiore, avvicinandosi al filtro ideale. Nei filtri analogici, la complessità è proporzionale al numero di elementi di accumulo di energia nel circuito (condensatori). Ogni dispositivo di accumulo di energia indipendente porta a un ordine aggiuntivo di un polinomio nel denominatore dell'equazione della funzione di trasferimento.  Filtri analogici: Ordine La complessità di un filtro è equivalente al numero di poli (radici dell’equazione del denominatore) della funzione di trasferimento. Ad esempio, un filtro del secondo ordine o a due poli ha una funzione di trasferimento con una polinomiale del secondo ordine nel denominatore e conterrebbe due elementi di accumulo di energia indipendenti.  Filtri analogici: Nitidezza inziale Sia la pendenza che la nitidezza iniziale aumentano con l'ordine del filtro (numero di poli), ma l'aumento dell'ordine del filtro aumenta anche la complessità e quindi il costo del filtro.  Filtri analogici: Nitidezza inziale È possibile aumentare l'acutezza iniziale delle caratteristiche di attenuazione del filtro senza aumentare l'ordine del filtro, se si è disposti ad accettare qualche irregolarità o ondulazione nella banda passante. La Figura mostra due filtri passa-basso del quarto ordine che hanno la stessa frequenza di taglio, ma differiscono per l'acutezza iniziale dell'attenuazione. Conversione ADC  Conversione ADC L'ultimo elemento analogico del tipico sistema di misura è l'interfaccia tra il mondo analogico e quello digitale: l'ADC. Come suggerisce il nome, questo componente elettronico converte una tensione analogica in un numero digitale equivalente. Nel processo di conversione ADC, una forma d'onda analogica o continua, x(t), viene convertita in una forma d'onda discreta, x[n], una funzione di numeri reali definiti come numeri interi in punti discreti del tempo. Questi numeri sono chiamati campioni e i punti discreti nel tempo sono di solito presi a intervalli regolari denominati intervallo di campionamento, Ts. L'intervallo di campionamento può anche essere definito da una frequenza denominata frequenza di campionamento  Conversione ADC La digitalizzazione di un segnale continuo (in alto a sinistra), richiede di tagliare il segnale sia nel tempo che nell'ampiezza (a destra). Il risultato è una serie di numeri discreti (quadrati) che approssimano il segnale originale. Il segnale digitalizzato risultante, in basso a sinistra, consiste in una fs=1Hz serie di valori numerici discreti campionati a intervalli di tempo discreti. In questo esempio, x[n] = 2, 3, 2, 1, 0, 2, 3, 5, e 8. Campionamento.m MATLAB Esempio in Matlab Generare un'onda sinusoidale discreta a 2 Hz. Frequenza di campionamento di 100Hz. Campionamento.m MATLAB Esempio in Matlab clear; close all; Fs = 100; % Sample Rate (Hz) Ts=1/Fs; % Intervallo di campionamento (s) TT = 1; % Total time = 1 sec f = 2; % Frequency of the sin = 2 Hz %% MAIN t = Ts:Ts:TT; % Vettore dei tempi x = sin(2*pi*f*t); % Generate desired sine wave plot(t,x,'.k'); % Plot sine wave as discrete points xlabel('Time (sec)'); % and label ylabel('x(t)');  Quantizzazione La suddivisione dell'ampiezza del segnale in livelli discreti, viene definita quantizzazione. Il valore numerico equivalente del segnale quantizzato può solo approssimare il livello del segnale analogico e il grado di approssimazione (risoluzione) dipende dal numero di valori diversi utilizzati per rappresentare il segnale. Poiché i segnali digitali sono rappresentati come numeri binari, la lunghezza dei bit dei numeri binari utilizzati determina il livello di quantizzazione. La tensione minima che può essere risolta e la dimensione della fetta di ampiezza sono note come livello di quantizzazione, q. La dimensione della fetta in volt è l'intervallo di tensione del convertitore diviso per il numero di valori discreti disponibili, supponendo che tutti i bit siano convertiti accuratamente. Se il convertitore produce un numero binario con b bit precisi, il numero di valori non nulli è (2b - 1), dove b è il numero di bit nell'uscita binaria. Se l'intervallo di tensione del convertitore varia tra 0 e VMAX volt, la dimensione del passo di quantizzazione, q, in volt, è data da:  Quantizzazione Se avessimo un ADC a 12 bit, e il range fosse 0-10V, la risoluzione equivalente dell’ADC al bit meno significativo (LSB) sarebbe: Solitamente i sistemi di acquisizione di biosegnali hanno una risoluzione di 12 o massimo 16 bit, sufficiente per coprire la gamma dinamica dei segnali di interesse (es. amplificatori audio possono arrivare a 24 bit)  Quantizzazione L'effetto della quantizzazione può essere visto come un rumore aggiunto al segnale originale. Questo rumore dipende dal livello di quantizzazione, q. Se esiste un numero sufficiente di livelli di quantizzazione, l'errore di quantizzazione può essere modellato come un rumore bianco additivo indipendente con una distribuzione uniforme tra ±q/2 e media zero.  Campionamento La suddivisione del segnale in punti discreti nel tempo viene definita campionamento. Tale suddivisione campiona la forma d'onda continua, x(t), in punti discreti nel tempo, 1Ts , 2Ts , 3Ts ,...., nTs , dove Ts è l'intervallo di campionamento. Poiché lo scopo del campionamento è produrre una copia accettabile della forma d'onda originale, la questione critica è: quanto bene questa copia rappresenta l'originale? In altre parole, è possibile ricostruire l'originale dalla copia digitalizzata? Se sì, la copia è chiaramente adeguata. La risposta a questa domanda dipende dalla frequenza di campionamento della forma d'onda analogica rispetto alle frequenze che contiene  Campionamento ll teorema di campionamento di Shannon afferma che qualsiasi forma d'onda sinusoidale può essere ricostruita in modo univoco a condizione che venga campionata almeno due volte in un periodo. (Cioè, la frequenza di campionamento, fs , deve essere > 2 fsinusoide). In altre parole, per specificare in modo univoco una sinusoide sono necessari solo due campioni equidistanti, che possono essere prelevati in qualsiasi punto del ciclo.  Campionamento La Figura mostra un'onda sinusoidale definita da due punti che coprono un periodo di tempo leggermente inferiore a 1 ciclo dell'onda sinusoidale. Secondo il teorema di Shannon, non esistono altre sinusoidi di frequenza inferiore che possano passare per questi due punti. Pertanto, questi due punti definiscono in modo univoco un'onda sinusoidale di frequenza minima.  Campionamento Esistono però molte sinusoidi di frequenza superiore che sono definite da questi due punti, ad esempio tutti i punti a frequenza multipla dell'originale. Non si può escludere dunque che questi due punti rappresentino anche tutte una serie di altre sinusoidi (idealmente infinite, a frequenza superiore all’originale)  Campionamento La Figura presenta un esempio di spettro del segnale prima (a) e dopo (b) il campionamento. Un singolo punto dello spettro rappresenta una singola forma d'onda sinusoidale; quindi, prima del campionamento, il segnale è una miscela di sette sinusoidi alle frequenze [1-7] Hz. Dopo il campionamento, lo spettro originale è ancora presente, ma in più si trova in molti altri punti dello spettro (sovrarmoniche)  Campionamento … … Le frequenze sono generate dalla matematica che descrive l'operazione di campionamento. Sebbene siano costrutti matematici, hanno un'influenza sul segnale campionato perché generano l'onda sinusoidale aggiunta riflessa a sinistra di fs (e a sinistra dei multipli di fs).  Campionamento … … Posto che il segnale campionato non è lo stesso dell'originale, la domanda cruciale è: è possibile recuperare il segnale originale dalla sua versione campionata? Banalmente, se possiamo ricostruire lo spettro originale non campionato dallo spettro campionato, allora il segnale digitale è una rappresentazione adeguata dell'originale.  Campionamento La Figura mostra un ingrandimento di un solo segmento dello spettro, il periodo compreso tra 0 e fs Hz. Confrontando questo spettro con lo spettro del segnale originale, si nota che i due spettri sono uguali per la prima metà dello spettro fino a fs /2 e la seconda metà è solo l'immagine speculare (frequenze negative teoriche). Sembrerebbe che si possa ottenere uno spettro di frequenza identico all'originale se si eliminano tutte le frequenze al di sopra di fs/2. In effetti, possiamo eliminare le frequenze al di sopra di fs/2 filtrando il segnale stesso. Finché riusciamo a tornare allo spettro originale, i dati campionati sono un utile riflesso dei dati originali. La frequenza fs/2 è definita: frequenza di Nyquist  Campionamento Abbiamo 4 sinusoidi con frequenze 100, 200, 300, 400 Hz, e campionando il segnale a 1000Hz, riusciremmo, ignorando le frequenze maggiori di fs/2 a ricostruire il segnale originale (figura a). Al contrario, se ci fossero frequenze superiori alla frequenza di Nyquist (es, 650 e 850 Hz), queste sarebbero riflesse nella metà inferiore dello spettro. Questo effetto è definito aliasing  Campionamento Questa strategia di ignorare tutte le frequenze superiori alla frequenza di Nyquist (fs /2) può essere utilizzata solo se il segnale originale non ha componenti spettrali a o superiori a fs /2. Teorema del campionamento di Shannon: Il segnale originale può essere recuperato da un segnale campionato a condizione che la frequenza di campionamento sia più del doppio della frequenza massima contenuta nell'originale.  Campionamento /2 /2 F Nyquist F Nyquist La frequenza di campionamento deve dunque essere sufficiente da coprire le informazioni spettrali di interesse, secondo il teorema di Shannon. Chiaramente aumentando troppo la frequenza di campionamento di rischia di appesantire il sistema informativo, acquisendo un numero ridondante di dati. Campionamento2.m MATLAB Esempio in Matlab clear; close all; Fs=2000; Ts = 1/Fs; % Intervallo di campionamento TT = 1; % Total time = 1 sec f1 = 100; % Frequency = 100 Hz f2 = 700; % Frequency = 700 Hz A1 = 1; % Ampiezza Sin 1 A2 = 2; % Ampiezza Sin 1 Campionamento2.m MATLAB Esempio in Matlab clear; close all; Fs=2000; Ts = 1/Fs; % Intervallo di campionamento TT = 1; % Total time = 1 sec f1 = 100; % Frequency = 100 Hz f2 = 700; % Frequency = 700 Hz %% MAIN t = Ts:Ts:TT; % Vettore dei tempi x = sin(2*pi*f1*t)+sin(2*pi*f2*t); % Generate desired sine wave plot(t,x); % Plot sine wave as discrete points xlabel('Time (sec)'); % and label ylabel('x(t)'); Y=fft(x); figure; plot(abs(Y))  Elaborazione Dati e Segnali Biomedici II A.A. 2023-24  Dr. Pietro Aricò, PhD [email protected]  Filtri Digitali Elaborazione Dati e Segnali Biomedici II A.A. 2023-24  Pietro Aricò, PhD [email protected] 04/10/2022 Summary Media mobile Grand average Trasformata Z Filtri FIR Filtri IIR Filtro di Wiener Esempi in Matlab  Filtri digitali Definizione. Un filtro digitale reale Tn è una qualsiasi funzione che restituisce un segnale a valori reali (output) per ogni segnale in ingresso (input). Si può esprimere la relazione tra input e output di un filtro digitale con la notazione con x(·) il segnale in input e y(n) il segnale in output al tempo n.  Filtri digitali Un filtro per essere lineare deve soddisfare due proprietà: ✓ Ridimensionamento: l’ampiezza dell’output deve essere proporzionale all’ampiezza dell’input. ✓ Sovrapposizione: se si sommano due segnali e si filtrano, l’output del filtro sarà lo stesso che si ottiene filtrando singolarmente i due segnali e sommandoli successivamente.  Filtri digitali Un filtro tempo-invariante è un filtro tale che l’output generato da un segnale ritardato è uguale all’output generato dal segnale originale, ritardato della stessa quantità. In particolare, Un filtro digitale Tn è tempo-invariante se, per ogni segnale x, si ha: dove ogni operatore di shift dell’N-esimo campionamento è definito come: cioè l’operatore ritarda il segnale in input di N campioni movAver.m MATLAB Media mobile Consideriamo la precedente equazione. La media di tre campioni successivi genera un nuovo segnale con componenti ad alta frequenza ridotte. dove x[n] è il segnale originale e y[n] è il nuovo segnale. Si può notare che questo processo di media riduce l'influenza dei campioni «fuori scala», o «outliers» e attenua le transizioni brusche.  Media mobile La media mobile è un processo che viene spesso descritto come un filtro passa basso. I filtri vengono utilizzati per migliorare l'SNR di un segnale e sono spesso pensati in termini di rimodellamento dello spettro di un segnale. Quindi, mentre il filtro della media mobile è un processo nel dominio del tempo, i suoi effetti sono solitamente concettualizzati nel dominio della frequenza. La maggior parte del rumore è a banda larga (il rumore più ampio è il rumore bianco con uno spettro piatto), mentre la maggior parte dei segnali è a banda stretta. Di conseguenza, dovrebbe essere sempre possibile migliorare l'SNR rimodellando opportunamente lo spettro di una forma d'onda mediante un filtro. noise_freq.m Rumore bianco MATLAB Media mobile La media mobile è un processo che viene spesso descritto come un filtro passa basso. I filtri vengono utilizzati per migliorare l'SNR di un segnale e sono spesso pensati in termini di rimodellamento dello spettro di un segnale. Quindi, mentre il filtro della media mobile è un processo nel dominio del tempo, i suoi effetti sono solitamente concettualizzati nel dominio della frequenza. La maggior parte del rumore è a banda larga (il EEG rumore più ampio è il rumore bianco con uno spettro piatto), mentre la maggior parte dei segnali è a banda stretta. Di conseguenza, dovrebbe essere sempre possibile migliorare l'SNR rimodellando opportunamente lo spettro di una forma d'onda mediante un filtro.  Riduzione del rumore attraverso grand average Quando si effettuano più misure, vengono generati più valori o segnali. Se queste misurazioni vengono combinate o sommate tra loro, il valore o il segnale combinato ha una media che è la media delle singole medie calcolate su ogni segnale. Lo stesso vale per la varianza: le varianze si sommano e la varianza media della misura combinata è la media delle singole varianze: Quando si sommano più segnali, si può dimostrare che la deviazione standard del rumore si riduce di un fattore pari a 1/ 𝑁 , dove N è il numero di misure che vengono mediate. In altre parole, la media delle misure provenienti da sensori diversi, o la media di più misure provenienti dalla stessa fonte, riduce la deviazione standard della variabilità della misura della radice quadrata del numero di medie.  Riduzione del rumore attraverso grand average Il grand average è una tecnica di elaborazione del segnale semplice ma potente per ridurre il rumore quando sono possibili osservazioni multiple del segnale. Tali osservazioni multiple possono provenire da più sensori, ma in molte applicazioni biomediche le osservazioni multiple provengono da risposte ripetute allo stesso stimolo. Ad esempio, questa tecnica è enormemente usata per l’analisi dei potenziali evocati, in cui il rapporto segnale rumore del singolo potenziale è talmente basso da non riuscire a visualizzare il potenziali stesso, immerso all’interno del rumore (segnale EEG di base). P300ex.m MATLAB Esempio Matlab: ERPs and grand average Sia data una matrice con 1152 ERPs su stimoli target, calcolare un grand average costituito da 16 medie.  Riduzione del rumore attraverso grand average Il grand average è una tecnica di elaborazione del segnale semplice ma potente per ridurre il rumore quando sono possibili osservazioni multiple del segnale. Tali osservazioni multiple possono provenire da più sensori, ma in molte applicazioni biomediche le osservazioni multiple provengono da risposte ripetute allo stesso stimolo. Ad esempio, questa tecnica è enormemente usata per l’analisi dei potenziali evocati, in cui il rapporto segnale rumore del singolo potenziale è talmente basso da non riuscire a visualizzare il potenziali stesso, immerso all’interno del rumore (segnale EEG di base).  Equazione alle differenze finite Una delle rappresentazioni di un filtro è l’equazione alle differenze finite. Essa calcola un output in base agli input passati e presenti al tempo n. Possiamo definirla nella maniera più generale possibile per un filtro LTI (lineare e tempo invariante, causale): con x segnale in input, y segnale in output, e le costanti ai; bi sono coefficienti. L’output finale è dovuto, oltre all’input, anche all’output degli istanti precedenti (ad esempio y(n − 1)), l’uso di questi fattori è detto feedback e ogni filtro avente almeno un elemento di feedback è detto ricorsivo. In particolare i coefficienti di questi elementi di feedback ai sono detti feedback coefficients, mentre gli altri bi sono detti feedforward coefficients. Deduciamo quindi che un filtro è ricorsivo se e solo se ai ≠ 0 per qualche i > 0. I filtri ricorsivi sono anche chiamati a infinite impulse response o IIR, se invece non abbiamo elementi di feedback (ai = 0 ∀i > 0), il filtro è detto finite impulse response o FIR.  Risposta all’impulso e stabilità Oltre che tramite l’equazione alle differenze finite, ogni filtro LTI può essere rappresentato in base alla sua risposta ad uno specifico segnale, detto impulso. Questa risposta è chiamata appunto risposta all’impulso. Definizione. Il segnale impulso è indicato con δ(n) e definito come: Definizione. Indichiamo la risposta di un filtro a δ(n) con h(n), ed è definita come: Definizione. Un filtro LTI è detto stabile se la risposta all’impulso h(n) tende a 0 per n che va all’infinito. Poichè solo i coefficienti di feedback possono causare instabilità, ogni filtro di ordine finito non ricorsivo è stabile.  Convoluzione Deriviamo ora quella che è detta rappresentazione di convoluzione per un filtro LTI, nel modo più generale possibile. Innanzitutto rappresentiamo un segnale x(·) come una combinazione lineare di impulsi spostati: dove ∗ è l’operatore di convoluzione. Quindi ogni campione può essere visto come un impulso di una certa ampiezza ad un certo istante, ognuno dei quali produce una risposta all’impulso. Se il nostro filtro oltre ad essere lineare è anche invariante nel tempo possiamo scrivere h(n; i) = h(n; −i), di conseguenza la formule precedente diventa Quindi l’output del filtro y è la convoluzione dell’input x con la risposta all’impulso h  Convoluzione e risposta all'impulso L'ingresso, x(t), e l'uscita, y(t), sono collegati tramite convoluzione attraverso la funzione h(t). Poiché sono tutte funzioni del tempo, la convoluzione è un'operazione nel dominio del tempo. Nel dominio digitale, le tre funzioni temporali diventano vettori discreti campionati: x[n], y[n] e h[n].  Convoluzione e risposta all'impulso La funzione h(t) è nota come risposta all'impulso. Come dice il nome, è la risposta del sistema a un ingresso impulsivo standard. Un ingresso impulsivo (definito anche δ(t)) è un impulso molto breve con un'area unitaria.  Convoluzione e risposta all'impulso Se si conosce la risposta del sistema a un impulso, è possibile determinare la sua risposta a qualsiasi ingresso semplicemente dividendo l'ingresso in una sequenza di impulsi. Ciascuna fetta di tempo genererà la propria piccola risposta all'impulso. L'ampiezza e la posizione di questa risposta all'impulso sono determinate dall'ampiezza e dalla posizione del segmento di segnale di ingresso associato. Se la sovrapposizione e l'invarianza temporale sono valide, l'uscita può essere determinata sommando (o integrando per funzioni continue) le risposte all'impulso di tutti i segmenti del segnale di ingresso.  Convoluzione e risposta all'impulso Ogni segmento di un segnale di ingresso può essere rappresentato da un impulso e ciascuno di essi produce una propria risposta all'impulso. L'ampiezza e la posizione della risposta all'impulso sono determinate dall'ampiezza e dalla posizione del segmento di ingresso associato.  Convoluzione e risposta all'impulso Quando si implementa la convoluzione, è necessario invertire il segnale di ingresso a causa del modo in cui i segnali sono graficizzati. Il lato sinistro del segnale è in realtà il lato a tempo basso; è la prima porzione che entra nel sistema. Pertanto, il segmento più a sinistra del segnale di ingresso produce la prima risposta all'impulso e il segnale di ingresso procede attraverso il sistema all'indietro rispetto al modo in cui è tracciato.  Convoluzione e risposta all'impulso Il processo di convoluzione viene espresso matematicamente come: dove n è la variabile che fa scorrere l'ingresso invertito, x[-k], attraverso la risposta all'impulso, h[k], e K è la lunghezza della funzione più breve, di solito h[k].  Funzione di trasferimento La funzione di trasferimento H(z) è una rappresentazione algebrica di un filtro LTI nel dominio della frequenza. Essa è definita come Y (z)/X(z), dove Y (z) è la trasformata z del segnale di output y(n), mentre X(z) è la trasformata z dell’input x(n). Teorema dello spostamento. Il ritardo di ∆ campioni nel dominio del tempo corrisponde a una moltiplicazione per z−∆ nel dominio della frequenza, dove con l’operatore Z indichiamo la trasformata z.  Funzione di trasferimento La funzione di trasferimento H(z) è una rappresentazione algebrica di un filtro LTI nel dominio della frequenza. Essa è definita come Y (z)/X(z), dove Y (z) è la trasformata z del segnale di output y(n), mentre X(z) è la trasformata z dell’input x(n). Teorema della convoluzione. Per ogni coppia di segnali causali x e y, la convoluzione nel dominio del tempo diventa moltiplicazione nel dominio z. Cioè:  Funzione di trasferimento La funzione di trasferimento H(z) è una rappresentazione algebrica di un filtro LTI nel dominio della frequenza. Essa è definita come Y (z)/X(z), dove Y (z) è la trasformata z del segnale di output y(n), mentre X(z) è la trasformata z dell’input x(n). Ora ricordando che l’output y di un LTI è dato dalla convoluzione dell’input x con la risposta all’impulso h possiamo usare i teoremi appena dimostrati effettuando la trasformata z di y ottenendo con H(z) trasformata di h. Otteniamo cosi la funzione di trasferimento:  Analisi poli-zeri Ogni filtro digitale può essere descritto dai suoi poli e zeri, inoltre essi possono dare un’idea della risposta di un filtro. Come abbiamo visto si parte dalla funzione di trasferimento, se prendiamo un filtro LTI ricorsivo, applicando la trasformata Z all’equazione alle differenze finite avremo: Fattorizzando il primo coefficiente del numeratore b0, e chiamandolo g Fattorizzando numeratore e denominatore: Ora i numeri q1;…; qM sono le radici o gli zeri della funzione. Se z è uguale a uno di questi valori la funzione di trasferimento va a 0. Allo stesso modo i numeri p; …; pN sono detti poli della funzione, quando z assume uno di questi valori l’ampiezza della funzione di trasferimento va all’infinito.  Analisi poli-zeri Tramite l’analisi di poli e zeri è possibile determinare anche la stabilità di un filtro. Infatti possiamo dire che un filtro è stabile se e solo se tutti i suoi poli sono all’interno della circonferenza unitaria nel piano z. Dimostrazione: consideriamo una risposta all’impulso causale nella forma, la cui trasformata z è: (serie di Taylor) Quindi, l’ultima uguaglianza è vera solo per , cioè se. La funzione di trasferimento ha dunque un solo polo a , ed esiste solo per. Consideriamo ora cosa accade se R >1. Il polo di H(z) va all’esterno della circonferenza unitaria, e la risposta all’impulso aumenta esponenzialmente, quindi la definizione di stabilità data in precedenza viene violata. Visto che la trasformata sappiamo che esiste solo se , deduciamo che R ≥ 1 implica che la trasformata z non esiste più sulla circonferenza unitaria, e la risposta in frequenza non è più definita.  Analisi poli-zeri Abbiamo quindi provato che un filtro a un polo è stabile se e solo se il polo è all’interno della circonferenza unitaria. Per quanto riguarda una funzione di trasferimento qualsiasi, la sua espansione in fratti semplici mostra che il comportamento vicino ad ogni polo è quello di un filtro polare avente solo quel polo. Quindi possiamo estendere la definizione di stabilità alla presenza di tutti i poli all’interno della circonferenza.  Funzione di trasferimento La maggior parte delle funzioni di trasferimento saranno composte da polinomi di z sia al numeratore che al denominatore: dove le b[k] sono coefficienti costanti del numeratore e le a[l] sono coefficienti del denominatore. Dalla funzione di trasferimento digitale, H[z], è possibile determinare l'uscita per qualsiasi ingresso, dove Y-1[z] è la trasformata Z inversa.  Risposta in frequenza La risposta in frequenza di un filtro LTI può essere definita come lo spettro del segnale di output diviso lo spettro del segnale di input. Una proprietà della trasformata z è che, se calcolata sulla circonferenza unitaria z=ejωT, si ottiene lo spettro (FFT(x)). Per dimostrarlo basta sostituire z=ejωT nella definizione della trasformata, per ottenere Essendo z complesso può essere espresso da un’ampiezza r = |z| e una fase θ = , è possibile quindi definire la risposta in ampiezza e la risposta in fase del filtro. La risposta in ampiezza a valori reali specifica il gain del filtro per ogni frequenza La risposta in fase Θ a valori reali da lo spostamento di fase, che ogni componente sinusoidale subisce per ogni frequenza.  Risposta in frequenza La risposta in ampiezza a valori reali specifica il gain del filtro per ogni frequenza La risposta in fase Θ a valori reali da lo spostamento di fase, che ogni componente sinusoidale subisce per ogni frequenza. Risposta in ampiezza e in fase di un filtro ellittico  Elaborazione Dati e Segnali Biomedici II A.A. 2023-24  Dr. Pietro Aricò, PhD [email protected]  Filtri Digitali Elaborazione Dati e Segnali Biomedici II A.A. 2023-24  Pietro Aricò, PhD [email protected] 11/10/2023  Progettazione filtro La progettazione di un filtro digitale consiste nel trovare i coefficienti a[l] e b[k] che forniscono lo shape spettrale desiderato. Questo processo di progettazione è agevolato da routine MATLAB che generano i coefficienti a e b che producono la risposta in frequenza desiderata. Nel dominio Z, gli spettri di frequenza della funzione di trasferimento H(z) possono essere ottenuti sostituendo la frequenza reale, ejω, alla frequenza complessa, z. dove FT indica la trasformata di Fourier. Come per tutte le trasformate di Fourier, la frequenza effettiva può essere ottenuta dal numero armonico m dopo averla moltiplicata per fs /N o 1/(NTs) hfilters.m MATLAB Esempio Matlab Sistema con la seguente funzione di trasferimento: fs = 1 kHz, N=512. Calcolare le risposta del filtro ad un impulso.  Esempio Matlab Sistema con la seguente funzione di trasferimento: fs = 1 kHz, N=512. Calcolare le risposta del filtro ad un impulso. Notiamo che a=1; a=-0.2; a=0.8; b=0.2; b=0.5. Applichiamo l’equazione Convertiamo le frequenze in Hz, moltiplicando m * fs/N Esempio Matlab  Filtri FIR e IIR I filtri aventi risposta all’impulso finita (FIR) o infinita (IIR) sono particolari sistemi lineari e tempo-invarianti, e possono essere implementati e simulati su macchine digitali. Offrono numerosi vantaggi sui filtri analogici, come ad esempio la riprogrammabilità del software senza cambiare hardware, oppure la modifica in tempo reale dei coefficienti, rendendoli filtri adattivi. Nelle applicazioni pratiche si preferisce in genere un filtro FIR, se è richiesta una fase lineare, se la distorsione della fase è invece tollerabile si preferisce l’uso dei filtri IIR, che richiedono meno parametri quindi sono meno complessi ed occupano meno spazio in memoria. Poiché i filtri a fase lineare hanno un ritardo di gruppo (group delay, definito come derivata della fase rispetto alla pulsazione) costante, tutte le componenti di frequenza hanno un uguale ritardo nel tempo. Ciò significa che non vi è distorsione dovuta alla risposta in frequenza; in molte applicazioni, questo ritardo di gruppo costante è vantaggioso. Per contro, un filtro con una "fase non lineare", ha un ritardo di gruppo che varia con la frequenza, che porta come risultato ad una distorsione di fase.  Filtri con risposta all’impulso finita (FIR) Definizione. Un sistema LTI è detto FIR se la risposta h(n) all’impulso unitario è finita, nel senso che h(n) = 0 per n < 0 e per n ≥ K per un opportuno K > 0. I filtri FIR sono filtri a media mobile, generalmente con pesi non uniformi. Nei filtri FIR, i pesi del filtro (coefficienti) definiscono la risposta all'impulso, quindi ovviamente la risposta all'impulso sarà finita (stabilità del filtro). I filtri FIR vengono implementati convolvendo i pesi con il segnale di ingresso, il che equivale matematicamente a fare una media ponderata del segnale di ingresso (equazione alle differenze finite per un filtro FIR). b[k] definisce i pesi, e K il numero di pesi, o la lunghezza del filtro, x[n] rappresenta il segnale in ingresso, e y[n] quello in uscita.  Filtri con risposta all’impulso finita (FIR) I filtri FIR possono essere implementati in MATLAB usando la routine conv, o meglio ancora attraverso la funzione filter. La funzione di trasferimento di un filtro FIR, rappresenta un caso particolare di trasformata Z, con i coefficienti a=0 (resta solo a0=1). j j  Filtri con risposta all’impulso finita (FIR) Nel filtro FIR a media mobile definito dall'equazione, ogni campione in uscita è la media degli ultimi N punti in ingresso. Un filtro di questo tipo, che utilizza solo dati attuali e passati, viene definito filtro causale. Questi filtri seguono il principio di causa-effetto, in quanto l'uscita è una funzione degli ingressi passati e presenti. Tutti i sistemi in tempo reale sono causali, poiché non hanno accesso al futuro. Non hanno altra scelta che operare sui valori attuali e passati di un segnale. Se invece si stanno effettuando delle analisi offline (ad esempio i dati sono memorizzati su un pc), è possibile impiegare filtri (definiti non causali) che usano valori anche futuri del segnale. FIR_blink.m MATLAB Filtri con risposta all’impulso finita (FIR) Filtro Causale Filtro Non Causale I filtri causali soffrono di un ritardo Causale K= 100 Non Causale  Filtri con risposta all’impulso finita (FIR) Per progettare un filtro, si parte dalla caratteristica di frequenza desiderata. A questo punto basta fare la trasformata di Fourier inversa. Per un filtro FIR ci basta solo lo spettro di magnitudo, perché i filtri FIR hanno una fase lineare, che è determinata in modo univoco dallo spettro di magnitudo.  Filtri con risposta all’impulso finita (FIR) Sebbene possano esserci circostanze in cui è necessaria una caratteristica di frequenza specifica, di solito si vuole semplicemente separare una gamma di frequenze desiderate da tutto il resto e si vuole che tale separazione sia il più netta possibile: finestra rettangolare Lo spettro di magnitudo di un filtro a finestra rettangolare è un filtro ideale. Viene mostrata anche la porzione di frequenza negativa di questo filtro perché è necessaria per il calcolo della trasformata di Fourier complessa inversa. L'asse delle frequenze è espresso in radianti al secondo ed è normalizzato al valore massimo di 1,0  Filtri con risposta all’impulso finita (FIR) k -ωc ωc L (numero di coefficienti del filtro) deve essere un numero dispari per rendere b[k] simmetrico rispetto a k = 0. Quando k = 0, il denominatore va a zero e il valore effettivo di b può essere ottenuto applicando i limiti e notando che (sin(x)/x)→ 1 quando x → 0.  Filtri con risposta all’impulso finita (FIR) La risposta impulsiva ad una finestra rettangolare (filtro ideale) è quindi una sinc → sin(x)/x. 0.05 0.2 La risposta all'impulso di un filtro a finestra rettangolare per L = 65 coefficienti (k varia tra ±32). Le frequenze di taglio sono indicate rispetto alla frequenza di campionamento, fs , come è comune per le frequenze dei filtri digitali.  Filtri con risposta all’impulso finita (FIR) Poiché stiamo progettando filtri FIR con un numero finito di coefficienti del filtro (cioè L è finito), abbiamo bisogno di una risposta all'impulso che sia anch'essa finita. Sfortunatamente, la soluzione all'Equazione è una risposta all'impulso teoricamente infinita: non va a zero per valori finiti di L. In pratica, dobbiamo semplicemente troncare l'equazione della risposta all'impulso a qualche valore finito L. Il troncamento ha due effetti negativi: o il filtro non ha più un taglio infinitamente netto(roll-out) o si producono oscillazioni nello spettro del filtro (ripple). FIR_rect.m MATLAB Esempio Matlab Trovare lo spettro di magnitudine del filtro a finestra rettangolare per due diverse lunghezze: L=18, L=66. Fc (frequenza di taglio)=300Hz*; Fs (Frequenza campionamento)=1000Hz. Fenomeno di Gibbs Poiché gli artefatti di Gibbs sono dovuti al troncamento di una *For digital filters, the funzione infinita, potremmo ridurli cutoff frequencies se la risposta all'impulso tendesse must lie between 0 and a zero dolcemente anziché 1, where 1 troncata bruscamente. corresponds to the Nyquist rate. In questa implementazione normalizziamo rispetto alla fs, perchè consideriamo tutta la h[n], anche negative.  Filtri con risposta all’impulso finita (FIR) Esistono diverse funzioni finestra, ma le due più popolari e più utili per le risposte all'impulso FIR sono la finestra di Hamming e la finestra di Blackman-Harris. w = hamming(N); w = blackmanharris(N); N rappresenta la lunghezza della finestra, che si fa corrispondere alla lunghezza del filtro (L). Tipicamente  Filtri con risposta all’impulso finita (FIR) Esistono diverse funzioni finestra, ma le due più popolari e più utili per le risposte all'impulso FIR sono la finestra di Hamming e la finestra di Blackman-Harris. Il roll-off del filtro non è ancora quello di un filtro ideale, ma diventa più ripido con l'aumentare della lunghezza del filtro. L’aumento della lunghezza del filtro aumenta il tempo di calcolo (e per i filtri causali il anche il ritardo) necessario per applicare il filtro a un set di dati, un tipico compromesso ingegneristico  Elaborazione Dati e Segnali Biomedici II A.A. 2023-24  Dr. Pietro Aricò, PhD [email protected]  Convoluzione e risposta all'impulso L'ingresso, x(t), e l'uscita, y(t), sono collegati tramite convoluzione attraverso la funzione h(t). Poiché sono tutte funzioni del tempo, la convoluzione è un'operazione nel dominio del tempo. Nel dominio digitale, le tre funzioni temporali diventano vettori discreti campionati: x[n], y[n] e h[n].  Convoluzione e risposta all'impulso La funzione h(t) è nota come risposta all'impulso. Come dice il nome, è la risposta del sistema a un ingresso impulsivo standard. Un ingresso impulsivo (definito anche δ(t)) è un impulso molto breve con un'area unitaria.  Convoluzione e risposta all'impulso Se si conosce la risposta del sistema a un impulso, è possibile determinare la sua risposta a qualsiasi ingresso semplicemente dividendo l'ingresso in una sequenza di impulsi. Ciascuna fetta di tempo genererà la propria piccola risposta all'impulso. L'ampiezza e la posizione di questa risposta all'impulso sono determinate dall'ampiezza e dalla posizione del segmento di segnale di ingresso associato. Se la sovrapposizione e l'invarianza temporale sono valide, l'uscita può essere determinata sommando (o integrando per funzioni continue) le risposte all'impulso di tutti i segmenti del segnale di ingresso.  Convoluzione e risposta all'impulso Ogni segmento di un segnale di ingresso può essere rappresentato da un impulso e ciascuno di essi produce una propria risposta all'impulso. L'ampiezza e la posizione della risposta all'impulso sono determinate dall'ampiezza e dalla posizione del segmento di ingresso associato.  Convoluzione e risposta all'impulso Quando si implementa la convoluzione, è necessario invertire il segnale di ingresso a causa del modo in cui i segnali sono graficizzati. Il lato sinistro del segnale è in realtà il lato a tempo basso; è la prima porzione che entra nel sistema. Pertanto, il segmento più a sinistra del segnale di ingresso produce la prima risposta all'impulso e il segnale di ingresso procede attraverso il sistema all'indietro rispetto al modo in cui è tracciato.  Convoluzione e risposta all'impulso Il processo di convoluzione viene espresso matematicamente come: dove n è la variabile che fa scorrere l'ingresso invertito, x[-k], attraverso la risposta all'impulso, h[k], e K è la lunghezza della funzione più breve, di solito h[k].  Funzione di trasferimento La funzione di trasferimento H(z) è una rappresentazione algebrica di un filtro LTI nel dominio della frequenza. Essa è definita come Y (z)/X(z), dove Y (z) è la trasformata z del segnale di output y(n), mentre X(z) è la trasformata z dell’input x(n). Teorema dello spostamento. Il ritardo di ∆ campioni nel dominio del tempo corrisponde a una moltiplicazione per z−∆ nel dominio della frequenza, dove con l’operatore Z indichiamo la trasformata z.  Funzione di trasferimento La funzione di trasferimento H(z) è una rappresentazione algebrica di un filtro LTI nel dominio della frequenza. Essa è definita come Y (z)/X(z), dove Y (z) è la trasformata z del segnale di output y(n), mentre X(z) è la trasformata z dell’input x(n). Teorema della convoluzione. Per ogni coppia di segnali causali x e y, la convoluzione nel dominio del tempo diventa moltiplicazione nel dominio z. Cioè:  Funzione di trasferimento La funzione di trasferimento H(z) è una rappresentazione algebrica di un filtro LTI nel dominio della frequenza. Essa è definita come Y (z)/X(z), dove Y (z) è la trasformata z del segnale di output y(n), mentre X(z) è la trasformata z dell’input x(n). Ora ricordando che l’output y di un LTI è dato dalla convoluzione dell’input x con la risposta all’impulso h possiamo usare i teoremi appena dimostrati effettuando la trasformata z di y ottenendo con H(z) trasformata di h. Otteniamo cosi la funzione di trasferimento:  Analisi poli-zeri Ogni filtro digitale può essere descritto dai suoi poli e zeri, inoltre essi possono dare un’idea della risposta di un filtro. Come abbiamo visto si parte dalla funzione di trasferimento, se prendiamo un filtro LTI ricorsivo, applicando la trasformata Z all’equazione alle differenze finite avremo: Fattorizzando il primo coefficiente del numeratore b0, e chiamandolo g Fattorizzando numeratore e denominatore: Ora i numeri q1;…; qM sono le radici o gli zeri della funzione. Se z è uguale a uno di questi valori la funzione di trasferimento va a 0. Allo stesso modo i numeri p; …; pN sono detti poli della funzione, quando z assume uno di questi valori l’ampiezza della funzione di trasferimento va all’infinito.  Analisi poli-zeri Tramite l’analisi di poli e zeri è possibile determinare anche la stabilità di un filtro. Infatti possiamo dire che un filtro è stabile se e solo se tutti i suoi poli sono all’interno della circonferenza unitaria nel piano z. Dimostrazione: consideriamo una risposta all’impulso causale nella forma, la cui trasformata z è: (serie di Taylor) Quindi, l’ultima uguaglianza è vera solo per , cioè se. La funzione di trasferimento ha dunque un solo polo a , ed esiste solo per. Consideriamo ora cosa accade se R >1. Il polo di H(z) va all’esterno della circonferenza unitaria, e la risposta all’impulso aumenta esponenzialmente, quindi la definizione di stabilità data in precedenza viene violata. Visto che la trasformata sappiamo che esiste solo se , deduciamo che R ≥ 1 implica che la trasformata z non esiste più sulla circonferenza unitaria, e la risposta in frequenza non è più definita.  Analisi poli-zeri Abbiamo quindi provato che un filtro a un polo è stabile se e solo se il polo è all’interno della circonferenza unitaria. Per quanto riguarda una funzione di trasferimento qualsiasi, la sua espansione in fratti semplici mostra che il comportamento vicino ad ogni polo è quello di un filtro polare avente solo quel polo. Quindi possiamo estendere la definizione di stabilità alla presenza di tutti i poli all’interno della circonferenza.  Funzione di trasferimento La maggior parte delle funzioni di trasferimento saranno composte da polinomi di z sia al numeratore che al denominatore: dove le b[k] sono coefficienti costanti del numeratore e le a[l] sono coefficienti del denominatore. Dalla funzione di trasferimento digitale, H[z], è possibile determinare l'uscita per qualsiasi ingresso, dove Y-1[z] è la trasformata Z inversa.  Risposta in frequenza La risposta in frequenza di un filtro LTI può essere definita come lo spettro del segnale di output diviso lo spettro del segnale di input. Una proprietà della trasformata z è che, se calcolata sulla circonferenza unitaria z=ejωT, si ottiene lo spettro (FFT(x)). Per dimostrarlo basta sostituire z=ejωT nella definizione della trasformata, per ottenere Essendo z complesso può essere espresso da un’ampiezza r = |z| e una fase θ = , è possibile quindi definire la risposta in ampiezza e la risposta in fase del filtro. La risposta in ampiezza a valori reali specifica il gain del filtro per ogni frequenza La risposta in fase Θ a valori reali da lo spostamento di fase, che ogni componente sinusoidale subisce per ogni frequenza.  Risposta in frequenza La risposta in ampiezza a valori reali specifica il gain del filtro per ogni frequenza La risposta in fase Θ a valori reali da lo spostamento di fase, che ogni componente sinusoidale subisce per ogni frequenza. Risposta in ampiezza e in fase di un filtro ellittico  Elaborazione Dati e Segnali Biomedici II A.A. 2023-24  Dr. Pietro Aricò, PhD [email protected]  Filtri Digitali Elaborazione Dati e Segnali Biomedici II A.A. 2023-24  Pietro Aricò, PhD [email protected] 11/10/2023  Progettazione filtro La progettazione di un filtro digitale consiste nel trovare i coefficienti a[l] e b[k] che forniscono lo shape spettrale desiderato. Questo processo di progettazione è agevolato da routine MATLAB che generano i coefficienti a e b che producono la risposta in frequenza desiderata. Nel dominio Z, gli spettri di frequenza della funzione di trasferimento H(z) possono essere ottenuti sostituendo la frequenza reale, ejω, alla frequenza complessa, z. dove FT indica la trasformata di Fourier. Come per tutte le trasformate di Fourier, la frequenza effettiva può essere ottenuta dal numero armonico m dopo averla moltiplicata per fs /N o 1/(NTs) hfilters.m MATLAB Esempio Matlab Sistema con la seguente funzione di trasferimento: fs = 1 kHz, N=512. Calcolare le risposta del filtro ad un impulso.  Esempio Matlab Sistema con la seguente funzione di trasferimento: fs = 1 kHz, N=512. Calcolare le risposta del filtro ad un impulso. Notiamo che a=1; a=-0.2; a=0.8; b=0.2; b=0.5. Applichiamo l’equazione Convertiamo le frequenze in Hz, moltiplicando m * fs/N Esempio Matlab  Filtri FIR e IIR I filtri aventi risposta all’impulso finita (FIR) o infinita (IIR) sono particolari sistemi lineari e tempo-invarianti, e possono essere implementati e simulati su macchine digitali. Offrono numerosi vantaggi sui filtri analogici, come ad esempio la riprogrammabilità del software senza cambiare hardware, oppure la modifica in tempo reale dei coefficienti, rendendoli filtri adattivi. Nelle applicazioni pratiche si preferisce in genere un filtro FIR, se è richiesta una fase lineare, se la distorsione della fase è invece tollerabile si preferisce l’uso dei filtri IIR, che richiedono meno parametri quindi sono meno complessi ed occupano meno spazio in memoria. Poiché i filtri a fase lineare hanno un ritardo di gruppo (group delay, definito come derivata della fase rispetto alla pulsazione) costante, tutte le componenti di frequenza hanno un uguale ritardo nel tempo. Ciò significa che non vi è distorsione dovuta alla risposta in frequenza; in molte applicazioni, questo ritardo di gruppo costante è vantaggioso. Per contro, un filtro con una "fase non lineare", ha un ritardo di gruppo che varia con la frequenza, che porta come risultato ad una distorsione di fase.  Filtri con risposta all’impulso finita (FIR) Definizione. Un sistema LTI è detto FIR se la risposta h(n) all’impulso unitario è finita, nel senso che h(n) = 0 per n < 0 e per n ≥ K per un opportuno K > 0. I filtri FIR sono filtri a media mobile, generalmente con pesi non uniformi. Nei filtri FIR, i pesi del filtro (coefficienti) definiscono la risposta all'impulso, quindi ovviamente la risposta all'impulso sarà finita (stabilità del filtro). I filtri FIR vengono implementati convolvendo i pesi con il segnale di ingresso, il che equivale matematicamente a fare una media ponderata del segnale di ingresso (equazione alle differenze finite per un filtro FIR). b[k] definisce i pesi, e K il numero di pesi, o la lunghezza del filtro, x[n] rappresenta il segnale in ingresso, e y[n] quello in uscita.  Filtri con risposta all’impulso finita (FIR) I filtri FIR possono essere implementati in MATLAB usando la routine conv, o meglio ancora attraverso la funzione filter. La funzione di trasferimento di un filtro FIR, rappresenta un caso particolare di trasformata Z, con i coefficienti a=0 (resta solo a0=1). j j  Filtri con risposta all’impulso finita (FIR) Nel filtro FIR a media mobile definito dall'equazione, ogni campione in uscita è la media degli ultimi N punti in ingresso. Un filtro di questo tipo, che utilizza solo dati attuali e passati, viene definito filtro causale. Questi filtri seguono il principio di causa-effetto, in quanto l'uscita è una funzione degli ingressi passati e presenti. Tutti i sistemi in tempo reale sono causali, poiché non hanno accesso al futuro. Non hanno altra scelta che operare sui valori attuali e passati di un segnale. Se invece si stanno effettuando delle analisi offline (ad esempio i dati sono memorizzati su un pc), è possibile impiegare filtri (definiti non causali) che usano valori anche futuri del segnale. FIR_blink.m MATLAB Filtri con risposta all’impulso finita (FIR) Filtro Causale Filtro Non Causale I filtri causali soffrono di un ritardo Causale K= 100 Non Causale  Filtri con risposta all’impulso finita (FIR) Per progettare un filtro, si parte dalla caratteristica di frequenza desiderata. A questo punto basta fare la trasformata di Fourier inversa. Per un filtro FIR ci basta solo lo spettro di magnitudo, perché i filtri FIR hanno una fase lineare, che è determinata in modo univoco dallo spettro di magnitudo.  Filtri con risposta all’impulso finita (FIR) Sebbene possano esserci circostanze in cui è necessaria una caratteristica di frequenza specifica, di solito si vuole semplicemente separare una gamma di frequenze desiderate da tutto il resto e si vuole che tale separazione sia il più netta possibile: finestra rettangolare Lo spettro di magnitudo di un filtro a finestra rettangolare è un filtro ideale. Viene mostrata anche la porzione di frequenza negativa di questo filtro perché è necessaria per il calcolo della trasformata di Fourier complessa inversa. L'asse delle frequenze è espresso in radianti al secondo ed è normalizzato al valore massimo di 1,0  Filtri con risposta all’impulso finita (FIR) k -ωc ωc L (numero di coefficienti del filtro) deve essere un numero dispari per rendere b[k] simmetrico rispetto a k = 0. Quando k = 0, il denominatore va a zero e il valore effettivo di b può essere ottenuto applicando i limiti e notando che (sin(x)/x)→ 1 quando x → 0.  Filtri con risposta all’impulso finita (FIR) La risposta impulsiva ad una finestra rettangolare (filtro ideale) è quindi una sinc → sin(x)/x. 0.05 0.2 La risposta all'impulso di un filtro a finestra rettangolare per L = 65 coefficienti (k varia tra ±32). Le frequenze di taglio sono indicate rispetto alla frequenza di campionamento, fs , come è comune per le frequenze dei filtri digitali.  Filtri con risposta all’impulso finita (FIR) Poiché stiamo progettando filtri FIR con un numero finito di coefficienti del filtro (cioè L è finito), abbiamo bisogno di una risposta all'impulso che sia anch'essa finita. Sfortunatamente, la soluzione all'Equazione è una risposta all'impulso teoricamente infinita: non va a zero per valori finiti di L. In pratica, dobbiamo semplicemente troncare l'equazione della risposta all'impulso a qualche valore finito L. Il troncamento ha due effetti negativi: o il filtro non ha più un taglio infinitamente netto(roll-out) o si producono oscillazioni nello spettro del filtro (ripple). FIR_rect.m MATLAB Esempio Matlab Trovare lo spettro di magnitudine del filtro a finestra rettangolare per due diverse lunghezze: L=18, L=66. Fc (frequenza di taglio)=300Hz*; Fs (Frequenza campionamento)=1000Hz. Fenomeno di Gibbs Poiché gli artefatti di Gibbs sono dovuti al troncamento di una *For digital filters, the funzione infinita, potremmo ridurli cutoff frequencies se la risposta all'impulso tendesse must lie between 0 and a zero dolcemente anziché 1, where 1 troncata bruscamente. corresponds to the Nyquist rate. In questa implementazione normalizziamo rispetto alla fs, perchè consideriamo tutta la h[n], anche negative.  Filtri con risposta all’impulso finita (FIR) Esistono diverse funzioni finestra, ma le due più popolari e più utili per le risposte all'impulso FIR sono la finestra di Hamming e la finestra di Blackman-Harris. w = hamming(N); w = blackmanharris(N); N rappresenta la lunghezza della finestra, che si fa corrispondere alla lunghezza del filtro (L). Tipicamente  Filtri con risposta all’impulso finita (FIR) Esistono diverse funzioni finestra, ma le due più popolari e più utili per le risposte all'impulso FIR sono la finestra di Hamming e la finestra di Blackman-Harris. Il roll-off del filtro non è ancora quello di un filtro ideale, ma diventa più ripido con l'aumentare della lunghezza del filtro. L’aumento della lunghezza del filtro aumenta il tempo di calcolo (e per i filtri causali il anche il ritardo) necessario per applicare il filtro a un set di dati, un tipico compromesso ingegneristico  Elaborazione Dati e Segnali Biomedici II A.A. 2023-24  Dr. Pietro Aricò, PhD [email protected]  Filtri Digitali Elaborazione Dati e Segnali Biomedici II A.A. 2023-24  Pietro Aricò, PhD [email protected] 16/10/2023 FIR_blinkLowPass.m MATLAB Esempio Matlab Applicare un filtro a finestra rettangolare passa-basso al segnale EEG contenente blink oculari di 30 secondi, che si trova come variabile dataBlink all’interno del file DataBlink.mat. Il segnale viene campionato a 256 Hz. Utilizzare una frequenza di taglio di 10 Hz e una lunghezza di filtro di 65. Utilizzare la finestra di BlackmanHarris per troncare la risposta all'impulso del filtro. Tracciare il segnale originale e quello filtrato.  Esempio Matlab Applicare un filtro a finestra rettangolare passa-basso al segnale EEG contenente blink oculari di 30 secondi, che si trova come variabile dataBlink all’interno del file DataBlink.mat. Il segnale viene campionato a 256 Hz. Utilizzare una frequenza di taglio di 10 Hz e una lunghezza di filtro di 65. Utilizzare la finestra di BlackmanHarris per troncare la risposta all'impulso del filtro. Tracciare il segnale originale e quello filtrato. FIR_blinkBandPass.m MATLAB Esempio Matlab Applicare un filtro a finestra rettangolare passa-banda al segnale EEG contenente blink oculari di 30 secondi dove c’è una componente additiva continua accoppiata al segnale (200microV). Il segnale viene campionato a 256 Hz. Utilizzare una frequenza di taglio di 3 - 10 Hz e una lunghezza di filtro di 129. Utilizzare la finestra di BlackmanHarris per troncare la risposta all'impulso del filtro. Tracciare il segnale originale e quello filtrato.  Esempio Matlab Applicare un filtro a finestra rettangolare passa-banda al segnale EEG contenente blink oculari di 30 secondi dove c’è una componente additiva continua accoppiata al segnale (200microV). Il segnale viene campionato a 256 Hz. Utilizzare una frequenza di taglio di 3 - 10 Hz e una lunghezza di filtro di 129. Utilizzare la finestra di BlackmanHarris per troncare la risposta all'impulso del filtro. Tracciare il segnale originale e quello filtrato. FIR_noiseSin.m MATLAB Esempio Matlab Sinusoide a 50Hz, accoppiata ad un rumore uniforme con gain 100. Fs = 1000Hz; Applicare un filtro FIR con L=129; fs = 1000; % Sample frequency in Hz NoiseVect=100*(rand(1,1000)-0.5); sinvect=sin(2*pi*50*[1/fs:1/fs:1]); NoiseVect=NoiseVect+sinvect;  Esempio Matlab Sinusoide a 50Hz, accoppiata ad un rumore bianco con gain 100. Fs = 1000Hz; Applicare un filtro FIR con L=129; fs = 1000; % Sample frequency in Hz NoiseVect=100*(rand(1,1000)-0.5); sinvect=sin(2*pi*50*[1/fs:1/fs:1]); NoiseVect=NoiseVect+sinvect;  Filtri con risposta all’impulso infinita (IIR) Per aumentare la pendenza di attenuazione di un filtro, si possono applicare due o più filtri in sequenza. Nel dominio della frequenza, gli spettri di frequenza si moltiplicano. Questo porta a una pendenza di attenuazione più ripida nello spettro. Tuttavia, a parità di sforzo computazionale, è possibile utilizzare un singolo filtro di ordine doppio e ottenere una pendenza di attenuazione ancora più ripida 12° 24th Double 12th 12°+12°  Filtri con risposta all’impulso infinita (IIR) Un'alternativa interessante all'utilizzo di due filtri FIR in serie è quella di filtrare due volte gli stessi dati, utilizzando prima un filtro FIR a media mobile standard e poi reintroducendo l'uscita del filtro FIR in un altro filtro che contribuisce anch'esso all'uscita. In questo modo si ottiene un filtro ricorsivo a media mobile. Il filtro inferiore opera sull'uscita, ma solo sui valori passati dell'uscita che sono già stati determinati dal filtro superiore  Filtri con risposta all’impulso infinita (IIR) L'equazione di un filtro IIR può essere derivata da tale co