function [TDFKnM, TDFKn0, Enm]= param_tdfkM(Cnm, Initval)% function [TDFKnM, TDFKn0, Enm]= param_tdfkM(Cnm, Initval)%   This file is part of the TFPM toolbox v1.0 (c)%   michael.jachan@tuwien.ac.at and underlies the GPL.% % Computes the right outer Kamen root [Kamen, LINALG, vol98, pp263,% 1988] of a cyclic TV polynomial given through Cnm. If no initial% value is  given (Initval==[]), the TDFR with minimum magnitude% and phase is taken?. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%if(0)% TEST TEST TEST TEST TEST TEST TEST TEST TEST TEST%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%clear;tfpm;MAR  =   0;LAR  =   0;MMA  =   3;LMA  =   3;N    =  16;re_im= 'r';mo_no= 'n';tfpm_file_gen;%-------------Cml= Bml;Cnm= param_ml_to_nm(Cml, N);%Cnm= ones(N, 1)*BLTI;Initval= [];%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%end;% TEST TEST TEST TEST TEST TEST TEST TEST TEST TEST%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Dimensions:NM= size(Cnm);N= NM(1);M= NM(2)-1;% TV Amplification factorTDFKn0= [];% \parc{0}[0:N-1]for n= 0:N-1   cn0= Cnm(n+1, 1);% \parc{0}[n]   Cnm(n+1, :)= Cnm(n+1, :)/cn0;   TDFKn0= [TDFKn0; cn0];end;if(nargin==1)   Initval= [];end;if(~M)%% Multiplier has no roots!   TDFKnm= zeros(N, M);% \PORO{0:0}[0:N-1]   Enm= Cnm;else   if(M==1)%% Order-1 system has very simple roots!      TDFKnm= -Cnm(:, 2)./Cnm(:, 1);% \PORO{1:1}[0:N-1]= -\parc{1}[0:N-1]   else      if(isempty(Initval))% Systen is TI prior to n= 0, takes min-root in abs and angle!!%         Initval= sort(roots(Cnm(1, :)));Initval= Initval(end);%*ones(M-1, 1);         Initval= j;      end;% The recursion      if(M>2)         TDFKnM= Initval.*ones(M-1, 1);         Enm= [];         for n= 0:N-1            eMM= 0;            for MM= M-1:-1:1               eMM= [eMM; double(( eMM(M-MM)-Cnm(n+1, 2+MM) )/TDFKnM(n+M-MM))];            end;            Enm= [Enm; [1; flipud(eMM(2:end))].'];            TDFKnM= [TDFKnM; eMM(M)-Cnm(n+1, 2)];         end;         TDFKnM= TDFKnM(M:end);%         TDFKnm= [param_tdfk(Enm, Initval) TDFKnM];      else%% The final step, M==2         TDFKnM= [];	 Enm= [];         for n= 0:N-1            p1= Cnm(n+1, 3)/Initval;            p2= -p1-Cnm(n+1, 2);            Initval= p2;            TDFKnM= [TDFKnM; p2];	    Enm= [Enm; -p1*Cnm(n+1, 1)];         end;      end;%if(M>2)%% The recursion   end;end;%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%if(0)% TEST TEST TEST TEST TEST TEST TEST TEST TEST TEST%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%Cml= param_rand(3, 1, N, 'r', 'n');Cnm= param_ml_to_nm(Cml, N);Initval= [];M= MMA;[TDFKnM, Enm]= param_tdfkM(Cnm, []);Initval= TDFKnM(N-M+1:N)TDFKnM[TDFKnm, TDFKn0]= param_tdfk(Cnm, []);figure(1);clf;tfpm_plot_roots(TDFKnm, max(1, max(max(abs(real(TDFKnm))))), ...			              max(1, max(max(abs(imag(TDFKnm))))), 'x')figure(2);clf;plot(real(TDFKnm8))mesh(abs(TDFKnm2(N+1:end, :) - TDFKnm8(7*N+1:end, :)))%[TDFKnm2, TDFKn02]= param_tdfk([Cnm;Cnm]);[TDFKnm8, TDFKn08]= param_tdfk([Cnm;Cnm;Cnm;Cnm;Cnm;Cnm;Cnm;Cnm]);[TDFKnm, TDFKn0]= param_tdfk(Cnm);clftfpm_plot_roots(TDFKnm, max(max(abs(real(TDFKnm)))), max(max(abs(imag(TDFKnm)))), 'o')%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%end;% TEST TEST TEST TEST TEST TEST TEST TEST TEST TEST%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%