tf2par.m

这是数字信号处理——matlab释义与实现一书的源代码

function [C,B,A] = tf2par(b,a);
%《数字信号处理教程——MATLAB释义与实现》
% 直接型到并联型的转换子程序
% 电子工业出版社出版  陈怀琛编著 2004年9月
%
% [C,B,A] = tf2par(b,a)
% --------------------------------------
%  C = 当length(b) >= length(a)时的多项式直通部分
%  B = 包含各bk的K乘2维实系数矩阵
%  A = 包含各ak的K乘3维实系数矩阵
%  b = 直接型的分子多项式系数
%  a = 直接型的分母多项式系数
% 
M = length(b); N = length(a);
[r1,p1,C] = residuez(b,a);		% 先求系统的单根p1,对应的留数r1及直接项C
p = cplxpair(p1,1e-9);	    % 用配对函数cplxpair由p1找共扼复根p，1e-9为容差
I = cplxcomp(p1,p);	        % 找p1变为p时的排序变化，MATLAB无此子程序，列于后面
r = r1(I);					% 让r1的排序同样变化成为r,以保持与极点对应
% 变为二阶子系统
K = floor(N/2); B = zeros(K,2); A = zeros(K,3);	% 二阶子系统变量的初始化
if K*2 == N; 			%  N 为偶, A(z)的次数为奇,有一个因子是一阶的 
	for i=1:2:N-2
		pi = p(i:i+1,:); 		% 取出一对极点
		ri = r(i:i+1,:);		% 取出一对对应的留数
		[Bi,Ai] = residuez(ri,pi,[]);		% 二个极点留数转为二阶子系统分子分母系数
		B(fix((i+1)/2),:) = real(Bi);		% 取Bi的实部，放入系数矩阵B的相应行中
		A(fix((i+1)/2),:) = real(Ai); 		% 取Ai的实部，放入系数矩阵A的相应行中
	end
	[Bi,Ai] = residuez(r(N-1),p(N-1),[]);		% 处理实单根
	B(K,:) = [real(Bi) 0]; A(K,:) = [real(Ai) 0];
else					%  N 为奇, A(z)的次数为偶,所有因子都是二阶的
    for i=1:2:N-1
		pi = p(i:i+1,:); 		% 取出一对极点
		ri = r(i:i+1,:);		% 取出一对对应的留数
		[Bi,Ai] = residuez(ri,pi,[]); 		% 二个极点留数转为二阶子系统分子分母系数
		B(fix((i+1)/2),:) = real(Bi); 	% 取Bi的实部，放入系数矩阵B的相应行中
		A(fix((i+1)/2),:) = real(Ai); 	% 取Ai的实部，放入系数矩阵A的相应行中
	end     
end