1

みなさん、お疲れ様でした!

いくつかの既知のインパルス応答で畳み込みを観察することから元の信号を取得するという基本的な問題に取り組んでみました。

しかし、私が得た結果はどういうわけか完全に間違っており、おそらくここではさまざまな間違ったステップが混在しています。私はすでにここやdeveloppezのような他のサイトで同様のトピックを調べましたが、原因を突き止めることができませんでした。助けていただければ幸いです。

私の真の信号f [。]が時間1のインパルスであり、インパルス応答g [。]がガウスであるとしましょう。私はそれらの畳み込みh [。]を、で計算し、基本的に、これがf [。]であることを期待して、conv()を見つけたいと思います。ifft( fft[h]./fft[g] )

最初の問題は、 n + m-1要素conv()の配列を作成することです。ここで、 n、mは引数配列の長さです。したがって、実行するには、gの長さでsmthを実行する必要があります。それは私が間違って行動するかもしれない最初の疑わしい場所です(コードを参照してください)。それを行う正しい方法は何ですか?fft[h]./fft[g]

2番目の問題は、最初の真の信号とは非常に異なるものが得られることです。

3つ目の問題は、信号シフトに対処する方法がわからないことです。matlabでは、正の時間信号を操作する必要がありますが、たとえば、ガウスインパルス応答には時間の負の要素と時間の正の要素の両方があるため、ここで操作するには、「前方」にシフトする必要があります(ピークは右に移動します)、結果を「シフト解除」する必要がありますか?

ありがとう!

これが私のがらくたです:)

close all;

TrueSignal = zeros( 101, 1 ); % impulse in t = 1.
TrueSignal( 1 ) = 1;
ImpulseResp = normpdf(-1:0.02:1)/normpdf( 0 ); % 101 elements array

figure;
subplot( 2,2,1 );
title('True signal')
plot( TrueSignal );
subplot( 2,2,2 );
title('Impulse response')
plot( ImpulseResp );

Conv = conv( TrueSignal, ImpulseResp ); % produces 201 elements array.
subplot( 2,2,3 );
title('Convolution')
plot( Conv );

% Wrong? I need a 201 elements array to represent the impulse response.
ImpulseResp_sparse = normpdf( -1:0.01:1 )/normpdf( 0 );
FIR = fft( ImpulseResp_sparse )/201;

Inverse = ifft( fft( Conv )./FIR ); % UPD Added fft() according to one of comments, bad mistake, but still not preventing.

subplot( 2,2,4 );
title('What is that???')
plot( abs( Inverse ) ); % It's weird! With no abs(), result is even more weird! 
4

2 に答える 2

2
  1. 結果の長さの問題を回避するfilter(g, 1, f)代わりに、を使用することをお勧めします。conv(g, f)または、結果の配列を切り取ります。
  2. 1つのfftを作成するのを忘れたようです。Inverse = ifft( fft(Conv)./FIR );

これは私のために働きます:

   close all;

   TrueSignal = zeros( 101, 1 ); % impulse in t = 1.
   TrueSignal( 1 ) = 1;
   ImpulseResp = normpdf(-1:0.02:1)/normpdf( 0 ); % 101 elements array

   figure;
   subplot( 2,2,1 );
   title('True signal')
   plot( TrueSignal );
   subplot( 2,2,2 );
   title('Impulse response')
   plot( ImpulseResp );

   Conv = filter( TrueSignal, 1, ImpulseResp ); 
   subplot( 2,2,3 );
   title('Convolution')
   plot( Conv );

   fftConv = fft(Conv);

   FIR = fft( ImpulseResp );

   Inverse = ifft( fftConv./FIR );

   subplot( 2,2,4 );
   plot( abs( Inverse ) );
  1. 言われているように、FFTは信号が周期的であると想定しています。非周期信号を分析する場合は、「ウィンドウ化されたFFT」アルゴリズムを使用して品質を向上させる必要があります。これはほとんど同じですが、入力に特別なウィンドウ関数( Blackmanなど)を掛ける必要があります。パディングの追加も機能しますが、最も計算コストの高い方法です。
于 2013-03-11T22:49:00.100 に答える
2

fftfor convolutionを簡単に使用すると、循環畳み込みが発生しますが、必要なもの(および実行するものconv)は線形畳み込みです。したがって、このようなスキームをで実装するfftには、信号を長さまでゼロパッドする必要がありますm+n-1

convこれは、線形畳み込みの出力とfftベースの線形畳み込みの間の同等性を示す例です。

x=rand(4,1);y=rand(3,1); %sample data
out1=conv(x,y);          %output from conv()
X=fft(x,6);Y=fft(y,6);   %zero pad and compute fft
out2=ifft(X.*Y);         %output from fft based lin. conv.

out1あなたはそれをチェックすることができout2、同じです(FP精度の範囲内で)。

この方法で問題を再定式化すれば、うまくいくはずです。私はあなたがre:shiftsに何を求めているのか理解できませんが、あなたはとを調べたいと思うかもしれませfftshiftifftshift

于 2013-03-11T23:04:04.967 に答える