filtering in freq domain

Hi all, 

To begin with thanks to all those who have made this wonderful group.
I am a mechanical engineer working on hilbert transforms. 

I  have a time domain signal sampled at 2100pts/sec for 60secs and has
signals at 3 freqs. I do a fft, get into my freq domain and then
filter all freq except the freq of interest. 

My question is 
a) What is the best way to filter in freq domain. Iam trying
butterworth filter, and get [b,a] ouput from matlab butter function.
How do use this on my data.

If i  do manage to filter say a freq which has a signal then make my
postive freqs double and negative freq zero(computing the hilbert),
and then  do a ifft and
calculate the angle. I get a plot of my angle for signal, and
similarly i get one for noise, implying i filter a freq where there is
no signal content.

My question is
b) Should'nt i get different plots for signal and noise. My noise
should be random but my plot doesnt show my noise to be random and
shows some kind of distinct pattern. 

Any help is greatly appreciated

Regards

Sundar

snsundar@olemiss.edu (sundar) writes:

> Hi all, 
>
> To begin with thanks to all those who have made this wonderful group.
> I am a mechanical engineer working on hilbert transforms. 
>
> I  have a time domain signal sampled at 2100pts/sec for 60secs and has
> signals at 3 freqs. I do a fft, get into my freq domain and then
> filter all freq except the freq of interest. 
>
> My question is 
> a) What is the best way to filter in freq domain. Iam trying
> butterworth filter, and get [b,a] ouput from matlab butter function.
> How do use this on my data.

Use Matlab's filter() function. [b,a] are the coefficients of the transfer
function for the filter. b is a vector of the numerator coefficients and
a is a vector of the denominator coefficients.

Note that frequency-domain filtering is only useful for FIR filters, i.e.,
filters in which the a vector above is 1. Butterworth filters are not
FIR filters, so you won't be able to use frequency-domain filtering on them.

Note also that I'm using the term "frequency-domain filtering" in its
conventional sense, i.e., performing a linear filter in the frequency
domain. 
-- 
%  Randy Yates                  % "So now it's getting late,
%% Fuquay-Varina, NC            %    and those who hesitate
%%% 919-577-9882                %    got no one..."
%%%% <yates@ieee.org>           % 'Waterfall', *Face The Music*, ELO
http://home.earthlink.net/~yatescr

First, if you have real data, use a real to complex FFT and you don't have to 
worry about the Hilbert Transform.  

Second, you have a small number of points, you can do a single transform of 
all the data points.  If you still need filtered time series, you can window 
in the frequency domain and do an inverse FFT.  Windowing the frequency 
regions means that you multiply the desired region by values that reflect the 
desired frequency response there.  

Third, the Butterworth Filter is a time domain filtering technique.  You don't 
have to go to the frequency domain.

In article <2b5e8f90.0411071826.3ce88c26@posting.google.com>, 
snsundar@olemiss.edu (sundar) wrote:
>Hi all, 
>
>To begin with thank         stoallthosewhohavemadethiswonderfulgroup.
>I am a mechanical engineer working on hilbert transforms. 
>
>I  have a time domain signal sampled at 2100pts/sec for 60secs and has
>signals at 3 freqs. I do a fft, get into my freq domain and then
>filter all freq except the freq of interest. 
>
>My question is 
>a) What is the best way to filter in freq domain. Iam trying
>butterworth filter, and get [b,a] ouput from matlab butter function.
>How do use this on my data.
>
>If i  do manage to filter say a freq which has a signal then make my
>postive freqs double and negative freq zero(computing the hilbert),
>and then  do a ifft and
>calculate the angle. I get a plot of my angle for signal, and
>similarly i get one for noise, implying i filter a freq where there is
>no signal content.
>
>My question is
>b) Should'nt i get different plots for signal and noise. My noise
>should be random but my plot doesnt show my noise to be random and
>shows some kind of distinct pattern. 
>
>Any help is greatly appreciated
>
>Regards
>
>Sundar

Randy Yates <yates@ieee.org> wrote in message news:<actt9icg.fsf@ieee.org>...
> snsundar@olemiss.edu (sundar) writes:
> 
> > Hi all, 
> >
> > To begin with thanks to all those who have made this wonderful group.
> > I am a mechanical engineer working on hilbert transforms. 
> >
> > I  have a time domain signal sampled at 2100pts/sec for 60secs and has
> > signals at 3 freqs. I do a fft, get into my freq domain and then
> > filter all freq except the freq of interest. 
> >
> > My question is 
> > a) What is the best way to filter in freq domain. Iam trying
> > butterworth filter, and get [b,a] ouput from matlab butter function.
> > How do use this on my data.
> 
> Use Matlab's filter() function. [b,a] are the coefficients of the transfer
> function for the filter. b is a vector of the numerator coefficients and
> a is a vector of the denominator coefficients.
> 
> Note that frequency-domain filtering is only useful for FIR filters, i.e.,
> filters in which the a vector above is 1. Butterworth filters are not
> FIR filters, so you won't be able to use frequency-domain filtering on them.

Yes you are. With the filter coefficients grouped as a and b above, 
the frequency domain IIR filter is implemented as

X=fft(x,N);           % DFT of input signal
A=fft(a,N);
B=fft(b,N); 
Y= X.*B./A;
y= real(ifft(y));

where the FFT length N is large enough to contain both the input signal 
and startup and end transients. The only difficult part here is to 
find a suitable N. Which can be hard enough, though.

Rune

Hi!

> Yes you are. With the filter coefficients grouped as a and b above,
> the frequency domain IIR filter is implemented as
>
> X=fft(x,N);           % DFT of input signal
> A=fft(a,N);
> B=fft(b,N);
> Y= X.*B./A;
> y= real(ifft(y));
>
> where the FFT length N is large enough to contain both the input signal
> and startup and end transients. The only difficult part here is to
> find a suitable N. Which can be hard enough, though.

I find it difficult to understand how can overlapp-and-add work analytically
correct. But I find it even more difficult to imagine, how could a
continuous
filter work in this case. The recursive formula indicates that the impulse
response
can be very very long. (infinite?) and then you still have to glue together
the individual pieces.

Regards!
Atmapuri.

Randy Yates <yates@ieee.org> wrote in message news:<actt9icg.fsf@ieee.org>...
> snsundar@olemiss.edu (sundar) writes:
Dear Randy 

Thanks once again for your advice. I tried using the filter function,
but iam not sure if i got what is right. Now, ok i have a signal at
380 Hz, i use a butterworth band pass filter from 379.5 to 380.5. Get
[b.a], then use it on my freq domain data. Now shouldnt my fft plot
show only one peak at 380 hz.

I think iam getting lost in terminology (freq domain filtering etc. )
above mentioned is my objective what do you think is the best method
of doing it, i tried rectangular window and hanning, they dont seem to
help.

 Use Matlab's filter() function. [b,a] are the coefficients of the
transfer
> function for the filter. b is a vector of the numerator coefficients and
> a is a vector of the denominator coefficients.
> 
> Note that frequency-domain filtering is only useful for FIR filters, i.e.,
> filters in which the a vector above is 1. Butterworth filters are not
> FIR filters, so you won't be able to use frequency-domain filtering on them.
> 
> Note also that I'm using the term "frequency-domain filtering" in its
> conventional sense, i.e., performing a linear filter in the frequency
> domain.

Thanks once again.

snsundar@olemiss.edu (sundar) writes:

> Randy Yates <yates@ieee.org> wrote in message news:<actt9icg.fsf@ieee.org>...
> > snsundar@olemiss.edu (sundar) writes:
> Dear Randy 
> 
> Thanks once again for your advice. I tried using the filter function,
> but iam not sure if i got what is right. Now, ok i have a signal at
> 380 Hz, i use a butterworth band pass filter from 379.5 to 380.5. Get
> [b.a], then use it on my freq domain data. Now shouldnt my fft plot
> show only one peak at 380 hz.

Well, it depends on what you mean by "one peak". The windowing characteristic
of the FFT will smear a peak, but the maximum should still be at 380 Hz. Note
that depending on the sample rate and FFT length, 380 Hz may or may not correspond
to a specific FFT bin, so you may have more than one bin at the peak value. 

> I think iam getting lost in terminology (freq domain filtering etc. )
> above mentioned is my objective what do you think is the best method
> of doing it, 

For right now, while using Matlab, the way you're going (using filter())
is just fine. 

It sounds like what you need to do filter a signal (whether you use a 
time-domain filter or frequency-domain filter is an implementation issue
which you probably don't need to consider right now) and look at the 
resulting response in the frequency domain. So I think you're just saying
you want to use the FFT to view the spectrum of your filter's output, which
is a valid thing to do.

> i tried rectangular window and hanning, they dont seem to
> help.

sundar, can you post your Matlab code and your data to a web page somewhere?
That way we can see exactly what's going on.
-- 
Randy Yates
Sony Ericsson Mobile Communications
Research Triangle Park, NC, USA
randy.yates@sonyericsson.com, 919-472-1124

Hi Randy and all others, 

thank you for your help and timely replies. I am really desperate to
get this done asap, my graduation depends on this.

yes,i understand this, but my binwidht is one. 
> Well, it depends on what you mean by "one peak". The windowing characteristic
> of the FFT will smear a peak, but the maximum should still be at 380 Hz. Note that depending on the sample rate and FFT length, 380 Hz may or may not correspond
to a specific FFT bin, so you may have more than one bin at the peak
value.
or right now, while using Matlab, the way you're going (using
filter())

> sundar, can you post your Matlab code and your data to a web page somewhere?
That way we can see exactly what's going on.
Well i dont have a website..let me trying copying the code here...

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% initial parameters
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% sampling rate for 60 sec
fs=2100

%fft_size
fft_size=2100

%binwidht
bwidht=fs/fft_size

%desired frequency
freq=380

%size of filter to apply
df=1

% apply a butterworth filter across the desired freqs.

w1 = (freq-df/2)/fs
w2 = (freq+df/2)/fs

% this gives me w1 at 379.5 & w2 at 380.5

% order of the butterworth filter

n=2
Wn=[w1 w2]

[b,a]=butter(n,Wn);

y=filter(b,a,data); % data(2100,60)

ffty=fft(y); % at this step here i should be able to c my peak at
380hz but dont..

% hilbert tranform remove f<0, 2*f>0;

for i=1:1:fft_size
	if i<=fft_size/2
	nfftdata(i,:)=ffty(i,:).*2;
else
	nfftdata(i,:)=0;
end
end

hildata=ifft(nfftdata,fft_size);

ang=angle(hildata);

% hence this will give me 2100 angle points for my first second data.

I need to prove that my angle should be random for noise while it
should not be random for signals. The second part works fine with a
hann windowing but i get a pattern for noise which i dont think is
logical. while with this filter function i get the same angle plot for
signal and noise.

Thanks once again.

Regards

Sundar

snsundar@olemiss.edu (sundar) writes:
> [...]
> ffty=fft(y); % at this step here i should be able to c my peak at
> 380hz but dont..

I do. Are you typing, e.g., 

plot(abs(ffty(:,1)));

? (you have 60 columns so you could use 1 to 60 for the second index).
-- 
%  Randy Yates                  % "Watching all the days go by...    
%% Fuquay-Varina, NC            %  Who are you and who am I?"
%%% 919-577-9882                % 'Mission (A World Record)', 
%%%% <yates@ieee.org>           % *A New World Record*, ELO
http://home.earthlink.net/~yatescr

yes even i do , thank you so much. 

Randy, i have one more question, after i do a 
hilbert(y)
ang=angle(y);
plot(ang(:,1))

this plot that iam getting has a pattern for a signal and is not
random for noise. Can you tell me why, am i not getting a random phase
angle distribution for that one sec data for noise. I read some of the
old posts and think that i might have something to do with my
filtering, and since my filtering widht is very narrow, all my points
are sort of correlated and hence lose their random behavior and show a
pattern.

i tried 
angun=unwrap(ang);   % unwrap my phase
phi=angun(:,1)-2*pi*freq*(0:1/fs:1)'; % remove my carrier freq
plot(phi)
phi1=atan2(sin(phi1),cos(phi1));  % convert from -pi to pi
plot(phi1)

This again leads me nowhere.

Using the hilbert transform how can i differentiate between signal and
noise.

thanks so much once again...u could mail me and post the reply since
that way it wud  be faster for us to communicate.

Thanks so much

Regards

Sundar    snsundar@olemiss.edu