Main Content

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

istft

逆短時間フーリエ変換

説明

x = istft(s) は、s の逆短時間フーリエ変換 (ISTFT) を返します。

x = istft(s,fs) は、サンプルレート fs を使用して s の ISTFT を返します。

x = istft(s,ts) は、サンプル時間 ts を使用して ISTFT を返します。

x = istft(___,Name,Value) は、名前と値のペアの引数を使用して追加オプションを指定します。オプションには、FFT ウィンドウの長さとオーバーラップされたサンプルの数が含まれます。これらの引数を前の入力構文のいずれかに追加できます。

[x,t] = istft(___) は、ISTFT が評価される信号時間を返します。

すべて折りたたむ

1 kHz で 1 秒間サンプリングされた 3 つの異なるチャープで構成される 3 チャネルの信号を生成します。

  1. 最初のチャネルは凹二次チャープで構成され、瞬時周波数は t = 0 では 100 Hz であり、t = 1 秒では 300 Hz になります。初期位相は 45 度です。

  2. 2 番目のチャネルは凸二次チャープで構成され、瞬時周波数は t = 0 では 200 Hz であり、t = 1 秒では 600 Hz になります。

  3. 3 番目のチャネルは対数チャープで構成され、瞬時周波数は t = 0 では 300 Hz であり、t = 1 秒では 500 Hz になります。

長さが 256 の周期的ハミング ウィンドウと 15 サンプルのオーバーラップ長を使用して、マルチチャネル信号の STFT を計算します。

fs = 1e3;
t = 0:1/fs:1-1/fs;
x = [chirp(t,100,1,300,'quadratic',45,'concave');
      chirp(t,200,1,600,'quadratic',[],'convex');
      chirp(t,300,1,500,'logarithmic')]'; 
  
[S,F,T] = stft(x,fs,'Window',hamming(256,'periodic'),'OverlapLength',15);

最初と 2 番目のチャネルについて、元のバージョンと再構成されたバージョンをプロットします。

[ix,ti] = istft(S,fs,'Window',hamming(256,'periodic'),'OverlapLength',15);

plot(t,x(:,1)','LineWidth',1.5)
hold on
plot(ti,ix(:,1)','r--')
hold off
legend('Original Channel 1','Reconstructed Channel 1')

plot(t,x(:,2)','LineWidth',1.5)
hold on
plot(ti,ix(:,2)','r--')

legend('Original Channel 2','Reconstructed Channel 2')

フェーズ ボコーダーは周波数領域にオーディオを変換することによってタイム ストレッチとピッチ スケーリングを実行します。この図は、フェーズ ボコーダーの実装に関連する操作を示します。

フェーズ ボコーダーは、ホップ サイズ R1 の解析ウィンドウを持つ信号の STFT を実行してから、ホップ サイズ R2 の合成ウィンドウを持つ ISTFT を実行します。ボコーダーは WOLA 方法を利用します。信号のタイム ストレッチを行うために、解析ウィンドウは合成よりも大きい数のオーバーラップ サンプルを使用します。結果として、周波数成分が同じだけ残っているにもかかわらず、入力時点よりも出力時点でのサンプルが多くなります (NS,Out>NS,In)。高いサンプルレートで再生することでこの信号のピッチ スケーリングができます。つまり、元の持続時間のまま、より高いピッチで信号を作成します。

8192 Hz でサンプリングされたヘンデルの「ハレルヤ コーラス」の断章を含むオーディオ ファイルを読み込みます。

load handel

長さが 512 のルート-ハン ウィンドウを設計します。解析オーバーラップの長さを 192、合成オーバーラップの長さを 166 に設定します。

wlen = 512;
win = sqrt(hann(wlen,'periodic'));
noverlapA = 192;
noverlapS = 166;

オーバーラップ 192 の解析ウィンドウとオーバーラップ 166 の合成ウィンドウを使用して、フェーズ ボコーダーを実装します。

S = stft(y,Fs,'Window',win,'OverlapLength',noverlapA);
iy = istft(S,Fs,'Window',win,'OverlapLength',noverlapS);


%To hear, type soundsc(w,Fs), pause(10), soundsc(iw,Fs);

解析ウィンドウと合成ウィンドウが同じでオーバーラップの長さが変更される場合、調整が必要になるゲイン/損失が追加されます。これは、フェーズ ボコーダーを実装するための一般的な方法です。

ホップ率を計算して、それを使用して再構成後の信号のゲインを調整します。ホップ率を使用して、ピッチ シフトしたデータの周波数も計算します。

hopRatio = (wlen-noverlapS)/(wlen-noverlapA);
iyg = iy*hopRatio;
Fp = Fs*hopRatio;


%To hear, type soundsc(iwg,Fs), pause(15), soundsc(iwg,Fp);

元の信号と固定されたゲインでタイム ストレッチした信号をプロットします。

plot((0:length(iyg)-1)/Fs,iyg,(0:length(y)-1)/Fs,y)
xlabel('Time (s)')
xlim([0 (length(iyg)-1)/Fs])
legend('Time Stretched Signal with Fixed Gain','Original Signal','Location','best')

タイム ストレッチした信号とピッチ シフトした信号を同じプロットで比較します。

plot((0:length(iy)-1)/Fs,iy,(0:length(iy)-1)/Fp,iy)
xlabel('Time (s)')
xlim([0 (length(iyg)-1)/Fs])
legend('Time Stretched Signal','Pitch Shifted Signal','Location','best') 

ピッチ シフトするデータの効果をより理解するには、2 秒を超える次の周波数 Fs の正弦波について考えます。

t =  0:1/Fs:2;
x = sin(2*pi*10*t);

オーバーラップの長さが 192 の短時間フーリエ変換と 166 の逆短時間フーリエ変換を計算します。

Sx = stft(x,Fs,'Window',win,'OverlapLength',noverlapA);
ix = istft(Sx,Fs,'Window',win,'OverlapLength',noverlapS);

1 つのプロットに元の信号、別のプロットにタイム ストレッチした信号とピッチ シフトした信号をプロットします。

subplot(2,1,1)
plot((0:length(ix)-1)/Fs,ix,'LineWidth',2) 
xlabel('Time (s)')
ylabel('Signal Amplitude')
xlim([0 (length(ix)-1)/Fs])
legend('Time Stretched Signal') 



subplot(2,1,2)
hold on
plot((0:length(x)-1)/Fs,x)
plot((0:length(ix)-1)/Fp,ix,'--','LineWidth',2)
legend('Original Signal','Pitch Shifted Signal','Location','best')
hold off
xlabel('Time (s)')
ylabel('Signal Amplitude')
xlim([0 (length(ix)-1)/Fs])

周波数 1 kHz、持続時間 2 秒間の複素正弦波を生成します。

fs = 1e3;
ts = 0:1/fs:2-1/fs;

x = exp(2j*pi*100*cos(2*pi*2*ts));

長さが 100 の周期的ハン ウィンドウを設計し、オーバーラップ サンプルの数を 75 に設定します。COLA 準拠に対するウィンドウとオーバーラップの長さを確認します。

nwin = 100;
win = hann(nwin,'periodic');
noverlap = 75;

tf = iscola(win,noverlap)
tf = logical
   1

信号をゼロ パディングして、エッジの影響を除去します。切り捨てを避けるために、(length(xZero)-noverlap)(nwin-noverlap) が整数になるように入力信号をゼロでパディングします。FFT 長を 128 に設定します。複素信号の短時間フーリエ変換を計算します。

nPad = 100;
xZero = [zeros(1,nPad) x zeros(1,nPad)];
fftlen = 128;
s = stft(xZero,fs,'Window',win,'OverlapLength',noverlap,'FFTLength',fftlen);

逆短時間フーリエ変換を計算して、完全再構成のためにゼロを削除します。

[is,ti] = istft(s,fs,'Window',win,'OverlapLength',noverlap,'FFTLength',fftlen);
is(1:nPad) = [];
is(end-nPad+1:end) = [];
ti = ti(1:end-2*nPad);

元の信号と再構成後の信号の実数部をプロットします。信号の虚数部も完全に再構成されます。

plot(ts,real(x))
hold on
plot(ti,real(is),'--')
xlim([0 0.5])
xlabel('Time (s)')
ylabel('Amplitude (V)')
legend('Original Signal','Reconstructed Signal') 
hold off

2 kHz で 1 秒間サンプリングされた正弦波を生成します。

fs = 2e3;
t = 0:1/fs:1-1/fs;
x = 5*sin(2*pi*10*t);

長さが 120 の周期的ハミング ウィンドウを設計します。オーバーラップが 80 サンプルのウィンドウの COLA 制約を確認します。ウィンドウとオーバーラップの組み合わせは COLA 準拠です。

win = hamming(120,'periodic');
noverlap = 80;
tf = iscola(win,noverlap)
tf = logical
   1

FFT 長を 512 に設定します。短時間フーリエ変換を計算します。

fftlen = 512;
s = stft(x,fs,'Window',win,'OverlapLength',noverlap,'FFTLength',fftlen);

逆短時間フーリエ変換を計算します。

[X,T] = istft(s,fs,'Window',win,'OverlapLength',noverlap,'FFTLength',fftlen,'Method','ola','ConjugateSymmetric',true);

元の信号と再構成後の信号をプロットします。

plot(t,x,'b')
hold on
plot(T,X,'-.r')
xlabel('Time (s)')
ylabel('Amplitude (V)')
title('Original and Reconstructed Signal')
legend('Original Signal','Reconstructed Signal')
hold off

入力引数

すべて折りたたむ

短時間フーリエ変換。行列または 3 次元配列として指定します。単一チャネル信号の場合は、s を列方向に時間が増加し、行方向に周波数が増加する行列として指定します。マルチチャネル信号の場合は、s を 3 番目の次元がチャネルに対応する 3 次元配列として指定します。周波数ベクトルと時間ベクトルは stft の出力として取得されます。

メモ

istft を使用して s を反転させて、その結果を x と同じ長さになるようにする場合、(length(x)-noverlap)/(length(window)-noverlap) の値は整数でなければなりません。

データ型: double | single
複素数のサポート: あり

サンプルレート (ヘルツ単位)。正のスカラーとして指定します。

データ型: double | single

サンプル時間。duration スカラーで指定します。

例: seconds(1) は、連続する信号サンプル間の 1 秒間の時間差を表す duration スカラーです。

データ型: duration

名前と値のペアの引数

オプションの引数 Name,Value のコンマ区切りペアを指定します。Name は引数名で、Value は対応する値です。Name は引用符で囲まなければなりません。Name1,Value1,...,NameN,ValueN のように、複数の名前と値のペアの引数を、任意の順番で指定できます。

例: istft(s,'Window',win,'OverlapLength',50,'FFTLength',128) は、隣接するセグメント間で 50 サンプルがオーバーラップし 128 点 FFT をもつウィンドウ win を使用してデータをウィンドウ処理します。

ウィンドウ処理関数。'Window' とベクトルから構成されるコンマ区切りのペアとして指定します。ウィンドウを指定しない場合、またはウィンドウを空として指定する場合、関数は長さが 128 の周期的ハン ウィンドウを使用します。Window の長さは 2 以上でなければなりません。

利用可能なウィンドウについては、ウィンドウを参照してください。

例: hann(N+1)(1-cos(2*pi*(0:N)'/N))/2 は、いずれも長さ N + 1 のハン ウィンドウを指定します。

データ型: double | single

オーバーラップするサンプル数。OverlapLengthwindow の長さより小さな正の整数で構成されるコンマ区切りのペアとして指定します。'OverlapLength' を省略するか、空として指定する場合、ウィンドウの長さが 75% より少ない最大整数に設定されます。これは、既定のハン ウィンドウの 96 サンプルになります。

データ型: double | single

DFT 点の数。FFTLength と正の整数で構成されるコンマ区切りのペアとして指定します。完全な時間領域の再構成を行うには、stft で使用されるものに一致するように FFTLength を設定する必要があります。

データ型: double | single

オーバーラップ加算の方法。Method と以下のいずれかで構成されるコンマ区切りのペアとして指定します。

  • 'wola' — 重み付きオーバーラップ加算

  • 'ola' — オーバーラップ加算

元の信号の共役対称性。ConjugateSymmetric と、true または false で構成されるコンマ区切りのペアとして指定します。このオプションを true に設定する場合、istft は、入力 s が対称であると仮定します。そうでない場合は、対称性の仮定は作成されません。s が丸め誤差のために厳密に共役対称性でない場合、名前と値のペアの true への設定では STFT は共役対称として扱われます。s が共役対称である場合、逆変換の計算が速くなり、出力は実数になります。

周波数範囲。Centered と、true または false で構成されるコンマ区切りのペアとして指定します。このオプションを true に設定した場合、スペクトルは中央揃えになり、計算区間は -π ~ π です。そうでない場合、スペクトルの計算区間は 0 ~ 2π です。

入力時間次元。InputTimeDimension と、acrosscolumns または downrows から構成されるコンマ区切りのペアとして指定します。この値が downrows に設定されている場合、istfts の時間次元が行に沿い、周波数が列に沿うと見なします。この値が acrosscolumns に設定されている場合、関数 istfts の時間次元が列に沿い、周波数次元が行に沿うと見なします。

出力引数

すべて折りたたむ

時間領域の再構成後の信号。ベクトルまたは行列として返されます。

データ型: single | double

時点。ベクトルとして返されます。

  • サンプルレート fs が指定された場合、t には時間値 (秒単位) が含まれます。

  • 持続時間 ts が指定された場合、t は、時間形式が入力持続時間と同じになり、duration 配列になります。

  • 時間情報が指定されない場合、t はサンプル数を含みます。

データ型: double | single

詳細

すべて折りたたむ

逆短時間フーリエ変換

逆短時間フーリエ変換は、STFT の各 DFT ベクトルの IFFT を取り、逆変換した信号をオーバーラップ加算することで計算されます。ISTFT は次のように計算されます。

x(n)=1/21/2m=Xm(f)ej2πfndf=m=1/21/2Xm(f)ej2πfndf=m=xm(n)

ここで、R は連続する DFT 間のホップ サイズで、Xm は時間 mRxm(n)=x(n)g(nmR) 付近を中心としたウィンドウが適用されたデータの DFT です。逆 STFT は、"解析ウィンドウ" g(n) が元の信号のウィンドウ処理に使用されて c が定数の m=ga+1(nmR)=cn と同じ長さの元の信号の完全再構成です。次の図は、元の信号の再構成で使用されているステップを示しています。

定数オーバーラップ加算 (COLA) 制約

変更していないスペクトルを正しく再構成させるには、解析ウィンドウが COLA 制約を満たさなければなりません。通常、解析ウィンドウが条件 m=ga+1(nmR)=cn を満たす場合、ウィンドウは COLA 準拠であると見なされます。さらに、COLA 準拠は弱いまたは強いのいずれかで説明されます。

  • 弱い COLA 準拠とは、

    G(fk)=0,k=1,2,,R1,fkkR.

    となるフレーム レートの高調波で解析ウィンドウのフーリエ変換がゼロであることを意味します。

    エイリアスのキャンセルはスペクトル変更によって妨害されます。弱い COLA は、周波数領域のエイリアスのキャンセルに依存しています。したがって、信号に何らかのスペクトル変更が行われていなければ、弱い COLA 準拠ウィンドウを使用して完全な再構成が可能です。

  • 強い COLA 準拠の場合、ウィンドウのフーリエ変換は

    G(f)=0,f12R.

    となるフレーム レートによってダウンサンプリングしながら常に帯域制限しなければなりません。この方程式は、強い COLA 制約によって許可されるエイリアシングがないことを示します。さらに、強い COLA 準拠の場合、定数 c の値は、1 と等しくなければなりません。通常、短時間のスペクトルが何らかの方法で変更される場合、より強い COLA 準拠ウィンドウが推奨されます。

関数 iscola を使用して、弱い COLA 準拠を確認できます。COLA 準拠を確認するために使用される加算数は、ウィンドウの長さとホップ サイズによって決まります。通常、重み付きオーバーラップ加算 (WOLA) には、m=ga+1(nmR)=cna=1、オーバーラップ加算 (OLA) には、a=0 を使用するのが一般的です。既定では、istft は、オーバーラップ加算方法を実行する前に "合成ウィンドウ" を適用して、WOLA 方法を使用します。

通常、合成ウィンドウは解析ウィンドウと同じです。強い OLA ウィンドウの平方根を取ることで有用な WOLA ウィンドウを構築できます。すべての非負の OLA ウィンドウに対してこの方法を使用できます。たとえば、ルートーハン ウィンドウは WOLA ウィンドウの良い例です。

完全再構成

一般に、入力信号の STFT を計算して逆変換しても、完全再構成は行われません。ISTFT の出力を元の入力信号と可能な限り一致させるには、信号とウィンドウが以下の条件を満たす必要があります。

  • 入力サイズ — istft を使用して stft の出力を逆変換して、その結果を入力信号 x と同じ長さにする場合、k = (length(x)-noverlap)(length(window)-noverlap) の値は整数でなければなりません。

  • COLA 準拠 — COLA 準拠ウィンドウを使用して、信号の短時間フーリエ変換を変更していないと仮定します。

  • パディング — k の値が整数にならないような入力信号の長さの場合、短時間フーリエ変換を計算する前に信号をゼロ パディングします。信号を逆変換した後に、余分なゼロを削除します。

参照

[1] Crochiere, R. E. “A Weighted Overlap-Add Method of Short-Time Fourier Analysis/Synthesis.” IEEE Transactions on Acoustics, Speech and Signal Processing. Vol. 28, No. 1, Feb 1980, pp. 99–102.

[2] A. D. Gotzen, N. Bernardini and D. Arfib, “Traditional Implementations of a Phase-Vocoder: The Tricks of the Trade”, Proceedings of the COST G-6 Conference on Digital Audio Effects (DAFX-00), Verona, Italy, Dec 7-9, 2000.

[3] Griffin, D. W., and J. S. Lim. “Signal Estimation from Modified Short-Time Fourier Transform.” IEEE Transactions on Acoustics, Speech and Signal Processing. Vol. 32, No. 2, April 1984, pp. 236–243.

[4] Portnoff, M. R. “Time-Frequency Representation of Digital Signals and Systems Based on Short-Time Fourier analysis.” IEEE Transactions on Acoustics, Speech and Signal Processing. Vol. 28, No. 1, Feb 1980, pp. 55–69.

[5] Smith, J. O. Spectral Audio Signal Processing. https://ccrma.stanford.edu/~jos/sasp/, online book, 2011 edition, accessed Nov 2018.

[6] Sharpe, Bruce. Invertibility of Overlap-Add Processing. https://gauss256.github.io/blog/cola.html, accessed July 2019.

拡張機能

R2019a で導入