Reply by Steven G. Johnson●August 4, 20042004-08-04
steven_hyh@yahoo.com (Steve) wrote...
> They are both approximation to ideal hilbert transform. But they
> should generate the same result.
Sigh... Yes, they are both approximations to the ideal transform, but
they are not the *same* approximation, and they do *not* generate
*exactly* the same result.
Even in the example case you posted, if you zoom in on the plots you
can see clear differences between the results of the two methods (e.g.
the amplitudes differ by about 0.3%).
My question was, how do the errors compare, and are there
circumstances under which one approximation is known to be
significantly better than the other?
Cordially,
Steven G. Johnson
Reply by Steve●August 4, 20042004-08-04
> [ I had asked about computing Hilbert transforms via FIR filters ala
> Oppenheim & Schafer vs. the FFT approach ala Matlab, and whether the
> FFT approach's periodic boundaries were a problem. ]
>
> steven_hyh@yahoo.com (Steve) wrote...
> > This is a good observation. IMHO, both of them are correct but just
> > realistic or not.
>
> This answer doesn't make sense to me. The two approaches don't
> compute exactly the same thing -- the FFT and FIR make *different*
> approximations in estimating the Hilbert transform of an infinite
> signal from a finite-length sample. Depending upon the application,
> they are certainly not always equally "correct", with the difference
> being only one of implementation.
Ok, see if this example makes sense to you:
% Create a sinusoid real signal
t = 0 : 0.1 : 9.9;
Xr = cos(2*pi*t);
% the hilbert transform, in theory, we know it is a sinusoidal too.
% so, the first method, using matlab's hilbert implementation.
Matlab's method % uses FFT method as you desribed.
X = hilbert(Xr); % Note that X is Xr + jXi, not the real hilbert
sequence
% Now this is the ture hilbert transform of Xr
X_hilb = imag(X);
% Obviously, X_hilb is a sinusoidal function, which meets theory
% Now turn to FIR approach
B = remez(30,[0.1 0.9],[1 1],'Hilbert');
Y = filter(B,1,Xr);
% You need to take care of delay here. I skip this part.
% But you can plot to see what is going on
plot(Y)
hold on
plot(X_hilb)
% except some delay, they are the same
They are both approximation to ideal hilbert transform. But they
should generate the same result. If you want to prove the uniqueness
of hilbert transform, that is another story. As for the "correctness"
you keep mentioning, if you can provide some examples(codes or
equations) to show what you mean, that is bettern than words.
Reply by sundar●August 3, 20042004-08-03
Hi guys,
I want to thank you all for all your valuable suggestions. I have a
question, I read in this book that.....if i do
x(t)-----> fft-----> X(f)------->X(f)*exp(i*pi/2)---------->
ifft--------H[x(t)].
In other words, if i do a fft on x(t) and then phase shift it and then
do a ifft i am supposed to get the hilbert transform of x(t). Fine,
well and good, but in matlab after i do this above method and then
directly use the matlab hilbert function....hilbert(x(t))....i get
different values.
Have i understood something wrong....or whats happening.Any clue.
Secondly iam looking for a IEEE paper which talks about how hilbert
transform works in matlab. I had no luck finding it, does anyone have
any idea what iam talking about.
Please help. Any advise is heartly appreciated.
Regards
Sundar
Reply by Steven G. Johnson●August 3, 20042004-08-03
[ I had asked about computing Hilbert transforms via FIR filters ala
Oppenheim & Schafer vs. the FFT approach ala Matlab, and whether the
FFT approach's periodic boundaries were a problem. ]
steven_hyh@yahoo.com (Steve) wrote...
> This is a good observation. IMHO, both of them are correct but just
> realistic or not.
This answer doesn't make sense to me. The two approaches don't
compute exactly the same thing -- the FFT and FIR make *different*
approximations in estimating the Hilbert transform of an infinite
signal from a finite-length sample. Depending upon the application,
they are certainly not always equally "correct", with the difference
being only one of implementation.
Regards,
Steven G. Johnson
Reply by Steve●July 30, 20042004-07-30
This is a good observation. IMHO, both of them are correct but just
realistic or not.
The FIR approach is more conventional to LTI system concept: a filter
sitting there waiting for input signals and generating output. While
the FFT approach uses the analytic sequence properties and needs to
collect certain amount of data then apply FFT with parameters like
length, etc. It is easier to implement on paper but hard in practice.
As for the comparison of these two, as far as I know the 2-D HT are
normally using FFT approach. For 1-D, both.
stevenj@alum.mit.edu (Steven G. Johnson) wrote in message news:<27cfb406.0407280640.65cff90b@posting.google.com>...
> Stephan M. Bernsee <spam@dspdimension.com> wrote...
> > You can calculate the Hilbert transformed signal using the FFT though,
> > by zeroing out the negative frequencies and multiplying by 2 before
> > doing the inverse transform.
>
> I'm familiar with this approach (e.g. it is apparently implemented in
> the hilbert() function of Matlab's Signal Processing Toolbox,
> http://www-ccs.ucsd.edu/matlab/toolbox/signal/hilbert.html). You can
> equivalently gain an additional factor of two in speed by exploiting
> the (usual) case where the input data are real (essentially by
> swapping/sign-flipping the real and imaginary parts instead of
> zeroing).
>
> I have some code for this with FFTW which I thought about posting, but
> I was told that this is not really the "right" way to compute the
> Hilbert transform -- if the input data are non-periodic, to avoid
> artifacts you should really implement an FIR filter with +/- response
> for positive and negative frequencies, as described e.g. in Oppenheim
> & Schafer.
>
> I'm not a signal processing expert, so I'm curious if anyone has a
> reference to indicate whether the FFT-based method ala Matlab is
> actually useful for serious applications, and how it compares to the
> FIR approach.
>
> Cordially,
> Steven G. Johnson
Reply by Steven G. Johnson●July 28, 20042004-07-28
Stephan M. Bernsee <spam@dspdimension.com> wrote...
> You can calculate the Hilbert transformed signal using the FFT though,
> by zeroing out the negative frequencies and multiplying by 2 before
> doing the inverse transform.
I'm familiar with this approach (e.g. it is apparently implemented in
the hilbert() function of Matlab's Signal Processing Toolbox,
http://www-ccs.ucsd.edu/matlab/toolbox/signal/hilbert.html). You can
equivalently gain an additional factor of two in speed by exploiting
the (usual) case where the input data are real (essentially by
swapping/sign-flipping the real and imaginary parts instead of
zeroing).
I have some code for this with FFTW which I thought about posting, but
I was told that this is not really the "right" way to compute the
Hilbert transform -- if the input data are non-periodic, to avoid
artifacts you should really implement an FIR filter with +/- response
for positive and negative frequencies, as described e.g. in Oppenheim
& Schafer.
I'm not a signal processing expert, so I'm curious if anyone has a
reference to indicate whether the FFT-based method ala Matlab is
actually useful for serious applications, and how it compares to the
FIR approach.
Cordially,
Steven G. Johnson
Reply by Steve●July 28, 20042004-07-28
He is describing an exact frequency component algorithm by using
Hilbert transform(HT). Unfortunately, he didn't describe it correctly.
In frequency estimation, the conventional way to do that is on two
steps:
1. using whatever spectrum method or simply a N-point FFT to show
frequency peaks, this is called coarse estimation. The peak regions
are then filtered out for step 2 use.
2. exact frequency estimation
This step takes a lot of computation since the resolution should be
very high or FFT length very high. But the task is simple: just do a
frequency sweeping at left and right points of the peak value. With
HT, the speed can be increase since only the magnitude and phase of
the original sequence are needed. I skip the equations here, you can
find reference easily on the Web or the original poster may point out
the paper he is referring to.
The reason why HT is used is because a causal sequence has been used,
that means the analytical sequence's(defined as the FFT without the
negative frequency part) real part and imaginary part are related by
HT, ie., the phase and magnitude of the sequence is easily know after
HT to the original sequence.
Reply by sundar●July 27, 20042004-07-27
Hi jerry
Thank you for your point and apologies for my communication problem.
ok iam rewritting what i read....
"The inverse FFT was taken out of the signal both before and after the
hilbert transform was applied. A complex time series was then formed
in which the real part was the inverse FFT of the signal and the
imaginery part was the inverse part of the hilbert transformed signal
.The magnitude and phase were calculated directly from this complex
signal. "
If this sound more clear , i would appreciate any suggestions on this.
Regards
sundar
Jerry Avins <jya@ieee.org> wrote in message news:<410673e9$0$5633$61fed72c@news.rcn.com>...
> Isn't the IFFT of a single non-zero "bin" something you can write by
> inspection? Don't let details of the math isolate you from the reality
> of what you're trying to do.
>
> You're twisting yourself around to the point that, while you may know
> what you mean, you aren't communicating well. That shows very clearly
> when you write, "Both A and B are complex outputs, and A is the real
> while B is the imaginary ..."
>
> Jerry
> --
> Engineering is the art of making what you want from things you can get.
> �����������������������������������������������������������������������
> sundar wrote:
>
> > hi everyone, I know most of you by now must have gone mad with my
> > posts and endless silly queries. Well the bad news is that i have some
> > more.
> >
> > i read this in an report and am trying to follow this.
> >
> > I have a time domain signal, do a fft and get into freq domain. Then
> > take the freq bin of interest and zero out the rest of the
> > frequencies. Then i do my ifft for this filtered signal. lets say this
> > output is A. Then if i rotate my filtered signal , by pi/2 and then do
> > my ifft on it. we get say output B. Both A and B are complex outputs,
> > and A is the real while B is the imaginary part and can be used to
> > calculate the magnitude and phase.
> >
> > does this sound ok to you all.
> >
> > Waiting for some help.
> >
> > Regards
> >
> > Sundar
> >
> >
> >
> >
> >
> > Randy Yates <yates@ieee.org> wrote in message news:<y8l7iza7.fsf@ieee.org>...
> >
> >>Hi sundar,
> >>
> >>What is confusing me in your description is the nature of the three
> >>"tonals" (I assume you mean sinusoids). Are they separated in
> >>frequency by a known minimum? If so, then you can perform the
> >>operations I think you desire by first bandpass-filtering the input to
> >>isolate the one frequency of interest, then performing the Hilbert
> >>transform on the result to get y[n] = A[n]e^(j*theta), where
> >>theta = j*2*pi*f*n*T + phi). Computing (theta[n] - theta[n-1]) / T
> >>of y[n] gives the instantaneous frequency, A[n] is the envelope.
> >>
> >>--RY
> >>
> >>snsundar@olemiss.edu (sundar) writes:
> >>
> >>
> >>>hi Stephen,
> >>>
> >>>Iam sorry for the vagueness in my description, thats because honestly
> >>>iam little lost in this place, and trying hard to stay afloat.
> >>>
> >>>Ok i am trying to find how my one particular frequency changes with
> >>>time.Does that make sense....?..
> >>>
> >>>also i was trying to get a envelope of my time domain signal.
> >>>
> >>>regards.
> >>>
> >>>sundar
> >>>
> >>>
> >>>Stephan M. Bernsee <spam@dspdimension.com> wrote in message news:<2md87gFllv9hU1@uni-berlin.de>...
> >>>
> >>>>Hi Sundar,
> >>>>
> >>>>I have some more questions for you because I'm afraid your description
> >>>>is still a bit vague:
> >>>>
> >>>>On 2004-07-23 21:31:55 +0200, snsundar@olemiss.edu (sundar) said:
> >>>>
> >>>>
> >>>>>Iam now a little more clear as to where iam heading. ok so getting
> >>>>>back to where i was . Iam trying to find the instantaneous frequency
> >>>>>by using hilbert transforms. this is where i stand.
> >>>>
> >>>>The instantaneous frequency of what exactly are you trying to
> >>>>determine? Your time domain sequence? Some frequency range of that
> >>>>sequence? One particular frequency?
> >>>>
> >>>>And where is the connection to your time domain amplitude envelope?
Reply by Jerry Avins●July 27, 20042004-07-27
Isn't the IFFT of a single non-zero "bin" something you can write by
inspection? Don't let details of the math isolate you from the reality
of what you're trying to do.
You're twisting yourself around to the point that, while you may know
what you mean, you aren't communicating well. That shows very clearly
when you write, "Both A and B are complex outputs, and A is the real
while B is the imaginary ..."
Jerry
--
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������
sundar wrote:
> hi everyone, I know most of you by now must have gone mad with my
> posts and endless silly queries. Well the bad news is that i have some
> more.
>
> i read this in an report and am trying to follow this.
>
> I have a time domain signal, do a fft and get into freq domain. Then
> take the freq bin of interest and zero out the rest of the
> frequencies. Then i do my ifft for this filtered signal. lets say this
> output is A. Then if i rotate my filtered signal , by pi/2 and then do
> my ifft on it. we get say output B. Both A and B are complex outputs,
> and A is the real while B is the imaginary part and can be used to
> calculate the magnitude and phase.
>
> does this sound ok to you all.
>
> Waiting for some help.
>
> Regards
>
> Sundar
>
>
>
>
>
> Randy Yates <yates@ieee.org> wrote in message news:<y8l7iza7.fsf@ieee.org>...
>
>>Hi sundar,
>>
>>What is confusing me in your description is the nature of the three
>>"tonals" (I assume you mean sinusoids). Are they separated in
>>frequency by a known minimum? If so, then you can perform the
>>operations I think you desire by first bandpass-filtering the input to
>>isolate the one frequency of interest, then performing the Hilbert
>>transform on the result to get y[n] = A[n]e^(j*theta), where
>>theta = j*2*pi*f*n*T + phi). Computing (theta[n] - theta[n-1]) / T
>>of y[n] gives the instantaneous frequency, A[n] is the envelope.
>>
>>--RY
>>
>>snsundar@olemiss.edu (sundar) writes:
>>
>>
>>>hi Stephen,
>>>
>>>Iam sorry for the vagueness in my description, thats because honestly
>>>iam little lost in this place, and trying hard to stay afloat.
>>>
>>>Ok i am trying to find how my one particular frequency changes with
>>>time.Does that make sense....?..
>>>
>>>also i was trying to get a envelope of my time domain signal.
>>>
>>>regards.
>>>
>>>sundar
>>>
>>>
>>>Stephan M. Bernsee <spam@dspdimension.com> wrote in message news:<2md87gFllv9hU1@uni-berlin.de>...
>>>
>>>>Hi Sundar,
>>>>
>>>>I have some more questions for you because I'm afraid your description
>>>>is still a bit vague:
>>>>
>>>>On 2004-07-23 21:31:55 +0200, snsundar@olemiss.edu (sundar) said:
>>>>
>>>>
>>>>>Iam now a little more clear as to where iam heading. ok so getting
>>>>>back to where i was . Iam trying to find the instantaneous frequency
>>>>>by using hilbert transforms. this is where i stand.
>>>>
>>>>The instantaneous frequency of what exactly are you trying to
>>>>determine? Your time domain sequence? Some frequency range of that
>>>>sequence? One particular frequency?
>>>>
>>>>And where is the connection to your time domain amplitude envelope?
Reply by sundar●July 27, 20042004-07-27
hi everyone, I know most of you by now must have gone mad with my
posts and endless silly queries. Well the bad news is that i have some
more.
i read this in an report and am trying to follow this.
I have a time domain signal, do a fft and get into freq domain. Then
take the freq bin of interest and zero out the rest of the
frequencies. Then i do my ifft for this filtered signal. lets say this
output is A. Then if i rotate my filtered signal , by pi/2 and then do
my ifft on it. we get say output B. Both A and B are complex outputs,
and A is the real while B is the imaginary part and can be used to
calculate the magnitude and phase.
does this sound ok to you all.
Waiting for some help.
Regards
Sundar
Randy Yates <yates@ieee.org> wrote in message news:<y8l7iza7.fsf@ieee.org>...
> Hi sundar,
>
> What is confusing me in your description is the nature of the three
> "tonals" (I assume you mean sinusoids). Are they separated in
> frequency by a known minimum? If so, then you can perform the
> operations I think you desire by first bandpass-filtering the input to
> isolate the one frequency of interest, then performing the Hilbert
> transform on the result to get y[n] = A[n]e^(j*theta), where
> theta = j*2*pi*f*n*T + phi). Computing (theta[n] - theta[n-1]) / T
> of y[n] gives the instantaneous frequency, A[n] is the envelope.
>
> --RY
>
> snsundar@olemiss.edu (sundar) writes:
>
> > hi Stephen,
> >
> > Iam sorry for the vagueness in my description, thats because honestly
> > iam little lost in this place, and trying hard to stay afloat.
> >
> > Ok i am trying to find how my one particular frequency changes with
> > time.Does that make sense....?..
> >
> > also i was trying to get a envelope of my time domain signal.
> >
> > regards.
> >
> > sundar
> >
> >
> > Stephan M. Bernsee <spam@dspdimension.com> wrote in message news:<2md87gFllv9hU1@uni-berlin.de>...
> >> Hi Sundar,
> >>
> >> I have some more questions for you because I'm afraid your description
> >> is still a bit vague:
> >>
> >> On 2004-07-23 21:31:55 +0200, snsundar@olemiss.edu (sundar) said:
> >>
> >> > Iam now a little more clear as to where iam heading. ok so getting
> >> > back to where i was . Iam trying to find the instantaneous frequency
> >> > by using hilbert transforms. this is where i stand.
> >>
> >> The instantaneous frequency of what exactly are you trying to
> >> determine? Your time domain sequence? Some frequency range of that
> >> sequence? One particular frequency?
> >>
> >> And where is the connection to your time domain amplitude envelope?