IFFT of custom frequency array - confusion with the time domain sampling frequency

Hi all,

I am another guy getting stuck with the same old stuff! The relationship
between point in a frequency response and the sampling frequency/time
characteristics of the impulse response.

I have created a custom frequency array, representing a frequency response
from an algorithm, between 0 and 300Hz, with 0.2Hz resolution, ie 1500
points. I then need to create the impulse response from this to convolve
with a .wav file.

I am having trouble working out the time scale of the resulting impusle
response. Here is the code I have written so far:


fl=0.0001;       %Lower frequency
fu=300;          %Upper frequency
fi=0.2;          %Increment in frequency
f=fl:fi:fu;      %frequency array
points=size(f,2);

%(algorithm here produces a value for each f, resulting in the vector
'press')

%**********************************************************************
%Mirror Response then IFFT to create Impulse Repsonse
%**********************************************************************
press_mirror=cat(2,press,fliplr(press));
press_mirror( ((points/2)+1) : points*2 ) = conj( press_mirror(
((points/2)+1) : points*2 ));

time_mirror=ifft(press_mirror);		           %IFFT
inverse=real(time_mirror(1:size(time_mirror,2)/2));%take real of 1st half
inverse=inverse';                      		   %rotate
inverse=inverse/(max(abs(inverse)));   		   %normalise IR
%**********************************************************************
%**********************************************************************

fs=44100;
fs_lf=600;		%This is the bit I'm confused with

[p,q]=rat(fs/fs_lf,0.000001);    %resample to 44100Hz for .wav file
inverse=resample(inverse,p,q); 
wavwrite(inverse,44100,16,'out_impulse.wav');


Now I've come to realise that changing the fs_lf value here results in a
differing time scale. I originally had 1200 which is twice the highest
frequency I will be dealing with. I do not know how to get the correct time
information from the impulse response. Is the sampling frequency 1500. If
it is, could I only get 1 seconds worth of information as there are 1500
points in the impulse response ('inverse')?

Any help would be hugely appreciated.
Many thanks,
Matt

"matt_w" <mat.w@btinternet.com> wrote in message 
news:As6dnUTkfOf1JYTVnZ2dnUVZ_t2inZ2d@giganews.com...
> Hi all,
>
> I am another guy getting stuck with the same old stuff! The relationship
> between point in a frequency response and the sampling frequency/time
> characteristics of the impulse response.
>
> I have created a custom frequency array, representing a frequency response
> from an algorithm, between 0 and 300Hz, with 0.2Hz resolution, ie 1500
> points. I then need to create the impulse response from this to convolve
> with a .wav file.
>
> I am having trouble working out the time scale of the resulting impusle
> response. Here is the code I have written so far:
>
>
> fl=0.0001;       %Lower frequency
> fu=300;          %Upper frequency
> fi=0.2;          %Increment in frequency
> f=fl:fi:fu;      %frequency array
> points=size(f,2);
>
> %(algorithm here produces a value for each f, resulting in the vector
> 'press')
>
> %**********************************************************************
> %Mirror Response then IFFT to create Impulse Repsonse
> %**********************************************************************
> press_mirror=cat(2,press,fliplr(press));
> press_mirror( ((points/2)+1) : points*2 ) = conj( press_mirror(
> ((points/2)+1) : points*2 ));
>
> time_mirror=ifft(press_mirror);            %IFFT
> inverse=real(time_mirror(1:size(time_mirror,2)/2));%take real of 1st half
> inverse=inverse';                         %rotate
> inverse=inverse/(max(abs(inverse)));      %normalise IR
> %**********************************************************************
> %**********************************************************************
>
> fs=44100;
> fs_lf=600; %This is the bit I'm confused with
>
> [p,q]=rat(fs/fs_lf,0.000001);    %resample to 44100Hz for .wav file
> inverse=resample(inverse,p,q);
> wavwrite(inverse,44100,16,'out_impulse.wav');
>
>
> Now I've come to realise that changing the fs_lf value here results in a
> differing time scale. I originally had 1200 which is twice the highest
> frequency I will be dealing with. I do not know how to get the correct 
> time
> information from the impulse response. Is the sampling frequency 1500. If
> it is, could I only get 1 seconds worth of information as there are 1500
> points in the impulse response ('inverse')?
>

Matt,

The time span is the reciprocal of the frequency resolution.
So it would be the reciprocal of 0.2Hz or 5 seconds if everything else is 
correct.

However, you say the span in frequency is 0 to 300Hz to get the 1500 points 
at .2Hz resolution i.e. 300/1500=0.2.  But the sample rate must be at least 
600Hz here - which would require 3,000 points.

There must be samples from 0 to fs-0.2
I see code that says fs=44,100.
That would make 220,500 samples.
There are other issues with just the code it seems.... have you run it at 
all?

So, for starters you need to create a frequency array that contains both the 
positive (0 to 300Hz) part and the negative part which is -300Hz to 0. 
Normally the negative part is the upper half of the array so it goes from 
300Hz to fs and the value at fs is the same value as at zero - and is not 
included in the array.  The frequency array is *periodic* - so this is how 
it's represented for a single period in frequency.

Also, if the time domain transform is to be purely real then the frequency 
domain values need to be real and even around zero and fs or, to be more 
general, need to be:
- real and even
- imaginary and odd.

Fred 

> Matt,
>
> The time span is the reciprocal of the frequency resolution.
> So it would be the reciprocal of 0.2Hz or 5 seconds if everything else is
> correct.

Ok, that's good. Am I right in thinking I'll only take half of the
output of the IFFT (as it should be mirrored), meaning with a
frequency resoultion of 0.2, I can get 2.5 seconds of impulse
response?

>
> However, you say the span in frequency is 0 to 300Hz to get the 1500 points
> at .2Hz resolution i.e. 300/1500=0.2.  But the sample rate must be at least
> 600Hz here - which would require 3,000 points.

I do have 1500 points (300/0.2), but then I mirror this in the code
(using the fliplr then creating the conjugate values for the second
half), resulting in 3000 points. Is this correct? Giving me a sample
rate of 600Hz? But then I thought this equation was the case: df=fs/N
and N now equals 6000 after mirroring, df still equals 0.2, so fs=1200
- oh. hmmm.??!

>
> There must be samples from 0 to fs-0.2
> I see code that says fs=44,100.

fs in the code is just for resampling up to that frequency to create
a .wav file, sorry, the sampling frequency for all the frequency and
IFFT stuff in my code above is fs_lf (sorry for confusion)

> That would make 220,500 samples.
> There are other issues with just the code it seems.... have you run it at
> all?

I have run the code, it 'seems' to give me the desired result (when I
analyse as a .wav file) if sampling frequency (fs_lf) is 1200Hz. (from
the equation: df=fs/N )

>
> So, for starters you need to create a frequency array that contains both the
> positive (0 to 300Hz) part and the negative part which is -300Hz to 0.
> Normally the negative part is the upper half of the array so it goes from
> 300Hz to fs and the value at fs is the same value as at zero - and is not
> included in the array.  The frequency array is *periodic* - so this is how
> it's represented for a single period in frequency.

I thought I had done this with the mirroring of the frequency array
(press). Is this not the case?

>
> Also, if the time domain transform is to be purely real then the frequency
> domain values need to be real and even around zero and fs or, to be more
> general, need to be:
> - real and even
> - imaginary and odd.
>
> Fred

Finally, I have plotted the (real) output of the IFFT and it doesn't
come out symmetrical which I thought it did, hence removing the second
half in my code.

Thanks for your help. I really appreciate it.

Matt

Matt,

Answers interspersed below at "***"

Fred

<mattwankling-google@yahoo.co.uk> wrote in message 
news:bca9bcba-2366-42f9-8464-90e97e356e13@k13g2000hse.googlegroups.com...
>
>> Matt,
>>
>> The time span is the reciprocal of the frequency resolution.
>> So it would be the reciprocal of 0.2Hz or 5 seconds if everything else is
>> correct.
>
> Ok, that's good. Am I right in thinking I'll only take half of the
> output of the IFFT (as it should be mirrored), meaning with a
> frequency resoultion of 0.2, I can get 2.5 seconds of impulse
> response?
>

***No.  That would be closer to correct for the results of an FFT but not 
for an IFFT.  You should be very careful about throwing out what seems 
"redundant" because it really isn't.  That's because generally:

***The time series is real but not purely even or odd.
Thus the frequency series is complex with the real part even and the 
imaginary part odd around zero/fs and is zero at fs/2.

***If you start with a purely real frequency sequence that is also even then 
the IFFt will be real and even.  That's what might suggest you could "throw 
away" half of the result.  But if you do that then you are effectively 
increasing the frequency resolution by a factor of 2 to 0.4Hz *and* are 
changing the frequency response to one that's no longer real and even.  It's 
pretty easy to show that the magnitude response is also not at all the same.

***The time span is 5 seconds if the frequency resolution is 0.2Hz.

>>
>> However, you say the span in frequency is 0 to 300Hz to get the 1500 
>> points
>> at .2Hz resolution i.e. 300/1500=0.2.  But the sample rate must be at 
>> least
>> 600Hz here - which would require 3,000 points.
>
> I do have 1500 points (300/0.2), but then I mirror this in the code
> (using the fliplr then creating the conjugate values for the second
> half), resulting in 3000 points. Is this correct? Giving me a sample
> rate of 600Hz? But then I thought this equation was the case: df=fs/N
> and N now equals 6000 after mirroring, df still equals 0.2, so fs=1200
> - oh. hmmm.??!
>

***Let's get all these numbers straight:
Time span 5 seconds.
Frequency resolution 0.2Hz.
Bandwidth 300Hz.
Sample rate 600Hz.
N=600/0.2= 3,000

df=fs/N = 600/3000 = 0.2

***Mind you, "N" refers to the number of complex samples in general so you 
would have 2*N values for N samples.  Is that part of the confusion re: a 
factor of 2?

>>
>> There must be samples from 0 to fs-0.2
>> I see code that says fs=44,100.
>
> fs in the code is just for resampling up to that frequency to create
> a .wav file, sorry, the sampling frequency for all the frequency and
> IFFT stuff in my code above is fs_lf (sorry for confusion)

***I don't "get" fs-lf = fs .  lf has nothing to do with anything we've 
dealt with above.

>
>> That would make 220,500 samples.
>> There are other issues with just the code it seems.... have you run it at
>> all?
>
> I have run the code, it 'seems' to give me the desired result (when I
> analyse as a .wav file) if sampling frequency (fs_lf) is 1200Hz. (from
> the equation: df=fs/N )
>
>>
>> So, for starters you need to create a frequency array that contains both 
>> the
>> positive (0 to 300Hz) part and the negative part which is -300Hz to 0.
>> Normally the negative part is the upper half of the array so it goes from
>> 300Hz to fs and the value at fs is the same value as at zero - and is not
>> included in the array.  The frequency array is *periodic* - so this is 
>> how
>> it's represented for a single period in frequency.
>
> I thought I had done this with the mirroring of the frequency array
> (press). Is this not the case?

***If you did it correctly then presumably yes.  Just to be clear:
If you start with a purely real frequency response then it will mirror 
around zero - meaning that it mirrors left/right around fs/2.  If the array 
is of even length then there is a single value at fs/2 and if the array is 
of odd length then there are two equal values around fs/2 with fs/2 midway 
between those samples.  Like this:
Even: [1 2 3 4 5 6 5 4 3 2] length = 10 value at fs/2 = 6; value at zero and 
fs = 1.
Odd: [ 1 2 3 4 5 6 6 5 4 3 2] length = 11 value near fs/2 = 6; value at zero 
and fs = 1.
So neither of them fliplr exactly at fs/2 ....
>
>>
>> Also, if the time domain transform is to be purely real then the 
>> frequency
>> domain values need to be real and even around zero and fs or, to be more
>> general, need to be:
>> - real and even
>> - imaginary and odd.
>>
>> Fred
>
> Finally, I have plotted the (real) output of the IFFT and it doesn't
> come out symmetrical which I thought it did, hence removing the second
> half in my code.

***If the real part of the IFFT is not even then that means that the 
frequency samples weren't purely real.  Since there appears to be no reason 
why your frequency response would have imaginary values to start with, I'm 
going to assume that they are zero.  Then you have only one possibilities:
- frequency samples are (real and) even
If the frequency samples are (real and) even then the IFFT will be real and 
even.

***If the frequency samples are complex and properly formed so the real part 
is even and the imaginary part is odd, then the IFFT should be purely real 
and need not be even..

>
> Thanks for your help. I really appreciate it.
>
> Matt
> 

Fred Marshall wrote:
(snip)

>>Ok, that's good. Am I right in thinking I'll only take half of the
>>output of the IFFT (as it should be mirrored), meaning with a
>>frequency resoultion of 0.2, I can get 2.5 seconds of impulse
>>response?

> ***No.  That would be closer to correct for the results of an FFT but not 
> for an IFFT.  You should be very careful about throwing out what seems 
> "redundant" because it really isn't.  That's because generally:

> ***The time series is real but not purely even or odd.
> Thus the frequency series is complex with the real part even and the 
> imaginary part odd around zero/fs and is zero at fs/2.

You might be interested in the Hartley transform, and
corresponding FHT:

http://en.wikipedia.org/wiki/Discrete_Hartley_transform

where the transform of real data is real, of even data is
even, and of odd data is odd.  (And of imaginary data
is imaginary.)

The web page indicates that FHT is about the same number
of operations as FFT optimized for real data, but might
be faster in some cases.

-- glen