Hi, and welcome to module 4.3 of Digital
Signal Processing.
In this module we will look at the DFT
in practice by showing some analysis
example by showing
you how to label the DFT axis to relate
the DFT to physical quantities like real
world frequency.
We will talk about DFT synthesis, how we
move from the frequency domain to the time
domain.
And finally we will introduce
a close relative of DFT called the
discreet Fourier series, which
is just a different flavor of DFT applied
to periodic signals.
You remember that we tried to wet your
appetite with respect to Fourier
analysis in the beginning of module 4.1,
by showing you a mystery signal.
And by showing you that by changing the
basis in which the signal is represented,
we all of a sudden could see some
structure in the data.
So now we can go back to that example
and analyze it applying what we know about
DFT.
So we take the DFT of the mystery signal,
we plot the real and imaginary part of the
DFT,
and we see that the structure that appears
is
a couple of spikes at symmetric positions
in the spectrum.
Now, we already know that this
indicates the presence of a strong
sinusoidal component
and more importantly, a sinusoidal
component has a
frequency which is a multiple of the basic
frequency for the space that the signal
lives in.
Now, if you look in detail at what happens
in the spectrum, we see
that the peaks are for k equal to 64 and k
equal to 960.
We also see that the
peaks appear only in the real part of the
spectrum and we remember that
when this happens the underlying sinusoid,
the
sinusoid represented by those peaks, is a
cosine.
So from this simple visual inspection, we
can already
write our signal as such, it will be a
cosine component and we will have to
determine both
the frequency and the initial phase of
this cosine.
And there will be another part that we can
probably call
a noise component, in the sense that it
doesn't have any structure.
Since the imaginary part of the Fourier
transform doesn't exhibit any particular
peak, we can assume that the phase of the
cosine is zero.
And as far as the frequency is concerned,
we
know that the peak occurs at k equals 64.
We are in a space of 1,024 points, and
therefore omega will be equal to
2 pi over 1,024, the basic frequency for
the space, times 64.
If it blocked the two components
separately at this point, we see
that we have indeed a cosine that
oscillates 64 times in 1,024 points.
And we have an additional noise component
that doesn't really seem to have any
structure.
Let's now look at another signal, one that
we showed in the introduction to this
class.
And that is a time series that maps the
solar activity since the 1700s.
Astronomers have noticed that solar
activity could be
related to what they call a sunspot
number.
The details are not really important, but
fundamentally
a measure of how many solar spots are
at a given point in time on the face of
the sun.
So we have a data set that goes from 1749
to 2003,
that is equivalent to 2,900 months, solar
spots are computed every month.
And as we plot this data we see that there
is an oscillatory behavior.
We're interested in finding out if
there's any fundamental periodicity in
solar activity.
So we can take the Fourier transform.
The sunspot time series is a real signal.
Now, if we are interested just in the
magnitudes of the DFT
coefficients, you remember we need only
show the first half of them.
Even so, if you look at the first 1,500
coefficients you see that
after the 100th coefficient or so, their
magnitude becomes too small to be
relevant.
So in this plot, for clarity, we just
show the first 100 DFT coefficients in
magnitude.
And indeed we can see that there is not
just one peak but a series of peaks.
Now, what is the fundamental period?
What is the most dominant mode of this
time series?
The main peak happens at k equal 22.
There are 22 cycles over the entire data
set.
So the data set is 2900 months, so the
frequency of the solar spot activity is
2900 divided
by 22 which is approximately 11 years.
So solar activity seems to have an
inherent period of 11 years.
We can perform the same analysis with a
data set that records the daily
temperature over a total of 2920 days.
If we take the DFT of this time
series, we see there is actually a very
pronounced periodicity in the spectrum,
there is basically just one peak.
And if we plot the normalized DFT
coefficient, so
if we divide the DFT coefficients by the
length
of the temperature vector, we can actually
extract some
information about the temperature values
in this time series.
So remember, for instance, the DFT
coefficient for k equal to
0 is the known normalized average of all
the data points
because x of 0 is simply the sum of all
the points in
the time series for n that goes from 0 to
big N minus 1.
So if you plot the normalized
coefficients,
the coefficient in 0 will be the average
temperature for the time series, which in
this case happened to be 12.3 Celsius
degrees.
And then we remark that there is a peak in
the DFT for k equal to 8.
So if we sum up what we learned from the
DFT about the temperature signal, we
know the average value of the temperature
is
the 0th DFT coefficient normalized, so
12.3 degrees Celsius.
The main peak is at k equal to 8, for a
value of 6.4 degrees.
What that means is there are eight cycles
of the temperature signal over the entire
duration of
our data set, which is 2,920 days.
So the frequency is 2,900 divided by 8,
which happens to be 365 days.
So indeed the temperature has a yearly
periodicity, as we all know very well.
Now, the value of the DFT main peak is 6.4
degrees Celsius.
Now, if we have a sinusoid of the form A
cosine of omega
n, and we take the DFT of this signal, we
will have a peak in magnitude
at some, for some index k, that has a
value of A
over 2, remember the definition of the DFT
of a cosine function.
So from this we can say that the
temperature excursion around
the average is twice the value of the DFT
main peak.
And, so we can say,
that the yearly temperature at the point
of measurement for the
year was an average of 12.3 degrees plus
or minus 12.8.
So this is a second example in which, in
order to find the real world frequency of
a certain
DFT component, we take the total length of
the
signal and we divide by the location of
the peak.
This is actually a general method to label
the frequency
axis of a DFT plot and it will be very
useful in the future
when we start analyzing generic signals
that
have been sampled in a variety of
contexts.
So if you remember, in module 2.2, we
informally introduced
the notion of a clock for a digital signal
processing system.
And this is really equivalent to
associating a certain
time interval, Ts, measured in actual
physical seconds between successive
samples in the signal.
So, for instance for the solar spot
signal, we had a Ts of
one month, whereas in the temperature
signal, the Ts was equal to one day.
We will see that for audio signal Ts
becomes very small
because we will need to take at least
8,000 samples per second.
Now if we have a value for Ts, which is
determined
by the experimental setup, we can reason
like so, the fastest positive
frequency in a digital signal is omega
equal to pi.
So a sinosoid at that frequency needs two
samples
to do a full revolution, remember the unit
circle.
Something that moves at a speed of pi will
be here at instance n, here at n
plus one, here at n plus two, so two
samples to complete the full revolution.
Now, the clock Ts can
also be expressed as 1 over Fs, where Fs
is the frequency of the system.
This is a standard relationship between
period and frequency for any system.
So, if the real world period for the
fastest sinusoid in a digital system
is 2 times Ts measured in seconds, the
real world frequency for the fastest
sinusoid is Fs over 2.
So the maximum frequency,
once we've fixed the period between
samples, is
Fs over 2, or equivalently, 1 over 2 Ts.
So let's take a concrete example.
Suppose I give you an audio file that
records a train whistle
and I'm asking you to find out which notes
make up the whistle
[SOUND].
And of course your idea is to take the DFT
of this file and find out the peaks
and the Fourier transform, which will
probably correspond to
the fundamental frequencies of the notes
that are being played.
So in order to do so you need some extra
information and namely you
need to know the sampling frequency of
this file, namely the time between
successive samples.
So the sampling frequency is Fs
equal to 8000 hertz, which means that Ts
is 1 over 8000 seconds,
this is the time between samples, and also
the file contains 32768 samples.
We will see it later why we often choose
signal lengths that are a power of two.
So if you take the DFT of this file and
plot
its magnitude you see that there are three
peaks that probably
corresponds to the three notes being
played.
Now in order to label the axis remember
what we said, the
highest frequency for a digital system
corresponds to half the inherent sampling
rate.
And on the DFT vector, this point will
correspond to the midpoint in
the vector, namely the point for k is
equal to N over 2.
So here indeed we have the midpoint of the
DFT vector and
therefore, this will correspond to a
frequency of four kilohertz.
We also have three peaks to take place for
different values of k.
To find the real frequency there, we just
apply a
linear mapping from zero to four
kilohertz, and so the
frequency of any intermediate point will
be 4 kilohertz divided
by the total number of points multiplied
by the index k.
And with this rule, we find that the peaks
take place at 494 hertz, 559 hertz, 739
hertz.
And if we look up what these frequencies
correspond to in the standard
western scale, we find that the whistle is
a simple B minor chord.
Okay, so far we have looked at DFT as an
analysis tool, let's now spend
some time thinking about the synthesis
properties
of the DFT, how we build the signal
starting from DFT coefficients.
And to better understand how the DFT
reconstruction formula works, it's useful
to look at
the standard sinusoidal generator that
works at one
of the standard frequencies in the Fourier
basis.
So here you have a complex exponential
with a frequency which is a
multiple of 2 pi over big N and with an
initial phase, phi k.
The way this generator works is
by successively generating points on the
unit circle, starting at the phase
phi k and proceeding in increments of 2 pi
over big N.
We can draw this concisely with this
notation where here inside you
have the index of the Fourier frequency
that the generator will produce.
And the input parameters for the
generators are
of course the initial phase here, and
again
factor A of k, that will lead to A of k
times to
e to the j 2 pi over big N, kn plus phi k.
So the DFT reconstruction formula can be
drawn
graphically like a machine composed of big
N sinusoidal
generators, each one of which will be
initialized
with a gain factor and with a phase
factor.
All their outputs will be summed together
and
this will give us back the original
signal.
How do we initialize each one of these
blocks?
Well the amplitude will be the amplitude
of the corresponding Fourier analysis
coefficient normalized by
big N, and the phase will be
the phase of the corresponding Fourier
analysis coefficient.
So, lets look at an example.
Lets take a simple signal, a seven point
long signal that looks like a triangle,
and we perform a Fourier analysis
on this vector and we get seven Fourier
coefficients
that we list here in terms of amplitude
and phase.
Now we initialize the seven generators
with the As
and the phis that we found in the analysis
part.
And then we turn the crank seven times and
we
will sum the outputs together to obtain
the original signal.
So here is how it works.
In the first line,
we will show the output of each generator
in turn, starting
for k equal to 0 up to k equal to 6.
And here on the bottom line, we will show
the cumulative
sum of the outputs as k proceeds from 0 to
6.
So the first generator is the one for k
equal to 0,
so the frequency is 0, and the output will
be simply a constant.
As a matter of fact, the value of this
contes,
constant will be equal to A of 0.
And so the cumulative sum now has only one
term,
so the output at this point will look like
a constant.
For k equal to 1, we have now a sinusoid
that will have a certain initial phase and
a certain amplitude.
And we can see that as we add this to the
previous constant, here in
grey you see the previous sum, and here
the updated sum for k equal to 1,
the real part of the signal starts to look
a little bit like a triangle.
We also have a significant non-negligible
imaginary part which will have to
disappear in the end because remember, we
started from a real signal.
So as we proceed for k equal to 2, we see
that the coefficient
amplitude is very small and so this will
be just small adjustments in the sum.
It's very hard to see the difference
between this step and the previous step.
K equal to 3 is very similar, small
adjustments.
For k equal to 4, we're going in the,
counter-clockwise part of the spectrum,
so more adjustments still, that will start
to bring down the imaginary part to 0.
Here you see another small adjustment.
And finally k equal to 6, the last
coefficient will completely kill the
imaginary part and give us back exactly
the signal we were starting from.
So why do we think
about the DFT reconstruction formula as a
machine?
Because in reality these machines were
used in
the past for instance to predict the
tides.
The tide is a periodic phenomenon and once
you find the modes you can actually use
the experimental data to predict the
evolution of
the tide as a superposition of sinusoidal
components.
So this machine, that was originally
invented by
Lord Kelvin, is nothing but a mechanical
implementation
of the DFT reconstruction machine that we
just saw before.
So here's a question, what happens if you
turn
the crank of the machine more than big N
times?
Well it turns out that the output will
become periodic, in other words,
x of small n plus big N will be equal to x
of n.
This is actually quite apparent from the
structure
of the synthesis and analysis formulas of
the DFT.
Let's start from the synthesis formula.
So in
theory, the output index small n should go
from 0 to big N minus 1.
But you can see that it only appears here,
inside this complex exponential.
Now this guy is 2 pi periodic, so if you
push n over big N minus 1, this will
simply cycle over once again as if n was
looping over the 0 to n minus 1 range.
So in the end, you can actually safely
take n
to be from the set of integers and the
output
will be an N periodic signal in the time
domain.
Actually the same holds for the analysis
formula.
You can let the index k roam over the
entire
set of integers and since it appears only
in this
complex exponential here, again, it will
be as if k
was looping over the 0 to n minus 1 range.
So what happens
really, is that you can consider the
sequence of DFT
coefficients as an N periodic signal in
the frequency domain.
Now, when we talk explicitly about the
periodicity of the Fourier
transform, we prefer to use the term
discrete Fourier series or DFS for short.
So the DFS is just the DFT with the
periodicity explicit.
The DFS maps
an N-periodic signal onto an N-periodic
sequence of Fourier coefficients.
And the inverse DFS maps an N-periodic
sequence of Fourier coefficients onto an
N-periodic signal.
But mathematically the two are the same,
it's just a conceptual difference.
However, DFS helps understand, for
instance, why we said that the
natural shift for finite length sequences
is actually a periodic shift.
Let's consider an N-periodic sequence x
tilda of n.
For periodic sequences, shifts are
well-defined.
We can take x tilda of n minus big M
and it's just a normal shift of an
infinite two-sided sequence.
If we take the DFS of this periodic
sequence with a shift, we can work out
very easily that the DFT coefficient for
index k is equal to the DFT coefficient
for index k of the original sequence
without the shift.
So, x of k is simply the DFS of tilde x of
n
with a phase shift term, e to the minus j,
2 pi over n, big M times k.
Similarly, it's very easy to verify that
the
inverse discrete Fourier series of a
vector of
DFS coefficients multiplied by a phase
delay factor,
is simply the periodic sequence delayed by
M.
So this term, e to the minus j, 2 pi
over n, Mk is a delay factor in the
frequency domain.
So what happens for finite length signals?
We saw that in this case, a shift is not a
well-defined operation and we need to
embed the signal into an infinite length
signal, either
by building a finite support or a periodic
signal.
Now we can always build x tilda of n as x
of n modules big N.
Now mathematically, the DFS of this
periodic signal is equal to
the DFT of the original signal.
So the inverse DFT of the DFT
of the signal times this delay factor,
that we showed in the previous slide,
is numerically equal to the inverse DFS of
the same quantity.
But we know that the inverse DFS of this
will be a delayed periodic sequence.
And if we go back to the definition of how
we build this periodic
sequence this will be simply x of n minus
big N modulus big N.
What we have shown here is that if we
think of the Fourier representation of a
finite length signal, circular shifts are
the natural extension of the shift
operator to finite length signals because
the
underlying Fourier representation is the
same as
that what would use for a periodic symbol
built on the finite length signal.