%Prelab 1 Matlab
%This MATLAB file takes a set of frequency, peak to peak input/output, and
%phase measurements and matches them against a known transfer function to
%find the likely cutoff frequency, 1/RC. This frequency then can be used to
%find the unknown capacitance given the resistance is known.
%
%Written by Stephen Maloney
%
clear all; clc; close all;
%Input data from LABVIEW - given in the prelab
input = [120, 4.1202, 3.7137, -23.7380; 160, 4.0984, 3.6122, -27.6350; 200, ...
4.0939, 3.4766, -32.5860; 240, 4.0995, 3.3042, -37.1420; 280, 4.1080, ...
3.1847, -38.8100; 320, 4.0994, 3.0345, -41.8310; 360, 4.0994, 2.8817, ...
-45.8470; 400, 4.0833, 2.7716, -47.4650; 440, 4.0875, 2.5484, -52.3820];
%The frequency range we are interested in for the transfer function
w = 2*pi*input(:, 1)';
%The range to try for the cutoff frequency
wc = 2*pi*linspace(300, 400, 100);
%Build up the magnitude plot - output/input
mag = (input(:, 3)./input(:, 2))';
%Build up the phase plot
phase = input(:, 4)';
mError = zeros(1, length(wc));
pError = zeros(1, length(wc));
for idx = 1:length(wc),
%Calculate the frequency response given wc(idx) as the cutoff freq
trial = 1./(1+j*w./(wc(idx)));
%Find the mean square error for this generated response and the given
%response. Notice that dividing by n is not necessary as every trial
%would be divided by the same n
mError(idx) = sum((abs(trial) - mag).^2);
pError(idx) = sum((angle(trial)*180/pi - phase).^2);
%For debug - this will show the phase fit in progress
semilogx(w, angle(trial)*180/pi, w, phase);
pause(.01);
end
%Display the output
[minVal, minIdx] = min(mError);
disp('Magnitude Match Cutoff frequency (Hz):');
FcM = wc(minIdx)/(2*pi)
[minVal, minIdx] = min(pError);
disp('Phase Match Cutoff frequency (Hz):');
FcP = wc(minIdx)/(2*pi)
%If the resistor value is known, the cutoff frequency can be used to find
%the unknown capacitance; Tau = RC = 1/Fc, C = 1/R*fc
R = 1e3;
disp('Magnitude Fit Unkown Capacitor Value (F):');
CM = 1/(R*FcM)
disp('Phase Fit Unkown Capacitor Value (F):');
CM = 1/(R*FcP)