Phase shift diagram for low pass filter...

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)

[Jan]

Currently we know, that you are getting it half of what you want. But I do not understand, which half is still missing. Could you please elaborate this again with all details? Which units do you want to change?

Answer 1

[Jan]

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);