gld_score.m

这是盲信号的代码 都已经通过编译了 做这方面的同仁可以参考一下 我觉得蛮惯用的

function [phi, phi2]=gld_score(x,lambda);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% GLD_SCORE - Calculates the score function and its derivative
% for a distribution from the Generalized Lambda Distribution.
%
% Input:
%   x      = the realization matrix (sample values)
%   lambda = the distribution parameters vector 
%            [lambda(1),lambda(2),lambda(3),lambda(4)]
%
% Output:
%   phi   = the score function matrix corresponding to x  
%   phi2  = the derivative of the score function matrix corresponding to x  
%
% Copyright (c) Helsinki University of Technology,
% Signal Processing Laboratory,
% Jan Eriksson, Juha Karvanen, and Visa Koivunen.
%
% For details see the files readme.txt
% and gpl.txt provided within the same package.
%
% Last modified: 5.9.2000
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

if nargin<2, disp('Usege: gld_score(x,lambda)');
             error('Not enough input arguments.');  
end

cp=gld_cdf(x,lambda);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% The values of cp are forced to interval [10^(-tres),1-10^(-tres)]
% to ensure the numerical stability.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
tres=10;
cp(cp<10^(-tres))=10^(-tres);
cp(cp>1-10^(-tres))=1-10^(-tres);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% The functions are given by the equations (20) and (21) in
%   Eriksson, J., Karvanen, J., and Koivunen, V.:
%   "Source Distribution Adaptive Maximum Likelihood Estimation of ICA Model", 
%   Proceedings of Second International Workshop on
%   Independent Component Analysis and Blind Signal Separation,
%   Helsinki 2000, pp. 227--232
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
A=cp.^(lambda(3)-1);
B=(1-cp).^(lambda(4)-1);
L3=(1-lambda(3))*lambda(3);
L4=(1-lambda(4))*lambda(4);
denominator=lambda(3).*A+lambda(4).*B;
phi =lambda(2)*(A.*L3./cp-B.*L4./(1-cp))./denominator.^2;
phi2=2*phi.^2+(lambda(2)^2)*(A.*(lambda(3)-2).*L3./...
		(cp.^2)+B.*(lambda(4)-2).*L4./((1-cp).^2))./(denominator.^3);

% The end of the function %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%