Filtri Ideali e Reali

Introduzione

La teoria dei Segnali e un Circuito in regime sinusoidale si incontrano al bar, ma non si salutano.

Tra i sistemi fisici ed in particolare sotto i cosiddetti sistemi LTI (lineari tempo-invariati), è opportuno porre l’accento sui filtri ideali. Un filtro è un sistema LTI selettivo nel dominio delle frequenze.

I filtri ideali sono non causali e sono il riferimento a cui approssimare filtri reali, costituiti da circuiti elettrici (resistori, induttori e condensatori). Dunque, grazie alla teoria dei circuiti, è possibile realizzare e studiare fisicamente i filtri.

Le quattro categorie principali di filtri sono:

  1. LP: passa-basso (low-pass)
  2. HP: passa-alto (high-pass)
  3. BP: passa-banda (band-pass)
  4. BR: elimina-banda (band-rejection)

La banda è l’intervallo di frequenze in cui il segnale assume valori maggiori di zero. Si considera il solo semiasse positivo.

Decibel

Le scale logaritmiche aiutano nella visualizzazione di grandezze fisiche molto estese su uno (o più) assi. Dalla teoria del circuiti è noto che la potenza elettrica è data dal prodotto tra tensione e corrente: p=vip=vi.

Date due potenze P1P_1 e P2P_2, il rapporto di potenze R è:

R=P2P1=(V2V1)2=(I2I1)2R = \frac{P_2}{P_1} = \bigg(\frac{V_2}{V_1}\bigg)^2 = \bigg(\frac{I_2}{I_1}\bigg)^2

Il rapporto di potenze espresso in decibel è:

RdB=10log⁡10(P2P1)=10log⁡10(V2V1)2=20log⁡10(V2V1)R_{dB} = 10 \log_{10}\bigg(\frac{P_2}{P_1}\bigg) = 10 \log_{10}\bigg(\frac{V_2}{V_1}\bigg)^2 = 20 \log_{10}\bigg(\frac{V_2}{V_1}\bigg)

RdB rappresenta il valore di R espresso in decibel.

Il decibel (dB) è la decima parte del bel (simbolo B):

10[dB]=1[B]10 [dB] = 1 [B]

Il decibel viene adottato per evidenziare rapporti di potenza.

Se si vuole ottenere un valore assoluto, si sostituisce P1 con il riferimento di 1 milli-watt e si ottiene il dBm:

dBm=10log⁡10(P21  [mW])dB_{m} = 10 \log_{10} \bigg(\frac{P_2}{1\;[mW]}\bigg)

Il dBm può anche essere indicato come dBmW per esplicitare ancora meglio che si tratta di un rapporto di potenze.

A livello teorico, la definizione del dB per le tensioni (calcolate sulla stessa resistenza) è:

dB=20log⁡10(V2V1)dB = 20 \log_{10} \bigg(\frac{V_2}{V_1}\bigg)

Si ottiene, come nel caso precedente, una grandezza relativa. Se si vuole ottenere una grandezza in valore assoluto, allora si deve porre a denominatore 1 milli-volt:

dBmV=20log⁡10(V21  [mV])dB_{mV} = 20 \log_{10} \bigg(\frac{V_2}{1\;[mV]}\bigg)

Data la funzione X(f)X(f) in ingresso al filtro e Y(f)Y(f) in uscita dal filtro. H(f)H(f) è la risposta in frequenza del filtro ed è il rapporto tra grandezze fisiche che rappresentano delle ampiezze. La risposta in ampiezza di un filtro, ovvero il modulo di H(f)H(f), può essere scritta come:

A(f)=∣H(f)∣=∣Y(f)∣∣X(f)∣  ⟹  A(f)dB=20log⁡10A(f)A(f) = |H(f)| = \frac{|Y(f)|}{|X(f)|} \implies A(f)_{dB} = 20 \log_{10} A(f)
R=V22/V12R=V_2^2 / V_1^2V2/V1V_2/V_1RdBR_{dB}
101010\sqrt{10}1010
10010010102020
10001000101010\sqrt{10}3030
111100
0.10.10.1\sqrt{0.1}−10-10
0.010.010.10.1−20-20
222\sqrt{2}33
1/21/21/21/\sqrt{2}−3-3

Dalla tabella soprastante si nota che la banda a −3dB-3 dB del segnale H(f)H(f) si ottiene A(f)A(f).

Grafici di Filtri Fondamentali

In questo articolo sono presenti simulazione effettuate con Simulink e codice scritto in MATLAB.

%%%
% sample usage:
% frequency_response(10, 2, 1, 4)
%
% parameters:
% #1: horizontal plots limit (both positive and negative) 
% #2: band value
% #3: lower bound for band-pass and band-rejection plots
% #4: upper bound for band-pass and band-rejection plots
%%%

function frequency_response(f_lim, B, fL, fH)

    % fprintf('»»» number of arguments: %d «««\n', nargin);

    if nargin ~= 4
        fprintf(['»»» USAGE: ideal_filters(f_lim, B, fL, fH)\n' ...
            '\tf_lim:\thorizontal plots limit (both positive and negative)\n' ...
            '\tB:\tband value\n' ...
            '\tfL:\tlower bound for band-pass and band-rejection plots\n' ...
            '\tfH:\tupper bound for band-pass and band-rejection plots\n' ...
            '«««\n']);
        return;
    end

    % take the absolute value of the first parameter and set as horizontal axes
    f_lim = abs(f_lim);
    f = linspace(-f_lim, f_lim, 1001);

    % define the functions
    H_lp = low_pass(f, B); % rect(f ./(2 * B));
    H_hp = high_pass(f, B);
    H_bp = band_pass(f, fL, fH);
    H_br = band_rejection(f, fL, fH);

    % ---

    figure('Name', 'Ideal Filters');
    plots = [H_lp; H_hp; H_bp; H_br];
    titles = ["Low Pass" "High Pass" "Band Pass" "Band Rejection"];

    % the use of i1 as loop indexes could cause some issuces since
    % it is mainly used to represent the immaginary unit
    for a = 1:4
        % fprintf('»»» i1 = %d', i1)
        subplot(2,2,a)
        plot(f, plots(a,:),'Color','blue','LineWidth',1.5)
        hold on
        plot(f, zeros(size(f)),zeros(size(f)),f,'Color','#808080','LineStyle','--','LineWidth',1.5)
        title(titles(a))
        hold off
        grid on
        ylim([-0.1, 1.1])
    end

end

% low_pass_plot.DataTipTemplate.DataTipRows(1).Label = "-B"
% datatip(lp_plot,-B,0);
% datatip(lp_plot,B,0);

function y = low_pass(x, B)
    y = rect(x ./ (2 * B));
    % den = 2 * B; % denominator
    % y = zeros(1, length(x));
    % y(abs(x./den) < 0.5) = 1;
    % y(abs(x./den) == 0.5) = 0.5;
end

function y = high_pass(f, B)
    y = 1 - low_pass(f, B);
    % den = 2 * B; % denominator
    % y = zeros(1, length(x));
    % y(abs(x./den) > 0.5) = 1;
    % y(abs(x./den) == 0.5) = 0.5;
end

function y = band_pass(f, fL, fH)
    B = fH - fL;
    f0 = (fH + fL) / 2;

    rect1 = rect((f-f0) ./ B);
    rect2 = rect((f+f0) ./ B);
    y = rect1 + rect2;
end

function y = band_rejection(f, fL, fH)
    y = 1 - band_pass(f, fL, fH);
end

function y = rect(x)
    y = zeros(1, length(x));
    y(abs(x) < 0.5) = 1;
    y(abs(x) == 0.5) = 0.5;
end

Le risposte in frequenza dei quattro filtri ideali più noti sono:

risposte in frequenza

%%%
% sample usage:
% impulse_response(10, 2, 1, 4)
%
% parameters:
% #1: horizontal plots limit (both positive and negative) 
% #2: band value
% #3: lower bound for band-pass and band-rejection plots
% #4: upper bound for band-pass and band-rejection plots
%%%

function impulse_response(t_lim, B, fL, fH)
    if nargin ~= 4
        fprintf(['»»» USAGE: impulse_response(t_lim, B, fL, fH)\n' ...
            '\tt_lim:\thorizontal plots limit (both positive and negative)\n' ...
            '\tB:\tband value\n' ...
            '\tfL:\tlower bound for band-pass and band-rejection plots\n' ...
            '\tfH:\tupper bound for band-pass and band-rejection plots\n' ...
            '«««\n']);
        return;
    end
    
    % take the absolute value of the first parameter and set as horizontal axes
    t_lim = abs(t_lim);
    t = linspace(-t_lim, t_lim, 10001);
    B2 = (2 * B)+0.2;

    % define the functions
    h_lp = impulsive_low_pass(t, B);
    h_hp = impulsive_high_pass(t, B);
    h_bp = impulsive_band_pass(t, fL, fH);
    h_br = impulsive_band_rejection(t, fL, fH);

    % ---

    figure('Name', 'Impulsive Responses');
    plots = [h_lp; h_hp; h_bp; h_br];
    titles = ["Low Pass" "High Pass" "Band Pass" "Band Rejection"];

    % create separate plots for each impulse response
    % the outer loop determinates which figure to plot (all togheter or not)
    % the inner loop draws the subplots
    % the use of i1 and j1 as loop indexes could cause some issuces since
    % they are mainly used to represent the immaginary unit
    for a = 1:2
        % fprintf('»»» \n i1 = %d', i1)

        for b = 1:4
            % fprintf('»»» \n j1 = %d', j1)

            if a == 1
                subplot(2,2,b)
            elseif a == 2
                figure('Name', titles(b))
            end
            
            plot(t, plots(b,:), 'Color', 'blue','LineWidth',1.5)
            hold on
            plot(t, zeros(size(t)), zeros(size(t)), t, 'Color', '#808080', 'LineStyle','--','LineWidth',1.5)
            title(titles(b))
            hold off
            grid on
            ylim([-B2, B2])
        end
    end
end

function y = impulsive_low_pass(x, B)
    y = 2 * B * sinc(2*B .* x);

    % true_h_lp = zeros(1,length(t));
    % two_b = B * 2;
    % true_h_lp(t==0)=1;
    % true_h_lp(t~=0)= two_b * sin(two_b * pi*t(t~=0))./(two_b * pi*t(t~=0));
end

function y = impulsive_high_pass(x, B)
    y = dirac(x) - 2 * B * sinc(2*B .* x);
end

function y = impulsive_band_pass(x, fL, fH)
    B = fH - fL;
    f0 = (fH + fL) / 2;

    y = 2 * B * sinc(B .* x) .* cos(2 * pi * f0 .* x);
end

function y = impulsive_band_rejection(x, fL, fH)
    y = dirac(x) - impulsive_band_pass(x, fL, fH);
end

function y = sinc(x)
    y = zeros(1, length(x));
    y(x==0) = 1;
    y(x~=0) = sin(pi*x(x~=0))./(pi*x(x~=0));
end

Le risposte impulsive (nel dominio temporale) dei quattro filtri ideali più noti sono:

risposte impulsive

Passa-Basso Ideale

Simbologia: due sinusoidi con quella più in alto barrata. Le frequenze più alte sono tagliate. La sigla LP in pedice significa low pass.

Il parametro B è la banda del filtro passa-basso.

HLP(f)=rect(f2B);hLP(t)=2B⋅sinc(2Bt)H_{\rm LP}(f) = \text{rect}\bigg(\frac{f}{2B}\bigg) \quad ; \quad h_{\rm LP}(t) = 2B \cdot \text{sinc}(2Bt)

La banda passante è compresa tra 00 e B. La banda oscura si estende da B a ∞\infty.

LP risposta in frequenza

LP risposta impulsiva

Passa-Basso Reale RC

Un filtro passa-basso reale è un circuito formato da un resistore con valore di resistenza R in serie ad un condensatore con capacità elettrica C. Lo schema elettrico potrebbe trarre in inganno facendo sembrare i componenti in parallelo. Il segnale in ingresso x(t)x(t) è la tensione ai capi dei due componenti in serie, mentre il segnale di uscita y(t)y(t) è la tensione ai capi del solo condensatore. La corrente i(t)i(t) scorre nel circuito (idealmente dall’ingresso all’uscita).

Passa-Basso Reale RC

Il 1° metodo di risoluzione consiste nel considerare il sistema LTI:

y(t)=T[x(t)]y(t) = \mathcal{T} \Big[ x(t)\Big]

Al circuito è associato il seguente sistema di equazioni:

{i(t)=Cddty(t)x(t)=y(t)+Ri(t)\sis{ i(t)=C \frac{d}{dt}y(t)\\ \\ x(t)=y(t)+R i(t) }

La risposta impulsiva del sistema passa-basso reale è la risposta all’impulso Delta di Dirac. Inoltre, si ricorda che la Delta di Dirac è la derivata della funzione gradino.

→u(t)ddt→δ(t)T→h(t)\Large \xrightarrow{u(t)} \boxed{\frac{d}{dt}} \xrightarrow{\delta(t)} \boxed{\mathcal{T}} \xrightarrow{h(t)}

Un sistema LTI a cascata è invertibile, dunque:

→u(t)T→z(t)ddt→h(t)\Large \xrightarrow{u(t)} \boxed{\mathcal{T}} \xrightarrow{z(t)} \boxed{\frac{d}{dt}} \xrightarrow{h(t)}

con z(t)≠δ(t)z(t)\neq \delta (t).

Allora z(t)=T[u(t)]z(t)=\mathcal{T}[u(t)] rappresenta il transitorio di carica del condensatore:

u(t)={1t≥00t<0  ⟹  z(t)={1−et/RCt≥00t<0u(t)=\sis{1 \quad t \geq 0 \\ 0\quad t < 0} \implies z(t)=\sis{ 1-e^{t/RC}\quad &t \geq 0 \\ 0 &t < 0 }

La risposta impulsiva è la derivata del transitorio di carica:

h(t)=ddtz(t)={1RC⋅e−t/RCt≥00t<0\Large h(t)=\frac{d}{dt}z(t) = \sis{ \frac{1}{RC}\cdot e^{-t/RC}\quad &t \geq 0 \\ 0 &t < 0 }

Si nota la notevole somiglianza con la definizione di esponenziale monolatero:

monoExp(t)=e−t/T⋅u(t)={e−t/Tt>00t≤0\text{monoExp}(t) = e^{-t/T} \cdot u(t)= \sis{e^{-t/T} \quad & t > 0\\ 0 & t \leq 0 }

La risposta in frequenza è la trasformata di Fourier della risposta impulsiva:

H(f)=∫−∞+∞x(t)e−j2πftdt==∫0+∞1RCe−t/(RC)e−j2πftdt==1RC∫0+∞e−t(1/(RC)+j2πf)dt==1RClim⁡c→+∞∫0ce−t(1/(RC)+j2πf)dt==1RCRC1+j2πfRC==11+j2πfRC\eq{ H(f) &= \int_{-\infty}^{+\infty}x(t)e^{-j2\pi ft}dt =\\ &= \int_{0}^{+\infty} \frac{1}{RC} e^{-t/(RC)} e^{-j2\pi ft}dt =\\ &= \frac{1}{RC}\int_{0}^{+\infty} e^{-t(1/(RC)+j2\pi f)}dt =\\ &= \frac{1}{RC}\lim_{c \to +\infty}\int_{0}^c e^{-t(1/(RC)+j2\pi f)}dt=\\ &= \frac{1}{RC}\frac{RC}{1+j2\pi fRC} =\\ &= \frac{1}{1+j2\pi fRC} }

La risposta in ampiezza è il modulo della risposta in frequenza:

A(f)=11+(2πfRC)2A(f) = \frac{1}{\sqrt{1+(2\pi f RC)^2}}

La risposta in fase è l’argomento della risposta in frequenza:

Θ(f)=arg⁡(H(f))=−arctan(2πfRC)\Theta(f) = \arg \bigg(H(f)\bigg) = -\text{arctan}(2\pi f RC)

Il 2° metodo di risoluzione consiste nell’utilizzo del metodo dei fasori, il quale prevede il passaggio ai numeri complessi per rappresentare circuiti elettrici che operano in regime sinusoidale. L’analisi circuitale in regime sinusoidale prevede che ad ogni elemento sia associato un fasore, un valore di impedenza e la corrispondente ammettenza.

L’articolo Analisi in Regime Sinusoidale tratta in modo molto approfondito questo argomento. Si può usare la tabella che segue per ottenere l’impedenza relativa ad un componente:

ElementoTempoFasoreImpedenzaAmmettenza
R: Resistorev=Riv=RiV=RIV=RIZ=RZ = RY=1/RY=1/R
L: Induttorev=Ldidtv = L \frac{di}{dt}V=jωLIV = j \omega LIZ=jωLZ = j \omega LY=1/(jωL)Y=1/(j \omega L)
C: Condensatorei=Cdvdti = C \frac{dv}{dt}V=I/(jωC)V = I/(j \omega C)Z=1/(jωC)Z = 1/(j \omega C)Y=jωCY = j \omega C

Da notare che la pulsazione è ω=2πf\omega = 2\pi f. Nelle formule che seguono, si esplicita la ff in quanto variabile indipendente dei grafici H(f)H(f), A(f)A(f) e Θ(f)\Theta(f).

Si applica un partitore tra i valori di impedenza della resistenza e del condensatore. Dato il fasore VX della funzione in ingresso x(t)x(t). Il fasore VY della tensione in uscita y(t)y(t) è dato da:

VY=VX⋅ZCZR+ZCV_Y = V_X \cdot \frac{Z_C}{Z_R + Z_C}

La risposta in frequenza è il rapporto tra il fasore della tensione in uscita e quello della tensione in ingresso:

H(f)=VYVX=ZCZR+ZC=11+j2πfRCH(f) = \frac{V_Y}{V_X} = \frac{Z_C}{Z_R + Z_C} = \frac{1}{1+j2\pi fRC}

La risposta in ampiezza è il modulo della risposta in frequenza:

A(f)=11+(2πfRC)2A(f) = \frac{1}{\sqrt{1+(2\pi f RC)^2}}

La risposta in fase è l’argomento della risposta in frequenza:

Θ(f)=arg⁡(H(f))=−arctan(2πfRC)\Theta(f) = \arg \bigg(H(f)\bigg) = -\text{arctan}(2\pi f RC)
%%%
% sample usage:
% real_low_pass_RC(1, 0.1)
%
% parameters:
% #1: resistor value in Ohm
% #2: capacitor value in Farad
%%%

function real_low_pass_RC(R, C)
	f = linspace(-40, 40, 10001);
    f_size = zeros(size(f));
    w = 2*pi.*f;
    Z_R = R;
    Z_C = 1./(1j*w*C);
    Hf = Z_C ./ (Z_R + Z_C);
    
    % amplitude response: A(f) = |H(f)|
    Af = 1 ./ sqrt(1 + (2 * pi * R * C .* f) .^ 2);

    % the filter is a low pass: f0 = 0;
    Af0 = 1 / sqrt(2);
    
    % B_(-3) --> band value at -3 dB
    B3 = 1 / (2 * pi * R * C);

    figure(Name="Real Low Pass");
    tiledlayout(1, 2)
    nexttile 
    plot(f, Af, 'LineWidth', 2.5, 'Color', 'blue')
    hold on
    plot(f, f_size, f_size, f, 'Color', '#808080', 'LineStyle','--')
    hold on
    yline(Af0, 'Color', 'red', 'LineStyle','-.')
    hold on
    xline(B3, 'Color', 'red', 'LineStyle','-.')
    hold off
    grid on
    ylim([-0.1, 1.1])
    xlim([-20 20])
    ax = gca;
    chart = ax.Children(5);
    datatip(chart,B3,Af0);
    chart = ax.Children(4);
    datatip(chart,B3,0);
    chart = ax.Children(3);
    datatip(chart,0,Af0,"Location","southwest");
    title('Real LP - Amplitude Response', FontSize=20)

    % ---

    % risposta in fase
    nexttile
    Fasef = angle(Hf);
    plot(f, Fasef, 'LineWidth', 2.5, 'Color', 'red')
    hold on
    plot(f, f_size, f_size, f, 'Color', '#808080', 'LineStyle','--')
    hold off
    grid on
    ylim([-2 2])
    xlim([-20 20])
    title('Real LP - Phase Response', FontSize=20)
end

Real LP

Passa-Basso Reale LR

Un filtro passa-basso reale, a differenza di quanto visto sopra, può anche essere formato ad un induttore con valore di induttanza L in serie ad un resistore con valore di resistenza R. Il segnale in ingresso x(t)x(t) è la tensione ai capi dei due componenti in serie, mentre il segnale di uscita y(t)y(t) è la tensione ai capi del solo resistore.

Passa-Basso Reale LR

Utilizzando il metodo dei fasori, si ottiene il fasore della tensione in uscita:

VY=VX⋅ZRZR+ZLV_Y = V_X \cdot \frac{Z_R}{Z_R + Z_L}

Nel dominio delle frequenze si può scrivere:

Y(f)=X(f)⋅ZRZR+ZL=X(f)⋅RR+j2πfLY(f) = X(f) \cdot \frac{Z_R}{Z_R + Z_L} = X(f) \cdot \frac{R}{R + j 2 \pi fL}

La risposta in frequenza è il rapporto tra il fasore della tensione in uscita e quello della tensione in ingresso:

H(f)=VYVX=Y(f)X(f)=ZRZR+ZL==RR+j2πfL=11+(j2πfL)/R==11+j2πfLR\eq{ H(f)&=\frac{V_Y}{V_X}=\frac{Y(f)}{X(f)}=\frac{Z_R}{Z_R + Z_L} =\\ &=\frac{R}{R+j2\pi fL}=\frac{1}{1+(j2\pi fL)/R}=\\ &=\frac{1}{1+j2\pi f \frac{L}{R}} }

La risposta in ampiezza è il modulo della risposta in frequenza:

A(f)=∣H(f)∣=R2R2+(2πfL)2=11+(2πfLR)2A(f) = |H(f)| = \frac{R^2}{\sqrt{R^2 + (2\pi fL)^2}} = \frac{1}{\sqrt{1+\Big(2\pi f \frac{L}{R}\Big)^2}}

Studiando la risposta in ampiezza al limite, si ricava che:

∣H(f)∣→f→01;∣H(f)∣→f→∞0|H(f)| \xrightarrow{f \to 0} 1 \quad ; \quad |H(f)| \xrightarrow{f \to \infty} 0\\

La risposta in ampiezza trova il proprio massimo per f0=0f_0 = 0, quindi si tratta di un filtro passa-basso.

f−3f_{-3} si ottiene grazie alla relazione che segue:

∣H(f−3)∣∣H(f0)∣=12\frac{|H(f_{-3})|}{|H(f_0)|} = \frac{1}{\sqrt{2}}
∣H(f−3)∣=12∣H(f0)∣  ⟹  11+(2πf−3LR)2=12  ⟹  1+(2πf−3LR)2=2  ⟹  f−32(2πLR)2=1  ⟹  f−32=R24π2L2  ⟹  f−3=±R2πL\eq{ |H(f_{-3})| &= \frac{1}{\sqrt{2}} |H(f_0)| \\ \implies & \frac{1}{\sqrt{1+\Big(2\pi f_{-3} \frac{L}{R}\Big)^2}}=\frac{1}{\sqrt{2}}\\ \implies & 1+\Big(2\pi f_{-3} \frac{L}{R}\Big)^2 = 2\\ \implies & f_{-3}^2\Big(2\pi\frac{L}{R}\Big)^2 = 1 \\ \implies & f_{-3}^2 = \frac{R^2}{4\pi^2L^2} \\ \implies & f_{-3}= \pm \frac{R}{2\pi L} \\ }

La banda a −3dB-3dB del filtro è il valore positivo (tra i due possibili) di f−3f_{-3}:

B−3dB≡f−3=R2πLB_{-3 dB} \equiv f_{-3}= \frac{R}{2\pi L}

La risposta in fase è l’argomento della risposta in frequenza:

Θ(f)=arg⁡(H(f))=−arctan(2πfLR)\Theta(f) = \arg \bigg(H(f)\bigg) = -\text{arctan}\bigg(\frac{2\pi f L}{R}\bigg)
%%%
% sample usage:
% real_low_pass_LR(1, 12.57)
%
% parameters:
% #1: inductor value in Henry
% #2: resistor value in Ohm
%%%

function real_low_pass_LR(L, R)
    f = linspace(-40, 40, 10001);
    f_size = zeros(size(f));
    % pulsazione (omega)
    w = 2*pi.*f;
    Z_R = R;
    Z_L = 1j*w*L;
    Hf = Z_R ./ (Z_R + Z_L);
    % amplitude response: A(f) = |H(f)|
    Af = 1 ./ sqrt(1 + (2 * pi * (L/R) .* f) .^ 2);
    Af0 = 1 / sqrt(2);
    
    % B_(-3) --> band value at -3 dB
    B3 = R / (2 * pi * L);
    fprintf("»»» valore di banda B = %f\n", B3);

    figure(Name="Real Low Pass - LR");
    tiledlayout(1, 2)
    nexttile 
    plot(f, Af, 'LineWidth', 2.5, 'Color', 'blue')
    % ...
    % risposta in fase
    nexttile
    Fasef = angle(Hf);
    plot(f, Fasef, 'LineWidth', 2.5, 'Color', 'red')
    % ...
end

Real Low Pass - LR

Calcolo dell’uscita dato l’ingresso

→x(t)Low  Pass→y(t)\xrightarrow{x(t)} \boxed{\rm Low \; Pass} \xrightarrow{y(t)}

Dato un ingresso nella forma:

x(t)=αcos⁡(2πf0t)x(t) = \alpha \cos(2\pi f_0 t)

N.B. La f0f_0 a cui si fa riferimento in questo paragrafo non è quella del paragrafo precedente.

Allora l’uscita è data da:

y(t)=α∣H(f0)∣cos⁡(2πf0t+arg⁡(H(f0)))y(t) = \alpha \Big|H(f_0)\Big| \cos\Big(2\pi f_0 t+\arg\big(H(f_0)\big)\Big)

Esempio:

Sia x(t)=8cos⁡(2πf0t)x(t) = 8 \cos(2\pi f_0 t) con 2πf0=103s−12\pi f_0=10^3 s^{-1} e L/R=3⋅10−3sL/R = \sqrt{3}\cdot 10^{-3}s.

Allora α=8\alpha = 8:

2πf0=103  ⟹  f0=1032πLR=31032\pi f_0=10^3 \implies f_0 = \frac{10^3}{2\pi}\\ \frac{L}{R} = \frac{\sqrt{3}}{10^3}

La risposta in ampiezza è:

∣H(f0)∣=∣H(1032π)∣==11+(2π1032π3103)2==11+(3)2=14=12\eq{ \Big|H(f_0)\Big| &= \bigg|H\bigg(\frac{10^3}{2\pi}\bigg)\bigg|=\\ &= \frac{1}{\sqrt{1+\Big(2\pi \frac{10^3}{2\pi} \frac{\sqrt{3}}{10^3}\Big)^2}}=\\ &= \frac{1}{\sqrt{1+(\sqrt{3})^2}} = \frac{1}{\sqrt{4}}=\frac{1}{2} }

La risposta in fase è:

Θ(f0)=−arctan(2πf0LR)==−arctan(2π1032π3103)==−arctan(3)==−π3\eq{ \Theta(f_0) &= -\text{arctan}\bigg(\frac{2\pi f_0 L}{R}\bigg) =\\ &= -\text{arctan}\bigg(2\pi \frac{10^3}{2\pi} \frac{\sqrt{3}}{10^3}\bigg) =\\ &= -\text{arctan}\Big(\sqrt{3}\Big) =\\ &= -\frac{\pi}{3}\\ }

Allora l’uscita è data da:

y(t)=α∣H(f0)∣cos⁡(2πf0t+arg⁡(H(f0)))==812cos⁡(2π1032πt−π3)==4cos⁡(103t−π3)\eq{ y(t) &= \alpha \Big|H(f_0)\Big| \cos\Big(2\pi f_0 t+\arg\big(H(f_0)\big)\Big)=\\ &= 8 \frac{1}{2}\cos\bigg(2\pi \frac{10^3}{2\pi} t-\frac{\pi}{3}\bigg)=\\ &= 4\cos\bigg(10^3 t-\frac{\pi}{3}\bigg) }

Passa-Alto Ideale

Simbologia: due sinusoidi con quella più in basso barrata. Le frequenze più basse sono tagliate. La sigla HP in pedice significa high pass.

HHP(f)=1−rect(f2B);hHP(t)=δ(t)−2B⋅sinc(2Bt)H_{\rm HP}(f) = 1-rect\bigg(\frac{f}{2B}\bigg) \quad ; \quad h_{\rm HP}(t) = \delta(t) -2B \cdot sinc(2Bt)

La banda passante si estende all’infinito. La banda oscura è compresa tra 00 e B.

HP risposta in frequenza

HP risposta impulsiva

Passa-Alto Reale CR

Un filtro passa-alto reale è un circuito formato da un condensatore con capacità elettrica C in serie ad un resistore con valore di resistenza R. Lo schema elettrico potrebbe trarre in inganno facendo sembrare i componenti in parallelo. Il segnale in ingresso x(t)x(t) è la tensione ai capi dei due componenti in serie, mentre il segnale di uscita y(t)y(t) è la tensione ai capi del solo resistore. La corrente i(t)i(t) scorre nel circuito (idealmente dall’ingresso all’uscita).

Passa-alto Reale

Utilizzando il metodo dei fasori come spiegato per il passa basso, si ottiene il fasore della tensione in uscita. In questo caso, nella formula del partitore, l’impedenza del resistore e l’impedenza del condensatore sono invertite.

VY=VX⋅ZRZR+ZCV_Y = V_X \cdot \frac{Z_R}{Z_R + Z_C}

La risposta in frequenza è il rapporto tra il fasore della tensione in uscita e quello della tensione in ingresso:

H(f)=VYVX=ZRZR+ZC=11−j2πfRCH(f) = \frac{V_Y}{V_X} = \frac{Z_R}{Z_R + Z_C} = \frac{1}{1-j2\pi fRC}

La risposta in ampiezza è il modulo della risposta in frequenza:

A(f)=∣H(f)∣=11+1(2πfRC)2\Large A(f) = |H(f)| = \frac{1}{\sqrt{1+\frac{1}{(2\pi f RC)^2}}}

Studiando la risposta in ampiezza per f→0f \to 0 e poi per f→∞f \to \infty si ottengono i valori:

f^0=0  ⟹  ∣X(f^0)∣=0f~0=∞  ⟹  ∣X(f~0)∣=1\eq{ \widehat{f}_0 &=0\implies |X(\widehat{f}_0)| = 0\\ \widetilde{f}_0 &=\infty\implies |X(\widetilde{f}_0)| = 1\\ }

Allora si ottiene:

max⁡∣X(f)∣=∣X(f0)∣=1  ⟹  f0=f~0=∞\max|X(f)| = |X(f_0)| = 1\implies f_0 = \widetilde{f}_0 = \infty

La risposta in fase è l’argomento della risposta in frequenza:

Θ(f)=arg⁡(H(f))=arctan(12πfRC)\Theta(f) = \arg \bigg(H(f)\bigg) = \text{arctan}\bigg(\frac{1}{2\pi fRC}\bigg)
%%%
% sample usage:
% real_high_pass(1, 0.1)
%
% parameters:
% #1: resistor value in Ohm
% #2: capacitor value in Farad
%%%

function real_high_pass(R, C)
    if nargin ~= 2
        fprintf(['»»» USAGE: real_high_pass(R, C)\n' ...
            '\tR:\tresistor value [Ohm]\n' ...
            '\tC:\tcapacitor value [Farad]\n' ...
            '«««\n']);
        return;
    end
    
    f = linspace(-40, 40, 10001);
    f_size = zeros(size(f));
    w = 2*pi.*f; % pulsazione (omega)
    Z_R = R;
    Z_C = 1./(1j*w*C);
    Hf = Z_R ./ (Z_R + Z_C);

    % amplitude response: A(f) = |H(f)|
    Af = 1 ./ sqrt(1 +(1./(2 * pi * R * C .* f) .^ 2));

    Af0 = 1 / sqrt(2);
    
    % B_(-3) --> band value at -3 dB
    B3 = 1 / (2 * pi * R * C);
    fprintf("»»» valore di banda B = %f\n", B3);
    figure(Name="Real High Pass");
    tiledlayout(1, 2)
    nexttile 
    plot(f, Af, 'LineWidth', 2.5, 'Color', 'blue')
    hold on
    plot(f, f_size, f_size, f, 'Color', '#808080', 'LineStyle','--')
    hold on
    yline(Af0, 'Color', 'red', 'LineStyle','-.')
    hold on
    xline(B3, 'Color', 'red', 'LineStyle','-.')
    hold off
    grid on
    ylim([-0.1, 1.1])
    xlim([-20 20])
    ax = gca;
    chart = ax.Children(5);
    datatip(chart,B3,Af0);
    chart = ax.Children(4);
    datatip(chart,B3,0);
    chart = ax.Children(3);
    datatip(chart,0,Af0,"Location","southwest");
    title('Real HP - Amplitude Response', FontSize=20)

    % ---

    % risposta in fase
    nexttile
    Fasef = angle(Hf);
    plot(f, Fasef, 'LineWidth', 2.5, 'Color', 'red')
    hold on
    plot(f, f_size, f_size, f, 'Color', '#808080', 'LineStyle','--')
    hold off
    grid on
    
    ylim([-2 2])
    xlim([-20 20])
    title('Real HP - Phase Response', FontSize=20)
end

Real HP

Passa-Alto RL

Un filtro passa-alto reale di tipo RL è un circuito formato da un resistore con valore di resistenza R in serie ad un induttore con valore di induttanza L. Lo schema elettrico potrebbe trarre in inganno facendo sembrare i componenti in parallelo. Il segnale in ingresso x(t)x(t) è la tensione ai capi dei due componenti in serie, mentre il segnale di uscita y(t)y(t) è la tensione ai capi del solo induttore. La corrente i(t)i(t) scorre nel circuito (idealmente dall’ingresso all’uscita).

Passa-alto Reale

Utilizzando il metodo dei fasori come spiegato per il passa basso, si ottiene il fasore della tensione in uscita.

VY=VX⋅ZLZR+ZLV_Y = V_X \cdot \frac{Z_L}{Z_R + Z_L}

La risposta in frequenza è il rapporto tra il fasore della tensione in uscita e quello della tensione in ingresso:

H(f)=VYVX=ZLZR+ZL=j2πfLR+j2πfL=11+R/(j2πfL)H(f) = \frac{V_Y}{V_X} = \frac{Z_L}{Z_R + Z_L} = \frac{j2\pi fL}{R+j2\pi fL} = \frac{1}{1+R/(j2\pi fL)}

La risposta in ampiezza è il modulo della risposta in frequenza:

A(f)=11+(R2πfL)2\Large A(f) = \frac{1}{\sqrt{1+\Big(\frac{R}{2\pi fL}\Big)^2}}

La risposta in fase è l’argomento della risposta in frequenza:

Θ(f)=arg⁡(H(f))=arctan(R2πfL)\Theta(f) = \arg \bigg(H(f)\bigg) = \text{arctan}\bigg(\frac{R}{2\pi fL}\bigg)
%%%
% sample usage:
% real_high_pass_RL(50, 1)
%
% parameters:
% #1: resistor value in Ohm
% #2: inductor value in Henry
%%%

function real_high_pass_RL(R, L)
    if nargin ~= 2
        fprintf(['»»» USAGE: real_high_pass_RL(R, L)\n' ...
            '\tR:\tresistor value [Ohm]\n' ...
            '\tL:\tinductor value [Henry]\n' ...
            '«««\n']);
        return;
    end
    
    f = linspace(-40, 40, 10001);
    f_size = zeros(size(f));
    w = 2*pi.*f; % pulsazione (omega)
    Z_R = R;
    Z_L = 1j*w*L;
    Hf = Z_L ./ (Z_R + Z_L);

    % amplitude response: A(f) = |H(f)|
    Af = 1./(sqrt(1+(R./(2*pi*L .*f)).^2));

    Af0 = 1 / sqrt(2);
    
    % B_(-3) --> band value at -3 dB
    B3 = R / (2 * pi * L);
    fprintf("»»» valore di banda B = %f\n", B3);
    figure(Name="Real High Pass");
    tiledlayout(1, 2)
    nexttile 
    plot(f, Af, 'LineWidth', 2.5, 'Color', 'blue')
    hold on
    plot(f, f_size, f_size, f, 'Color', '#808080', 'LineStyle','--')
    hold on
    yline(Af0, 'Color', 'red', 'LineStyle','-.')
    hold on
    xline(B3, 'Color', 'red', 'LineStyle','-.')
    hold off
    grid on
    ylim([-0.1, 1.1])
    xlim([-20 20])
    ax = gca;
    chart = ax.Children(5);
    datatip(chart,B3,Af0);
    chart = ax.Children(4);
    datatip(chart,B3,0);
    chart = ax.Children(3);
    datatip(chart,0,Af0,"Location","southwest");
    title('Real HP - Amplitude Response', FontSize=20)

    % ---

    % risposta in fase
    nexttile
    Fasef = angle(Hf);
    plot(f, Fasef, 'LineWidth', 2.5, 'Color', 'red')
    hold on
    plot(f, f_size, f_size, f, 'Color', '#808080', 'LineStyle','--')
    hold off
    grid on
    
    ylim([-2 2])
    xlim([-20 20])
    title('Real HP - Phase Response', FontSize=20)
end

Real HP - RL

Passa-banda Ideale

Simbologia: tre sinusoidi con quelle esterne barrate. Le frequenze tagliate sono quelle più basse di fLf_L (low) e più alte di fHf_H (high). La sigla BP in pedice significa band pass.

La banda passante è compresa tra fLf_L e fHf_H, si può quindi affermare che B=fH−fLB = f_H-f_L.

La frequenza di centro-banda è:

f0=fL+fH2f_0 = \frac{f_L + f_H}{2}

Il fattore di selettività (o qualità) è: Q=f0/BQ = f_0/B

HBP=rect(f−f0B)+rect(f+f0B)hBP=2B⋅sinc(Bt)⋅cos⁡(2πf0t)\eq{ H_{\rm BP} &= rect\bigg(\frac{f-f_0}{B}\bigg) + rect\bigg(\frac{f+f_0}{B}\bigg) \\ h_{\rm BP} &= 2B \cdot sinc(Bt) \cdot \cos(2\pi f_0 t) }

BP risposta in frequenza

BP risposta impulsiva

Passa-Banda Reale RLC

Un bipolo è una rete elettrica caratterizzata da due soli terminali (o poli) di accesso. Può essere costituito da uno o più elementi circuitali connessi fra loro.

Un filtro passa-banda reale è un circuito formato da un resistore con valore di resistenza R in serie ad un bipolo costituito da un condensatore con capacità elettrica C in parallelo ad un induttore con induttanza L. Il segnale in ingresso x(t)x(t) è la tensione ai capi dei due componenti in serie (resistore e bipolo), mentre il segnale di uscita y(t)y(t) è la tensione ai capi del solo bipolo (condensatore in parallelo ad induttore).

Passa-banda Reale

Il metodo dei fasori prevede il passaggio a numeri complessi per rappresentare circuiti elettrici che operano in regime sinusoidale. Si ottengono i fasori:

ZR=R;ZL=j2πfL;ZC=1j2πfCZ_R = R \quad ; \quad Z_L = j 2\pi f L \quad ; \quad Z_C = \frac{1}{j 2\pi f C}

Si calcola ZP, ovvero l’impedenza in parallelo tra ZL e ZC:

ZP=ZL∥ZC=ZL⋅ZCZL+ZC==j2πfL−(2πf)2LC+1\eq{ Z_P &= Z_L \parallel Z_C = \frac{Z_L \cdot Z_C}{Z_L + Z_C} = \\ &= \frac{j 2\pi f L}{-(2\pi f)^2 LC + 1} }

Dato il fasore VX della funzione in ingresso x(t)x(t). Il fasore VY della tensione in uscita y(t)y(t) è dato da:

VY=VX⋅ZPZR+ZPV_Y = V_X \cdot \frac{Z_P}{Z_R + Z_P}

La risposta in frequenza è il rapporto tra il fasore della tensione in uscita e quello della tensione in ingresso:

H(f)=VYVX=ZPZR+ZP=j2πfLR[1−(2πf)2LC]+j2πfLH(f) = \frac{V_Y}{V_X} = \frac{Z_P}{Z_R + Z_P} = \frac{j 2\pi f L}{R[1-(2\pi f)^2 LC ] + j 2\pi f L }

La risposta in ampiezza è il modulo della risposta in frequenza:

A(f)=∣H(f)∣=11+[1−(2πf)2LC]2(2πf)2L2⋅R2\Large A(f) = |H(f)| = \frac{1}{\sqrt{1+\frac{[1-(2\pi f)^2 LC]^2}{(2\pi f)^2 L^2}\cdot R^2}}

La frequenza massima della risposta in ampiezza è:

fmax=max⁡[A(f)]  ⟹  1−(2πf)2LC=0  ⟹  fmax=±12πLC\eq{ f_{\rm max}&=\max[A(f)] \implies 1-(2\pi f)^2 LC = 0 \\ \implies f_{\rm max} &= \pm \frac{1}{2\pi \sqrt{LC}} }

La frequenza di centro banda coincide con la frequenza massima, nel semi-piano positivo:

f0=12πLCf_0 = \frac{1}{2\pi \sqrt{LC}}

La frequenza massima (sull’asse orizzontale) interseca il grafico della risposta in ampiezza per A(fmax)=1A(f_{\rm max})=1. Le frequenze fLf_L ed fHf_H intersecano il grafico per A(f)=1/2A(f)=1/\sqrt{2}:

A(f−3)A(f0)=12\frac{A(f_{-3})}{A(f_0)} = \frac{1}{\sqrt{2}}

Bisogna fare qualche utile precisazione:

f0≡fmax  ⟹  A(f0)=A(fmax)=1  ⟹  A(f−3)1=12  ⟹  A(f−3)=12\eq{ f_0 \equiv f_{\rm max} \implies & A(f_0)=A(f_{\rm max})=1\\ \implies & \frac{A(f_{-3})}{1} = \frac{1}{\sqrt{2}}\\ \implies & A(f_{-3}) = \frac{1}{\sqrt{2}}\\ }

La banda B−3=−3[dB]B_{-3}=-3[dB] è data da:

f−3L≡fL[dB]  ;  f−3H≡fH[dB]B−3=f−3H−f−3L\eq{ & f_{-3L} \equiv f_L[dB] \;;\; f_{-3H} \equiv f_H[dB]\\ & B_{-3} = f_{-3H}-f_{-3L}\\ }
%%%
% sample usage:
% real_band_pass(5, 0.4, 0.01)
%
% parameters:
% #1: resistor value in Ohm
% #2: inductor value in Henry
% #3: capacitor value in Farad
%%%

function real_band_pass(R, L, C)
    f = linspace(-40, 40, 10001);
    f_size = zeros(size(f));
    w = 2*pi.*f;
    % calcolo di H(f)
    Z_R = R;
    Z_L = 1j*w*L;
    Z_C = 1./(1j*w*C);
    Z_P = (Z_L .* Z_C)./(Z_L + Z_C);

    Hf = Z_P ./ (Z_R + Z_P);

    % amplitude response: A(f) = |H(f)|
    Af = abs(Hf); % amplitude response: A(f) = |H(f)|
    % 1 ./ sqrt(1 + (2 * pi * R * C .* f) .^ 2);

    % center frequency: frequenza di centro banda
    f_0 = 1 / (2*pi*sqrt(L*C)); % x-axis
    % quality factor: Il fattore di selettività (o qualità)
    Q = 1/R * sqrt(L/C);
    % band value
    B = f_0 ./ Q;
    f_low = f_0 - (B/2);
    f_high = f_0 + (B/2);

    fprintf("»»» valore di banda B = %f\n", B);
    fprintf("»»» frequenza di centro banda f0 = %f\n", f_0);
    fprintf("»»» quality factor Q = %f\n", Q);
    fprintf("»»» f_low = %f\n", f_low);
    fprintf("»»» f_high = %f\n", f_high);
    figure('Name', 'Real Low Pass - Amplitude Response');
    plot(f, Af, 'LineWidth', 2.5, 'Color', 'blue')
    hold on
    % draw horizontal and vertical axes
    plot(f, f_size, f_size, f, 'Color', '#808080', 'LineStyle','--')
    % ...
    % risposta in fase
    Fasef = angle(Hf);
    figure('Name', 'Real Low Pass - Phase Response');
    plot(f, Fasef, 'LineWidth', 2.5, 'Color', 'red')
    % ...
end

real band pass - A(f)

real band pass - fase

Elimina-Banda

Simbologia: tre sinusoidi con quella interna barrata. Le frequenze tagliate sono quelle più alte di fLf_L (low) e più basse di fHf_H (high). La sigla BR in pedice significa band rejection.

Il valore di banda B e di f0f_0 sono uguali al filtro passa-banda.

HBR=1−HBP=1−rect(f−f0B)−rect(f+f0B)hBR=δ(t)−[2B⋅sinc(Bt)⋅cos⁡(2πf0t)]\eq{ H_{\rm BR}&=1-H_{\rm BP}=1-rect\bigg(\frac{f-f_0}{B}\bigg)-rect\bigg(\frac{f+f_0}{B}\bigg) \\ h_{\rm BR}&=\delta(t)-\Big[2B \cdot sinc(Bt) \cdot \cos(2\pi f_0 t)\Big] }

La banda passante è compresa tra 00 e fLf_L e si estende da fHf_H a ∞\infty.

BR risposta in frequenza

BR risposta impulsiva

Calcolo della Banda

I filtri ideali sono definiti utilizzando il valore B di banda passante. Un filtro reale cerca di approssimare il comportamento di un filtro ideale ed è definito tramite circuito composto da resistenze, condensatori e induttori.

I filtri reali sono costruiti per ottenere un certo valore di banda, che diventa quindi il punto di arrivo dello sviluppo.

Passa-Basso RC e Passa-Alto CR

Il valore B di banda passante di un filtro passa-alto/basso reale è legato ai valori di resistenza R e capacità elettrica C.

La costante di tempo del filtro è: T=RCT=RC.

Grazie alla natura del filtro, la banda coincide con la frequenza di taglio (cutoff):

B=f0=12πT=12πRCB = f_0 = \frac{1}{2\pi T}=\frac{1}{2\pi RC}

Il calcolo della banda è uguale per entrambi i filtri. La differenza tra passa-alto e passa-basso sta nelle frequenze che possono passare: sopra o sotto il valore di banda.

N.B. La frequenza di taglio viene indicata con la c o lo zero in pedice: fc≡f0f_c \equiv f_0.

Passa-Basso LR

Il valore B di banda passante di un filtro passa-basso LR (serie) è legato ai valore di induttanza L e resistenza R. La costante di tempo del filtro è:

T=LRT = \frac{L}{R}

La banda coincide con la frequenza di taglio (cutoff):

B=fc=12πT=R2πLB = f_c = \frac{1}{2\pi T}=\frac{R}{2\pi L}

(Esisterà sicuramente un passa-alto del tipo RL. La sua banda si calcolerà sicuramente nel medesimo modo).

Passa-Banda RLC

In un filtro passa-banda RLC, le correlazioni tra i valori di resistenza R, induttanza L e capacità C forniscono i parametri caratteristici del filtro.

La frequenza di centro banda è:

f0=12πLCf_0 = \frac{1}{2\pi \sqrt{LC}}

Il fattore di selettività (o qualità) è:

Q=1RLCQ = \frac{1}{R} \sqrt{\frac{L}{C}}

La banda passante è:

B=f0QB = \frac{f_0}{Q}

Le frequenze di cutoff inferiore e superiore sono:

fL=f0−B2;fH=f0+B2f_L = f_0 - \frac{B}{2} \quad ; \quad f_H = f_0 + \frac{B}{2}

Il coefficiente di smorzamento è:

ζ=12Q=R2CL\zeta = \frac{1}{2Q} = \frac{R}{2} \sqrt{\frac{C}{L}}

La frequenza angolare (natural frequency) è:

ω0=1LC\omega_0 = \frac{1}{\sqrt{LC}}

La costante di tempo si ottiene grazie al coefficiente di smorzamento e la frequenza angolare:

T=1ω01−ζ2T = \frac{1}{\omega_0 \sqrt{1-\zeta^2}}

Applicazioni pratiche

I filtri studiati nella prima parte di questo articolo trovano importanti applicazioni “reali” nell’editing di immagini e nella modifica di file audio.

Filtri applicati a file audio

Lo script filtro_audio, che verrà illustrato alla fine del paragrafo, permette di applicare i quattro principali tipi di filtro ad una traccia audio fornita dall’utente. Lo script è stato opportunamente modificato, rispetto alla versione illustrata in classe, per permettere all’utente di selezionare i seguenti parametri:

  • filename: file audio di ingresso,
  • filter_type: tipo di filtro da applicare al file di ingresso,
  • B: dimensione della banda in Hz,
  • f0: frequenza centrale in Hz.

Si prenda come esempio di file audio un arpeggio, una serie di accordi suonati con una chitarra classica o la canzone Scar Tissue dei Red Hot Chili Peppers.

Passando come filter_type il valore bp o band-pass, lo script filtro_audio filtra le frequenze al di fuori della banda B centrate nella frequenza f0f_0. Si può verificare che, diminuendo progressivamente il valore di B, l’audio in uscita risulterà sempre più selettivo attorno alla nota musicale relativa alla frequenza f0 passata come parametro. Mantenendo invece B costante e facendo variare f0, grazie al grafico delle trasformate dei segnali, è possibile apprezzare lo scostamento della banda passante.

Scar Tissue - passa-basso - A=5 Scar Tissue - passa-alto - A=5 Scar Tissue - passa-banda - A=5 Scar Tissue - elimina-banda - A=5

%% Utilizzo:
% filtro_audio("ScarTissue.wav", 'bp', 3.0e3, 3.0e3)
% filename ~ nome del file da leggere (path, URL, ...)
% filter_type ~ una stringa tra le seguenti:
% 'bp','band-pass', 'br','band-rejection', 'lp','low-pass', 'hp', 'high-pass'
% B ~ dimensione di banda [Hz]
% f0 ~ frequenza centrale [Hz]
function filtro_audio(filename, filter_type, B, f0)
% --- Lettura campioni --- %
info = audioinfo(filename);      % informazioni sul file audio
freqC = 44.1e3;                  % [Hz]: frequenza di campionamento
tempoC = 1 / freqC;              % [s]: tempo di campionamento
durata = info.Duration;          % [s]: durata del segnale audio
numCampioni = freqC * durata;    % 
%numCampioni=info.TotalSamples;  % numero totale di campioni
% funziona solo se il file è campionato a 44.1e3, altrimenti dà errore

[xstereo,fc] = audioread(filename,[1,numCampioni]);
x = (xstereo(:,1))'; % solo canale sinistro
% la frequenza di campionamento dipende dal file sorgente, es. CD: 44.1 kHz
fprintf('»»» frequenza di campionamento = %d [Hz]\n', fc);

T = 50/B; % durata risposta impulsiva [s]
tempoFiltro=0:tempoC:T;
tempoNorm = tempoFiltro-T/2;
% --- applica un filtro in base al parametro filter_type --- %
switch filter_type
    case {'bp','band-pass'} % --- Filtro passa-banda ideale --- %
        % più opzioni per uno stesso blocco
        figure(Name="Filtro passa banda")
        g=2*B*sinc(B*tempoNorm).*rectpuls(tempoNorm/T).*cos(2*pi*f0*tempoNorm);
    case {'br','band-rejection'} % --- Filtro elimina-banda ideale --- %
        figure(Name="Filtro elimina banda")
        g=2*B*sinc(B*tempoNorm).*rectpuls(tempoNorm/T).*(cos(2*pi*f0*tempoNorm).^2);
    case {'lp','low-pass'}  % --- Filtro passa-basso --- % 
        figure(Name="Filtro passa basso")
        g=2*B*sinc(2*B*tempoNorm).*rectpuls(tempoNorm/T);
    case {'hp', 'high-pass'}% --- Filtro passa-alto --- % 
        figure(Name="Filtro passa alto")
        g=2*B*sinc(B*tempoNorm).*rectpuls(tempoNorm/T).*cos(2*pi*f0*tempoNorm);
    otherwise
        warning('Unexpected filter type.')
        return;
end

% --- Uscita filtrata --- %
w = conv(g,x) * tempoC; % convoluzione tra g ed x
w = w(length(g):length(w));

tempo = 0:tempoC:durata-tempoC;
tempoW = tempo(1:length(w))+T/2;

durataTransitorio = 0.0250; % [s]
A = 5; % amplificazione: rende più visibile il segnale filtrato (è un valore empirico)

tiledlayout(1, 2); nexttile;
plot((tempo-durataTransitorio)*1e3, x, 'Color', 'cyan', 'LineWidth', 2.5);
hold on;
plot((tempoW-durataTransitorio)*1e3, A*w, 'Color', 'black', 'LineWidth', 1.5);
grid on;
xlabel('Tempo (ms)', FontSize=14);
ylabel('Segnali temporali', FontSize=14);
legend('x(t)', 'A\cdot{w(t)}', FontSize=12);
axis([2000 2200 -1*max(abs(w)) max(abs(w)) ]);
titoloGrafico=sprintf('Filtraggio %s con banda B = %.2f [Hz] e frequenza centrale %.2f [Hz]', filter_type, B, f0);
title(titoloGrafico,'FontSize',16);
hold off;

% --- Calcolo della trasformata di Fourier dell'ingresso --- %
lunghezzaFft = 2^nextpow2(length(x));
X = fft(x,lunghezzaFft)*tempoC;
X = [X(lunghezzaFft/2+1:lunghezzaFft) X(1:lunghezzaFft/2+1)];
frequenza = freqC*linspace(-0.5,0.5,lunghezzaFft+1);

% --- Calcolo della trasformata di Fourier dell'uscita --- %
%lunghezzaFft=2^nextpow2(length(w));
W = fft(w,lunghezzaFft)*tempoC;
W = [W(lunghezzaFft/2+1:lunghezzaFft) W(1:lunghezzaFft/2+1)];
%frequenza=freqC*linspace(-0.5,0.5,lunghezzaFft+1);

nexttile; % aggiunge una "piastrella" alla figura
plot(frequenza/1e3,abs(X), 'Color', 'cyan', 'LineWidth', 1.5);
hold on;
plot(frequenza/1e3,abs(W), 'Color', 'black', 'LineWidth', 1.5);
grid on;
xlabel('Frequenza (kHz)', FontSize=14);
ylabel('Spettro di ampiezza', FontSize=14);
legend('|X(f)|', '|W(f)|', FontSize=12);
axis([0 8 0 1.2*max(abs(W))]);
title('Trasformate di Fourier dei segnali convoluti', 'FontSize',16)
hold off;

file = split(filename,".");
name = string(file(1));
out_filename = sprintf('./output/%s_output_%s_b%.0f_f%.0f.wav',name, filter_type, B, f0);
audiowrite(out_filename,[w.'*0.99/max(abs(w)),w.'*0.99/max(abs(w))],freqC);

end

Clipping

Lo script clip_audio prende tutto ciò che viene fatto in filtro_audio ed aggiunge un blocco non lineare a cui il segnale (relativo al file audio) viene dato in pasto a seguito del filtraggio. Questo blocco ha un andamento lineare finché |y|<yM. Superato questo limite, il valore di z satura a yM o –yM.

Lo scopo del blocco è quello di tagliare il segnale in ingresso, aggiungendovi dei disturbi udibili nel file di output generato al termine del programma. I disturbi sono inversamente proporzionali al valore di yM. Le trasformazioni non lineari aggiungono inoltre dei contributi in frequenza e contrastano quindi il lavoro dei filtri.

Per semplicità, clip_audio permette di applicare solamente un filtro passa-basso.

Dato il seguente blocco:

% --- Blocco non lineare --- %
switch block_number
    case 1 % --- blocco #1 --- %
        yM=0.12;      % ampiezza della saturazione
        z=y;          % copia non distorta del segnale di ingresso
        z(y>yM)=yM;   % clipping dei valori maggiori di yM
        z(y<-yM)=-yM; % clipping dei valori minori di -yM
    case 2
        yM=0.10;      % ampiezza della saturazione
        z=y;          % copia non distorta del segnale di ingresso
        z(y>0)=yM;    % clipping dei valori maggiori di 0
        z(y<0)=-yM;   % clipping dei valori minori di 0
    case 3
        yM=0.10;      % ampiezza della saturazione
        z=y;          % copia non distorta del segnale di ingresso
        z(y>0)=yM*(1-exp(-y(y>0)./yM)); % clipping dei valori maggiori di 0
        z(y<0)=-yM*(1-exp(y(y<0)./yM));% clipping dei valori minori di 0
    otherwise
        warning('Unexpected block number.')
        return;
end

Fornita in ingresso una sinusoide di ampiezza 0.2 e periodo uguale alla dimensione della banda, è interessante notare il confronto del segnale prima e dopo l’applicazione del blocco non lineare.

sin blocco #1 sin blocco #2 sin blocco #3

Applicando il 3° blocco distorsivo all’introduzione della canzone Scar Tissue dei Red Hot Chili Peppers, si ottengono i seguenti grafici.

Scar Tissue cliping - grafico 1 Scar Tissue cliping - grafico 2 Scar Tissue cliping - grafico 3

Si può chiaramente apprezzare il disturbo dovuto al blocco non-lineare. I grafici mostrano in modo chiaro la distorsione, ma ascoltare la traccia audio generata dal comando audiowrite è sicuramente il modo migliore per comprendere l’entità della distorsione.

%% Utilizzo:
% clip_audio('Scar_Tissue.mp3', 3, 1.0e3)
% filename: nome del file da leggere (path, URL, ...)
% block_number: numero del blocco non-lineare da applicare
% B: dimensione di banda [Hz]
function clip_audio(filename, block_number, B)

% --- Lettura campioni --- %
freqC = 44.1e3;        % [Hz]: frequenza di campionamento = 44.1 kHz
tempoC = 1/freqC;      % [s] tempo di campionamento
durataTransitorio=0.1; % [s]
durata = 10.0;         % [s]
numeroCampioni = durata * freqC; % numero totale di campioni
inizioCampioni = durataTransitorio*freqC;
[xstereo,~]=audioread(filename,[inizioCampioni+1,inizioCampioni+1+numeroCampioni]);
x = (xstereo(:,1))'; % solo canale sinistro
tempo=0:tempoC:durata;
% x = 0.2 .* sin(B.*tempo);

% --- Filtro passa-basso --- %
T = 50/B; % durata risposta impulsiva [s]
tempoFiltro=0:tempoC:T;
tempoNorm = tempoFiltro-T/2;
h=2*B*sinc(2*B*tempoNorm).*rectpuls(tempoNorm/T);

% --- Filtraggio --- %
y=conv(x,h)*tempoC; % uscita filtrata
y=y(length(h):length(y));
tempoY=tempo(1:length(y))+T/2;
tempoY_plot=tempoY-durataTransitorio;

figure;
set(gcf,'defaultaxesfontname','Courier New')
plot(tempo-durataTransitorio, x, 'Color', 'cyan', 'LineWidth', 1.5);
hold on;
plot(tempoY_plot, y, 'Color', 'black', 'LineWidth', 1.5);
grid on;
xlabel('Tempo (s)','FontSize',12);
ylabel('Segnali temporali','FontSize',12)
legend('x(t)', 'y(t)');
axis([0 0.4 -0.8001 0.8001]);

% --- Blocco non lineare --- %
switch block_number
    case 1 % --- blocco #1 --- %
        yM=0.12;      % ampiezza della saturazione
        z=y;          % copia non distorta del segnale di ingresso
        z(y>yM)=yM;   % clipping dei valori maggiori di yM
        z(y<-yM)=-yM; % clipping dei valori minori di -yM
    case 2 % --- blocco #2 --- %
        yM=0.10;      % ampiezza della saturazione
        z=y;          % copia non distorta del segnale di ingresso
        z(y>0)=yM;    % clipping dei valori maggiori di 0
        z(y<0)=-yM;   % clipping dei valori minori di 0
    case 3 % --- blocco #3 --- %
        yM=0.10;      % ampiezza della saturazione
        z=y;          % copia non distorta del segnale di ingresso
        z(y>0)=yM*(1-exp(-y(y>0)./yM)); % clipping dei valori maggiori di 0
        z(y<0)=-yM*(1-exp(y(y<0)./yM));% clipping dei valori minori di 0
    otherwise
        warning('Unexpected block number.')
        return;
end

figure;
plot(tempoY_plot, y, 'Color', 'cyan', 'LineWidth', 1.5); hold on;
plot(tempoY_plot, z, 'Color', 'black', 'LineWidth', 1.5); grid on;
xlabel('Tempo (s)','FontSize',12);
ylabel('Segnali temporali','FontSize',12);
legend('y(t)', 'z(t)');
axis([0 0.4 -0.2001 0.2001]);

% --- Calcolo della trasformata di Fourier dell'ingresso --- %
lunghezzaFft=2^nextpow2(length(y));
Y=fft(y,lunghezzaFft)*tempoC;
Y=[Y(lunghezzaFft/2+1:lunghezzaFft) Y(1:lunghezzaFft/2+1)];
frequenza=freqC*linspace(-0.5,0.5,lunghezzaFft+1);

% --- Calcolo della trasformata di Fourier dell'uscita --- %
%lunghezzaFft=2^nextpow2(length(z));
Z=fft(z,lunghezzaFft)*tempoC;
Z=[Z(lunghezzaFft/2+1:lunghezzaFft) Z(1:lunghezzaFft/2+1)];
%frequenza=freqC*linspace(-0.5,0.5,lunghezzaFft+1);

figure;
plot(frequenza/1e3,20*log10(abs(Y)./max(abs(Y))), 'Color', 'cyan', 'LineWidth', 1.5);
hold on;
plot(frequenza/1e3,20*log10(abs(Z)./max(abs(Y))), 'Color', 'black', 'LineWidth', 1.5);
plot(frequenza/1e3,20*log10(abs(Y)./max(abs(Y))), 'Color', 'cyan', 'LineWidth', 0.5);
grid on;
xlabel('Frequenza (kHz)','FontSize',12);
ylabel('Spettro di ampiezza (dB)','FontSize',12);
legend('|Y(f)|', '|Z(f)|');
axis([0 8 -120 0]);

% scrittura del file audio a cui è stato applicato il filtro
% N.B. se il filename è un path relativo e inizia con uno o più punti,
% questo modo per estrapolare il nome senza estensione non funziona.
% Dunque 'filename' deve essere un path assoluto, a meno che il file non si
% trovi nella cartella nella quale è situato lo script
file = split(filename,".");
name = string(file(1));
out_filter= sprintf('./output/%s_output_filter_b%.0f.wav',name,B);
out_dist = sprintf('./output/%s_output_dist_b%.0f.wav',name,B);

audiowrite(out_filter,[y.',y.'],freqC);
audiowrite(out_dist,[z.',z.'],freqC);

Filtri applicati ad immagini

MATLAB è in grado di leggere il contenuto di un immagine mediante il comando imread, ad esempio:

x=imread('unife', '.tif');
% oppure, con singolo parametro:
x=imread('unife.tif');

Ad x viene assegnata la matrice bidimensionale contenente il segnale di intensità luminosa dell’immagine unife.tif presente nel path corrente. Le componenti a bassa frequenza descrivono le caratteristiche generali dell’immagini (ai quali l’occhio umano è più sensibile). Le componenti ad alta frequenza descrivono invece i dettagli dell’immagine.

Il comando imgaussfilt applica un filtraggio passa-basso di tipo Gaussiano. Si ha la possibilità di passare il parametro opzionale deviazione standard σ\sigma (sigma):

sigma = 20; % default: 0.5
y1 = imgaussfilt(x, sigma);

Si può applicare un filtro passa-alto sottraendo al segnale originale il segnale ricavato dal filtraggio passa-basso.

z1=x-y1;

Si può ulteriormente filtrare l’immagine sottraendo il segnale filtrato ad un quadrato totalmente bianco.

quadratoBianco=255*ones(size(x),'uint8');
y2=quadratoBianco-z1;

Lo script che segue applica questi concetti per generare e mostrare le immagini filtrate, a partire da un’immagine di partenza fornita come parametro, è il seguente:

%% Utilizzo:
% filtro_lp_hp_immagine("Unife.tif")
function filtro_lp_hp_immagine(filename)

% estrapolazione di nome ed estensione dal parametro in ingresso
file = split(filename,".");
name = file(1);
ext = file(2);

x=imread(name, ext); % lettura immagine

figure;
imshow(x); % effettua il grafico dell'immagine
title('Immagine originale');

%% filtraggio passa-basso di tipo Gaussiano
sigma = 20;
y1 = imgaussfilt(x, sigma);

figure
imshow(y1);
title('Immagine elaborata con il filtro passa-basso');

% genera un'immagine a partire dal vettore bidimensionale in ingresso
% il nome di tale immagine è una concatenazione di stringhe
imwrite(y1,strcat(name,'_lp.', ext),ext);

%% Filtraggio passa-alto
z1=x-y1;

figure
imshow(z1);
title('Immagine elaborata con il filtro passa-alto');
imwrite(z1,strcat(name,'_hp.', ext),ext);

% si può filtrare l'immagine ottenuta sottraendo il segnale filtrato ad un quadrato totalmente bianco 
quadratoBianco=255*ones(size(x),'uint8');
y2=quadratoBianco-z1;

figure
imshow(y2);
title('passa-alto (bianco)');
imwrite(y2,strcat(name,'_white_hp.', ext),ext);

A partire dall’immagine del Dipartimento di Ingegneria di Unife, sono illustrati i risultati dello script precedente:

Inge originaleInge lp
Inge hpInge hp bianco

Applicazioni della funzione fspecial

Tra gli innumerevoli modi per applicare filtri alle immagini, è bene parlare della funzione fspecial, la quale permette di applicare un filtro avente una determinata risposta impulsiva h passata come parametro. Per applicare un filtro con fspecial bisogna fare uso di imfilter, ad esempio:

# H = fspecial('nome-filtro', altri-parametri);

H = fspecial('motion',20,45);
B = imfilter(I,H);
imshow(B);

Estrazione delle componenti RGB

È possibile applicare altre tipologie di filtraggio, tra le quali l’estrazione delle componenti RGB (red, green, blue).

function rgb_extractor(filename)
% estrapolazione di nome ed estensione dal parametro in ingresso
% N.B. se il filename è un path relativo e inizia con uno o più punti,
% questo modo per estrapolare il nome senza estensione non funziona.
% Dunque 'filename' deve essere un path assoluto, a meno che il file non si
% trovi nella cartella nella quale è situato lo script
file = split(filename,".");
name = file(1);
ext = file(2);

my_image=imread(filename);

% componenti rossa, verde e blu in termini di intensità luminosa
red_brigh=my_image(:,:,1);
green_brigh=my_image(:,:,2);
blue_brigh=my_image(:,:,3);
black=zeros(size(my_image,1),size(my_image,2));
white=my_image(:,:,1:1:1);

% componenti rossa, verde e blu in termini di colore
red_component=cat(3,red_brigh,black,black);
green_component=cat(3,black,green_brigh,black);
blue_component=cat(3,black,black,blue_brigh);
just_black=cat(3,black,black,black);
just_white=cat(3,white,white,white);

figure
imshow(my_image);

figure
imshow(red_component);
imwrite(red_component,strcat('./output/', name,'_red.', ext),ext);

figure
imshow(green_component);
imwrite(green_component,strcat('./output/', name,'_green.', ext),ext);

figure
imshow(blue_component);
imwrite(blue_component,strcat('./output/', name,'_blue.', ext),ext);

figure
imshow(just_black);
figure
imshow(just_white);
end

La tabella che segue illustra i risultati ottenuti con lo script appena illustrato:

originalerosso
verdeblu