DSPRelated.com
Forums

interesting Filtering problem

Started by Abhinav Kashyap October 11, 2004
I am trying to compare two data sets with equal no. of points (it has
several peaks in it).I am using correlation technique .First I filter
the noise present in signal.
For that I first FFT the data, multiply low pass filter on real part
of FFT and than take IFFT.
I do same procedure for other data set.
Than I do cross  correlation and auto correlations.

I am not sure is this right approach or not. End result what i am
getting is like a triangle around y-axis.
Has anybody used this technique earlier....any references?
I am Using Matlab for programming.

Thanks in advance
Abhinav Kashyap wrote:

> I am trying to compare two data sets with equal no. of points (it has > several peaks in it).I am using correlation technique .First I filter > the noise present in signal. > For that I first FFT the data, multiply low pass filter on real part > of FFT and than take IFFT. > I do same procedure for other data set. > Than I do cross correlation and auto correlations. > > I am not sure is this right approach or not. End result what i am > getting is like a triangle around y-axis. > Has anybody used this technique earlier....any references? > I am Using Matlab for programming. > > Thanks in advance
If the FFT has an imaginary part, you can't ignore it. You need it to put the signal back together. Lowpass both parts. There may be other, more efficient ways to filter your data, but that might not matter to you. Jerry -- Engineering is the art of making what you want from things you can get. �����������������������������������������������������������������������
"Jerry Avins" <jya@ieee.org> wrote in message
news:2t0ep6F1qdjg1U1@uni-berlin.de...
> Abhinav Kashyap wrote: > > > I am trying to compare two data sets with equal no. of points (it has > > several peaks in it).I am using correlation technique .First I filter > > the noise present in signal. > > For that I first FFT the data, multiply low pass filter on real part > > of FFT and than take IFFT. > > I do same procedure for other data set. > > Than I do cross correlation and auto correlations. > > > > I am not sure is this right approach or not. End result what i am > > getting is like a triangle around y-axis. > > Has anybody used this technique earlier....any references? > > I am Using Matlab for programming. > > > > Thanks in advance > > If the FFT has an imaginary part, you can't ignore it. You need it to > put the signal back together. Lowpass both parts.
When you lowpass your real part, are you multiplying point by point? If so, then I suspect you are trying to apply a filter mask (amplitude vs freq descriptor) to your freq domain data. However, your real (or imaginary) output from the FFT does *not* by itself describe your signal as 'amplitude vs freq' - this means, you cannot apply your filter in this manner (a simple way to state this is that your units don't match). You have to convert your FFT output (complex data) into an amplitude vs freq trace, apply your filter, convert back to complex data and then back to time domain (IFFT).....OR create your filter mask as complex data in the freq domain and apply to your FFT output. Cheers Bhaskar
> > There may be other, more efficient ways to filter your data, but that > might not matter to you. > > Jerry > -- > Engineering is the art of making what you want from things you can get. > &#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;
Bhaskar Thiagarajan wrote:

> When you lowpass your real part, are you multiplying point by point? If so, > then I suspect you are trying to apply a filter mask (amplitude vs freq > descriptor) to your freq domain data. However, your real (or imaginary) > output from the FFT does *not* by itself describe your signal as 'amplitude > vs freq' - this means, you cannot apply your filter in this manner (a simple > way to state this is that your units don't match). > You have to convert your FFT output (complex data) into an amplitude vs freq > trace, apply your filter, convert back to complex data and then back to time > domain (IFFT).....OR create your filter mask as complex data in the freq > domain and apply to your FFT output.
All that's necessary is to apply the filter mask to the real and imaginary parts separately. That will scale the magnitudes and preserve the phases. The probably faux pas is making the transition too abrupt, with the consequent ringing that causes. Jerry -- Engineering is the art of making what you want from things you can get. &#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;
"Jerry Avins" <jya@ieee.org> wrote in message
news:2t0mhlF1puk5mU1@uni-berlin.de...
> Bhaskar Thiagarajan wrote: > > > When you lowpass your real part, are you multiplying point by point? If
so,
> > then I suspect you are trying to apply a filter mask (amplitude vs freq > > descriptor) to your freq domain data. However, your real (or imaginary) > > output from the FFT does *not* by itself describe your signal as
'amplitude
> > vs freq' - this means, you cannot apply your filter in this manner (a
simple
> > way to state this is that your units don't match). > > You have to convert your FFT output (complex data) into an amplitude vs
freq
> > trace, apply your filter, convert back to complex data and then back to
time
> > domain (IFFT).....OR create your filter mask as complex data in the freq > > domain and apply to your FFT output. > > All that's necessary is to apply the filter mask to the real and > imaginary parts separately. That will scale the magnitudes and preserve > the phases. The probably faux pas is making the transition too abrupt, > with the consequent ringing that causes. > > Jerry > -- > Engineering is the art of making what you want from things you can get. > &#4294967295;&#4294967295;&#4294967295;&#4294967295;
&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295; As Jerry Avins says , if you make sure that you scale the real and imaginary parts of your FFT equally according to the magitude of your transfer function then you preserve the relative phase of all the components. Two more points : 1/ If you are using matlab you may want to fftshift your fft data before you apply this weighting (or fftshift the weights from natural increasing frequency order if that's what they are generated in). 2/ I remember seeing several threads in this group over the course of time which related to correlation and convolution using FFT so you may find that you can get a more efficient arrangement than doing the IFFT then cross-correlation. Lastly , why are you doing autocorrelation if you are comparing the two sample sets? Do you autocorrelate each one and then look for similarities by eye? I'm just curious. Best of Luck - Mike
akashya@uncc.edu (Abhinav Kashyap) wrote in message news:<c56359ad.0410111340.40c0c7c9@posting.google.com>...
> I am trying to compare two data sets with equal no. of points (it has > several peaks in it).I am using correlation technique .First I filter > the noise present in signal. > For that I first FFT the data, multiply low pass filter on real part > of FFT and than take IFFT. > I do same procedure for other data set. > Than I do cross correlation and auto correlations. > > I am not sure is this right approach or not. End result what i am > getting is like a triangle around y-axis.
Well, the frequency-domain filter operation in matlab goes as %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% H=fft(h,N); X=fft(x,N); Y=X.*H; % Multiply *complex* DFT of x with H y=real(ifft(Y)); % Filtered data should be real-valued %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% where h is the (time domain) impulse response of the filter (equivalently, H is the frequency-domain transfer function of the filter), x is the time-domain data, y is the time-domain filtered data, and N is an integer which must meet the requirement N >= length(x)+length(h)-1; Note that N is the "significant length" if the impulse response in time domain, i.e. the number of "non-vanishing" coefficients in the time-domain impulse response. If you know what you're doing, this is a very efficient way, computationally speaking, of implementing a digital filter in off-line applications. There are hazards, though, that you need to be aware of. You haven't said much about where you got the LP filter function from, or how it looks. I assume (rightly or wrongly) that you basically set all DFT coefficients above a certain bin to zero and leave the others alone. As others already have pointed out, you may see ringing in time domain if the transition band in your filter is too narrow. I would recommend that you design a filter by some of the "proper" methods and use this filter function. The results will be better than by the quick'n dirty method of setting spectrum coefficients to zero.
> Has anybody used this technique earlier....any references?
Setting spectrum coefficients to zero is a common way of implementing quick'n dirty filters. As always with quick'n dirty methods, it works but there are problems. One of the main problems here is that you need to avoid circular convolution, hence the restriction on N mentioned above. The other problem is that narrow transition bands introduce artifacts in the filtered time-domain data, which is why you would be better off designing a filter by the proper methods. So if you are serious about this job, it would be well worth the while to take the time to learn how to design filters with matlab. If you have the Sig Proc toolbox available, you have everything you need.
> I am Using Matlab for programming. > > Thanks in advance
Rune
"Rune Allnor" <allnor@tele.ntnu.no> wrote in message
news:f56893ae.0410112207.44e9f071@posting.google.com...
> akashya@uncc.edu (Abhinav Kashyap) wrote in message
news:<c56359ad.0410111340.40c0c7c9@posting.google.com>...
> > I am trying to compare two data sets with equal no. of points (it has > > several peaks in it).I am using correlation technique .First I filter > > the noise present in signal. > > For that I first FFT the data, multiply low pass filter on real part > > of FFT and than take IFFT. > > I do same procedure for other data set. > > Than I do cross correlation and auto correlations. > > > > I am not sure is this right approach or not. End result what i am > > getting is like a triangle around y-axis. > > Well, the frequency-domain filter operation in matlab goes as > > %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% > H=fft(h,N); > X=fft(x,N); > > Y=X.*H; % Multiply *complex* DFT of x with H > y=real(ifft(Y)); % Filtered data should be real-valued > %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% > > where h is the (time domain) impulse response of the filter > (equivalently, H is the frequency-domain transfer function of the > filter), x is the time-domain data, y is the time-domain filtered > data, and N is an integer which must meet the requirement > > N >= length(x)+length(h)-1; > > Note that N is the "significant length" if the impulse response in > time domain, i.e. the number of "non-vanishing" coefficients in the > time-domain impulse response. If you know what you're doing, this is > a very efficient way, computationally speaking, of implementing a > digital filter in off-line applications. > > There are hazards, though, that you need to be aware of. You haven't > said much about where you got the LP filter function from, or how it > looks. I assume (rightly or wrongly) that you basically set all DFT > coefficients above a certain bin to zero and leave the others alone. > > As others already have pointed out, you may see ringing in time domain > if the transition band in your filter is too narrow. I would recommend > that you design a filter by some of the "proper" methods and use this > filter function. The results will be better than by the quick'n dirty > method of setting spectrum coefficients to zero. > > > Has anybody used this technique earlier....any references? > > Setting spectrum coefficients to zero is a common way of implementing > quick'n dirty filters. As always with quick'n dirty methods, it works > but there are problems. One of the main problems here is that you need > to avoid circular convolution, hence the restriction on N mentioned > above. The other problem is that narrow transition bands introduce > artifacts in the filtered time-domain data, which is why you would be > better off designing a filter by the proper methods. So if you are > serious about this job, it would be well worth the while to take the > time to learn how to design filters with matlab. If you have the Sig > Proc toolbox available, you have everything you need. > > > I am Using Matlab for programming. > > > > Thanks in advance > > Rune
Just one more thing : elsewhere in this group there are a number of threads dealing with cross-correlation using FFTs and if you look them up you may well find that you can use your filtered FFT values without having to do the IFFTs before you cross-correlate. You could probably save yourself some run-time if you want to take the trouble. Best of luck - Mike
akashya@uncc.edu (Abhinav Kashyap) wrote in message news:<c56359ad.0410111340.40c0c7c9@posting.google.com>...
> I am trying to compare two data sets with equal no. of points (it has > several peaks in it).I am using correlation technique .First I filter > the noise present in signal. > For that I first FFT the data, multiply low pass filter on real part > of FFT and than take IFFT. > I do same procedure for other data set. > Than I do cross correlation and auto correlations. > > I am not sure is this right approach or not. End result what i am > getting is like a triangle around y-axis. > Has anybody used this technique earlier....any references? > I am Using Matlab for programming. > > Thanks in advance
You may also save some time and effort by exploiting the Weiner-Kintchine theorem which states that the autocorrelation of a time function is mathematically equivalent to the inverse of the Power (i.e. magnitude squared) spectrum. I suspect that the triangle you are seeing at the origin is due to the autocorrelation of noise. Essentially, when you look at an autocorrelation function, the dimension of the x axis is correlation length so as you get closer to zero, you are observing effects that don't correlate with themselves over time i.e. noise and non periodic events. You can experiment with Matlab by generating a few cycles of a sine wave with added noise. The ACF of this will be a cos wave with the same amplitude with a 'blip' at zero with amplitude roughly equivalent to he noise rms.