MATLAB:Transfer Functions

Example

Analytical

Assume you have a transfer function \(\mathbb{H}(s)\) for a low-pass \(RC\) circuit where you have a function for the ratio of the Laplace transform of the voltage across the capacitor to the voltage across the series combination of a capacitor and a resistor:

\( \begin{align} \mathbb{H}(s)=\frac{\frac{1}{sC}}{R+\frac{1}{sC}}=\frac{1}{(RC)s+1} \end{align} \)

There are a few different ways to examine the magnitude and phase content of the Fourier version of this transfer function:

\( \begin{align} \mathbb{H}(j\omega)=\frac{1}{(RC)j\omega+1} \end{align} \)

over a range of frequencies. This document will show how to use MATLAB's tf function to set up and analyze the magnitude and phase of the transfer function of circuit.

%% Circuit constants
R = 10000;
C = 22e-9;

%% Set up transfer function
%  Create "s" as a transfer function for use later
s = tf([1 0], [1]);
%  Use s to generate transfer function for circuit
H = 1 / (s * R * C + 1);
%  Generate list of frequencies; must use angular frequencies in bode command
F = logspace(1, 5, 1000);
Omega = 2 * pi * F;
%  Use bode command to analyze transfer function
[HMag, HPhase, HOmega] = bode(H, Omega);
HMag   = squeeze(HMag);
HPhase = squeeze(HPhase);

%% Make plot
figure(1); clf
%  Magnitude plot on top
subplot(2,1,1)
semilogx(HOmega, 20*log10(HMag), 'k-')
xlabel('\omega, rad/s')
ylabel('|H|, dB')
%  Phase plot on bottom
subplot(2,1,2)
semilogx(HOmega, HPhase, 'k-')
xlabel('\omega, rad/s')
ylabel('\angle H, rad')

Theoretical

If you have a data set and want to find the theoretical transfer function between two variables in the set, you can have MATLAB come up with a transfer function estimate using the tfestimate command. In the example below, assume you are trying to find an estimate for the transfer function:

\( \begin{align} \mathbb{H}_{\mbox{est}}(j\omega)&=\frac{\mathbb{V}_{out}(j\omega)}{\mathbb{V}_{in}(j\omega)} \end{align} \)

to generate a magnitude plot and a phase plot of an experimentally determined transfer function. MATLAB's tfestimate will produce a numerical estimate of the magnitude and phase of a transfer function given an input signal, an output signal, and possibly other information. The specific form of this command is:

[EstH, EstF] = tfestimate(Vin, Vout, [], [], [], Fs)
EstMag   = abs(EstH);
EstPhase = angle(EstH);
EstOmega = EstF*2*pi;

where:

• Vin - vector containing the input voltage values
• Vout - vector containing the output voltage values
• Fs - the sampling frequency for the voltage values, in Hz
• EstH - vector containing complex numbers that contain the amplitudes and phases of the estimate of the transfer function
• EstF - vector containing the corresponding frequencies, in Hz, for the magnitudes and phases stored in EstH
• EstMag - magnitude of the transfer function estimate
• EstPhase - phase of the transfer function estimate
• EstOmega - corresponding angular frequencies for the estimate values

The square brackets in the third through fifth arguments are placeholders for parameters whose default values are fine for this experiment. Essentially, this command will determine the frequency content of the input and of the output using subsections of the data; it will then compute the magnitude ratio and phase difference between the input and the output and provide an approximation of the transfer function at particular frequencies contained in EstF.

One issue that comes into play is you must make sure the input signal has energy at as many frequencies as possible to give tfestimate values to work with. Using a single-frequency cosine as an input, for example, might lead to a disaster - if that frequency happens to exactly hit one of the frequencies tfestimate is using - since there would only be one input frequency on which to base an estimate for the response to all frequencies. This would be similar to estimating the acoustics of a concert hall by hitting a single tuning fork.

It is important to note, however, that what MATLAB is really doing is assuming that your input is one period's worth of some periodic input - meaning \(T\) is \({t}_{\mbox{end}}-{t}_{\mbox{start}}\) and \(\omega_0=2\pi/T\). Furthermore, the tfestimate command breaks the signal up into several windows, further reducing the possibility that a problem of this nature will occur. If, however, the input signal does not have any energy in the higher frequency ranges, the tfestimate program will not be able to compensate. You should therefore make sure that whatever signal gets used contains a wide range of frequencies. Ideally, the input signal would be an impulse function since that contains all frequencies at equal amounts. Barring that - and the Laws of Nature do a pretty nice job of barring that - generating a noise signal or a sweeping frequency signal with energy in the bands of interest works well.

To make a plot, you can use code similar to:

%% Make plot
figure(1); clf
%  Magnitude plot on top
subplot(2, 1, 1)
semilogx(EstOmega, 20*log10(EstMag), 'k-')
xlabel('\omega, rad/s')
ylabel('|H|, dB')
%  Phase plot on bottom
subplot(2,1,2)
semilogx(EstOmega, EstPhase, 'k-')
xlabel('\omega, rad/s')
ylabel('\angle H, rad')

Comparison Between Theoretical and Experimental

Take the plotting examples and combine the semilogx arguments. Use different line styles. Add a legend. Comment to self about how easy that was!