butterworthncs(n, freq)

function [denpoles, gain] = butterworthncs(n, freq)

% The Butterworth Filter transfer function
%  Use zero-pole-gain block in Simulnik with the kth pole given by the
%  formula:
%                                            (1/2n)
%  Poles with real part <0 among    p  = i(-1)      2*pi*freq
%                                    k
%    Where:
%          freq  = the Butterworth filter design frequency in Hz.
%          n     = order of the desired filter
%          i     = sqrt(-1)

%  denroots are roots of -1, denpoles are left half plane roots only:
denroots = (i*roots([1 zeros(1,2*n-1) 1]))*2*pi*freq; 
denpoles = denroots(find(real(denroots)<=0));         

%  If there is a real pole, zero out the imaginary part.  There is only one
%  real pole (only when n is odd). Because of the numeric precision
%  for computing the roots, the imaginasry part will not be zero.

index           = find(abs(imag(denpoles))<1e-6); % Find 0 imag. part poles
denpoles(index) = real(denpoles(index));          % and make them real only

% To insure the steady state gain is 1, the Butterworth filter gain must be
% n times the magnitude of the Butterworth circle squared:

gain = (2*pi*freq)^n