scal2frq

音階から周波数

構文

frq = scal2frq(A,wname,delta)

frq = scal2frq(A,wname)

説明

frq = scal2frq(A,wname,delta) は、A で与えられた音階および wname とサンプリング周期 delta で指定されたウェーブレットに対応する疑似周波数を返します。出力 frq は実数値で、次元は A と同じになります。

frq = scal2frq(A,wname) は frq = scal2frq(A,wname,1) と等価です。

例

すべて折りたたむ

音階と疑似周波数

ライブスクリプトを開く

この例では、音階を 2 倍にしたときに疑似周波数がどのように変化するかを示します。

5 オクターブにわたってオクターブあたりの音が 10 となる音階のベクトルを作成します。

vpo = 10;
no = 5;
a0 = 2^(1/vpo);
ind = 0:vpo*no;
sc = a0.^ind;

音階の範囲が 5 オクターブをカバーしていることを確認します。

log2(max(sc)/min(sc))

ans = 5.0000

音階をプロットする場合、データカーソルを使用して、インデックス $n + 10$ の音階がインデックス $n$ の音階の 2 倍になることを確認できます。各オクターブをマークするように y 軸目盛りを設定します。

plot(ind,sc)
title('Scales')
xlabel('Index')
ylabel('Scale')
grid on
set(gca,'YTick',2.^(0:5))

音階を実数値 Morlet ウェーブレットの疑似周波数に変換します。まず、サンプリング周期を 1 と仮定します。オクターブあたり 10 音なので、表の 10 行目ごとに表示します。音階が 2 倍になるたびに疑似周波数が半分になっていることを確認します。

pf = scal2frq(sc,"morl");
T = [sc(:) pf(:)];
T = array2table(T, ...
    'VariableNames',{'Scale','Pseudo-Frequency'});
disp(T(1:10:end,:))

    Scale    Pseudo-Frequency
    _____    ________________

      1            0.8125    
      2           0.40625    
      4           0.20313    
      8           0.10156    
     16          0.050781    
     32          0.025391

100 Hz でデータがサンプリングされると仮定します。音階、対応する疑似周波数、および周期を含む table を作成します。表の 10 行目ごとに表示します。

Fs = 100;
DT = 1/Fs;
pf = scal2frq(sc,"morl",DT);
T = [sc(:)/Fs pf(:) 1./pf(:)];
T = array2table(T, ...
    'VariableNames',{'Scale','Pseudo-Frequency','Period'});
T(1:vpo:end,:)

ans=6×3 table
    Scale    Pseudo-Frequency     Period 
    _____    ________________    ________

    0.01           81.25         0.012308
    0.02          40.625         0.024615
    0.04          20.313         0.049231
    0.08          10.156         0.098462
    0.16          5.0781          0.19692
    0.32          2.5391          0.39385

scal2frq に $Δ t = \frac{1}{F s}$ の係数があることに注意してください。これは、音階から周波数への適切な変換を実現するために必要です。 $Δ t$ は生の音階を適切に調整するために必要です。たとえば、次のようにするとします。

f = scal2frq(1,'morl',0.01);

これは、マザー Morlet ウェーブレットを 0.01 倍に拡大した場合、ウェーブレットの中心周波数はどうなるかということです。言い換えれば、 $ψ (t)$ の代わりに $ψ (t / 0.01)$ を見た場合、中心周波数にどのような影響があるかということです。 $Δ t$ によって、音階の正しい調整係数が提供されます。

調整したサイズに音階をまず変換してから scal2frq を指定せずに $Δ t$ を使用することでも、同じ結果が得られます。

scadjusted = sc.*0.01;
pf2 = scal2frq(scadjusted,'morl');
max(pf-pf2)

ans = 0

周波数を使用した等高線図での CWT のプロット

ライブスクリプトを開く

この例では、Hz 単位の近似周波数を使用して CWT の等高線図を作成する方法を示します。

加法性ノイズのサポートが互いに素である 2 つの正弦波で構成される信号を作成します。1 kHz で信号がサンプリングされると仮定します。

Fs = 1000;
t = 0:1/Fs:1-1/Fs;
x = 1.5*cos(2*pi*100*t).*(t<0.25)+1.5*cos(2*pi*50*t).*(t>0.5 & t<=0.75);
x = x+0.05*randn(size(t));

入力信号の CWT を取得し、結果をプロットします。

[cfs,f] = cwt(x,Fs);
contour(t,f,abs(cfs).^2); 
axis tight;
grid on;
xlabel('Time');
ylabel('Approximate Frequency (Hz)');
title('CWT with Time vs Frequency');

入力引数

すべて折りたたむ

`A` — 音階
正の実数値のベクトル

音階。正の実数値のベクトルとして指定します。

`wname` — ウェーブレット
文字ベクトル | string スカラー

ウェーブレット。文字ベクトルまたは string スカラーとして指定します。詳細については、wavefun を参照してください。

`delta` — サンプリング周期
`1` (既定値) | 正の実数値スカラー

サンプリング周期。実数値のスカラーとして指定します。

例: pf = scal2frq([1:5],"db4",0.01)

詳細

すべて折りたたむ

疑似周波数

音階と周波数の関係については、おおよその答えしかありません。

ウェーブレット解析において、音階を周波数に関連付ける方法は、ウェーブレットの中心周波数 F_c を決定し、次の関係を使用することです。

$F_{a} = \frac{F_{c}}{a}$

ここで、

a は音階。
F_c はウェーブレットの中心周波数 (Hz)。
F_a は音階 a に対応する疑似周波数 (Hz)。

この考え方は、与えられたウェーブレットに周波数 F_c の純粋に周期的な信号を関連付けるというものです。ウェーブレット係数のフーリエ変換を最大化する周波数は F_c です。関数 centfrq は、指定されたウェーブレットの中心周波数を計算します。上述の関係から、音階は疑似周波数に反比例することがわかります。たとえば、音階が上がると、ウェーブレットはさらに広がり、疑似周波数は低くなります。

中心周波数とウェーブレットの対応の例を次の Figure に示します。

実数ウェーブレットと複素数ウェーブレットの中心周波数

このように、中心周波数ベースの近似 (赤) は、ウェーブレットの主要な振動 (青) を捉えています。中心周波数は、ウェーブレットの支配的な周波数を便利かつ簡単に特徴付けるものです。

参照

[1] Abry, P. Ondelettes et turbulence. Multirésolutions, algorithmes de décomposition, invariance d'échelles et signaux de pression. Diderot, Editeurs des sciences et des arts, Paris, 1997.