maxkurt.m

来自「这是盲信号的代码 都已经通过编译了 做这方面的同仁可以参考一下 我觉得蛮惯用的」· M 代码 · 共 30 行

function [y,num]=maxkurt(x)
%2005.12.21
[n T]=size(x);
y=x-mean(x,2)*ones(1,T);
y=inv(sqrtm(x*x'/T))*y;
encore=1;
num=1;
while encore, 
    encore=0;
    num=num+1;
    for p=1:n-1
        for q=p+1:n
            xi=y(p,:).*y(q,:);
            eta=y(p,:).*y(p,:)-y(q,:).*y(q,:);
            omega=atan2(4*(eta*xi'),eta*eta'-4*(xi*xi'));
            if abs(omega)>0.1/sqrt(T)
                encore=1;
                c=cos(omega/4);
                s=sin(omega/4);
                y([p q],:)=[c s;-s c]*y([p q],:);
            end
        end
    end
end
    return