www.pudn.com > SPICE.zip > aug_SPICE.m, change:2011-11-12,size:2192b


function  x = aug_SPICE(R,f,c,d,angleinterv) 
 
j = sqrt(-1); 
%R: augmented covariance matrix 
M0 = size(R,1); 
invR = pinv(R); 
%grid configuration 
theta = -90:angleinterv:90; 
N = length(theta); 
A = zeros(M0,N+M0); 
for  k = 1:N 
    A(:,k) = exp(-j*2*pi*(f/c)*(0:M0-1)'*d*sin((pi/180)*theta(k))); 
end 
A(:,N+1:N+M0) = eye(M0); 
%initialization: x_ini and sigma_ini 
x_ini = zeros(1,N); 
for k = 1:N 
    a = zeros(M0,1); 
    a(:) = A(:,k); 
    x_ini(k) = (a'*R*a)/(norm(a,2))^4; 
end 
[x0 inx0] = sort(x_ini); 
sigma_ini = 0; 
for k = 1:M0 
    a = zeros(M0,1); 
    a(:) = A(:,inx0(k)); 
    sigma_ini = sigma_ini+x_ini(inx0(k))*(norm(a,2))^2; 
end 
sigma_ini = sigma_ini/M0; 
%iteration 
x = zeros(1,N); 
W = zeros(1,N+M0); 
for k = 1:N+M0 
    a = zeros(M0,1); 
    a(:) = A(:,k); 
    W(k) = (a'*invR*a)/M0; 
end 
rr = sum(W(N+1:N+M0)); 
R_ini = A*diag(([x_ini sigma_ini*ones(1,M0)]).')*A'; 
rau = sqrt(rr)*sigma_ini*norm(pinv(R_ini)*R^(1/2),'fro'); 
for k = 1:N 
    a = zeros(M0,1); 
    a(:) = A(:,k); 
    rau = rau+sqrt(W(k))*x_ini(k)*norm(a'*pinv(R_ini)*R^(1/2),2); 
end 
for k = 1:N 
    a = zeros(M0,1); 
    a(:) = A(:,k); 
    x(k) = x_ini(k)*norm(a'*pinv(R_ini)*R^(1/2),2)/(sqrt(W(k))*rau); 
end 
sigma = sigma_ini*norm(pinv(R_ini)*R^(1/2),'fro')/(sqrt(rr)*rau); 
% while norm([x-x_ini sigma-sigma_ini],2)/norm([x sigma],2) >= 1e-10 
for p = 1:100 
    p 
    aa = norm([x-x_ini sigma-sigma_ini],2)/norm([x sigma],2) 
    x_ini(:) = x(:); 
    sigma_ini = sigma; 
    W = zeros(1,N+M0); 
    for k = 1:N+M0 
        a = zeros(M0,1); 
        a(:) = A(:,k); 
        W(k) = (a'*invR*a)/M0; 
    end 
    rr = sum(W(N+1:N+M0)); 
    R_ini = A*diag(([x_ini sigma_ini*ones(1,M0)]).')*A'; 
    rau = sqrt(rr)*sigma_ini*norm(pinv(R_ini)*R^(1/2),'fro'); 
    for k = 1:N 
        a = zeros(M0,1); 
        a(:) = A(:,k); 
        rau = rau+sqrt(W(k))*x_ini(k)*norm(a'*pinv(R_ini)*R^(1/2),2); 
    end 
    for k = 1:N 
        a = zeros(M0,1); 
        a(:) = A(:,k); 
        x(k) = x_ini(k)*norm(a'*pinv(R_ini)*R^(1/2),2)/(sqrt(W(k))*rau); 
    end 
    sigma = sigma_ini*norm(pinv(R_ini)*R^(1/2),'fro')/(sqrt(rr)*rau); 
end 
x = abs(x);