信道迭代.m

一个双线性滤波器的设计

% 信道容量C的迭代算法                              %
% 函数说明：                                       %
% [CC,Paa]=ChannelCap(P,k) 为信道容量函数        	%
% 变量说明：                                       %
% P：输入的正向转移概率矩阵，k：迭代计算精度       %
% CC：最佳信道容量，Paa：最佳输入概率向量          %
% Pa：初始输入概率向量，Pba：正向转移概率矩阵    	%
% Pb：输出概率矩阵                                	%
% C：初始信道容量， r：输入符号数，s：输出符号数   %
%************************************************  %
function[CC,Paa]=ChannelCap(P,K)


%初始化Pa
P=input('请输入转移概率矩阵P:\n');
[r,s]=size(P);            %求输入输出符号的个数分别为r、s
Pa=(1/r)*ones(1,r);       %设初始概率分布为均匀分布
Pba=P;                    %正向信道转移概率矩阵
k=input('请输入迭代计算精度k:\n');


%下面进行迭代计算
n=0;
C=0;                   %设置c和cc的处置，满足abs(CC-C)>=k即可，使进入循环执行
CC=1;
while abs(CC-C)>=k
    n=n+1;
%求Pb
Pb=zeros(1,s);
for j=1:s
  for i=1:r
        Pb(j)=Pb(j)+Pa(i)*Pba(i,j);
    end
end
%求Paa
a=zeros(1,r);
b=zeros(1,r);
  for i=1:r
     for j=1:s
        a(i)=a(i)+Pba(i,j)*log2((Pba(i,j)+eps)/(Pb(j)+eps));
     end
     b(i)=exp(a(i));
end
sum1=0;
for i=1:r
sum1=sum1+Pa(i)*(b(i)+eps);
end
sum2=0;
for i=1:r
    if sum2<b(i)
        sum2=b(i);
    end
end
for i=1:r
    Paa(i)=Pa(i)*b(i)/(sum1+eps);
end
%求信道容量
C=log(sum1);
CC=log(sum2);
Pa=Paa;
end



%打印输出结果
fprintf('最佳输入概率分布'),Paa
fprintf('信道容量CC'),CC
fprintf('迭代次数n'),n
fprintf('输入符号数r'),r
fprintf('输出符号数s'),s
fprintf('迭代精度k'),k