Reply by robert bristow-johnson●April 7, 20062006-04-07
having trouble sleeping, so...
krishna_sun82 wrote:
> Thank you all for your replies. I tried using both freqz and fftshift. In
> my case, there is not much difference in the end result. I don't really
> use arbitrary frequency values, but still freqz provides that flexibility
> over fftshift.
>
> Robert, can you please tell me why and in which case you would prefer
> fftshift over freqz? Is there any mathematical reason that I fail to
> comprehend?
it's not just fftshift() vs. freqz(). it's fftshift( fft( fftshift() )
) vs. freqz().
all fftshift() does is swap the first half with the second half of the
array. suppose your chunk of data of which spectrum you are interested
in was yanked out of a larger stream of samples (that is essentially
applying a window function, even if it's the rectangular window
function). also, suppose you used a "better" window function, like a
Hamming or Kaiser, so that the samples in the middle have more
amplitude than the samples at the two ends. now, where would you like
your "y-axis" to go through (where t=0 is), the middle of that chunk or
one of the ends? that's the reason to use fftshift on the data going
*in* to the FFT, so that the data in the second half of the input is
treated as the data left of the t=0 mark, the negative time data.
doing it on the output is so that the negative frequency output comes
before, to the left of, the positive frequency output.
>
> One more doubt: is there any way to find the inverse of the DTFT in
> MATLAB?
the inverse DTFT is an integral (of a continuous-frequency function).
unless you're doing it mathematically or symbolically, evaluating an
integral of a continuous function would entail an uncountably infinite
number of additions, something pretty hard to make a computer do. so
like the DFT is an evaluation of the DTFT at a finite number of equally
spaced points, the iDFT is an approximation of the iDTFT where the
integral is replaced by a Riemann summation.
so if you have a continuous-frequency function of w, X(e^jw) that is
the DTFT of an unknown discrete-time sequence, x[n], and you want to
numerically compute some approximation to x[n], then uniformly sample
X(e^jw) a gazillion (N) times from w=-pi to w=+pi, swap the halfs
(using fftshift) and inverse DFT that result and you'll have a
gazillion samples of x[n]. but since the real x[n] might never have
been periodic (which is why it is describable with the DTFT of it), the
output of the iDFT will be one period of a periodic sequence (of period
N) that should very well approximate your aperiodic x[n] for small
enough n. that is the consequence of sampling X(e^jw), since X(e^jw)
has more information in it than the samples of it: X[k] = X(e^jw) where
w = 2*pi*k/N .
> Also, DTFT output is actually an infinite sequence,
no, that's the INPUT to the DTFT. the output of it is a
continuous-frequency function from -pi to +pi.
> but in a computer we
> can get only finite points. To calculate the inverse DTFT, we have to
> integrate the scaled DTFT output over 2pi. But the computer output for
> DTFT contains finite points.
geez, i wish i would read through your post before typing. it might
save me some keystrokes.
> If there is no instant function in MATLAB,
> can I have summation in the place of integration?
yes, you have to. but what you're summing are samples of the DTFT
output (with something multipying them) which has the effect of
periodically repeating and overlap-adding the original infinite
sequence x[n] resulting in a periodic (or aliased) extension of x[n].
> What other factors do I
> need to consider to have a reasonably accurate output?
well your original x[n] that is an infinite sequence must decay to zero
rapidly enough in both directions (call those the "tails" of x[n]) so
that the DTFT summation converges to finite values of X(e^jw). when
you grab a finite segment of x[n] to pass to a DFT (which is what you
are using to approximate the DTFT), you are windowing it. to minimize
the effect of windowing, make sure the window is wide enough to include
most of the tails (so that what gets cut off is very small in
amplitude).
perhaps, just read the enormous back-and-forth we've been having in the
"Periodicity of the DFT" subthread we've had very recently. i don't
want to repeat myself all over again.
r b-j
Reply by krishna_sun82●April 7, 20062006-04-07
Thank you all for your replies. I tried using both freqz and fftshift. In
my case, there is not much difference in the end result. I don't really
use arbitrary frequency values, but still freqz provides that flexibility
over fftshift.
>i *still* recommend the two fftshift() operations.
>there is a way to do split screen plots, but i forgot it.
>
>r b-j
>
Robert, can you please tell me why and in which case you would prefer
fftshift over freqz? Is there any mathematical reason that I fail to
comprehend?
For split screen plots, you may use subplot(m, n, p) function.
One more doubt: is there any way to find the inverse of the DTFT in
MATLAB?
Also, DTFT output is actually an infinite sequence, but in a computer we
can get only finite points. To calculate the inverse DTFT, we have to
integrate the scaled DTFT output over 2pi. But the computer output for
DTFT contains finite points. If there is no instant function in MATLAB,
can I have summation in the place of integration? What other factors do I
need to consider to have a reasonably accurate output?
Reply by robert bristow-johnson●April 6, 20062006-04-06
Andor wrote:
>
> freqz() can evalute the DTFT on any arbitrary frequency points (not
> necessarily evenly spaced, as is the result of the FFT).
that's true, but it would still have to be evalutated at a finite
selection of points, though they need not be uniformly spaced.
i made a small screwup in the MATLAB pseudo-code in my previous post:
N = some_really_big_number_likely_a_power_of_two;
omega = linspace(-pi, pi(1-2/N), N);
x = some_really_big_array_of_size_N;
X = fftshift( fft( fftshift(x) ) );
plot(omega, abs(X)); % this is N points of the DFT magnitude in a
linear scale
pause:
plot(omega, angle(X)); % this is N points of the DFT wrapped phase in
a linear scale
i *still* recommend the two fftshift() operations.
there is a way to do split screen plots, but i forgot it.
r b-j
Reply by Andor●April 6, 20062006-04-06
robert bristow-johnson wrote:
...
> but the point is the same. since the output of the DTFT is a
> continuous-frequency function with infinite number of points (usually
> between w=-pi and w=+pi), and since computers have trouble dealing with
> any infinite quantity of any object, if you're gonna simulate, you have
> to decide on some finite number that is close enough. that would be
> the DFT "N".
freqz() can evalute the DTFT on any arbitrary frequency points (not
necessarily evenly spaced, as is the result of the FFT).
Reply by robert bristow-johnson●April 6, 20062006-04-06
Rick Lyons wrote:
> On Thu, 06 Apr 2006 00:53:27 -0500, "krishna_sun82"
> <krishna.sun@gmail.com> wrote:
>
> >How do I have to simulate the Discrete time Fourier Transform in MATLAB? Do
> >I have to make use of freqz function for it? Thanks in advance.
> >
> >- Krishna
>
> You simulate the DTFT in MATLAB by zero padding
> a time-domain sequence with many many zeros.
>
> Then perform the DFT on that extended-length sequence
> by whatever DFT method [fft() or freqz()] you wish.
he might not need to zero-pad if he can define the whole extended
length sequence meaningfully.
but the point is the same. since the output of the DTFT is a
continuous-frequency function with infinite number of points (usually
between w=-pi and w=+pi), and since computers have trouble dealing with
any infinite quantity of any object, if you're gonna simulate, you have
to decide on some finite number that is close enough. that would be
the DFT "N".
the other thing i would recommend is to fftshift() the input to the
fft() (to avoid that phase being off by +/- pi for odd X[k]) and to
also fftshift() the output of the fft() so that the negative
frequencies components come first. (if MATLAB had the ability to do
negative indices, a version of the DFT could be written so that none of
this swapping of the two halves would be necessary.)
N = some_really_big_number_likely_a_power_of_two;
omega = linspace(-pi, pi(1-2/N), N);
x = some_really_big_array_of_size_N;
X = fftshift( fft( fftshift(x) ) );
plot(omega, X); % this is N points of the DFT in a linear scale
r b-j
r b-j
Reply by Rick Lyons●April 6, 20062006-04-06
On Thu, 06 Apr 2006 00:53:27 -0500, "krishna_sun82"
<krishna.sun@gmail.com> wrote:
>How do I have to simulate the Discrete time Fourier Transform in MATLAB? Do
>I have to make use of freqz function for it? Thanks in advance.
>
>- Krishna
You simulate the DTFT in MATLAB by zero padding
a time-domain sequence with many many zeros.
Then perform the DFT on that extended-length sequence
by whatever DFT method [fft() or freqz()] you wish.
[-Rick-]
Reply by Andor●April 6, 20062006-04-06
krishna_sun82 wrote:
> How do I have to simulate the Discrete time Fourier Transform in MATLAB? Do
> I have to make use of freqz function for it? Thanks in advance.
For finite sequences, the DTFT is just the z-transform evaluated on the
unit circle. So yes, if x is your vector then you get the DTFT
evaluated at the frequencies F (in Hz) by
H = freqz(x,1,F,Fs),
where Fs is the sampling frequency.
Reply by krishna_sun82●April 6, 20062006-04-06
How do I have to simulate the Discrete time Fourier Transform in MATLAB? Do
I have to make use of freqz function for it? Thanks in advance.
- Krishna