thd

r = thd(x) は、実数値の正弦波信号 x の THD (全高調波歪み) を dBc 単位で返します。全高調波歪みは、基本周波数と最初の 5 つの高調波から、入力信号と同じ長さの修正ピリオドグラムを使用して求められます。この修正ピリオドグラムでは、β = 38 のカイザーウィンドウを使用します。

r = thd(x,fs,n) では、THD の計算で使用するサンプルレート fs および高調波の数 (基本波を含む) を指定します。

r = thd(pxx,f,"psd") では、入力 pxx を片側 PSD (パワースペクトル密度) 推定として指定します。f は pxx の PSD 推定に対応する周波数のベクトルです。

r = thd(pxx,f,n,"psd") では、THD の計算に使用する高調波の数 (基本波を含む) を指定します。

r = thd(sxx,f,rbw,"power") では、入力を片側パワースペクトルとして指定します。rbw はそれについて各パワー推定が積分される分解能帯域幅です。

r = thd(sxx,f,rbw,n,"power") では、THD の計算に使用する高調波の数 (基本波を含む) を指定します。

r = thd(___,"aliased") では、ナイキスト範囲内にエイリアシングされた基本波の高調波が報告されます。このオプションは、入力信号がアンダーサンプリングされた場合に使用します。このオプションを指定しない場合、または "omitaliases" に設定した場合、関数は、ナイキスト範囲外にある基本周波数の高調波を無視します。

[r,harmpow,harmfreq] = thd(___) は、高調波 (基本波を含む) のパワー (単位: dB) と周波数を返します。

出力引数を設定せずに thd(___) を使用すると、現在の Figure ウィンドウで、信号のスペクトルがプロットされ、高調波に注釈が付けられます。スペクトルの基本成分、高調波、DC レベルおよびノイズが、さまざまな色で描画されます。プロットの上に THD と表示されます。基本波と高調波にはラベルが付けられます。DC の項は測定対象にはならないため、ラベル付けされません。

例

2 つの高調波をもつ信号の THD の決定

この例では、基本波と 2 つの高調波で構成される信号の全高調波歪みを dBc 単位で陽的に求める方法を示します。この陽的な計算は、thd で返される結果と照合されます。

1 kHz でサンプリングされた信号を作成します。この信号は、振幅が 2 である 100 Hz の基本波と、周波数が 200 Hz および 300 Hz、振幅が 0.01 および 0.005 である 2 つの高調波で構成されています。全高調波歪みを陽的な計算と thd により求めます。

t = 0:0.001:1-0.001;
x = 2*cos(2*pi*100*t)+0.01*cos(2*pi*200*t)+0.005*cos(2*pi*300*t);
tharmdist = 10*log10((0.01^2+0.005^2)/2^2)

tharmdist = -45.0515

r = thd(x)

r = -45.0515

高調波の数の指定

1 kHz でサンプリングされた信号を作成します。この信号は振幅が 2 である 100 Hz の基本波と、周波数が 200 Hz、300 Hz、400 Hz、振幅が 0.01、0.005、0.0025 である 3 つの高調波で構成されています。

高調波の数を 3 に設定します。この中には基本波も含まれます。したがって、THD の計算には 100 Hz、200 Hz、300 Hz でのパワーが使用されます。

t = 0:0.001:1-0.001;
x = 2*cos(2*pi*100*t)+0.01*cos(2*pi*200*t)+ ...
    0.005*cos(2*pi*300*t)+0.0025*sin(2*pi*400*t);
r = thd(x,1000,3)

r = -45.0515

高調波の数を 3 に設定すると、THD の計算において 400 Hz でのパワーが無視されます。

高調波の数の指定 (PSD 入力)

信号のピリオドグラム PSD 推定を取得し、この PSD 推定を thd の入力として使用します。高調波の数を 3 に設定します。この中には基本波も含まれます。したがって、THD の計算には 100 Hz、200 Hz、300 Hz でのパワーが使用されます。

t = 0:0.001:1-0.001;
fs = 1000;
x = 2*cos(2*pi*100*t)+0.01*cos(2*pi*200*t)+ ...
    0.005*cos(2*pi*300*t)+0.0025*sin(2*pi*400*t);
[pxx,f] = periodogram(x,rectwin(length(x)),length(x),fs);
r = thd(pxx,f,3,'psd')

r = -45.0515

パワースペクトルからの THD

ハミングウィンドウで得たパワースペクトルとウィンドウの分解能帯域幅を入力して、THD を求めます。

10 kHz でサンプリングされた信号を作成します。この信号は振幅が 2 である 100 Hz の基本波と、周波数が 300 Hz、500 Hz、700 Hz、振幅が 0.01、0.005、0.0025 である 3 つの奇数番目の高調波で構成されています。高調波の数に 7 を指定します。THD を求めます。

fs = 10000;
t = 0:1/fs:1-1/fs;
x = 2*cos(2*pi*100*t)+0.01*cos(2*pi*300*t)+ ...
    0.005*cos(2*pi*500*t)+0.0025*sin(2*pi*700*t);
[sxx,f] = periodogram(x,hamming(length(x)),length(x),fs,'power');
rbw = enbw(hamming(length(x)),fs);
r = thd(sxx,f,rbw,7,'power')

r = -44.8396

エイリアシングされた高調波を使用した場合と使用しない場合の THD

入力が 2.1 kHz トーンの弱非線形増幅器の出力に類似した信号を生成します。信号は 10 kHz で 1 秒間サンプリングされます。信号のパワースペクトルを計算してプロットします。β = 38 のカイザーウィンドウを使用して計算します。

Fs = 10000;
f = 2100;

t = 0:1/Fs:1; 
x = tanh(sin(2*pi*f*t)+0.1) + 0.001*randn(1,length(t));

periodogram(x,kaiser(length(x),38),[],Fs,'power')

高調波は周波数 4.2 kHz、6.3 kHz、8.4 kHz、10.5 kHz、12.6 kHz、14.7 kHz でノイズから突出します。最初の周波数を除くすべての周波数は、ナイキスト周波数より大きくなります。高調波はそれぞれ 3.7 kHz、1.6 kHz、0.5 kHz、2.6 kHz、4.7 kHz にエイリアシングされます。

信号の全高調波歪みを計算します。thd は、既定で、エイリアシングされた高調波をノイズの一部として扱います。

thd(x,Fs,7);

計算を繰り返しますが、今度はエイリアシングされた高調波を信号の一部として扱います。

thd(x,Fs,7,'aliased');

高調波のパワーと対応する周波数

10 kHz でサンプリングされた信号を作成します。この信号は振幅が 2 である 100 Hz の基本波と、周波数が 300 Hz、500 Hz、700 Hz、振幅が 0.01、0.005、0.0025 である 3 つの奇数番目の高調波で構成されています。高調波の数に 7 を指定します。THD、高調波のパワー、対応する周波数を求めます。

fs = 10000;
t = 0:1/fs:1-1/fs;
x = 2*cos(2*pi*100*t)+0.01*cos(2*pi*300*t)+ ...
    0.005*cos(2*pi*500*t)+0.0025*sin(2*pi*700*t);
[r,harmpow,harmfreq] = thd(x,10000,7);
[harmfreq harmpow]

ans = 7×2

  100.0000    3.0103
  201.0000 -321.0983
  300.0000  -43.0103
  399.0000 -281.9259
  500.0000  -49.0309
  599.0000 -282.1066
  700.0000  -55.0515

偶数番目の高調波のパワーはおよそ $- 300$ dB で、 $1 0^{- 15}$ の振幅に相当します。

増幅信号の THD

50 kHz でサンプリングされる、周波数 2.5 kHz の正弦波を生成します。標準偏差 0.00005 のホワイトガウスノイズを信号に付加します。結果を弱非線形増幅器に通します。THD をプロットします。

fs = 5e4;
f0 = 2.5e3;
N = 1024;
t = (0:N-1)/fs;

ct = cos(2*pi*f0*t);
cd = ct + 0.00005*randn(size(ct));

amp = [1e-5 5e-6 -1e-3 6e-5 1 25e-3];
sgn = polyval(amp,cd);
thd(sgn,fs);

プロットには、比率の計算に使用されるスペクトルと、ノイズとして処理される領域が示されます。DC レベルは計算対象になりません。基本波と高調波にはラベルが付けられます。

入力引数

`x` — 実数値正弦波入力信号
ベクトル

実数値正弦波入力信号。行ベクトルまたは列ベクトルで指定されます。

例: cos(pi/4*(0:159))+cos(pi/2*(0:159))

データ型: single | double

`fs` — サンプルレート
正のスカラー

サンプルレート。正のスカラーで指定します。サンプルレートは単位時間あたりのサンプル数です。時間の単位が秒の場合、サンプルレートの単位は Hz です。

`n` — 高調波数
正の整数

高調波の数。正の整数として指定します。

`pxx` — 片側 PSD 推定
ベクトル

片側 PSD 推定。非負の実数値列ベクトルで指定されます。

パワースペクトル密度は、デシベル単位ではなく、線形単位で表さなければなりません。デシベル値をパワー値に変換するには、db2pow を使用します。

例: [pxx,f] = periodogram(cos(pi./[4;2]*(0:159))'+randn(160,2)) は、2π Hz でサンプリングされたノイズの多い 2 チャネル正弦波のピリオドグラム PSD 推定値とそれについて計算される周波数を指定します。

データ型: single | double

`f` — 巡回周波数
ベクトル

片側 PSD 推定 pxx に対応する巡回周波数。行ベクトルまたは列ベクトルで指定されます。f の最初の要素は 0 でなければなりません。

データ型: double | single

`sxx` — パワースペクトル
非負の実数値行ベクトルまたは列ベクトル

パワースペクトル。非負の実数値行ベクトルまたは列ベクトルで指定されます。

パワースペクトルは、デシベル単位ではなく線形単位で表さなければなりません。デシベル値をパワー値に変換するには、db2pow を使用します。

例: [sxx,w] = periodogram(cos(pi./[4;2]*(0:159))'+randn(160,2),'power') では、ホワイトガウスノイズに組み込まれる 2 チャネルの正弦波のピリオドグラムパワースペクトル推定と、計算に使う正規化周波数を指定します。

`rbw` — 分解能帯域幅
正のスカラー

分解能帯域幅。正のスカラーで指定されます。この分解能帯域幅は、離散フーリエ変換の周波数分解能と、ウィンドウの等価ノイズ帯域幅の積です。

出力引数

`r` — 全高調波歪み (dBc)
実数値のスカラー

全高調波歪み (dBc)。実数値のスカラーとして返されます。

`harmpow` — 高調波のパワー
実数値のスカラーまたはベクトル

高調波のパワー。dB 単位表記の実数値のスカラーまたはベクトルとして返されます。harmpow がスカラーであるかベクトルであるかは、入力引数 n として指定する高調波の数によって異なります。

`harmfreq` — 高調波の周波数
非負のスカラーまたはベクトル

高調波の周波数。非負のスカラーまたはベクトルとして返されます。harmfreq がスカラーであるかベクトルであるかは、入力引数 n として指定する高調波の数によって異なります。

詳細