Guess a system's transfer function

I have the step and ramp response of a system. How can I have MATLAB show me the transfer function of the system? I have searched 'System
Identification toolbox' but I haven't found anything useful till now.
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload
Zeinab Ghofrani wrote:

I don't know the MATLAB toolboxes, but getting an answer is pretty easy using least-mean squares estimation. Assuming that you want a discrete-time model you have a system with output y and input u:
y_k = a_y_ + ay_ ... a_0y_ + b_u_ ...
If you have K samples you make a 1 x K vector of all the y_k (call it Y) and a 2n x K matrix (call it P), where each row is
[y_ y_ ... y_ u_ u_ ... u_]
You want to find the gain vector B that satisfies Y = P*B, which you do by finding B = P\Y (if I have my slash in the right direction).
Now, the caveat is that I have given you 1/4 of one page out of any book on system identification. For a nice well-behaved system where you have a pretty good idea of the number of poles and zeros and you have given it input data that shows their characteristics this will work fine.
Otherwise go read the book. I'd recommend one, but what I know on the subject comes from Astrom's "Adaptive Control" and the IEEE's "The Control Handbook".
--

Tim Wescott
Wescott Design Services
  Click to see the full signature.
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload

...
Y)

book

have

It can be much easier, depending on what the output looks like. If you think a single-pole low-pass filter describes your data, that model is entirely described by the peak value and time constant, i.e., A[1-e^(-t/tau) --- there are many methods to determine the time constant.
If you need to use a higher order model, then things get more complex, but you can still find the closed-form solution for the step response, and work from there.
Shouldn't need to use adaptive control techniques unless your model is time-varying, or unless you need to constantly estimate your transfer function while simultaneously trying to control your system.
Scott
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload
Scott Seidman wrote:

It's not really an adaptive control technique, it's just a system ID technique often used by adaptive controllers.
Yes, if you have a pretty good idea that it's a first- or second-order system you can just go by the parameters.
--

Tim Wescott
Wescott Design Services
  Click to see the full signature.
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload

True enough. I've done just that to estimate system parameters in biological systems (no control involved), but those vary with time, so the method really shines. It's a PITA to start with a continuous time transfer function, do the bilinear transformation to discretize it, and recover your parameters from the system ID output. It's not hard, but a little tedious. I suppose you can just go for a generic low-order discrete time system, and then do the inverse bilinear transform to bring it back to the continuous world, but I don't recall that being the approach I took. Perhaps I made life hard on myself for no reason.
There are just much better tools available for estimating simple low order transfer functions for LTI systems, especially when the inputs are nice enough to give you closed form solutions.
If the systems are higher order, and the input can take any form (i.e., no closed-form solution) I'd use the optimization toolbox instead of the sysid toolbox. Take a fairly generic transfer function of the appropriate order, take it back to the differential equation, derive your jacobian (that's an optional step, though, but it makes the optimization converge faster), use the ode tools to solve for the output and jacobian, and pass the whole kit and kaboodle to the nonlinear optimizer. Damn--looking at that (and I've done it plenty of times), maybe it is better to just do a sysid. The nifty thing about the optimization solution, though, is that you can get confidence intervals on your id'ed parameters, calculate coefficients of variation using your variance/covariance matrix to give you a clue as to the precision of your estimates, and do a sensitivity analysis for sufficient excitation. Powerful stuff from the world of parameter estimation.
Another option might be to use transfer function estimation techniques from the world of spectral analysis, like the Welch Periodogram. Take the estimated function, and fit it in the frequency domain to any transfer function you like using simple curve fitting tools. There's even the bonus that you can throw out frequency points of low coherence, and just use your most robust gain and phase estimations.
BTW, when you use a sysid for a non time varying system, do you turn the forgetting factor off, or just let the algorithm have its way with the i/o?
Scott
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload
Scott Seidman wrote:

Don't take me for an expert on this! I just answered his question because no one else jumped in there!
If I _did_ use the above technique you'll note that it has no forgetting factor -- it's just a pure LMS fit to the coefficients. To date I've always been able to either model a system adequately from first principals or run a sine-wave sweep on it and do my tuning in the frequency domain.
I have a pretty jaundiced view of the highfalutin toolboxes -- it seems that you're always twiddling with things in the end anyway, because your model never matches reality. I've never had a system that either absolutely required careful characterization and hand tuning, or was just so non-critical that I could do it by guess and by gosh.
When I run into the right system I'm sure my attitude will change...
--

Tim Wescott
Wescott Design Services
  Click to see the full signature.
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload

I agree completely with Tim. I have tried quite a few of these techniques over the years and been thoroughly disappointed with the results. Even with systems that I know the transfer function of, the results obtained have been ambiguous at best or plain incorrect at worst. I think parameter ID is one of those things where theory and practice are quite widely separated and where theoretical assumptions are easily violated in practice. If there is significant measurement noise present, then things get even worse. I remember reading a paper where Kumpati Narendra at Yale stated that in 30 years of his experience with parameter estimation, convergence to the correct parameter values never happened even in simulation studies, even though the modelling error signal went to zero as required by the algorithms. This echoes what I have found.
For these reasons, I always try and do a frequency response determination rather than a parameter estimation of a practical plant and have had much more success. Using good engineering judgement and a few careful experiments seems to be the best way to get good results with model fitting.
Fred.
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload
fred snipped-for-privacy@hotmail.com (Fred Stevens) wrote in

Agreed to an extent-- its more of an art than a math. It starts out with your choice of a model. Then it goes to careful investigation of the condition number and the jacobian, then you calculate sensitivity of the output to each parameter as a function of time. You need to investigate the Coefficient of Variation of each estimated parameter.
A few simple rules of thumb seem to help. 1- Never estimate what you can measure 2- If you have too many parameters, things start to look silly in parameter space. If you're looking at more than five or six parameters, think long and hard before you embark on a parameter estimation 3- Always start from a variety of initial guesses fairly well distributed in parameter space. Many of these methods get stuck at local minima.
Some points we disagree on. If your error goes to zero, and you still don't have the "right" answer, then your model is overparametrized, or your input is not sufficient to stimulate the system. So, you're looking at a poor model choice, or a poor experimental design. Both situations are avoidable with care.
Also, sometimes the parameter estimation is unavoidable. Sure, dont estimate the parameters unless you need to know values of the parameters. A frequency response is fine, but once you're trying to fit it to a specific form of transfer function, you're either doing a curve fit or a parameter estimation (though some incorrectly use those terms interchangably). In fact, for much of my work, when I don't really care about the form of transfer function, I use random signals techniques to estimate the transfer function.
Last, but not least, the sysid tools are largely recursive least squares, as are the parameter estimation tools. They aren't that different.
Scott
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload
On Thu, 08 Jul 2004 16:54:28 -0700, Tim Wescott wrote:

<snip>
Maybe I got lost and maybe I'm late getting in on this thread, but how would you sweep a system with a motor fitted with a tach or position sensor? I understand sweeping an open loop for a PLL or amplifier, but I can't quite wrap my noodle around sweeping through a motor at this time.

<snip>
--
Best Regards,
Mike

Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload
Activ8 wrote:

You can sweep any system that will respond to a sine-wave input, and that won't be expensive to run a sweep on (running a 1 microhertz to 1 millihertz sweep on an oil refinery, for instance, may not be be cost effective).
You can run a sweep on a closed-loop system to get open loop response: you inject a sine signal at a summing junction that you've made for the purpose, then measure the amplitude and phase right _after_ the junction, and right before. Then you divide the output by the input for the value of the transfer function at that frequency.
You can buy an instrument (formerly control system analyzer, now a dynamic system analyzer) to do this work for you, or you can build it into your software.
--

Tim Wescott
Wescott Design Services
  Click to see the full signature.
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload

It also wouldn't be safe -- the operators would shoot you. Actually running a frequency sweep on a process system would yield very little of value since the system is highly nonlinear and all responses depend on signal amplitude not to mention factors such as component saturation. What is usually done is a small bump test near the setpoint.
Walter.
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload
Walter Driedger wrote:

I was thinking of some of your comments about sweeping when I wrote that.
Now that we've drawn you in, what do you do with your data once you've done your bump test?
--

Tim Wescott
Wescott Design Services
  Click to see the full signature.
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload
reply at bottom

running
since
amplitude
done
Good question. I think the original question was how to determine the transfer function of a process. I don't recall what the purpose of the transfer function was. If it was for controller tuning, what I am saying is basically to use Nicols Ziegler. For other purposes -- remember this, a real plant has many and complex components to its transfer function. For example, a temperature control loop has components ranging from the multiple lags of a temperature element to the sampling lag of the controller. This is followed by the control valve dynamics, the fuel gas flow dynamics, the entire dynamics of combustion, etc. etc. etc. It's bottomless. In reality, there will be one or two components that dominate the particular phenomenon of interest. A first order lag, a dead time, a process gain, that's usually all you need to produce a model good enough for most things. The basic bump test will provide all these. If there is a resonance present in the system a bump test will provide frequency and decay rate. Do you really need more?
Walter.
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload
On Mon, 19 Jul 2004 12:09:01 -0700, Tim Wescott wrote:

<snip>
I think my prob here was in thinking about cutting the loop and applying the sweep at any point in the loop. Naturally, if I apply a low (enough) level sweep to the motor input, it won't turn and we'll measure zero at the tach output.
Never mind :)
--
Best Regards,
Mike

Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload
Activ8 wrote:

Running a sweep assumes that the system behaves linearly -- in the case of a motor with friction you have a fairly severe nonlinearity that requires you to be careful about the amplitude of your sweep if you want "meaningful" results.
There's actually a technique for nonlinear control that says you run a sine-sweep of a nonlinear system & pretend that it's for a linear one, then design your controller. I don't believe that the results are guaranteed solid, but it gives you a pretty good idea of how to proceed, and it actually works in many cases, particularly if you sweep at a number of different amplitudes.
--

Tim Wescott
Wescott Design Services
  Click to see the full signature.
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload

I think the theoretical basis for this would be the describing function approach, wouldn't it?
Anyway, the sine sweep seems better for the real world because since the input energy is concentrated at one frequency at a time, one can operate the test at a lower drive overall level compared to other stimuli.
Incidently Hewlett Packard (or is it Agilent, I suppose so--god how I hate that name) has an old application note showing many different schemes for performing a open loop measurement on a closed loop system. It was an instructive read since there turns out to be many possible configurations and some are not so good.
dave y.
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload
dave y. wrote:

Yes. I was driving myself up the wall trying to remember that name when I wrote this post, then it came to me about two minutes after I hit the "send" button.

I have seen the driving function recommended with a step input (where you'd use N-Z or a more formal system ID/control design scheme). The rational is that if you're primarily controlling the system with a step input then following a step is a good thing.

I'd like to see that note -- do you have a link?

--

Tim Wescott
Wescott Design Services
  Click to see the full signature.
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload
On Sat, 24 Jul 2004 14:16:09 -0700, Tim Wescott
...

The first application note below is the one I was referring to. The second one is related. Agilent doesn't list these notes--maybe they didn't 'translate' them over from the HP format. Maybe your local Agilent sales office could produce a copy. If you need one bad let me know and I'll see what I can do.
"Control System Loop Gain Measurements" Application Note 243-5 Revised Hewlett Packard Part Number 5091-3809E Printed 3/92
"Control System Measurement Fundamentals Using Dynamic Signal Analyzers and Accessories" Application Note 243-6 Hewlett Packard part Number 5091-5886E Printed 10/92
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload
The system whose transfer functiuon I want to guess has a pole at zero and another unknown pole and a delay (exp(-tau*s))
so I want to find the pole and the time constant (tau). Which command should I use in MATLAB?
Thanks Zeinab Ghofrani
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload
snipped-for-privacy@gmail.com (Zeinab Ghofrani) wrote in

I would try to find the delay by doing a cross-correlation of the input and output, or simply measuring it as best I can by viewing the step response. The second pole is going to be highly related to the decay time of the exponential decay of the step response. The gain of the system can be directly measured.
So, given the trasfer function ATs/(sT+1) (ignoring the time delay for the moment), "A" will be equal to the maximum of the step response. "T" is just the time constant of the exponential decay, which is the time (referenced to the start of the step response) at which the output is equal to A/e.
You can work out the matlab implementation on your own.
Scott
Add pictures here
<% if( /^image/.test(type) ){ %>
<% } %>
<%-name%>
Add image file
Upload

Polytechforum.com is a website by engineers for engineers. It is not affiliated with any of manufacturers or vendors discussed here. All logos and trade names are the property of their respective owners.