Main Content

このページの翻訳は最新ではありません。ここをクリックして、英語の最新版を参照してください。

cfirpm

複素かつ非線形位相の等リップル FIR フィルターの設計

説明

b = cfirpm(n,f,fresp) では、関数 fresp によって返される望ましい周波数応答を最もよく近似する長さ n+1 の FIR フィルターが返されます。この関数は、関数ハンドル (@fresp) により呼び出されます。

b = cfirpm(n,f,fresp,w) では、w で指定された重みを使用して、各周波数帯域の近似に重み付けが行われます。

b = cfirpm(n,f,a) では、f の帯域エッジの振幅 a を指定します。この構文は、b = cfirpm(n,f,{@multiband,a}) と同じ結果を返します。

b = cfirpm(n,f,a,w) では、最適化中にオプションとして使用される、帯域ごとに 1 つの正加重のセットが適用されます。w を指定しない場合、関数は重みを 1 に設定します。

b = cfirpm(___,sym) では、設計のインパルス応答に対称的な制約を与えます。sym の指定に加え、前の構文の入力引数を任意に組み合わせて指定します。

b = cfirpm(___,debug) では、フィルター設計中に中間結果を表示または非表示にします。

b = cfirpm(___,lgrid) では、周波数グリッドの密度を制御します。

b = cfirpm(___,'skip_stage2') では、関数 cfirpm によって標準の firpm 誤差変換アルゴリズムで最適解が求められなかったと判定された場合にのみ実行される、第 2 段最適化アルゴリズムを無効にします。このアルゴリズムを無効にすると計算速度が向上する場合がありますが、精度が低下する可能性があります。既定の設定では、第 2 段最適化は有効になっています。

[b,delta] = cfirpm(___) では、最大リップル高 delta が返されます。

[b,delta,opt] = cfirpm(___) では、関数 cfirpm で計算されたオプションの結果が返されます。

すべて折りたたむ

31 タップの線形位相ローパス フィルターを設計します。その振幅応答と位相応答を表示します。

b = cfirpm(30,[-1 -0.5 -0.4 0.7 0.8 1],@lowpass);
fvtool(b,1,'OverlayedAnalysis','phase')

Figure Figure 1: Magnitude Response (dB) and Phase Response contains an axes object. The axes object with title Magnitude Response (dB) and Phase Response contains an object of type line.

exp(-jπfN/2+j4πf|f|) で近似的に設定される周波数応答をもつ、次数 22 の非線形位相オールパス FIR フィルターを設計します。ここで、f[-1,1] とします。

n = 22;                              % Filter order
f = [-1 1];                          % Frequency band edges
w = [1 1];                           % Weights for optimization
gf = linspace(-1,1,256);             % Grid of frequency points 
d = exp(-1i*pi*gf*n/2 + 1i*pi*pi*sign(gf).*gf.*gf*(4/pi));
                                     % Desired frequency response

cfirpm を使用して、この FIR フィルターを計算します。実際の振幅応答と近似の振幅応答を dB 単位で、位相応答を度単位でプロットします。

b = cfirpm(n,f,'allpass',w,'real');  % Approximation
freqz(b,1,256,'whole')

subplot(2,1,1)                       % Overlay response
hold on
plot(pi*(gf+1),20*log10(abs(fftshift(d))),'r--')

subplot(2,1,2)
hold on
plot(pi*(gf+1),unwrap(angle(fftshift(d)))*180/pi,'r--')
legend('Approximation','Desired','Location','SouthWest')

Figure contains 2 axes objects. Axes object 1 contains 2 objects of type line. These objects represent Approximation, Desired. Axes object 2 contains 2 objects of type line.

カスタム周波数応答関数 fresp を使用して次数 30 のローパス フィルターを設計します。関数 fresp のコードは、この例の終わりにあります。

[b,delta]= cfirpm(30,linspace(-1,1,32),@fresp);

FVTool を使用して、フィルターの振幅応答を可視化します。

fvtool(b,1)

Figure Figure 1: Magnitude Response (dB) contains an axes object. The axes object with title Magnitude Response (dB) contains an object of type line.

ユーザー定義の fresp 関数: ローパス フィルターの設計

関数 fresp を使用すると、ローパス フィルター、ハイパス フィルター、または微分器の設計を選択できます。フィルター次数 N と周波数配列 F を指定しなければなりません。周波数グリッド GF と重み W を指定しない場合、関数はこれらの値を自動的に決定します。

function [dh,dw] = fresp(N,F,GF,W)

W = [1;1]*(W(:).'); W = W(:);

type = 'lowpass';

mags = zeros(size(W));

switch type
    case 'lowpass'
        mags(10:end-10) = 1;
    case 'highpass'
        mags(1:10) = 1;
        mags(end-10:end) = 1;
    case 'differentiator'
        mags = abs(linspace(-pi,pi,length(mags)));
end

dh = interp1(F(:),mags,GF).*exp(-1j*pi*GF*N/2);
dw = interp1(F(:),W,GF);

end

入力引数

すべて折りたたむ

フィルター次数。正の実数スカラーとして指定します。

正規化周波数点。範囲 [–1, 1] の要素をもつ実数値のベクトルとして指定します。ここで、1 は正規化されたナイキスト周波数に対応します。周波数は増加する順でなければならず、f は偶数長であることが必要です。周波数帯域は、k が奇数の場合、f(k)f(k+1) の範囲です。奇数 k に対して f(k+1) ~ f(k+2) の区間は "遷移帯域"、すなわち最適化の "don't care" 領域です。

周波数応答。関数ハンドルとして指定します。詳細については、事前に定義された周波数応答関数ユーザー定義の周波数応答関数を参照してください。

f で指定した点での目的の振幅。ベクトルとして指定します。k が偶数の場合、点 f(k) と点 f(k+1) 間の周波数における希望する振幅は、点 (f(k),a(k)) と点 (f(k+1),a(k+1)) とを結ぶ線分です。

各周波数帯域の適合の調整に使用される重み付け。実数値ベクトルとして指定します。w の長さは f の長さの半分で、帯域あたり厳密に 1 つの重みに相当します。w を指定しない場合、関数は重みを 1 に設定します。

フィルター設計のインパルス応答に課される対称性の制約。次のいずれかの値を指定します。

  • 'none' — 対称性の制約を課しません。関数に負の帯域周波数を渡した場合、または fresp によって既定値が与えられない場合は、このオプションが既定値になります。

  • 'even' — 実数かつ偶数のインパルス応答に制限します。このオプションはハイパス、ローパス、オールパス、バンドパス、バンドストップ、逆 sinc およびマルチバンドの設計の場合の既定値です。

  • 'odd' — 実数かつ奇数のインパルス応答に制限します。このオプションは、ヒルベルト変換および微分器設計の場合の既定値です。

  • 'real' — 周波数応答に共役対称性の制約を課します。

'none' 以外の値を指定した場合、正の周波数上にのみ帯域エッジを指定しなければなりません (負の周波数領域は対称性によって埋められます)。sym を指定しない場合、関数は既定の設定を確認するために fresp を呼び出します。ユーザー定義の関数 fresp は、フィルター次数 n として 'defaults' を受け取れば、有効な sym オプションを返さなければなりません。

フィルター設計中の中間結果の表示。'off''trace''plots' または 'both' として指定します。

周波数グリッドの密度。整数の cell 配列として指定します。周波数グリッドには約 2^nextpow2(lgrid*n) 個の周波数点があります。

出力引数

すべて折りたたむ

フィルター係数。長さ n+1 の行ベクトルとして返されます。

最大リップル高。スカラーとして返されます。

関数 cfirpm によって計算されるオプションの結果。次のフィールドを含む構造体として返されます。

フィールド

説明

opt.fgrid

フィルター設計の最適化に使用される周波数グリッド ベクトル

opt.des

opt.fgrid の各点に対する望ましい周波数応答

opt.wt

opt.fgrid の各点の重み

opt.H

opt.fgrid の各点に対する実際の周波数応答

opt.error

opt.fgrid 内の各点の誤差

opt.iextr

極値周波数に対する opt.fgrid へのインデックスのベクトル

opt.fextr

極値周波数のベクトル

詳細

すべて折りたたむ

事前に定義された周波数応答関数

事前に定義された周波数応答関数 fresp は、この節のいくつかの一般的なフィルター設計で使用されます。カスタム関数 fresp の作成方法の詳細については、関数ハンドルの作成を参照してください。

すべての事前定義された周波数応答関数では、f が負の周波数を含まず、かつ d = 0 の場合、対称性を示すオプション sym には既定値の 'even' が使用されます。それ以外の場合、sym には既定値の 'none' が使用されます。詳細は、symを参照してください。すべての事前定義された周波数応答関数では、フィルターの応答がサンプリング区間の単位で n/2+d の群遅延をもつように d によって群遅延オフセットが指定されます。負の値では遅延が減少し、正の値では遅延が増大します。既定では、d = 0 です。

  • @lowpass, @highpass, @allpass, @bandpass, @bandstop

    これらの関数では共通の構文が使用されます。ここでは @lowpass の例を示します。

    b = cfirpm(n,f,@lowpass,...) および

    b = cfirpm(n,f,{@lowpass,d},...) では、線形位相 (n/2+d 遅延) フィルターを設計します。

    メモ

    @bandpass フィルターでは、周波数ベクトルの最初の要素がゼロ以下で、最後の要素はゼロ以上でなければなりません。

  • @multiband では、任意の帯域振幅をもつ線形位相周波数応答フィルターを設計します。

    b = cfirpm(n,f,{@multiband,a},...) および

    b = cfirpm(n,f,{@multiband,a,d},...) では、f の帯域エッジでの希望する振幅を含むベクトル a を指定します。k が偶数の場合、点 f(k) と点 f(k+1) 間の周波数における希望する振幅は、点 (f(k),a(k)) と点 (f(k+1),a(k+1)) とを結ぶ線分です。

  • @differentiator では、線形位相微分器を設計します。これらの設計では、ゼロ周波数は遷移帯域内でなければならず、帯域の重み付けは周波数に反比例して設定されます。

    b = cfirpm(n,f,{@differentiator,fs},...) および

    b = cfirpm(n,f,{@differentiator,fs,d},...) では、微分器応答の勾配を決定するのに使用されるサンプルレート fs を指定します。指定しなかった場合、fs には既定値の 1 が設定されます。

  • @hilbfilt では、線形位相ヒルベルト変換フィルターの応答を設計します。ヒルベルト設計では、ゼロ周波数は遷移帯域内になければなりません。

    b = cfirpm(n,f,@hilbfilt,...) および

    b = cfirpm(N,F,{@hilbfilt,d},...) では、線形位相 (n/2+d 遅延) ヒルベルト変換フィルターを設計します。

  • @invsinc では、線形位相逆 sinc フィルターの応答を設計します。

    b = cfirpm(n,f,{@invsinc,a},...) および

    b = cfirpm(n,f,{@invsinc,a,d},...) は、sinc(a*g) として計算された sinc 関数に対するゲイン a を指定します。ここで、g は範囲 [–1, 1] に正規化された最適化グリッド周波数を含むものとします。既定の設定では a = 1 です。群遅延オフセットは d です。フィルター応答はサンプル間隔の単位において n/2+d の群遅延をもちます。ここで、n はフィルターの次数です。負の値では遅延が減少し、正の値では遅延が増大します。既定の設定では d = 0 です。

ユーザー定義の周波数応答関数

fresp に対して事前定義された周波数応答関数の代わりに、ユーザー定義関数を使用することもできます。

関数 cfirpm は、この構文を使用してこのユーザー定義関数を呼び出します。

[dh,dw] = fresp(n,f,gf,w,p1,p2,...)

  • n はフィルター次数です。

  • f は、–1 ~ 1 の間で単調増加する周波数帯域エッジのベクトルです。ここで、1 はナイキスト周波数に相当します。

  • gf は、cfirpm によって指定された周波数帯域で線形に内挿されたグリッド点のベクトルです。入力 gf では、応答関数が評価される周波数グリッドが決定されます。関数 cfirpm は、opt 構造の fgrid フィールドでこのデータを返します。

  • w は、最適化中に使用される帯域ごとの正で実数の重みのベクトルです。w は、cfirpm を呼び出す際にはオプションです。この入力を指定しない場合、cfirpmfresp に渡す前にこれを重み 1 に設定します。

  • dhdw は、それぞれ望ましい複素周波数応答と帯域重みのベクトルで、グリッド gf の各周波数で評価されます。

  • p1,p2,... は、fresp に渡されるオプションのパラメーターです。

さらに、関数 cfirpm は既定の対称性 sym を決定するために、fresp の予備呼び出しを行います。cfirpm は、この構文を使用してこの呼び出しを行います。

sym = fresp('defaults',{n,f,[],w,p1,p2,...})
必要に応じて、適切な対称性の既定設定を決定するために、引数を使用できます。新しい周波数応答関数を生成する場合のテンプレートとして、ローカル関数 lowpass を使用できます。関数 lowpass を見つけるには、コマンド ラインで edit cfirpm と入力し、lowpass で関数コード cfirpm 内を検索してください。関数をコピーして、変更したり、名前を変更したり、パス内に保存することができます。

アルゴリズム

関数 cfirpm を使用すると、複素となる可能性のある FIR フィルターの設計において、任意の周波数領域の制約を指定できます。チェビシェフ (またはミニマックス) フィルター誤差が最適化され、等リップル FIR フィルター設計が作成されます。

複素数の場合、Remez 変換法の拡張型が実行されます。この変換法では、フィルターの等リップル特性が n+2 の極をもつように制限されている場合に、最適なフィルターが得られます。このフィルターが収束しない場合は、アルゴリズムが上昇-下降アルゴリズムに切り替わり、最適な解に収束するように処理が行われます。詳細については、参考文献を参照してください。

参照

[1] Demjanjov, V. F., and V. N. Malozemov. Introduction to Minimax. New York: John Wiley & Sons, 1974.

[2] Karam, L.J. Design of Complex Digital FIR Filters in the Chebyshev Sense. Ph.D. Thesis, Georgia Institute of Technology, March 1995.

[3] Karam, L.J., and J. H. McClellan. "Complex Chebyshev Approximation for FIR Filter Design." IEEE® Transactions on Circuits and Systems II: Analog and Digital Signal Processing 42, no. 3 (March 1995): 207–216.

拡張機能

バージョン履歴

R2006a より前に導入

参考

| | |