ドキュメンテーション

最新のリリースでは、このページがまだ翻訳されていません。 このページの最新版は英語でご覧になれます。

FIR ハーフバンド フィルターの設計

この例では、FIR ハーフバンド フィルターの設計法を示します。ハーフバンド フィルターは、マルチレート信号処理アプリケーションにおいて 2 の係数で内挿または間引きする際に幅広く使用されます。ハーフバンド フィルターは、ポリフェーズの形式で効率的に実装されます。その理由は、その係数の約半分がゼロに等しくなるためです。

ハーフバンド フィルターには 2 つの重要な特徴があります。つまり、通過帯域リップルと阻止帯域リップルは同じでなければならず、通過帯域エッジと阻止帯域エッジの周波数はハーフバンド周波数 Fs/4 (または正規化周波数の pi/2 ラジアン/サンプル) から等距離になければなりません。

ハーフバンド係数の取得

関数 firhalfband は、FIR ハーフバンド等リップル フィルターの係数を返します。簡単な例として、96 kHz でサンプリングされたデータと通過帯域周波数 22 kHzで処理を行うハーフバンド フィルターを考えます。

Fs  = 96e3;
Fp  = 22e3;
N   = 100;
NUM = firhalfband(N,Fp/(Fs/2));
hfvt = fvtool(NUM,'Fs',Fs,'Color','white');
hfvt.MagnitudeDisplay = 'Zero-phase';

応答にズーム インすることで、通過帯域と阻止帯域のピーク ツー ピーク リップルが同じであることを確認できます。また、Fs/4 (24 kHz) 点を中心とした対称性もあります。指定したとおりに通過帯域が最大 22 kHz まで拡張し、26 kHz で阻止帯域が始まります。インパルス応答を調べることで、係数が 1 つおきに 0 と等しいことも確認できます。これにより、係数 2 で内挿/間引きするためのフィルター実装を大変効率よく行えます。

hfvt.Analysis = 'impulse';

dsp.FIRHalfbandInterpolator および dsp.FIRHalfbandDecimator

関数 firhalfband には、この他にいくつかの設計オプションがあります。ただし、ほとんどの場合、dsp.FIRHalfbandInterpolator および dsp.FIRHalfbandDecimator で直接操作する方が適切です。これら 2 つの System object では、係数を設計するだけでなく、ポリフェーズの内挿/間引きを効率的に実装できます。これらは倍精度/単精度浮動小数点データだけでなく、固定小数点データのフィルター処理もサポートします。また、ARM® Cortex® M および ARM Cortex A に最適化されたコード生成に加え、C コードおよび HDL コード生成もサポートします。

HBi = dsp.FIRHalfbandInterpolator('SampleRate',Fs,...
    'Specification','Filter order and transition width','FilterOrder',N,...
    'TransitionWidth',4000);
fvtool(HBi,'Fs',2*Fs,'Color','white');

内挿を実行するには、関数 step を使用します。これはマルチレート フィルターであるため、サンプルレートが意味するものを定義することが重要です。この System object やその他すべての System object では、サンプルレートは入力信号のサンプルレートを指します。ただし、FVTool では、サンプルレートはフィルターが処理を行うレートとして定義されます。内挿の場合には、アップサンプリング後にフィルター処理を行います (概念上)。したがって、2 でアップサンプリングされるため、FVTool のサンプルレートを 2*Fs として指定する必要があります。

FrameSize = 256;
SA  = dsp.SpectrumAnalyzer('SampleRate',2*Fs,'SpectralAverages',5);
SW10 = dsp.SineWave('Frequency',10e3','SampleRate',Fs,...
    'SamplesPerFrame',FrameSize);
SW20 = dsp.SineWave('Frequency',20e3','SampleRate',Fs,...
    'SamplesPerFrame',FrameSize);
tic
while toc < 10
    x = step(SW10) + step(SW20) + 0.01.*randn(FrameSize,1); %  96 kHz
    y = step(HBi,x);                                         % 192 kHz
    step(SA,y);
end

スペクトルの複製は約 40 dB 減衰することに注意してください。これはハーフバンド フィルターによって生じるおおよその減衰量です。フィルターの群遅延を補正する際には、入力サンプルと内挿されたサンプルを重ねてプロットできます。フィルターの出力時に入力サンプルが変わっていないことに注意してください。これは、ハーフバンドのポリフェーズ分岐のいずれかが、入力サンプルを変化させない純粋な遅延分岐であるためです。

GD = 50;
n = 0:2:511;
stem(n(1:end-GD/2),x(1:end-GD/2),'k','filled')
hold on
nu = 0:511;
stem(nu(1:end-GD),y(GD+1:end))
legend('Input samples','Interpolated samples')

間引きの場合には、フィルター処理後にダウンサンプリング (概念上) するため、dsp.FIRHalfbandDecimator で指定したサンプルレートはフィルターのサンプルレートに対応します。そのため、間引きでは、FVTool で指定した Fs を係数で乗算する必要はありません。

FrameSize = 256;
FsIn = 2*Fs;
HBd = dsp.FIRHalfbandDecimator('SampleRate',FsIn,...
    'Specification','Filter order and transition width','FilterOrder',N,...
    'TransitionWidth',4000);
fvtool(HBd,'Fs',FsIn,'Color','white');
SA  = dsp.SpectrumAnalyzer('SampleRate',Fs,'SpectralAverages',5);
SW10 = dsp.SineWave('Frequency',10e3','SampleRate',Fs,...
    'SamplesPerFrame',FrameSize);
SW20 = dsp.SineWave('Frequency',20e3','SampleRate',Fs,...
    'SamplesPerFrame',FrameSize);
tic
while toc < 10
    x = step(SW10) + step(SW20) + 0.01.*randn(FrameSize,1); %  96 kHz
    y = step(HBi,x);                                         % 192 kHz
    xd = step(HBd,y);                                       %  96 kHz
    step(SA,xd);
end

フィルター係数の取得

関数 tf を使用すると、内挿/間引きからフィルター係数を抽出できます。

NUM = tf(HBi); % Or NUM = TF(HBd);

さまざまな設計仕様の使用

フィルター次数と遷移幅を指定する代わりに、所定の遷移幅と阻止帯域の減衰量を与えて最小次数のフィルターを設計できます。

Ast = 80; % 80 dB
HBi = dsp.FIRHalfbandInterpolator('SampleRate',Fs,...
    'Specification','Transition width and stopband attenuation',...
    'StopbandAttenuation',Ast,'TransitionWidth',4000);
fvtool(HBi,'Fs',2*Fs,'Color','white');

すべての内挿と同様に、絶対単位での通過帯域ゲインは内挿係数と等価です (ハーフバンドの場合は 2)。これは、6.02 dB の通過帯域ゲインに相当します。

フィルター次数と阻止帯域の減衰量を指定することも可能です。

HBd = dsp.FIRHalfbandDecimator('SampleRate',Fs,...
    'Specification','Filter order and stopband attenuation',...
    'StopbandAttenuation',Ast,'FilterOrder',N);
fvtool(HBd,'Fs',Fs,'Color','white');

内挿と異なり、間引きでは通過帯域に 1 (0 dB) のゲインがあります。

フィルター バンクでのハーフバンド フィルターの使用

ハーフバンドの内挿と間引きを使用すると、合成/解析フィルター バンクを効率的に実装できます。これまでに示したハーフバンド フィルターは、すべてローパス フィルターでした。1 つの加算器を追加すると、ローパス応答に加えてハイパス応答も得られ、2 つの応答をフィルター バンクの実装に使用できるようになります。

以下のコードは、直交ミラー フィルター (QMF) バンクをシミュレートします。1 kHz と 3 kHz の正弦波で構成される 8 kHz 信号が、ローパス/ハイパス ハーフバンド間引きを使用して 2 つの 4 kHz 信号に分割されます。ローパス信号は 1 kHz の正弦波を保持し、ハイパス信号は 3 kHz の正弦波 (ダウンサンプリング後に 1 kHz にエイリアシングされる) を保持します。次に、信号は合成フィルター バンクによって、ハーフバンド内挿を使用して統合されます。ハイパス分岐により、エイリアシングされた 1 kHz の正弦波がアップコンバートされ、3 kHz に戻ります。内挿された信号のサンプルレートは 8 kHz です。

Fs1 = 8000; % Units = Hz
Spec = 'Filter order and transition width';
Order = 52;
TW = 4.1e2; % Units = Hz

% Construct FIR Halfband Interpolator
HBi = dsp.FIRHalfbandInterpolator(...
    'Specification',Spec,...
    'FilterOrder',Order,...
    'TransitionWidth',TW,...
    'SampleRate',Fs1/2,...
    'FilterBankInputPort',true);

% Construct FIR Halfband Decimator
HBd = dsp.FIRHalfbandDecimator(...
    'Specification',Spec,...
    'FilterOrder',Order,...
    'TransitionWidth',TW,...
    'SampleRate',Fs1);

% Input
f1 = 1000;
f2 = 3000;
InputWave   = dsp.SineWave('Frequency',[f1,f2],'SampleRate',Fs1,...
    'SamplesPerFrame',1024,'Amplitude',[1 0.25]);

% Construct Spectrum Analyzer object to view the input and output
SA = dsp.SpectrumAnalyzer('SampleRate',Fs1,...
    'PlotAsTwoSidedSpectrum',false,'ShowLegend',true,'YLimits',...
    [-120 30],...
    'Title',...
    'Input Signal and Output Signal of Quadrature Mirror Filter');
SA.ChannelNames = {'Input','Output'};

tic
while toc < 10
    Input = sum(step(InputWave),2);
    NoisyInput = Input+(10^-5)*randn(1024,1);
    [Lowpass, Highpass] = step(HBd,NoisyInput);
    Output = step(HBi,Lowpass,Highpass);
    step(SA,[NoisyInput,Output]);
end

詳細設計オプション: さまざまな設計アルゴリズムの指定

これまでに示した設計は、すべて最適な等リップル設計でした。fdesign.interpolator および fdesign.decimator を使用すると、その他の設計アルゴリズムを利用できます。

Fs = 44.1e3;
N  = 90;
TW = 1000/Fs; % Transition width
d = fdesign.interpolator(2,'halfband','N,TW',N,TW);
HBeq = design(d,'equiripple','SystemObject',true); % Equiripple design
HBls = design(d,'firls','SystemObject',true);      % Least-squares design
HBka = design(d,'kaiserwin','SystemObject',true);  % Kaiser-window design

FVTool で設計を比較できます。さまざまな設計により、最小の阻止帯域の減衰量とより全体的な減衰量のトレードオフが可能になります。

hfvt = fvtool(HBeq,HBls,HBka,'Fs',2*Fs,'Color','white');
legend(hfvt, 'Equiripple design', 'Least-squares design', ...
    'Kaiser-window design')

阻止帯域の減衰の制御

代わりに次数と阻止帯域の減衰を指定することもできます。これにより、全体的な阻止帯域の減衰量と遷移幅とのトレードオフが可能になります。

Ast  = 60; % Minimum stopband attenuation
d = fdesign.interpolator(2,'halfband','N,Ast',N,Ast);
HBeq = design(d,'equiripple','SystemObject',true); % Equiripple design
HBka = design(d,'kaiserwin','SystemObject',true);  % Kaiser-window design
hfvt = fvtool(HBeq,HBka,'Fs',2*Fs,'Color','white');
legend(hfvt, 'Equiripple design', 'Kaiser-window design')

最小次数設計

設計仕様を満たすために必要な最小次数のフィルターを設計する場合には、等リップル設計に加えてカイザー ウィンドウ設計も使用できます。カイザー ウィンドウ設計のための実際の次数は、等リップル設計に必要な次数より大きくなりますが、代わりに全体的な阻止帯域の減衰量は改善されます。

Fs = 44.1e3;
TW = 1000/(Fs/2); % Transition width
Ast = 60;  % 60 dB minimum attenuation in the stopband
d = fdesign.decimator(2,'halfband','TW,Ast',TW,Ast);
HBeq = design(d,'equiripple','SystemObject',true); % Equiripple design
HBka = design(d,'kaiserwin','SystemObject',true);  % Kaiser-window design
hfvt = fvtool(HBeq,HBka,'Fs',Fs,'Color','white');
legend(hfvt, 'Equiripple Design', 'Kaiser-window design')

阻止帯域の減衰が増大する等リップル設計

カイザー ウィンドウ フィルターを設計する代わりに、変更した '等リップル' 設計を使用して阻止帯域の減衰量を増やすこともできます。

HBeq1 = design(d,'equiripple','StopbandShape','1/f','StopbandDecay',4,...
    'SystemObject',true);
HBeq2 = design(d,'equiripple','StopbandShape','linear', ...
    'StopbandDecay',53.333,'SystemObject',true);
hfvt = fvtool(HBeq1,HBeq2,'Fs',Fs,'Color','white');
legend(hfvt,'Stopband decaying as (1/f)^4','Stopband decaying linearly')

ハイパス ハーフバンド フィルター

ハイパス ハーフバンド フィルターは、ローパス ハーフバンド フィルターの 1 つおきの係数の符号を変更することで、ローパス ハーフバンド フィルターから取得できます。あるいは、'Type' プロパティを 'Highpass' に設定することで、ハイパス ハーフバンドを直接設計することもできます。

d  = fdesign.decimator(2,'halfband','Type','Highpass','TW,Ast',TW,Ast);
HBhp = design(d,'equiripple','StopbandShape','linear', ...
    'StopbandDecay',53.333,'SystemObject',true);
hfvt = fvtool(HBhp,HBeq2,'Fs',Fs,'Color','white');
legend(hfvt, 'Highpass halfband filter', 'Lowpass halfband filter')

この情報は役に立ちましたか?