tcf2.m

Software for design and tuninig of SISO and MIMO contol systems

function y=tcf2(g)
%Remove numerator and denomiator zeros at the origin
%[n,d,ts,td]=tfdata(g,'v');
[n,d]=tfdata(g,'v');
td=g.inputdelay;

ln=length(n);
ld=length(d);

kk=0; %counter for origin zeros of nummerator.
for i=ln:-1:0 % count from the end of polynomial.
   if n(i)==0
      kk=kk+1;
   else
      break
   end
end
nn=n(1:(ln-kk)); % numerator after remove origin zeros.

jj=0; %counter for origin zeros of denominator.
for i=ld:-1:0 % count from the end of polynomial.
   if d(i)== 0
      jj=jj+1;
   else
      break
   end
end
dd=d(1:(ld-jj));% denominator after remove origin zeros.
gtf=tf(nn,dd,'inputdelay',td); % form the new transfer function from 
							  % the origin zero free polynomials.
lz=kk-jj; % lz is the order of pure integrator which will be added later.

%TCF TIME CONSTANT FORM
%    tcf(tf) is used to change transfer function into the 
%    time constant form.

[z,p,k]=zpkdata(gtf,'v');

% Pre-treatment of the zeros and poles

[zr,zc]=size(z);
[pr,pc]=size(p);
for i=1:zr
   if real(z(i,1))~=0 & (abs(imag(z(i,1)))/abs(real(z(i,1))))<0.005
      z(i,1)=real(z(i,1));
   end
end
for i=1:pr
   if real(p(i,1))~=0 & (abs(imag(p(i,1)))/abs(real(p(i,1))))<0.005
      p(i,1)=real(p(i,1));
   end
end

%gain in time constant form-------------------

zerogain=1;
for i=1:zr
   zerogain=-z(i,1)*zerogain;
end
polegain=1;
for i=1:pr
   polegain=-p(i,1)*polegain;
end
gain=k*zerogain/polegain;
if imag(gain)<0.00000001
   gain=real(gain);
end
%gain

%numerator and denominator polynomier coefficients matrix 
%num_max den_max in time constant form---------------
zinit=zr;
pinit=pr;
znew=z;
pnew=p;
zb=[];
pb=[];
flag=0;
for i=1:zr
   for j=1:pinit
      if abs(pnew(j,1)-z(i,1))< 0.001
         for k=1:zinit-1
            if k < i-flag
               zb(k,1)=znew(k,1);
            end
            if k >= (i-flag)
               zb(k,1)=znew(k+1,1);
            end
         end
         
         for k=1:pinit-1
            if k < j
               pb(k,1)=pnew(k,1);
            end
            if (k >=j)
               pb(k,1)=pnew(k+1,1);
            end
         end   
         pnew=pb;
         znew=zb;
         pinit=size(pnew,1);
         pb=[];
         zb=[];
         flag=flag+1;
         break   
     end
   end
   zinit=size(znew,1);
end
z=znew;
p=pnew;
[zr,zc]=size(z);
[pr,pc]=size(p);


for i=1:zr
   zero_tau(i,1)=1/(-z(i,1));
   %   num=conv([zero_tau(i,1),1],num);
end
num_max=[];
den_max=[];
i=1;
flag=0;
while i<=zr
   if isreal(z(i,1))
         num_max(i-flag,1)=0;
         num_max(i-flag,2)=zero_tau(i,1);
         num_max(i-flag,3)=1;
         i=i+1;
    else
         num_max(i-flag,:)=conv([zero_tau(i,1),1],[zero_tau(i+1,1),1]);     
         i=i+2;
         flag=flag+1;
    end
 end 
 if lz > 0
    clz=ceil(lz/2);
    flz=floor(lz/2);
    for i=1:flz
       num_max=[num_max;1 0 0];
    end
    if (clz-flz)>0
       num_max=[num_max;0 1 0];
    end
 end

for i=1:pr
   pole_tau(i,1)=1/(-p(i,1));   
end
i=1;
flag=0;
while i<=pr
   if isreal(p(i,1))
         den_max(i-flag,1)=0;
         den_max(i-flag,2)=pole_tau(i,1);
         den_max(i-flag,3)=1;
         i=i+1;
    else
         den_max(i-flag,:)=conv([pole_tau(i,1),1],[pole_tau(i+1,1),1]);     
         i=i+2;
         flag=flag+1;
    end
 end 
 if lz < 0
        clz=ceil(abs(lz/2));
        flz=floor(abs(lz/2));
    for i=1:flz
       den_max=[den_max;1 0 0];
    end
    if (clz-flz)>0
       den_max=[den_max;0 1 0];
    end
 end

%num_max
%den_max
%Convert the num_max and den_max into polynomier form for output.
%----------------------------------------------------------------
%numstr=mtx2str(num_max,2);
%denstr=mtx2str(den_max,2);

if (zr==0) & (lz==0)
   numstr='';
else
   numstr=['(',poly2str(num_max(1,:),'s'),')'];
   for i=2:size(num_max,1)
      str_n=['(',poly2str(num_max(i,:),'s'),')'];
      numstr=[numstr,'*',str_n];
   end
end
denstr=['(',poly2str(den_max(1,:),'s'),')'];
for i=2:size(den_max,1);
   str_d=['(',poly2str(den_max(i,:),'s'),')'];
   denstr=[denstr,'*',str_d];
end
%Make use of the functions in IMCTUNE--

%Form the transfer function in time constant form----------------
%Time Constant From:
%gain*numerator
%--------------exp(-td*s)
%Denomenator
%
disp(' ');
disp('Time Constant Form:');
disp(' ');
if (zr==0) & (lz==0)
   numstr=[num2str(gain)];
   disp(numstr);
else
   numstr=[num2str(gain),'*',numstr];
   disp(numstr);
end
disp('--------------------------------------------------------');
disp(denstr);
disp(' ');
strdelay=['input delay: ',num2str(td)];
disp(strdelay);