《MATLAB輔助現(xiàn)代工程數(shù)字信號處理》課件第3章

第3章數(shù)字濾波器的結(jié)構(gòu)與分析3.1數(shù)字濾波器及其實(shí)現(xiàn)3.2IIR數(shù)字濾波器的基本結(jié)構(gòu)3.3FIR數(shù)字濾波器的基本結(jié)構(gòu)3.4數(shù)字濾波器的格型結(jié)構(gòu)3.5數(shù)字濾波器的MATLAB時(shí)頻分析3.6數(shù)字濾波器的量化效應(yīng)3.7小結(jié)

3.1數(shù)字濾波器及其實(shí)現(xiàn)

3.1.1數(shù)字濾波器概述

數(shù)字濾波器是數(shù)字信號處理的基礎(chǔ)部分，與模擬濾波器相比，數(shù)字濾波器具有精度高、可靠性高、靈活性高、便于大規(guī)模集成和多維過濾等優(yōu)點(diǎn)，已廣泛應(yīng)用于現(xiàn)代各類工程領(lǐng)域。

數(shù)字濾波器是具有一定傳輸特性的數(shù)字信號處理裝置，它的輸入和輸出都是離散數(shù)字信號，它借助于數(shù)字器件和一定的數(shù)值計(jì)算方法，對輸入信號進(jìn)行處理，改變輸入信號，進(jìn)而去掉信號中的無用成分而保留有用成分。如果在數(shù)字處理系統(tǒng)前、后分別加上A/D轉(zhuǎn)換器和D/A轉(zhuǎn)換器，就可以處理模擬信號。數(shù)字濾波器的輸入輸出是一個(gè)時(shí)間序列。設(shè)H(z)為數(shù)字濾波器的系統(tǒng)函數(shù)，h(n)為其相應(yīng)的脈沖序列，則在時(shí)域內(nèi)有

y(n)=x(n)*h(n)

(3.1)

在z域內(nèi)有

Y(z)=H(z)X(z)

(3.2)

式中，X(z)和Y(z)分別為輸入x(n)和輸出y(n)的Z變換。

在頻域內(nèi)有

Y(jω)=H(jω)X(jω)

(3.3)

式中，H(jω)為數(shù)字濾波器的頻率特性;X(jω)和Y(jω)分別為輸入x(n)和輸出y(n)的頻譜。3.1.2數(shù)字濾波器的分類與實(shí)現(xiàn)

數(shù)字濾波器有多種分類方法，按功能，數(shù)字濾波器可分為低通濾波器、高通濾波器、帶通濾波器和帶阻濾波器。按脈沖響應(yīng)的長度，數(shù)字濾波器可分為無限脈沖響應(yīng)(IIR)濾波器和有限脈沖響應(yīng)(FIR)濾波器。

IIR濾波器的差分方程為

(3.4)IIR濾波器的系統(tǒng)函數(shù)為(3.5)

FIR濾波器的差分方程為(3.6)FIR濾波器的系統(tǒng)函數(shù)為(3.7)

IIR濾波器在結(jié)構(gòu)上存在輸出到輸入的反饋;FIR濾波器在結(jié)構(gòu)上不存在輸出到輸入的反饋，信號流圖中不存在環(huán)路。數(shù)字濾波器既可利用專用處理器(如DSP)實(shí)現(xiàn)，也可直接利用計(jì)算機(jī)和通用軟件編程實(shí)現(xiàn)，其具體實(shí)現(xiàn)過程如圖3.1所示。圖3.1數(shù)字濾波器的實(shí)現(xiàn)一個(gè)輸出序列是其過去N點(diǎn)的線性組合加上當(dāng)前輸入序列與過去N點(diǎn)輸入序列的線性組合。輸出y(n)除了與當(dāng)前的輸入x(n)有關(guān)外，同時(shí)還與過去的輸入和過去的輸出有關(guān)，系統(tǒng)是帶有記憶的。

對于上面的算式，可以化成不同的計(jì)算形式，如直接計(jì)算、分解為多個(gè)有理函數(shù)相加、分解為多個(gè)有理函數(shù)相乘等，不同的計(jì)算形式表現(xiàn)出不同的計(jì)算結(jié)構(gòu)，而不同的計(jì)算結(jié)構(gòu)可能會帶來不同的效果，或者是實(shí)現(xiàn)簡單、編程方便，或者是計(jì)算精度較高等。由于數(shù)字信號是通過采樣和轉(zhuǎn)換得到的，而轉(zhuǎn)換的位數(shù)是有限的，所以存在著量化誤差。此外，計(jì)算機(jī)中數(shù)的表示也總是有限的，經(jīng)此表示的濾波器的系數(shù)同樣存在著量化誤差，故在計(jì)算過程中因有限字長也會造成誤差。

量化誤差主要有三種:A/D變換量化效應(yīng)，系數(shù)的量化效應(yīng)，數(shù)字運(yùn)算的有限字長效應(yīng)。3.1.3數(shù)字濾波器的運(yùn)算結(jié)構(gòu)

數(shù)字濾波器無論采用硬件實(shí)現(xiàn)還是采用軟件實(shí)現(xiàn)，首先都應(yīng)確定數(shù)字濾波器的運(yùn)算結(jié)構(gòu)。運(yùn)算結(jié)構(gòu)可以用方框圖表示，也可用信號流圖來表示。為了簡單起見，通常用信號流圖來表示其運(yùn)算結(jié)構(gòu)。信號流圖表示濾波器有加法、乘法及延遲三種基本運(yùn)算單元。

在信號流圖中，只有輸出支路的節(jié)點(diǎn)稱為輸入節(jié)點(diǎn)或源點(diǎn);只有輸入支路的節(jié)點(diǎn)稱為輸出節(jié)點(diǎn)或阱點(diǎn);既有輸入支路又有輸出支路的節(jié)點(diǎn)叫做混合節(jié)點(diǎn)。通路是指從源點(diǎn)到阱點(diǎn)之間沿著箭頭方向的連續(xù)的一串支路，通路的增益是該通路上各支路增益的乘積?；芈肥侵笍囊粋€(gè)節(jié)點(diǎn)出發(fā)沿著支路箭頭方向到達(dá)同一個(gè)節(jié)點(diǎn)的閉合通路，它象征著系統(tǒng)中的反饋回路。組成回路的所有支路增益的乘積通常稱為回路增益。在信號流圖中，系統(tǒng)函數(shù)按梅遜(Mason)公式計(jì)算:(3.8)式中，Tk為從輸入節(jié)點(diǎn)(源點(diǎn))到輸出節(jié)點(diǎn)(阱點(diǎn))的第k條前向通路增益;Δ為流圖的特征式，即(3.9)式中，為所有不同回路增益之和;為每兩個(gè)互不接觸回路增益之和;Δk是不接觸第k條前向通路的特征式余因子。信號流圖的轉(zhuǎn)置定理:對于單個(gè)輸入、單個(gè)輸出的系統(tǒng)，通過反轉(zhuǎn)網(wǎng)絡(luò)中的全部支路的方向，并且將其輸入和輸出互換，得出的流圖具有與原始流圖相同的系統(tǒng)函數(shù)。信號流圖轉(zhuǎn)置可以轉(zhuǎn)變運(yùn)算結(jié)構(gòu)，并可以驗(yàn)證所計(jì)算的流圖系統(tǒng)函數(shù)的正確與否。

運(yùn)算結(jié)構(gòu)對濾波器的實(shí)現(xiàn)很重要，尤其對于一些定點(diǎn)運(yùn)算的處理機(jī)，結(jié)構(gòu)的不同將會影響系統(tǒng)的精度、誤差、穩(wěn)定性、經(jīng)濟(jì)性以及運(yùn)算速度等許多重要的性能。對于IIR數(shù)字濾波器與FIR數(shù)字濾波器，它們在結(jié)構(gòu)上各有自己不同的特點(diǎn)，下面分別加以討論。3.1.4數(shù)字濾波器的MATLAB實(shí)現(xiàn)

MATLAB信號工具箱提供了很多內(nèi)部函數(shù)來設(shè)計(jì)數(shù)字濾波器，具體函數(shù)將在以后各章中逐一介紹，這里僅舉幾例來介紹數(shù)字濾波的MATLAB實(shí)現(xiàn)。

【例3.1】

用MATLAB的內(nèi)部函數(shù)分別設(shè)計(jì)一個(gè)巴特沃斯和切比雪夫數(shù)字濾波器。

MATLAB程序如下:

%MATLABPROGRAM3-1

clc;

fp=2;

%2kHz

fs=10;

%10kHz

Rp=3;

%3dB

Rs=60;

%60dB

Wp=fp/fp;

Ws=fs/fp;

[n，Wn]=buttord(Wp，Ws，Rp，Rs，‘s′);%巴特沃斯濾波器階數(shù)選擇

[z，p，k]=buttap(n);

%巴特沃斯濾波器設(shè)計(jì)

figure(1);

[b，a]=zp2tf(z，p，k);

[h，w]=freqz(b，a，512，′whole′);

plot(w，abs(h));

[n1，Wn1]=cheb1ord(Wp，Ws，Rp，Rs，′s′);

%切比雪夫Ⅰ型濾波器階數(shù)選擇

[z1，p1，k1]=cheb1ap(n，Rp);

%切比雪夫Ⅰ型濾波器設(shè)計(jì)

figure(2);

[b，a]=zp2tf(z1，p1，k1);

[h，w]=freqz(b，a，512，′whole′);

plot(w，abs(h));

程序運(yùn)行結(jié)果如圖3.2所示。圖3.2兩種數(shù)字濾波器的幅頻曲線

【例3.2】

設(shè)離散信號為x(n)=s(n)+d(n)，其中，原始無損信號s(n)=2n(0.8)n，d(n)為隨機(jī)噪聲信號，利用MATLAB設(shè)計(jì)一個(gè)滑動平均濾波器。

MATLAB程序如下:

%MATLABPROGRAM3-2

R=50;

d=rand(R，1)－0.5;

m=0:1:R－1;

s=2*m.*(0.9.^m);

x=s+d′;

plot(m，d，′r－′，m，s，′b*′，m，x，′m:′);

xlabel(′時(shí)間序號n′);ylabel(′振幅′);

legend(′r－′，′d[n]′，′b*′，′s[n]′，′m:′，′x[n]′);

pause

M=input(′輸入樣本數(shù)=′);

b=ones(M，1)/M;

y=filter(b，1，x);

plot(m，s，′r－′，m，y，′b*′);

legend(′r－′，′s[n]′，′b*′，′y[n]′);

xlabel(′時(shí)間序號n′);ylabel(′振幅′);

程序運(yùn)行過程中，輸入樣本數(shù)50，程序運(yùn)行結(jié)果如圖3.3所示。圖3.3滑動平均濾波輸出的幅值曲線

3.2IIR數(shù)字濾波器的基本結(jié)構(gòu)

IIR數(shù)字濾波器的結(jié)構(gòu)存在反饋環(huán)路，具有遞歸型結(jié)構(gòu)特點(diǎn)。對于同一系統(tǒng)函數(shù)，IIR數(shù)字濾波器有各種不同的結(jié)構(gòu)形式。由式(3.4)和式(3.5)表示的差分方程或系統(tǒng)函數(shù)可知，其基本結(jié)構(gòu)有直接型、級聯(lián)型和并聯(lián)型三種。

3.2.1直接型結(jié)構(gòu)

直接型結(jié)構(gòu)是直接由IIR濾波器的差分方程所得到的網(wǎng)絡(luò)結(jié)構(gòu)，包括直接Ⅰ型結(jié)構(gòu)、直接Ⅱ型結(jié)構(gòu)以及轉(zhuǎn)置結(jié)構(gòu)。

1.直接Ⅰ型結(jié)構(gòu)

設(shè)，則令有由此可知，H1(z)實(shí)現(xiàn)了系統(tǒng)的零點(diǎn)，H2(z)實(shí)現(xiàn)了系統(tǒng)的極點(diǎn)。H(z)由這兩部分級聯(lián)構(gòu)成。因此，可得到如圖3.4所示的IIR濾波器的直接Ⅰ型結(jié)構(gòu)。圖3.4IIR濾波器的直接Ⅰ型結(jié)構(gòu)直接Ⅰ型結(jié)構(gòu)具有以下缺點(diǎn):

(1)需要2N個(gè)延遲器(z－1)，延遲器太多。

(2)系數(shù)ai、bi對濾波器性能的控制不直接，對零、極點(diǎn)的控制難，一個(gè)ai、bi的改變會影響系統(tǒng)的零點(diǎn)或極點(diǎn)分布。

(3)對字長變化敏感，即對ai、bi的準(zhǔn)確度要求嚴(yán)格。

(4)該結(jié)構(gòu)極易不穩(wěn)定，階數(shù)高時(shí)，上述影響更大。

2.直接Ⅱ型結(jié)構(gòu)

直接Ⅰ型結(jié)構(gòu)可看做兩個(gè)獨(dú)立的網(wǎng)絡(luò)H1(z)和H2(z)兩部分串接構(gòu)成的總的系統(tǒng)函數(shù)。由線性系統(tǒng)函數(shù)的不變性，交換兩個(gè)網(wǎng)絡(luò)次序，可得到如圖3.5所示的直接Ⅱ型結(jié)構(gòu)。直接Ⅱ型結(jié)構(gòu)的延遲線減少一半，為N個(gè)，可節(jié)省寄存器或存儲單元。圖3.5IIR濾波器的直接Ⅱ型結(jié)構(gòu)

3.直接型轉(zhuǎn)置結(jié)構(gòu)

如果將原網(wǎng)絡(luò)中所有支路的方向加以反轉(zhuǎn)，并將輸入和輸出相互交換，則網(wǎng)絡(luò)的系統(tǒng)函數(shù)不會改變，這就是轉(zhuǎn)置定理。由此可得到如圖3.6所示的直接型轉(zhuǎn)置結(jié)構(gòu)。

IIR數(shù)字濾波器的直接型結(jié)構(gòu)簡單直觀，但存在直接Ⅰ型的缺點(diǎn)。對于三階以上的IIR濾波器，幾乎都不采用直接型結(jié)構(gòu)，而是采用級聯(lián)型、并聯(lián)型等其他形式的結(jié)構(gòu)，把高階變成不同組合的低階系統(tǒng)來實(shí)現(xiàn)。

MATLAB信號處理工具箱提供的函數(shù)filter可實(shí)現(xiàn)利用IIR直接型濾波器結(jié)構(gòu)計(jì)算濾波器對輸入的響應(yīng)。

【例3.3】

求巴特沃斯濾波器的單位脈沖響應(yīng)。圖3.6IIR濾波器的轉(zhuǎn)置結(jié)構(gòu)(Ⅰ型)MATLAB程序如下:

%MATLABPROGRAM3-3

x=[1，zeros(1，120)];

[b，a]=butter(12，500)/1000;

y=filter(b，a，x);

stem(y);

grid;

程序運(yùn)行結(jié)果如圖3.7所示。圖3.7巴特沃斯濾波器的單位脈沖響應(yīng)3.2.2級聯(lián)型結(jié)構(gòu)

級聯(lián)型將濾波器系統(tǒng)函數(shù)H(z)的分子和分母分解為一階和二階實(shí)系數(shù)因子之積的形式，即把它的分子、分母都表達(dá)為因子形式:

(3.11)由于系數(shù)ai、bi都是實(shí)數(shù)，零、極點(diǎn)為實(shí)根或共軛復(fù)根，故可將相互共軛的零點(diǎn)(極點(diǎn))合并起來，形成一個(gè)實(shí)系數(shù)的二階多項(xiàng)式:(3.12)式中，gi、pi為實(shí)根;βi、αi為復(fù)根，且N1+2N2=N，M1+2M2=N。

為了簡化級聯(lián)形式，將實(shí)系數(shù)的兩個(gè)一階因子組合成二階因子，則式(3.12)可寫成二階因子的形式:(3.13)式中，M表示(N+1)/2中的最大整數(shù)。Hi(z)稱為二階基本環(huán)節(jié)，可用直接Ⅱ型實(shí)現(xiàn)，H(z)由各Hi(z)級聯(lián)而成。級聯(lián)型結(jié)構(gòu)如圖3.8所示。圖3.8級聯(lián)型IIR結(jié)構(gòu)級聯(lián)型結(jié)構(gòu)的優(yōu)點(diǎn)如下：

(1)簡化硬件實(shí)現(xiàn)，只用一個(gè)二階環(huán)節(jié)，通過變換系數(shù)進(jìn)行時(shí)分復(fù)用就可實(shí)現(xiàn)整個(gè)系統(tǒng);可流水線操作。

(2)零、極點(diǎn)可單獨(dú)控制，且每一個(gè)基本環(huán)節(jié)系數(shù)變化只影響該子系統(tǒng)的零、極點(diǎn)。調(diào)整系數(shù)β1i、β2i可單獨(dú)調(diào)整第i對零點(diǎn)，而不影響其他零、極點(diǎn)。調(diào)整系數(shù)α1i、α2i可單獨(dú)調(diào)整第i對極點(diǎn)，而不影響其他零、極點(diǎn)。

(3)各二階環(huán)節(jié)零、極點(diǎn)的搭配可互換位置，優(yōu)化組合以減小運(yùn)算誤差;對系數(shù)變化的敏感度小，受有限字長的影響比直接型低。同時(shí)，級聯(lián)型結(jié)構(gòu)還存在著二階環(huán)節(jié)電平難控制的缺點(diǎn):電平大，易溢出;電平小則使信噪比減小。

在MATLAB中，信號處理工具箱定義了SOS模型來表示級聯(lián)型結(jié)構(gòu)，若已知系統(tǒng)函數(shù)H(z)，則可以借助內(nèi)部函數(shù)tf2sos將濾波器的直接型結(jié)構(gòu)轉(zhuǎn)換為級聯(lián)型結(jié)構(gòu)。

【例3.4】

濾波器的系統(tǒng)函數(shù)為求該系統(tǒng)的級聯(lián)結(jié)構(gòu)形式和系統(tǒng)的單位脈沖響應(yīng)。MATLAB程序如下:

%MATLABPROGRAM3-4

num=[1，－5，10，－32，19];%系統(tǒng)函數(shù)的分子系數(shù)

den=[15，13，5，－6，－2];

%系統(tǒng)函數(shù)的分母系數(shù)

[z，p，k]=tf2zp(num，den);

%求系統(tǒng)函數(shù)的零、極點(diǎn)和增益

sos=zp2sos(z，p，k);

%直接型結(jié)構(gòu)轉(zhuǎn)換為級聯(lián)型結(jié)構(gòu)

disp(′sos=′);disp(sos);

%顯示級聯(lián)型結(jié)構(gòu)

num1=conv(sos(1，1:3)，sos(2，1:3));

den1=conv(sos(1，4:6)，sos(2，4:6));

x=[1，zeros(1，40)];

n=[0:length(x)－1];

y=filter(num，den，x);

%直接型結(jié)構(gòu)系統(tǒng)的單位脈沖響應(yīng)

y1=filter(num1，den1，x);

%級聯(lián)型結(jié)構(gòu)系統(tǒng)的單位脈沖響應(yīng)

subplot(2，2，1);

plot(n，y);

title(′ImpulseResponse′);

legend(′Original′);

subplot(2，2，2);

plot(n，y1);

title(′ImpulseResponse′);

legend(′Cascadeform′);

程序執(zhí)行結(jié)果為

sos=

0.0667－0.32530.19531.0000－0.2644－0.1668

1.0000－0.11996.48521.00001.13110.7992系統(tǒng)的SOS級聯(lián)型結(jié)構(gòu)為濾波器的單位脈沖響應(yīng)如圖3.9所示。從圖中可以看出，直接型結(jié)構(gòu)和級聯(lián)型結(jié)構(gòu)的單位脈沖響應(yīng)相同，這說明結(jié)構(gòu)形式的改變不會改變?yōu)V波器的時(shí)域特性。圖3.9濾波器不同結(jié)構(gòu)的單位脈沖響應(yīng)3.2.3并聯(lián)型結(jié)構(gòu)

將系統(tǒng)函數(shù)展開成部分分式之和，可用并聯(lián)方式構(gòu)成濾波器:式中，di是極點(diǎn)，Ai為部分分式展開的系數(shù)，常數(shù)γ0=bN/aN。式(3.14)所蘊(yùn)含的結(jié)構(gòu)如圖3.10所示，它包含了一組單極點(diǎn)濾波器。圖3.10IIR系統(tǒng)的并聯(lián)結(jié)構(gòu)為避免復(fù)數(shù)乘法，可以組合復(fù)共軛極點(diǎn)對形成雙極點(diǎn)子系統(tǒng)。另外，可將實(shí)值極點(diǎn)對形成雙極點(diǎn)子系統(tǒng)，得(3.15)式中，L為(N+1)/2的整數(shù)部分。當(dāng)N為奇數(shù)時(shí)，其中一個(gè)Hk(z)是一個(gè)單極點(diǎn)系統(tǒng)。并聯(lián)型結(jié)構(gòu)的特點(diǎn)如下:

(1)系統(tǒng)實(shí)現(xiàn)簡單，只需一個(gè)二階環(huán)節(jié)，系統(tǒng)通過改變輸入系數(shù)即可完成。

(2)極點(diǎn)位置可單獨(dú)調(diào)整，但不能像級聯(lián)型那樣單獨(dú)調(diào)整零點(diǎn)的位置，因?yàn)椴⒙?lián)型各子系統(tǒng)的零點(diǎn)并非整個(gè)系統(tǒng)函數(shù)的零點(diǎn)。當(dāng)需要準(zhǔn)確的傳輸零點(diǎn)時(shí)，級聯(lián)型最合適。

(3)由于基本環(huán)節(jié)并聯(lián)，故可同時(shí)對輸入信號進(jìn)行運(yùn)算，因此該結(jié)構(gòu)運(yùn)算速度快。

(4)各并聯(lián)二階網(wǎng)絡(luò)的誤差互不影響，運(yùn)算誤差最小，對字長要求低。

在MATLAB中，信號處理工具箱沒有直接提供該結(jié)構(gòu)的生成信號，但借助函數(shù)residuez可以獲得IIR濾波器并聯(lián)型結(jié)構(gòu)的參數(shù):先用residuez求出H(z)的部分分式展開式，然后把分式兩兩合并成式(3.15)的形式。若分式項(xiàng)為偶數(shù)，則并聯(lián)型結(jié)構(gòu)全為二階環(huán)節(jié);若分式項(xiàng)為奇數(shù)，則并聯(lián)型結(jié)構(gòu)除二階環(huán)節(jié)外，還包括一階環(huán)節(jié)。

【例3.5】

濾波器的系統(tǒng)函數(shù)為

求該系統(tǒng)的級聯(lián)結(jié)構(gòu)形式和系統(tǒng)的單位脈沖響應(yīng)。

MATLAB程序如下:

%MATLABPROGRAM3-5

clc;

num=[1，-3，11，-28，19];

den=[16，13，2，-5，-1];

[r，p，c]=residuez(num，den);

n1=length(r);

ifrem(n1，2)==0N=n1;

elseN=n1-1;

end

num1=zeros(N/2，2);

den1=zeros(N/2，2);

fori=1:N/2

num1(i，:)=r(i*2-1:i*2)′;

den1(i，:)=p(i*2-1:i*2)′;

end

fori=1:N/2

den2(i，:)=conv(poly(den1(i，1))，poly(den1(i，2)));

num2(i，:)=num1(i，1)*poly(den1(i，2))+num1(i，2)*poly(den1(i，1));

end

ifrem(n1，2)~=0;

nn1=n1-N/2;

fori=1:nn1

den2(i+N/2，:)=[poly(p(N+i))0];

num2(i+N/2，:)=[r(N+i)0];

end

end

disp(′并聯(lián)結(jié)構(gòu)的整體形式′)

sys1=filt(num2(1，:)，den2(1，:))

sys2=filt(num2(2，:)，den2(2，:))程序執(zhí)行結(jié)果為

并聯(lián)結(jié)構(gòu)的整體形式:

Transferfunction:

－8.124－4.612z^－1

----------------------------

1+1.136z^－1+0.5978z^－2

Samplingtime:unspecified

Transferfunction:

27.19－13.67z^－1

----------------------------

1－0.324z^－1－0.1046z^－2

Samplingtime:unspecified

c=

－19

3.3FIR數(shù)字濾波器的基本結(jié)構(gòu)

FIR數(shù)字濾波器主要是非遞歸結(jié)構(gòu)，無反饋，但在頻率采樣結(jié)構(gòu)等某些結(jié)構(gòu)中也包含有反饋的遞歸部分。FIR的基本結(jié)構(gòu)包括直接型、級聯(lián)型、線性相位型以及頻率采樣型。

3.3.1直接型結(jié)構(gòu)

根據(jù)FIR濾波器的差分方程式(3.6)和系統(tǒng)函數(shù)式(3.7)，M階FIR濾波器的直接型結(jié)構(gòu)如圖3.11所示。該結(jié)構(gòu)由M+1個(gè)乘法器、M個(gè)延遲器以及M個(gè)加法器組成。

直接型也叫橫向?yàn)V波器、卷積型結(jié)構(gòu),FIR濾波器的直接型結(jié)構(gòu)與IIR相似，可以用MATLAB的內(nèi)部函數(shù)filter來實(shí)現(xiàn)。圖3.11FIR濾波器的直接型結(jié)構(gòu)3.3.2級聯(lián)型結(jié)構(gòu)

當(dāng)需要控制濾波器的傳輸零點(diǎn)時(shí)，可將系統(tǒng)函數(shù)H(z)分解為若干個(gè)一階實(shí)系數(shù)因子和二階實(shí)系數(shù)因子相乘的形式:(3.16)

由上式可得到FIR濾波器的級聯(lián)型結(jié)構(gòu)，如圖3.12所示。圖3.12FIR濾波器的級聯(lián)型結(jié)構(gòu)級聯(lián)型結(jié)構(gòu)每一個(gè)一階因子控制一個(gè)實(shí)數(shù)零點(diǎn)，每一個(gè)二階因子控制一對共軛零點(diǎn)。調(diào)整零點(diǎn)位置比直接型方便，但是它所需要的系數(shù)a比直接型h(n)多，因而需要的乘法器也多。

FIR濾波器的級聯(lián)型結(jié)構(gòu)可用MATLAB編程實(shí)現(xiàn)，具體與IIR濾波器類似。3.3.3線性相位型結(jié)構(gòu)

若FIR濾波器的單位脈沖響應(yīng)h(n)為實(shí)數(shù)，且滿足h(n)=±h(N－1－n)，則稱該濾波器具有線性相位。FIR濾波器的重要特點(diǎn)是，可設(shè)計(jì)成具有嚴(yán)格線性相位的濾波器，此時(shí)h(n)滿足偶對稱或奇對稱條件，但無論h(n)為偶對稱還是奇對稱，都有以下性質(zhì):

(1)N為偶數(shù)，系統(tǒng)函數(shù)為

(3.17)(2)N為奇數(shù)，系統(tǒng)函數(shù)為(3.18)由以上兩式可得到線性相位FIR濾波器的線性相位型結(jié)構(gòu)，如圖3.13所示。圖3.13FIR濾波器的線性相位型結(jié)構(gòu)線性相位型結(jié)構(gòu)本質(zhì)上是直接型的，但乘法次數(shù)較少。當(dāng)N為偶數(shù)時(shí)，乘法次數(shù)減為N/2;當(dāng)N為奇數(shù)時(shí)，乘法次數(shù)減為(N+1)/2。對于線性相位FIR濾波器來說，如果零點(diǎn)在單位圓上，則零點(diǎn)以共軛對出現(xiàn);如果零點(diǎn)在實(shí)軸上，則零點(diǎn)互為倒數(shù)出現(xiàn);如果零點(diǎn)既在單位圓上又在實(shí)軸上，則此時(shí)只有一個(gè)零點(diǎn)。

FIR濾波器的線性相位型結(jié)構(gòu)也可用MATLAB編程實(shí)現(xiàn)，具體與IIR濾波器的直接型結(jié)構(gòu)幾乎相同。3.3.4頻率采樣型結(jié)構(gòu)

若在頻率域進(jìn)行等間隔采樣，則相應(yīng)的時(shí)域信號會以采樣點(diǎn)數(shù)為周期進(jìn)行周期性延拓。如果在頻率域的采樣點(diǎn)數(shù)N大于等于原序列的長度M，則不會引起信號失真。

若h(n)是長為N的序列，對系統(tǒng)函數(shù)H(z)在單位圓上作N等份采樣，這個(gè)采樣值就是h(n)的離散傅里葉變換H(k)。原序列的Z變換H(z)與頻域采樣值H(k)滿足下面的關(guān)系式:

(3.19)式(3.19)為FIR濾波器提供了一種稱為頻率采樣的結(jié)構(gòu)。利用頻率采樣的內(nèi)插公式，可得(3.20)式中，

由此可知，H(z)由梳狀濾波器Hc(z)和N個(gè)一階網(wǎng)絡(luò)Hk(z)兩部分級聯(lián)而成。第一部分H(z)是一個(gè)由N節(jié)延時(shí)器組成的梳狀濾波器，它在單位圓上有N個(gè)等份的零點(diǎn):(3.21)其頻率響應(yīng)為(3.22)其幅度響應(yīng)為(3.23)第二部分由N個(gè)一階網(wǎng)絡(luò)Hk(z)并聯(lián)而成，每個(gè)一階網(wǎng)絡(luò)在單位圓上有一個(gè)極點(diǎn):(3.24)該網(wǎng)絡(luò)在處的頻率響應(yīng)為∞，是一個(gè)諧振頻率為的諧振器。這些并聯(lián)諧振器的極點(diǎn)正好各自抵消一個(gè)梳狀濾波器的零點(diǎn)，從而使這個(gè)頻率點(diǎn)的響應(yīng)等于H(k)。兩部分級聯(lián)后，就得到FIR濾波器的頻率采樣型結(jié)構(gòu)，如圖3.14所示。圖3.14FIR濾波器的頻率采樣型結(jié)構(gòu)這一結(jié)構(gòu)的最大特點(diǎn)是它的系數(shù)H(k)直接就是濾波器在ω=處的響應(yīng)。因此，控制濾波器的響應(yīng)很直接。只要h(n)的長度N相同，則對于任何頻響形狀，其梳狀濾波器部分和N個(gè)一階網(wǎng)絡(luò)部分的結(jié)構(gòu)就完全相同，只是各支路增益H(k)不同。這樣，相同部分便于標(biāo)準(zhǔn)化和模塊化。

然而，上述頻率采樣型結(jié)構(gòu)亦有如下兩個(gè)缺點(diǎn):

(1)所有的系數(shù)WN－k和H(k)都是復(fù)數(shù)，而硬件乘法器要完成復(fù)數(shù)乘法運(yùn)算，需采用更加復(fù)雜的硬件結(jié)構(gòu)。

(2)所有諧振器的極點(diǎn)都在單位圓上，考慮到系數(shù)量化的影響，有些極點(diǎn)實(shí)際上不能與梳狀濾波器的零點(diǎn)相抵消，故使系統(tǒng)的穩(wěn)定性變差。為了克服上述兩個(gè)缺點(diǎn)，對頻率采樣型結(jié)構(gòu)需進(jìn)行如下修正:

(1)將單位圓上所有的零、極點(diǎn)向單位圓內(nèi)收縮，收縮到半徑r略小于1的圓上。當(dāng)由于某種原因，零、極點(diǎn)不能抵消時(shí)，極點(diǎn)位置可仍在單位圓內(nèi)，保持了系統(tǒng)穩(wěn)定。H(z)變?yōu)?3.25)

(2)共軛根合并。將一對復(fù)數(shù)一階網(wǎng)絡(luò)合并成一個(gè)實(shí)系數(shù)的二階網(wǎng)絡(luò)，從而將復(fù)數(shù)乘法運(yùn)算變成了實(shí)數(shù)運(yùn)算，這些共軛根在圓周上對稱。因h(m)是實(shí)數(shù)，故其DFT也是圓周共軛對稱。即

H(N－k)=H*(k)

(3.26)因此，將第k及第(N－k)個(gè)諧振器合并為一個(gè)二階網(wǎng)絡(luò)Hk(z):(3.27)式中，該二階網(wǎng)絡(luò)是一個(gè)有限Q值的諧振器，諧振頻率。修正后的頻率采樣型結(jié)構(gòu)如圖3.15所示。圖3.15修正后的頻率采樣型結(jié)構(gòu)除以上共軛極點(diǎn)外，還有實(shí)數(shù)極點(diǎn)，分如下兩種情況:

(1)當(dāng)N為偶數(shù)時(shí)，有兩個(gè)實(shí)數(shù)極點(diǎn)z=±r，對應(yīng)有兩個(gè)一階網(wǎng)絡(luò)，系統(tǒng)函數(shù)為

(3.28)式中，對應(yīng)的兩個(gè)一階網(wǎng)絡(luò)分別為

(2)當(dāng)N為奇數(shù)時(shí)，只有一個(gè)實(shí)數(shù)極點(diǎn)z=r，對應(yīng)有一個(gè)一階網(wǎng)絡(luò)，系統(tǒng)函數(shù)為(3.29)式中，對應(yīng)的一階網(wǎng)絡(luò)為

頻率采樣型結(jié)構(gòu)的選頻性好，適于窄帶濾波，大部分H(k)為0，且只有較少的二階網(wǎng)絡(luò);不同的FIR濾波器，若長度相同，則可通過改變系數(shù)用同一個(gè)網(wǎng)絡(luò)實(shí)現(xiàn)，因此復(fù)用性好。但頻率采樣型的結(jié)構(gòu)復(fù)雜，因而采用的存儲器也多。頻率采樣型結(jié)構(gòu)適合于任何FIR系統(tǒng)函數(shù)。

3.4數(shù)字濾波器的格型結(jié)構(gòu)

格型(Lattice)結(jié)構(gòu)是數(shù)字濾波器的另一種實(shí)現(xiàn)結(jié)構(gòu)，它不僅適用于IIR濾波器，也適用于FIR濾波器。格型濾波器廣泛應(yīng)用于語音處理和自適應(yīng)濾波器實(shí)現(xiàn)中。

格型結(jié)構(gòu)可分為全零點(diǎn)(AZ)濾波器、全極點(diǎn)(AP)濾波器和零極點(diǎn)(AZAP)濾波器三類。

3.4.1全零點(diǎn)濾波器的格型結(jié)構(gòu)

全零點(diǎn)濾波器的系統(tǒng)函數(shù)為

(3.30)

全零點(diǎn)濾波器的格型結(jié)構(gòu)如圖3.16所示。圖3.16全零點(diǎn)濾波器的格型結(jié)構(gòu)根據(jù)系統(tǒng)函數(shù)，由高階系數(shù)遞推各低階反射系數(shù)Kp，即(3.31)

MATLAB信號處理工具箱提供了函數(shù)poly2rc來實(shí)現(xiàn)從FIR濾波器直接型結(jié)構(gòu)到格型結(jié)構(gòu)的計(jì)算，即在已知FIR濾波器的全零點(diǎn)多項(xiàng)式模型的條件下，計(jì)算FIR濾波器格型結(jié)構(gòu)的反射系數(shù)Kp。函數(shù)調(diào)用格式為

K=poly2rc(b)式中，b為FIR濾波器多項(xiàng)式系數(shù)實(shí)向量，b(1)不能為0;K為反射系數(shù)。

MATLAB還提供了FIR格型結(jié)構(gòu)反射系數(shù)Kp求濾波器多項(xiàng)式系數(shù)的函數(shù)rc2poly，其調(diào)用格式為

b=rc2poly(K)

【例3.6】

已知FIR濾波器的差分方程為

y(n)=x(n)+0.5984x(n－1)+0.8945x(n－2)－0.0561x(n－3)+0.0132x(x－5)

求該系統(tǒng)的格型結(jié)構(gòu)。

MATLAB程序如下:

%MATLABPROGRAM3-6

b=[1，0.5984，0.8945，－0.0561，0，0.0132];

k=poly2rc(b);

disp(′k=′);disp(k′);

程序運(yùn)行結(jié)果為

k=

0.33830.9441－0.0632－0.00790.01323.4.2全極點(diǎn)濾波器的格型結(jié)構(gòu)

全極點(diǎn)濾波器的系統(tǒng)函數(shù)為(3.32)

全極點(diǎn)濾波器的格型結(jié)構(gòu)如圖3.17所示。圖3.17全極點(diǎn)濾波器的格型結(jié)構(gòu)圖中，Kp為反射系數(shù)。

MATLAB信號處理工具箱提供了函數(shù)tf2latc，以實(shí)現(xiàn)從系統(tǒng)函數(shù)分母多項(xiàng)式a求反射系數(shù)K。該函數(shù)也可用于生成IIR傳輸函數(shù)的格型結(jié)構(gòu)。3.4.3零極點(diǎn)濾波器的格型結(jié)構(gòu)

零極點(diǎn)濾波器的系統(tǒng)函數(shù)為(3.33)

含極點(diǎn)和零點(diǎn)濾波器的格型結(jié)構(gòu)如圖3.18所示。圖3.18零極點(diǎn)濾波器的格型結(jié)構(gòu)圖中的方框?yàn)閳D3.19所示的基本格型單元。圖3.19零極點(diǎn)濾波器的基本格型單元在零極點(diǎn)濾波器的格型結(jié)構(gòu)中，Ki為濾波器格型系數(shù)，即反射系數(shù)；ci為濾波器梯形系數(shù)。在此種情況下，格型結(jié)構(gòu)中K、C參數(shù)可按下列步驟進(jìn)行確定。

(1)利用AZ系統(tǒng)反射系數(shù)Kp的遞推公式，遞推出K參數(shù)；

(2)確定cp

：cp=bp；

(3)遞推求出C參數(shù)。(3.34)

MATLAB信號處理工具箱提供了函數(shù)tf2latc，以實(shí)現(xiàn)一般IIR濾波器系統(tǒng)函數(shù)H(z)求濾波器的格型結(jié)構(gòu)參數(shù)K和C。該函數(shù)的調(diào)用格式為

[K，C]=tf2latc(num，den)

K=tf2latc(1，den)

[K，C]=tf2latc(1，den)

K=tf2latc(num)

其中，num和den分別是以z的降冪排列的系統(tǒng)函數(shù)的分子多項(xiàng)式及分母多項(xiàng)式的系數(shù)向量;K為濾波器格型系數(shù)(反射系數(shù));C為濾波器梯形系數(shù)。

第一種調(diào)用格式可獲得零極點(diǎn)格型結(jié)構(gòu)濾波器;第二、三種可獲得全極點(diǎn)格型結(jié)構(gòu)濾波器;第四種調(diào)用格式可獲得全零點(diǎn)格型結(jié)構(gòu)濾波器。因此，該函數(shù)也適用于FIR濾波器。

MATLAB信號處理工具箱還提供了由格型結(jié)構(gòu)濾波器求濾波器系統(tǒng)函數(shù)的latc2tf函數(shù)。該函數(shù)的調(diào)用格式為

[num，den]=latc2tf(K，C)

[num，den]=latc2tf(K，′iir′)

num=latc2tf(K，′fir′)

num=latc2tf(K)

其中，iir表示K是全極點(diǎn)IIR格型濾波器的反射系數(shù);fir表示K是全零點(diǎn)FIR格型濾波器的反射系數(shù)。

【例3.7】

求下面系統(tǒng)的零極點(diǎn)IIR濾波器的格型結(jié)構(gòu)參數(shù)，系統(tǒng)函數(shù)為

MATLAB程序如下:

%MATLABPROGRAM3-7

num=[1，3，5，-2];

den=[1，0.665，0.875，－.468];

[K，C]=tf2latc(num，den);

disp(′格型結(jié)構(gòu)反射系數(shù)K＝′);

disp(K′);

disp(′格型結(jié)構(gòu)梯形系數(shù)C＝′);

disp(C′);

程序運(yùn)行結(jié)果為

格型結(jié)構(gòu)反射系數(shù)K＝ 0.5462

1.5189－0.4680

格型結(jié)構(gòu)梯形系數(shù)C＝ －7.3881－3.9591

6.3300－2.0000

3.5數(shù)字濾波器的MATLAB時(shí)頻分析

3.5.1時(shí)域分析

在MATLAB中，求解離散時(shí)間系統(tǒng)時(shí)域響應(yīng)的方法，同樣適用于數(shù)字濾波器。這里重點(diǎn)介紹常用的濾波器時(shí)間響應(yīng)分析工具函數(shù)filter、filtic、fftfilt及impz。

1.函數(shù)filter

函數(shù)filter用于實(shí)現(xiàn)IIR和FIR濾波器對數(shù)據(jù)濾波，并計(jì)算濾波器對輸入的響應(yīng)。具體用法參見第2章相關(guān)內(nèi)容。該函數(shù)可用于全極點(diǎn)濾波器、全零點(diǎn)濾波器以及零極點(diǎn)濾波器。

2.函數(shù)filtic

函數(shù)filtic用于從濾波器過去值y和x求濾波器狀態(tài)的初始值。調(diào)用格式為

z=filtic(num，den，y，x)

z=filtic(num，den，y)

其中，y=[y(－1)，y(－2)，…，y(－na)]為輸出y的過去值向量，x=[x(－1)，x(－2)，…，x(－nb)]為輸入x的過去值向量，nb=length(num)－1，na=length(den)－1;num，den分別為濾波器分子、分母多項(xiàng)式系數(shù)向量;z為濾波器初始狀態(tài)。

函數(shù)filter設(shè)立狀態(tài)初始值zi和zf，可將特別長的信號“切斷”來分段濾波處理，再把它們銜接起來而不影響全信號的濾波效果。

【例3.8】

設(shè)濾波器系統(tǒng)函數(shù)為對隨機(jī)長信號濾波，試?yán)肕ATLAB編程驗(yàn)證分段濾波的效果。

MATLAB程序如下:

%MATLABPROGRAM3-8

num=[1，0.9，1.2];

den=[0.8，0.3，－0.6];

x=rand(60，1);

y=filter(num，den，x);

x1=x(1:15);

x2=x(16:30);

x3=x(31:45);

x4=x(46:60);

[y1，Zf1]=filter(num，den，x1);Zi2=Zf1;

[y2，Zf2]=filter(num，den，x2，Zi2);Zi3=Zf2;

[y3，Zf3]=filter(num，den，x3，Zi3);Zi4=Zf3;

[y4，Zf4]=filter(num，den，x4，Zi4);

yd=[y1;y2;y3;y4];

n=[0:length(y)-1];

plot(n，y，′r－′，n，yd，′b*′);

legend(′r－′，′y[n]′，′b*′，′yd[n]′);

xlabel(′n′);title(′兩種濾波曲線′);

程序運(yùn)行結(jié)果如圖3.20所示。圖3.20信號分段濾波輸出對比圖由圖可知，信號分段濾波和全程濾波完全吻合。分段濾波必須采用函數(shù)filter。

3.函數(shù)fftfilt

函數(shù)fftfilt基于FFT重疊相加算法來實(shí)現(xiàn)對數(shù)據(jù)的濾波，該函數(shù)只適用于FIR濾波器。該函數(shù)的調(diào)用格式為

y=fftfilt(num，x)

y=ftfilt(num，x，n)

其中，num為濾波器的系數(shù)向量;x為輸入序列;n為FFT長度，缺省時(shí)，選擇最佳的FFT長度;y為濾波器輸出。

【例3.9】FIR低通濾波器的截止頻率為250Hz，采樣頻率fs為1kHz，對信號x(t)進(jìn)行濾波，x(t)=2sin(2πf1t)+

3cos(2πf2t+π/3)，f1為50Hz，f2為300Hz，求濾波輸出。

MATLAB程序如下:

%MATLABPROGRAM3-9

clc;

N=1000;

Fs=1000;

n=[0:N－1];

t=n/Fs;

x=2*sin(2*pi*50*t)+3*cos(2*pi*300*t+pi/3);

num=fir1(40，250/500);%產(chǎn)生截止頻率為250Hz的FIR低通濾波器

yfft=fftfilt(num，x，256);%fft重疊相加算法濾波

n1=[81:161];

t1=t(n1);

x1=x(n1);

y1=yfft(n1);

subplot(2，2，1);

plot(t1，x1);title(′原始信號x′);

subplot(2，2，2);

plot(t1，y1);title(′濾波輸出信號y′);

grid

程序運(yùn)行結(jié)果如圖3.21所示，通過低頻濾波器后的信號只包含50Hz頻率的正弦波。圖3.21fftfilt實(shí)現(xiàn)低頻濾波

4.函數(shù)impz

函數(shù)impz用于產(chǎn)生數(shù)字濾波器的脈沖響應(yīng)。其調(diào)用格式為

[h，t]=impz(num，den)

[h，t]=impz(num，den，n)

[h，t]=impz(num，den，n，F(xiàn)s)

impz(num，den)

其中，num和den分別為系統(tǒng)函數(shù)的多項(xiàng)式分子系數(shù)和分母系數(shù);n為采樣點(diǎn)數(shù)，缺省時(shí)自動選擇;Fs為采樣頻率，默認(rèn)為1;h為濾波器單位脈沖響應(yīng)，缺省時(shí)，繪制濾波器脈沖響應(yīng)曲線;t為h對應(yīng)的時(shí)間向量。

【例3.10】

設(shè)計(jì)一個(gè)8階的巴特沃斯帶通濾波器，通帶為120~280Hz，采樣頻率為1kHz，繪制濾波器的單位脈沖響應(yīng)曲線。

MATLAB程序如下:

%MATLABPROGRAM3-10

N=101;

Fs=1000;

[num，den]=butter(8，2*[120250]/Fs);

Impz(num，den，N，F(xiàn)s);

title(′ImpulseResponse′);

grid

程序運(yùn)行結(jié)果如圖3.22所示。圖3.22脈沖響應(yīng)曲線3.5.2頻域分析

頻域分析主要包括幅頻響應(yīng)、相頻響應(yīng)和零極點(diǎn)位置等。

1.濾波器的頻率響應(yīng)

MATLAB提供了函數(shù)freqs、freqz等工具來分析濾波器的頻率響應(yīng)。

1)函數(shù)freqs

函數(shù)freqs用于求模擬濾波器的頻率響應(yīng)。其調(diào)用格式為

h=freqs(b，a，w)

[h，w]=freqs(b，a)

[h，w]=freqs(b，a，n)

freqs(b，a)

其中，b，a分別為模擬濾波器傳遞函數(shù)的分子和分母系數(shù)向量;n為頻率點(diǎn)數(shù)，缺省值為200;h為頻率響應(yīng)，以復(fù)數(shù)表示，當(dāng)輸出缺省時(shí)，繪制模擬濾波器的幅頻響應(yīng)和相頻響應(yīng)曲線;w為頻率向量，以實(shí)數(shù)表示。

MATLAB信號處理模塊還提供了兩個(gè)函數(shù)abs(x)和angle(x)，以求取頻率響應(yīng)H(ejω)的幅頻響應(yīng)|H(ejω)|和相頻響應(yīng)∠H(ejω)，單位為弧度。

2)函數(shù)freqz

函數(shù)freqz用于求數(shù)字濾波器的頻率響應(yīng)。此函數(shù)已在第1章介紹了部分用法，其他調(diào)用格式為

[h，w]=freqz(num，den，n，F(xiàn)s)

[h，w]=freqz(num，den，n，′whole′，F(xiàn)s)

h=freqz(num，den，f，F(xiàn)s)

freqz(num，den)

其中，num和den分別為多項(xiàng)式分子系數(shù)和分母系數(shù);n為復(fù)頻率響應(yīng)的計(jì)算點(diǎn)數(shù)，默認(rèn)為512;Fs為采樣頻率;h為復(fù)頻率響應(yīng);w為n點(diǎn)頻率矢量(弧度);′whole′表示返回的w值為整個(gè)單位圓上的N點(diǎn)等間距的頻率矢量;當(dāng)函數(shù)freqz的返回值缺省時(shí)，繪制頻率響應(yīng)曲線。

3)函數(shù)unwrap

函數(shù)freqz輸出的頻率向量ω在0～2π之間。為了獲得一個(gè)濾波器真正的相頻特性曲線，需對相位角ω進(jìn)行修正。函數(shù)unwrap用于展開ω，其調(diào)用格式為

P=unwrap(w)

【例3.11】設(shè)計(jì)一個(gè)20階的FIR低通濾波器，截止頻率為0.5，求濾波器的頻率響應(yīng)。

MATLAB程序如下:

%MATLABPROGRAM3-11

clc;

num=fir1(20，0.4);%產(chǎn)生截止頻率為0.5的FIR低通濾波器

[h，w]=freqz(num，1，512);

magH=20*log10(abs(h));

phaH1=unwrap(angle(h));

phaH=180*phaH1/pi;

wnyq=w/pi;

subplot(2，1，1);

plot(wnyq，magH);

xlabel(′頻率f′);ylabel(′幅值(dB)′);grid;

subplot(2，1，2);

plot(wnyq，phaH);

xlabel(′頻率f′);ylabel(′相位(degree)′);grid;

程序運(yùn)行結(jié)果如圖3.23所示。圖3.23數(shù)字濾波器的頻率響應(yīng)

2.零極點(diǎn)圖

濾波器的零、極點(diǎn)位置決定了濾波器的穩(wěn)定性和性能，MATLAB信號處理模塊提供有繪制離散時(shí)間系統(tǒng)零極點(diǎn)圖的工具函數(shù)zplane，其調(diào)用格式為

zplane(z，p)

zplane(num，den)

其中，z和p分別為系統(tǒng)的零、極點(diǎn)向量;num和den分別為系統(tǒng)函數(shù)的分子和分母系數(shù)向量。

3.6數(shù)字濾波器的量化效應(yīng)

3.6.1二進(jìn)制數(shù)的表示

1.定點(diǎn)表示

整個(gè)運(yùn)算中，若小數(shù)點(diǎn)在數(shù)碼中的位置固定不變，則稱為定點(diǎn)制。定點(diǎn)制總把數(shù)限制在±1之間，最高位為符號位，0為正，1為負(fù)，小數(shù)點(diǎn)緊跟在符號位之后。若數(shù)的本身只有小數(shù)部分，則稱為“尾數(shù)”;定點(diǎn)數(shù)作加減法時(shí)結(jié)果可能會超出±1，這稱為“溢出”;乘法運(yùn)算不溢出，但字長要增加一倍。為保證字長不變，乘法后，一般要對增加的尾數(shù)作截尾或舍入處理，這將帶來誤差。定點(diǎn)數(shù)的缺點(diǎn)是動態(tài)范圍小，有溢出。定點(diǎn)數(shù)的表示分為原碼、反碼及補(bǔ)碼三種。

若設(shè)(b+1)位碼定點(diǎn)數(shù)為β0β1β2…βb，則

(1)原碼表示:

(3.35)

(2)反碼表示:正數(shù)同原碼，負(fù)數(shù)則將原碼中的尾數(shù)按位求反。即(3.36)

(3)補(bǔ)碼表示:正數(shù)同原碼，負(fù)數(shù)則將原碼中的尾數(shù)求反加1。即(3.37)

2.浮點(diǎn)表示

浮點(diǎn)表示為(3.38)式中，M是尾數(shù)，c為階數(shù)。浮點(diǎn)制運(yùn)算的優(yōu)點(diǎn)是動態(tài)范圍大，一般不溢出，但存在當(dāng)執(zhí)行相乘或相加運(yùn)算時(shí)要對尾數(shù)作量化處理的缺點(diǎn)。一般情況下，浮點(diǎn)數(shù)都用較長的字長，精度較高，因此討論誤差影響主要針對的是定點(diǎn)制。3.6.2定點(diǎn)制的量化誤差

定點(diǎn)制中的乘法，運(yùn)算完畢后會使字長增加。假定定點(diǎn)數(shù)是b位字長,則運(yùn)算后將增長到b1位,需對尾數(shù)作量化處理使b1位字長降低到b位。量化處理包括截尾處理和舍入處理。

1.截尾處理

截尾處理是保留b位，拋棄余下的尾數(shù)。用[·]T表示截尾處理。量化寬度或量化階q=2－b，代表b位字長可表示的最小數(shù)。

1)正數(shù)

一個(gè)b1位的正數(shù)的截尾誤差為

(3.39)－q≤ET≤0

2)負(fù)數(shù)(β0=1)

負(fù)數(shù)的三種碼表示方式不同，所以誤差也不同。

原碼的截尾誤差為(3.40)0≤ET≤q

補(bǔ)碼的截尾誤差為(3.41)

反碼的截尾誤差為(3.42)

補(bǔ)碼的截尾誤差均是負(fù)值，原碼、反碼的截尾誤差取決于數(shù)的正負(fù)，正數(shù)時(shí)為負(fù)，負(fù)數(shù)時(shí)為正。

【例3.12】

試用MATLAB編寫函數(shù)a2dT，其功能是將原碼十進(jìn)制數(shù)轉(zhuǎn)化為二進(jìn)制數(shù)，對于這個(gè)二進(jìn)制數(shù)，用給定位數(shù)進(jìn)行截尾處理得到小數(shù)部分，并將結(jié)果表示為十進(jìn)制。

MATLAB程序如下:

%MATLABPROGRAM3-12

functionbeq=a2dT(d，n)

%產(chǎn)生一個(gè)十進(jìn)制向量d的二進(jìn)制表示的十進(jìn)制等效beq，其中n為截尾處理位數(shù)

m=1;d1=abs(d);

whilefix(d1)>0

d1=abs(d)/(10^m);

m=m+1;

end

beq=0;

fork=[1:n]

beq=fix(d1*2)/(2^k)+beq;

d1=(d1*2)-fix(d1*2);

end

beq=sign(d).*beq*10^(m-1);

2.舍入處理

舍入處理就是按最接近的值取b位碼，通過b+1位上加1后作截尾處理實(shí)現(xiàn)。該法即為通常的四舍五入法，按最接近的數(shù)取量化，所以不論正數(shù)、負(fù)數(shù)，還是原碼、補(bǔ)碼、反碼，誤差總是在±q/2之間，以[x]R表示對x作舍入處理。舍入處理的誤差比截尾處理的誤差小，所以對信號進(jìn)行量化時(shí)多用舍入處理。

【例3.13】

試用MATLAB編寫函數(shù)a2dR，其功能是將原碼十進(jìn)制數(shù)轉(zhuǎn)化為二進(jìn)制數(shù)，對于這個(gè)二進(jìn)制數(shù)，用給定位數(shù)進(jìn)行舍入處理得到小數(shù)部分，并將結(jié)果表示為十進(jìn)制。

MATLAB程序如下:

%MATLABPROGRAM3-13

functionbeq=a2dR(d，n)

%產(chǎn)生一個(gè)十進(jìn)制向量d的二進(jìn)制表示的十進(jìn)制等效beq，其中n為舍入處理位數(shù)

m=1;d1=abs(d);

whilefix(d1)>0

d1=abs(d)/(10^m);

m=m+1;

end

beq=0;d1=d1+2^(－n－1);

fork=[1:n]

beq=fix(d1*2)/(2^k)+beq;

d1=(d1*2)－fix(d1*2);

end

beq=sign(d).*beq*10^(m－1);3.6.3A/D變換的量化效應(yīng)

A/D變換器分為采樣和數(shù)字編碼兩部分。數(shù)字編碼對采樣序列作舍入或截尾處理，得到有限字長的數(shù)字信號。本節(jié)將討論這一過程中的量化效應(yīng)。

1.量化噪聲

對一個(gè)采樣數(shù)據(jù)x(n)作截尾和舍入處理，則截尾量化誤差為

(3.43)

舍入量化誤差為(3.44)式(3.43)和式(3.44)給出了量化誤差的范圍，但要精確知道誤差的大小卻很困難。一般總是通過分析量化噪聲的統(tǒng)計(jì)特性來描述量化誤差。對量化誤差e(n)的統(tǒng)計(jì)特性作如下假定:

(1)e(n)是平穩(wěn)隨機(jī)序列，具有均勻等概率分布;

(2)e(n)與信號x(n)不相關(guān);

(3)e(n)序列中的任意兩個(gè)值不相關(guān)，即為白噪聲。

由上述假定知，量化誤差是一個(gè)與信號序列完全不相關(guān)的白噪聲序列，稱為量化噪聲。

截尾量化噪聲的均值和方差分別為

(3.45)舍入量化噪聲的均值和方差分別為(3.46)由此可見，量化噪聲的方差與A/D變換的字長直接有關(guān)，字長越長，量化噪聲越小。

2.量化信噪比

量化信噪比定義為

(3.47)式中，σx2為輸入信號的方差，也即輸入信號的功率;σ2e為噪聲的方差，也即量化噪聲功率。量化信噪比用對數(shù)表示為(3.48)由此可知，字長每增加1位，量化信噪比增加6個(gè)分貝;信號能量越大，量化信噪比越高。因信號本身有一定的信噪比，故單純提高量化信噪比并無意義。

【例3.14】

已知噪聲在－1至1之間均勻分布，求b=8位和b=12位時(shí)A/D的SNR。

解由于噪聲均勻分布，故均值和方差分別為當(dāng)b=8位時(shí)，SNR=54dB;當(dāng)b=12位時(shí)，SNR=78dB。

3.量化噪聲通過線性系統(tǒng)

為了單獨(dú)分析量化噪聲通過系統(tǒng)后的影響，可將系統(tǒng)近似看做完全理想的(即具有無限精度的線性系統(tǒng))。在輸入端線性相加的噪聲，在系統(tǒng)的輸出端也是線性相加的。系統(tǒng)的輸出為

(3.49)

輸出噪聲為

ef(n)=e(n)*h(n)

(3.50)若e(n)為舍入噪聲，則輸出噪聲的方差為由于e(n)為白噪聲，其各變量間互不相關(guān)，故

由帕斯維爾定理和留數(shù)定理得(3.52)則歸一化輸出噪聲方差為(3.53)

【例3.15】

一個(gè)8位A/D變換器(b=7)，其輸出作為IIR濾波器的輸入，求濾波器輸出端的量化噪聲功率和歸一化噪聲功率。已知IIR濾波器的系統(tǒng)函數(shù)為

解由于A/D的量化效應(yīng)，故濾波器輸入端的噪聲功率為

濾波器的輸出噪聲功率為其積分值等于單位圓內(nèi)所有極點(diǎn)留數(shù)的和。單位圓內(nèi)有一個(gè)極點(diǎn)z=0.999，故用MATLAB編寫的求取歸一化輸出噪聲方差程序如下:

%MATLABPROGRAM3-15

num=1;

den=[1，－0.999];

x=1;

order=max(length(num)，length(den))－1;

si=[zeros(1，order)];

nvar=0;k=1;

whilek>0.0000001

[y，sf]=filter(num，den，x，si);

si=sf;k=abs(y)*abs(y);

nvar=nvar+k;

x=0;

end

disp(′歸一化輸出噪聲方差:′);disp(nvar);程序運(yùn)行結(jié)果為

輸出噪聲方差:500.25013.6.4有限字長運(yùn)算對數(shù)字濾波器的影響

數(shù)字濾波器的實(shí)現(xiàn)涉及到相乘和求和兩種運(yùn)算。在定點(diǎn)制運(yùn)算中，每一次乘法運(yùn)算之后都要作一次舍入(截尾)處理，因此，引入了非線性。采用統(tǒng)計(jì)分析的方法，將舍入誤差作為獨(dú)立噪聲e(n)疊加在信號上。對舍入噪聲e(n)作如A/D變換量化誤差相同的假設(shè)。根據(jù)假設(shè)，整個(gè)系統(tǒng)就可作為線性系統(tǒng)處理。每一個(gè)噪聲可用第1章所講的線性離散系統(tǒng)的理論求出其輸出噪聲，所有輸出噪聲經(jīng)線性疊加后即可得到總的噪聲輸出。

1.IIR的有限字長效應(yīng)

由于e(n)是疊加在輸入端的，故由e(n)造成的輸出誤差為

ef=e(n)*h(n)

輸出噪聲方差

(3.54)

【例3.16】

已知一個(gè)二階IIR低通數(shù)字濾波器，其系統(tǒng)函數(shù)為要求采用定點(diǎn)制運(yùn)算，尾數(shù)作舍入處理，分別計(jì)算其直接型、級聯(lián)型和并聯(lián)型三種結(jié)構(gòu)的舍入誤差。

解

(1)直接型:圖3.24直接型結(jié)構(gòu)流圖

直接型結(jié)構(gòu)流圖如圖3.24所示。圖中，e0(n)、e1(n)、e2(n)分別為系數(shù)0.04、1.7、－0.72相乘后引入的舍入噪聲。采用線性疊加方法，可看出輸出噪聲ef(n)是這三個(gè)舍入噪聲通過網(wǎng)絡(luò)形成的，因此

ef(n)=(e0(n)+e1(n)+e2(n))*h0(n)

其中，h0(n)是H0(z)的單位脈沖響應(yīng)。

輸出噪聲的方差為

將σ2e=q2/12和B(z)代入，利用留數(shù)定理可得

σ2f=22.4q2

(2)級聯(lián)型:

將H(z)分解為圖3.25級聯(lián)型結(jié)構(gòu)流圖

級聯(lián)型結(jié)構(gòu)流圖如圖3.25所示。由圖中可見，噪聲e0(n)、e1(n)通過網(wǎng)絡(luò)H1(z)，噪聲e2(n)只通過網(wǎng)絡(luò)H2(z)。

則有其中，h1(n)和h2(n)分別是H1(z)和H2(z)的單位脈沖響應(yīng)，因此將B1(z)=1－0.9z－1，B2(z)=1－0.8z－1，代入，得

σ2f=15.2q2

(3)并聯(lián)型:

將H(z)分解為部分分式圖3.26并聯(lián)型結(jié)構(gòu)流圖

并聯(lián)型結(jié)構(gòu)流圖如圖3.26所示。并聯(lián)型結(jié)構(gòu)有4個(gè)系數(shù)，故有4個(gè)舍入噪聲，其中，(e0(n)+e1(n))只通過網(wǎng)絡(luò)1/B1(z)，(e2(n)+e3(n))通過網(wǎng)絡(luò)1/B2(z)。

輸出噪聲方差為

代入B1(z)和B2(z)及σe2值，得

σf2=1.34q2通過比較三種結(jié)構(gòu)的誤差大小可知，直接型誤差最大，并聯(lián)型誤差最小。這是由于直接型結(jié)構(gòu)的所有舍入誤差都經(jīng)過全部網(wǎng)絡(luò)的反饋環(huán)節(jié)，反饋過程中誤差積累，故輸出誤差很大;級聯(lián)型結(jié)構(gòu)的每個(gè)舍入誤差只通過其后面的反饋環(huán)節(jié)，而不通過它前面的反饋環(huán)節(jié)，故誤差小于直接型;并聯(lián)型的每個(gè)并聯(lián)網(wǎng)絡(luò)的舍入誤差只通過本身的反饋環(huán)節(jié)，與其他并聯(lián)網(wǎng)絡(luò)無關(guān)，