ifft returns NaN when plotting the impulse response function

조회 수: 28 (최근 30일)
粤轩 杨
粤轩 杨 2023년 4월 10일
댓글: 粤轩 杨 2023년 4월 10일
Hello, I seem to be having issues using MATLAB's fft and ifft functions.
An impulse response function h(t) has the following formula: inj(t) * h(t) = AIF(t). We know that the graphs of inj(t) and AIF(t) are as followed. I wrote the following code to do the deconvolution but h(t) in the output graph is zero.
I realized that the h returns from ifft is NaN but I don't know how to correct the code. Thanks!
load("AIF_1.mat");
load("inj_1.mat");
inj_1 = inj(201:2400);
inj1_FFT = fft(inj_1);
AIF_1 = AIF(201:2400);
AIF1_FFT = fft(AIF_1);
h_FFT = AIF1_FFT ./ inj1_FFT;
h_FFT(isnan(h_FFT)==1) = 0;
h = ifft(h_FFT);
X = zeros(1,200)
ht = [X,h];
plot(time,real(ht));title('h(t)');

채택된 답변

Paul
Paul 2023년 4월 10일
Hi 粤轩 杨,
It looks like h_FFT also has a few values that are inf, in addition to the NaNs
load("AIF_1.mat");
load("inj_1.mat");
inj_1 = inj(201:2400);
inj1_FFT = fft(inj_1);
AIF_1 = AIF(201:2400);
AIF1_FFT = fft(AIF_1);
h_FFT = AIF1_FFT ./ inj1_FFT;
sum(isnan(h_FFT))
ans = 3
sum(isinf(h_FFT))
ans = 2
% quick correction to get things to run, not sure if this is really what
% should be done
inj1_FFT(inj1_FFT == 0) = eps;
h_FFT = AIF1_FFT ./ inj1_FFT;
sum(isnan(h_FFT))
ans = 0
sum(isinf(h_FFT))
ans = 0
h = ifft(h_FFT);
X = zeros(1,200);
ht = [X,h];
plot(time,real(ht));title('h(t)');
  댓글 수: 1
粤轩 杨
粤轩 杨 2023년 4월 10일
Hi Paul,
The output graph is a little strange. But you finging of inf values in h_FFT made me realize that I didn't intercept the inj correctly.
After changing to
inj_1 = inj(200:2400)
the output graph is good!
Thank you for your inspiration!

댓글을 달려면 로그인하십시오.

추가 답변 (1개)

粤轩 杨
粤轩 杨 2023년 4월 10일
There are two related attachments for your information. Thanks!
  댓글 수: 1
粤轩 杨
粤轩 杨 2023년 4월 10일
Just as mentioned above, the final code should be
load("AIF_1.mat");
load("inj_1.mat");
inj_1 = inj(200:2400);
inj1_FFT = fft(inj_1);
AIF_1 = AIF(200:2400);
AIF1_FFT = fft(AIF_1);
h_FFT = AIF1_FFT ./ inj1_FFT;
h = ifft(h_FFT);
X = zeros(1,199);
ht = [X,h];
plot(time,real(ht));title('h(t)');
And the output graph looks like

댓글을 달려면 로그인하십시오.

태그

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by