Using IFFT for obtaining time response of measured freq response

Hello, I have numerical frequency response data (G(s=j2*pi*f), f) for a system. I want to obtain the impulse response in time domain. I should be able to obtain this using the IFFT function. But the scaling is not clear to me. Here is the code I am using. What am I missing? Thanks. Amit.

    clc; close all;
    fmax        = 1e6; 
    L           = 1024;  
    fdelta      = 2*fmax/(L-1)
    fgrid       = [-fmax:fdelta:fmax]; %fmax = fs/2 => deltaT = 1/(2*fmax)
    tau         = 1/(2*pi*5e4)
       %  Y(s) = 1/(1+s.tau) Given Freq Domain data
    Y           = 1./(1+j*2*pi*fgrid*tau);  
    y_t         = ifft((Y));
    Tdelta      = 1/(2*fmax);
    t           = [0:Tdelta:Tdelta*(length(y_t)-1)];

     % compare with known time domain function
    figure; plot(t, y_t, 'bx', t, exp(-t/tau), 'r'); 
    grid; axis([0 1e-4 0 1])

Answer 1

[Rick Rosson]

Please try:

   Fs = 2e6;
   L = 1024;
   dF = Fs/L;
   f = (-Fs/2:dF:Fs/2-dF)';
   s = j*2*pi*f;

   Fc = 50e3;
   alpha = 2*pi*Fc;
   tau = 1/alpha;
   Y = alpha./(s+alpha);

   dt = 1/Fs;
   t = dt*(0:L-1)';
   y = L*ifft(ifftshift(Y));

   x = exp(-t/tau);

HTH.

Rick

[Amit]

Rick,

Thanks for the reply.

I tried the code you suggested. I still get a scaling error.

>> [(x(1:10)) abs(y(1:10))']

ans =

1.0000e+000 77.8669e+000

854.6360e-003 152.1368e+000

730.4027e-003 109.5818e+000

624.2284e-003 105.7690e+000

533.4881e-003 81.7647e+000

455.9381e-003 76.5834e+000

389.6611e-003 59.9684e+000

333.0184e-003 55.8890e+000

284.6095e-003 43.7460e+000

243.2376e-003 40.9324e+000

>> [(x(1:10))./abs(y(1:10))']

ans =

12.8424e-003

5.6176e-003

6.6654e-003

5.9018e-003

6.5247e-003

5.9535e-003

6.4978e-003

5.9586e-003

6.5060e-003

5.9424e-003

>> 2*pi

ans =

6.2832e+000

>> 2*pi/1024

ans =

6.1359e-003

Looks like there should be scaling factor of 2*pi instead of L=1024.

However, increasing the sampling frequency changes the scaling.

With fs=200e6.

>> [(x(1:10))./abs(y(1:10))']

828.7821e-003

463.9740e-003

517.4152e-003

484.3603e-003

507.4654e-003

489.3110e-003

504.0767e-003

491.5127e-003

502.3784e-003

492.7540e-003

I am really not sure what is happening. I wish there were a matlab function to take care of all this and just give me a time domain impulse response!

Amit.

[LI QZH]

Y = alpha./(s+alpha); i think this is not the right Laplace transform

y=1./(s+alpha) is the right form

am i right ? thanks

Answer 2

[Rick Rosson]

Do you have access to either the Control Systems Toolbox or the Signal Processing Toolbox? If so, which one (or both)?

Answer 3

[Rick Rosson]

Please check this related answer:

Scaling the FFT and the IFFT

HTH.

Answer 4

[Rick Rosson]

I modified the code I posted earlier to correct the scaling factor:

   Fs = 2e6;
   L = 1024;
   dF = Fs/L;
   f = (-Fs/2:dF:Fs/2-dF)';
   s = j*2*pi*f;

   Fc = 50e3;
   alpha = 2*pi*Fc;
   tau = 1/alpha;
   Y = alpha./(s+alpha);

   dt = 1/Fs;
   t = dt*(0:L-1)';
   G = 2*pi;
   y = G*abs(ifft(ifftshift(Y)));

   x = exp(-t/tau);

   figure;
   plot(t,x,t,y);

HTH.

Rick

[Greg Heath]

FFTing and IFFTing is tricky business because you have to be perfectly clear about original assumptions. For example, time domain signals contain N measurements at a sampling frequency

Fs. If the signal is real the corresponding nonnegative frequency spectrum measurements have a spacing df = Fs/N and

have either the range

f = df*[0:N/2], Neven

fmax = Fs/2

L = length(f) = N/2+1

with L odd when N is even

or

f = df*[0:(N-1)/2], N odd

fmax = Fs/2 - df/2

L = (N+1)/2

with L even when N is odd

Now since you seem to have L = 1024 measurements, my conclusion is that

N = 2*L + 1 = 2049 is odd

and

Fs = 2*N*fmax/(N-1) = 2.001e6

Hope this helps.

Greg