Determining The Transfer Function Of A System Using FFT Methods

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

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

The resonant peak lets you compute the damping ratio, z. (Click here to go to the point in the lessons where the resonant peak-damping 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 one-dimensional 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 log-log 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:

G(s) = Y(s)/U(s).

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

G(jw) = Y(jw)/U(jw).

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:

G(jw) = jwY(jw).

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 s-plane. 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

Call that time Dt.

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:N-2]'/(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 mid-point 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 input-output 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 non-zer 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 non-zero 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.