MATLAB: Phase shift diagram for low pass filter…

low-pass filterphase shift plot

Hello All, Well I am designing a low pass filter. I was able to make a lot of progress however fell short at the phase shift diagram, if any one knows if or what I'm doing wrong it would help a lot. I'm getting "half" of what I want and not quite right on the units. it looks like the image under "Phase Response" on this website: http://www.play-hookey.com/ac_theory/filters/lo_pass_filters.html
Here is the code so far, (Just a note it takes about 30 seconds, any advice on making it faster will help too!):
clc
clear all
format long
R = 50.6;
C = 1.4*(10^-6);
F = 1/(2*pi*R*C);
Fc = [0:1:100000];
w = 1/(R*C);
wc = (2*pi*Fc)';
Vin = 1;
for n = 1:length(wc)
    Vout(n) = Vin*(1/sqrt(1+(wc(n)^2)*(R^2)*(C^2)));
end
loglog(Fc',Vout')
grid on
axis([0 10^5 10^-2 10])
title('Cut Off Frequency Chart')
ylabel('Attnuation (V_{out} / V_{in})')
xlabel('Frequency (Hz)')
hold on
Fc_actual = [1 5 10 20 30 150 200 300 1000 2000 2151 2500 2750 3100 4000 5200 7000 10000 15000 20000 32500];
Vout_actual = [1.101 1.088 1.016 1.132 1.157 1.149 1.137 1.112 0.938 0.730 0.707 0.643 0.605 0.554 0.447 0.35 0.25 0.16 0.08 0.04 0.01];
loglog(Fc_actual,Vout_actual,'r.')
loglog(F,0.707,'g*')
legend('Calculated','Actual','-3 dB point')
yL = get(gca,'YLim');
line([F F],yL,'Color','g');
xL = get(gca,'XLim');
line(xL,[0.707 0.707],'Color','g');
hold off
% this is the code for the phase shift
x = wc;
for m = 1:length(wc)
    y(m) = 1/(-atand(w/wc(m)));
end
loglog(x,y)

Best Answer

About the speed:
1. Omit clear all because it is a waste of time to remove all loaded functions from the memory.
2. Pre-allocate:
Vout = zeros(1, length(wc));  % Pre-allocate!!!
for n = 1:length(wc)
  Vout(n) = Vin / sqrt(1+(wc(n)^2)*(R^2)*(C^2));
end
...
y = zeros(1, length(wc));  % pre-allocate!!!
for m = 1:length(wc)
  y(m) = 1/(-atand(w/wc(m)));
end
3. Or use the vectorized notation instead of the FOR loop:
Vout = Vin ./ sqrt(1 + wc .^ 2 *(R ^ 2 * C ^ 2));
...
y = -1 ./ atand(w ./ wc);
Related Question