See Second Order RLC Filters (pdf file).
Most of the experiments in this lab use the "LC box" shown in the following figure.
The schematic of the box is
E1. LPF Step Response Measurement. Use the LC box (with L = 15 mH and C = 6.8 uF) to implement the following second order RLC LPF.
Use the function generator (rectangular signal with frequency approximately 100 Hz and amplitude approximately 1 Vpp) and the oscilloscope to display the step response of the circuit. Choose the resistor R such that the system is critically damped. Write down the value of R. Repeat the step response measurement with R = 0 and estimate the value of the damping ratio zeta from the overshoot in the step response. Compare your experimental values of R in the first case and zeta in the second case with what you would expect from a (nearly) ideal RLC LPF. You should find quite a bit of discrepancy. Explain where this discrepancy comes from and how it imposes a fundamental limit on what can be measured in the lab using the "step response". Hint: Remember Thévenin and the difference between ideal and real sources.
E2. LPF Frequency Response and Signal Processing in Matlab. An alternative to measuring the step response of a circuit with a function generator (that is designed to match a 50 ohm load) is to use a sinusoidal sweep and measure the magnitude and the phase of the frequency response. To get started, use the same circuit as in E1 with R = 10 ohms and make a recording using the soundcard probe and Audacity. Use the function generator in sweep mode to generate a sine signal that sweeps logarithmically from 20 Hz to 20 kHz in about 4 seconds. Use a signal level of about 500 mVpp and make sure the recording has sufficient level but is not distorted. Select one sweep and save it, e.g., as lpf01.wav.
Use the following Matlab script file to display the signals at the input (left channel) and ouptut (right channel) of the LPF. To become familiar with the basic features of Matlab, look at the Introduction to Matlab.
| Frequency Response Script V1.0: f_resp10.m |
|---|
%f_resp10 Frequency Response, V1.0
% Displays vS(t) and vO(t) of Frequency Response wav-file
% 2-13-07, P. Mathys
[y Fs] = wavread('lpf01.wav'); %Read recorded wav-file
%containing one sweep
vS = y(:,1)'; %Input vS(t) is left channel
vO = y(:,2)'; %Output vO(t) is right channel
t = [0:length(y)-1]/Fs; %Generate time axis (in sec)
subplot(2,1,1) %Selects plot in upper half of
%plot window
plot(t,vS,'-g') %Plot vS(t)
grid
title('Second Order LPF, Sinusoidal Sweep 20 Hz ... 20 kHz')
ylabel('Input v_S(t)')
subplot(2,1,2) %Selects plot in lower half of
%plot window
plot(t,vO,'-r') %Plot vO(t)
grid
ylabel('Output v_O(t)')
xlabel('t [sec]')
figure(gcf) %Brings current figure to foreground
|
Using the sample data in lpf01.wav the following graphs were obtained.
Since the input signal does not have constant amplitude (due to the fact that the input impedance of the RLC LPF is frequency dependent), it is necessary to look at the ratio AO/AS, where AS is the amplitude of the sine at the input of the LPF and AO is the amplitude of the sine at the output of the LPF, to obtain the magnitude of the frequency response. Thus, the first signal processing task that needs to be done in Matlab is to find the amplitude of a sinusoidal waveform.
The next figure shows how this can be made into a Matlab function (called avgT and saved as avgT.m).
The Matlab script file f_resp11.m shown below uses the avgT function to compute and plot the amplitudes AS and AO.
| Frequency Response Script V1.0: f_resp11.m |
|---|
%f_resp11 Frequency Response, V1.1
% Displays vS(t) and vO(t) of Frequency Response wav-file
% and computes amplitudes AS and AO
% 2-13-07, P. Mathys
[y Fs] = wavread('lpf01.wav'); %Read recorded wav-file
%containing one sweep
T = 0.02; %Averaging time constant
vS = y(:,1)'; %Input vS(t) is left channel
AS = sqrt(2*avgT(vS.^2,Fs,T)); %Amplitude of vS(t)
vO = y(:,2)'; %Output vO(t) is right channel
AO = sqrt(2*avgT(vO.^2,Fs,T)); %Amplitude of vO(t)
t = [0:length(y)-1]/Fs; %Generate time axis (in sec)
subplot(2,1,1) %Selects plot in upper half of
%plot window
plot(t,vS,':g',t,AS,'-b') %Plot vS(t)
grid
title('Second Order LPF, Sinusoidal Sweep 20 Hz ... 20 kHz')
ylabel('Input v_S(t), amplitude A_S')
subplot(2,1,2) %Selects plot in lower half of
%plot window
plot(t,vO,':r',t,AO,'-b') %Plot vO(t)
grid
ylabel('Output v_O(t), amplitude A_O')
xlabel('t [sec]')
figure(gcf) %Brings current figure to foreground
|
Using again the sample data in lpf01.wav the graphs shown below were obtained.
For the choice of T (20 ms for the above graph) a compromise needs to be made. If T is increased, the ripple that is quite visible in the graph at the lower frequencies is reduced. But this would introduce a larger delay at higher frequencies which shifts the amplitude to the right, some of which is also visible in the above graph.
The magnitude of the frequency response can now be easily obtained using |H| = AO/AS. To display it in dB compute 20 log10(AO/AS) an then plot the result as follows in Matlab:
magHdB = 20*log10(AO./AS); %Magnitude of H in dB subplot(2,1,1) %Selects plot in upper half of plot window plot(t,magHdB,'-b') %Plot |H| in dB
The result for the data in lpf01.wav is shown below.
This looks quite reasonable, except for the high frequency end where noise and/or distortion may limit the amount of attenuation that can be measured. One major shortcoming, however, still exists, and that is that the x-axis is still labeled in time and not in frequency.
Extract Frequency from Sinusoid. To extract the frequency from the measured signal vS(t), make use of the fact that the derivative of cos(2*pi*f*t) with respect to t is -2*pi*f*sin(2*pi*f*t), i.e., it is another sinusoid, but with amplitude proportional to f. Thus, the following Matlab statements can be used to extract f from vS(t) (in vS):
AS = sqrt(2*avgT(vS.^2,Fs,T)); %Amplitude of vS(t)
vScos = vS./AS; %Normalized cosine from vS(t)
dvScosdt = Fs*[0 diff(vScos)]; %d[vScos]/dt
%Derivative of normalized cosine
ffest = sqrt(2*avgT(dvScosdt.^2,Fs,T/2))/(2*pi);
%Estimated frequency in Hz
In principle it is now possible to obtain a plot of the magnitude of the frequency response using the Matlab commands
plot(ffest,magHdB,'-b') grid
which results in the following figure:
However, since the measurement was made using a logarithmic sweep and since frequency responses are easier to interpret when f is on a logarithmic scale, the Matlab command semilogx (type help semilogx at the Matlab prompt to see more details) should be used instead of the plot command. Also, there are some transient effects at the beginning of the ffest vector and some unreliable data for frequencies at the upper end which can be eliminated from the graph by dropping the first T seconds from ffest and limiting the display to the range from 10 Hz to 10 kHz. Here are the Matlab commands to achieve this:
ixdisp = [round(T*Fs):length(ffest)];
%Indexes for display (skip first T seconds)
semilogx(ffest(ixdisp),magHdB(ixdisp),'-b')
%Plot |H| in dB versus f on log10 scale
xlim([10 10000]) %Display only 10Hz to 10 kHz range
grid
The resulting graph for the data in lpf01.wav is shown below.
A Matlab function, called shift90 (save as shift90.m before using in Matlab), that implements the Hilbert transform (or -90 deg phase shift) is given below.
| -90 deg Phase Shift (Hilbert Transform) Function: shift90.m |
|---|
function y = shift90(x,Fs,ord)
%shift90 -90 deg phase shift (Hilbert transform filter)
% >>>>> y = shift90(x,Fs,ord) <<<<<
% where
% y: Filter output y(t), sampling rate Fs
% x: Filter input x(t), sampling rate Fs
% Fs: Sampling rate for x(t), y(t)
% ord: Filter order (must be even integer)
% 2-14-07, P. Mathys
x = reshape(x,1,length(x)); %Make x row vector
ord2 = floor(ord/2); %Half of filter order
t = [-ord2:ord2]/Fs; %Time axis for h(t)
h = zeros(size(t)); %Initialize h(t)
ix = find(t~=0);
h(ix) = 1./(pi*t(ix)); %Hilbert transform response
h(ix) = h(ix)-cos(pi*Fs*t(ix))./(pi*t(ix));
%H(f) limited to |f|<=Fs/2
y = filter(h,1,[x zeros(1,ord2)])/Fs; %Filter output y(t)
y = y(ord2+1:end); %Filter delay compensation
|
Use the Matlab script file given below to test the shift90 function. In particular, look at the filter order ord that is needed to convert a cosine to a sine (watch for both the correct phase shift and the correct amplitude, both should be within about 10% of the theoretical values) at the lower end (20 Hz) and the higher end (20 kHz) of the range of interest for the measurements made with the soundcard probes. Note that the filter order should be an even number and may have to be chosen quite high, on the order of several hundred to a few thousand.
Now all the bits and pieces to compute the phase shift phi between the measured input signal vS(t) and the measured output signal vO(t) are available and it remains to assemble them in the right order for Matlab to execute. Assume that vScos (= cos wt) and vO (= vO(t)=AOcos(wt+phi)) are available and then use the following Matlab statements to obtain phi:
ord = 1000; %Filter order vSsin = shift90(vScos,Fs,ord); %Shift vScos by -90 deg => vSsin a = avgT(2*vO.*vScos,Fs,T); %a = AO*cos(phi) b = avgT(2*vO.*vSsin,Fs,T); %b = -AO*sin(phi) phi = atan2(-b,a); %Phase phi of vO(t), modulo 2*pi %Use the unwrap command to remove jumps in phi by 360 deg %For BPF and HPF it may be better to leave out the unwrap command phi = unwrap(phi); %Phase unwrapped phi = (180/pi)*phi; %Phase in degreees
Note that phi is computed using the command atan2 which resolves the angle phi over the full range of ±pi or ±180 degrees (as opposed to atan which has only one input argument and a range of ±pi/2 or ±90 degrees). In some cases (especially for LPFs and BSFs) it is useful to use the unwrap command to remove jumps in the phase by ±2pi, thereby attempting to make the phase a continuous function.
Using
ixdisp = [round(T*Fs):length(ffest)];
%Indexes for display (skip first T seconds)
subplot(2,1,2) %Selects plot in lower half of plot window
semilogx(ffest(ixdisp),phi(ixdisp),'-r')
%Plot angle(H) in deg versus f on log10 scale
xlim([10 10000]) %Display only 10Hz to 10kHz range
grid
together with the magnitude of the frequency response computed earlier results in the following final frequency response display for the data in lpf01.wav:
The complete Matlab script file that was used to generate these graphs is shown below.
Use the Matlab code developed here to display the frequency response (both magnitude and phase) of your 2'nd order LPF filter. Compare the results from the Matlab graphs with frequency response measurements obtained using the oscilloscope. Point out any differences that you see, try to explain where the differences come from, and make suggestions how the differences could be minimized.
E3. Design and Measurement of BPF. Use the LC box (with L = 15 mH and C = 6.8 uF) to implement the following second order RLC BPF.
Choose R such that you obtain a BPF with with a center frequency of approximately 500 Hz and a half-power (or -3dB) bandwidth of 100 Hz. Measure and plot the magnitude and the phase of the frequency response and adjust R as necessary. Compare the performance of this BPF to the one you built in the previous lab from two 1'st order circuits.
E4. Design and Measurement of HPF or BSF. Choose to design either a HPF (highpass filter) or a BSF (bandstop filter) from the LC box. Use zeta approximately equal to 0.1 and measure the important parameters and the performance of the filter you built.
Lab worksheet in PDF format: lws03.pdf
Note: Each student needs to turn in a lab worksheet. The raw measured data for each student in a lab team will be the same, but the conclusions drawn from it are individual to each partner in the team.
©2003-2007, P. Mathys. Last revised: 2-27-07, PM.