Using
FFT Methods
Why Use FFT Methods?
Goals For This
Lesson
What Is Involved
In Using The FFT For Identification?
Using Impulse Response
Data?
Step Response
Response
To General Inputs
You are at: Basic Concepts
 System Models  FFT Methods
Click
here to return to the Table of Contents
Why
Use FFT Methods?
FFT methods can be used to obtain transfer function models for systems.

In a linear system, the
frequency response function can be obtained by using sinusoidal inputs.
However, some inputs  including random signals  can be expressed as a
sum or integral of sinusoids. Finite length signals can be expressed
as a sum of sinusoids, for example. These signals can then be used
to compute the frequency response function, and once you have the frequency
response function you can use Bode' plots to determine a numerical transfer
function.

Using more than one time
response will let you check things like linearity of the system or whether
the system varies in time, however.

Getting time response
data into a computer for analysis is easily automated.
Goals
For This Lesson
What are you trying to do in this lesson?

Given a
linear system for which you need the transfer function, and for which you
can get the impulse response or the step response of the system

Be able
to use knowledge of the input and the measured response to get the transfer
function of the system using FFT methods.
What
Is Involved In Using The FFT For Identification?
Using the FFT to get a transfer function of a system is not really overly
complex.

There must be an input
to the system which you can observe and record  preferably in a computer
file.

There must be an output
from the system which you can observe and record  again, preferably in
a computer file.

The input and output data
has to be able to be read into an analysis program, like Mathcad or Matlab,
that can take the Fast Fourier Transform  the FFT  of both input and
output data records.
We
will start with a relatively clean set of data.

The data is shown below.
Assume that it is the response of a system to a unit pulse input  a pulse
of very short duration with a unit area.

A basic tenet of linear
analysis is that the unit impulse response of a system and the transfer
function of the system are a Laplace transform pair, i.e. the transfer
function is the Laplace transform of the unit impulse response.

The implication is that
we can ge the transfer function by getting the Laplace transform of the
unit impulse response.

We will assume that we
have data, and it will never be as "clean" as the data we use in this example.
Our goal is to get the transform of the data that we have, and the path
we take to get the transform is circuitous because we need to allow for
noisy data and not having a nice mathematical function that we can apply
our analysis skills to.

This response would then
be very close to being the impulse response of the system.

The algorithm to compute
the transfer function starts with taking the FFT of this data set.
The complete algorithm is given next.
The
algorithm is as follows.

Load the input data and
output data into the mathematical analysis program.

FFT the input data and
output data.

Divide the output FFT
by the input FFT and take the magnitude. Since the input for our
example is a unit impulse, the input FFT is 1.0.

Plot the result as a Bode'
plot.

Treat the resulting Bode'
plot as a frequency response  which it really is  and use frequency response
methods to fit a transfer function to the calculated Bode' plot.
Looking
at the first step, you should realize that there are different programming
environments in which you can do this. You can load data into:

Mathcad

Matlab

A Spreadsheet
The point is that you have a number of options,
and there are several programming environments in which you can load data,
and FFT the data you load. The code to load a set of data into a
Mathcad workspace is what we have here.

READPRN reads data from
a flat file and the argument for the function is the complete file name
including the complete path.

Next, you need to FFT
the data. Mathcad supplies a number of FFT functions.

CFFT will FFT data for
cases where the number of data points is not a power of 2. (Powers
of 2 will make the computation faster!)
If we FFT the data in
the plot above, and plot the magnitude of the FFT, we get the plot below.
The data plot
is repeated at the left.

The resultant plot doesn't
look particularly informative  and it isn't! We will fix that soon,
but there are some points to note on this plot.
This
plot has several typical characteristics of FFT magnitude plots.

Detail is lost when linear
scales are used.

The magnitude plot is
symmetric around the center frequency. There is no additional information
in the higher frequency portion of the plot and we
will not use the data above the middle of the plot.

This plot is plotted agains
the index, i, of the data array. We may need to correct scales to
account for the actual sampling times.
We
can change the plot to a Bode' plot by doing the following.

Use logarithmic scales
on both the vertical and horizontal axes.

Plot the magnitude of
the FFTed data.

Display the plot only
for frequencies below half the maximum frequency.
Here's
the plot that results. This clearly looks like a Bode' magnitude
plot.
There are some questions you need to answer
about this system.
Problem
P1

The low frequency asymptote
lets you compute the DC gain.

What is the DC Gain for
the example?
Next
you can find the damping ratio for this system.
Problem
P1

The resonant peak lets
you compute the damping ratio, z.
(Click here to
go to the point in the lessons where the resonant peakdamping ratio relationship
is discussed.)

How high is the resonant
peak?

The resonant peak in this
system is 3.3 since gain at the resonant peak is .033 and the DC gain is
.01. Now, you need to convert that to a damping ratio.

Here's the expression
that relates resonant peak gain to damping ratio.

Gain at Resonant Peak
= G_{dc}/2z

Now, what is the damping
ratio?
Finally,
we need to determine the resonant frequency. First, we need to clarify
the horizontal axis and what is being plotted there. The plot we
have been using is plotted with the array index along the horizontal axis.
In other words, our data array was onedimensional and was just an array
of the measured data points.

The total time interval
for the data record, T, determines the frequency spacing in the FFT plot.

For the ith point on the
plot, the frequency for that point is i/T.

Now, with some data we
can determine the resonant frequency.
We're
ready to compute the resonant frequency.

First determine the point
at which the peak occurs.

Determine the index, i,
for that point?
Now, we need
some data.

We will assume that the
total time for the data was 1.024 seconds, and there were 1024 points in
the data record, with one millisecond between points.

What is the frequency
for an index, i = 16? Answer is (16/1024)/(.001) (One millisecond)
SET UP PROBLEM.
So
far we haven't had any problems with this technique. Here's a review
of the technique and some questions.

Assuming that you have
the impulse response, you can get a Bode' plot that is equivalent to what
you would get taking a frequency response, by FFTing the impulse response,
and plotting it on loglog scales or db vs log f scales.

We've demonstrated the
technique for a second order system with complex poles, but it seems clear
 although you should check this  that we can do the same for data for
other systems, e.g. systems with several real poles, or general combinations
of real and complex poles and zeroes.

The one real question
is "What if we don't have the impulse response?", or in particular, "What
if we measure a step response  a widely used test signal  and all we
have is a step response measurement?".
Using
Step Response Data
In this section we are going to use step response data to get a Bode' plot
equivalent to evaluate system transfer functions. Note the following.

If the input to a system
is u(t), with a transform U(s), and the output of the system is y(t), with
a transform Y(s), the transfer function is:

We can evaluate the transfer
function by getting a plot of G(jw).
That should work out to be:

So, we should be able
to transform y(t) and u(t) numerically, and divide them and plot the result.

That's the theory, but
it only works if both functions actually have transforms. A tempting
line of reasoning is the following.

If the input to a system
is a step, then the transform of the input is 1/s.

We can substitute s =
jwin
1/s, and then compute the transfer function as:

Plot the Bode' plot and
we've got the transfer function, right?

Wrong!  but for a subtle
reason.
The
line of reasoning above is tempting but wrong  as so many things in life
seem to be.

It's true that the transform
of a step is 1/s.

It's not true that you
can substitute s = jw
in 1/s. Here's the subtle thing! The transform of a step 
1/s  is only defined for s in the right half of the splane. It's not
defined on the imaginary axis. That seems like such a small technicality.

However, if you try to
FFT a step  a sequence consisting entirely of 1's  you will find
out something interesting. The FFT will not drop off with frequency.
It will be a constant for zero frequency, and zero for all other frequencies!
A step only contains one frequency  zero frequency.

That makes it really hard
to divide the output by the input at all frequencies because the input
FFT is zero at most frequencies.
Are
there any conclusions we can draw from this?

We can't manipulate the
transform of the step response  doing calculations in the frequency domain
 in a way to get the transfer function!

That doesn't mean we can't
do something in the time domain.

In the time domain, the
step response is the integral of the impulse response. On the other hand,
the impulse response is the derivative of the step response  and we can
differentiate numerically!
To
differentiate data numerically, the simplest approach is the following.

Form the difference between
succesive samples of the data.

Diff_{i}
= Data_{i+1}  Data_{i}

The time between successive
data points should be the same for all pairs of consecutive data points

Approximate the derivative
with the expression:

Derivative_{i}
= Diff_{i} / Dt.

The derivative is really
the mathematical limit of this expression as Dt
> 0. Differentiating numerically in this fashion is based on the basic
idea of a derivative, but there are some precautions.

It's not clear how small
Dt
needs to be to make the computation acurate.

As Dt
gets smaller and smaller, the graininess of the data  determined by the
number of bits in the numerical computer representation will eventually
limit how accurate the difference can be computed.
Now, we can look at how to do the numerical differentiation in Mathcad
and Matlab.
Using
Step Response Data In Matlab
Here's a Matlab workspace.
cd a
load 'a:\StepData.txt'
time = StepData(:,1);
DT = time(2)  time(1);
Response = StepData(:,2);
plot (time,Response)
Response1 = Response;
Response1([1],:) =
[ ];
Response2 = Response;
N = length(time);
Response2([N],:) =
[ ];
ImpulseResponse =
(Response1  Response2)/DT;
time([N],:) = [ ];
plot (time,ImpulseResponse)
FreqResponse = abs(fft(ImpulseResponse)*DT);
Frequency = [0:N2]'/(N*DT);
loglog(Frequency,FreqResponse)

This workspace loads a
data file.

Then, time is extracted
from the first column of the data and Dt is computed.

Next, the step response
data is extracted and the step response is stored in Response1.

Now, the first element
is extracted from Response1, and the last element is extracted from a copy
named Response2.

The impulse response is
computed using the difference expression.

Then, take the last element
from the time array.

Plot the impulse response
vs time.

FFT the impulse response,
and normalize the result, and plot the Bode' plot.
After you have the Bode' plot, then you need to be able to extract the
transfer function from that data. Click
here if you want to review getting a transfer function from frequency
response data.
You should also observe the
following precautions.

Be sure that you have the response to a unit impulse.

If your data is the response to an impulse of magnitude
5, then you need to divide the data by 5 before doing this analysis.

Remember that the data is really only valid up to
the center point of the frequency data you get after FFTing, and that you
should be careful in your analysis to throw away the data beyond the midpoint
of the FFT.
Using
Step Response Data In Mathcad
Here's a Mathcad workspace. (Actually, a picture of a Mathcad workspace.
It's not live, obviously.)
Here are the steps in the Mathcad example workspace.

Data is generated in this
workspace. We did that just to show the technique. In general,
you would read data from a file.

Plot the data.

Compute the derivative
 the impulse response.

FFT the impulse response,
normalize the result, and plot the Bode' plot.
Using
General Response Data
You may now realize that you can use a wide variety of inputoutput data
and calculate a transfer function with FFT techniques. The algorithm
is simple.

There must be an input
to the system which you can observe, record in a file, load into an analysis
program and compute the FFT. Caution  the input should not contain
just one frequency. The input must contain all frequencies.
You will have trouble if you try to divide by zero at any frequency.

What's true for the input
must also be true for the output, so load it also and be sure all frequency
components are nonzer magnitude.

The algorithm is simple
 FFT the input and output, and divide the output FFT by the input FFT.

Plot the quotient as a
Bode' plot.

The only caveat is that
the input FFT must be nonzero at all frequencies for which you do the
computation. That's the problem that kept the step input from working
with this technique. It won't be a good idea to use sinusoidal inputs
for this reason.
Problems
Links
To Related Lessons
Lessons on Getting Models
From Data
Send
us your comments on these lessons.