function [tz2,tn2,sz2,sn2,Ntz2,Nsz2,A_all2]=kh_get_pure_terms_svdcheck(tz,tn,sz,sn,Ntz,Nsz,rp)
%
% Script to calculate the proportion of variance explained by each of the
% three pure subduction processes (Aw, Af, Ah). Outputs the
% tz,sz,tn,sn,ntz,nsz terms that are explained by the combination of all 3
% pure processes
%
% All input terms are density levels*profiles (ie 40*m)
%
% Written by Kieran Helm (21-12-2007)
%
% Includes SVD check (as on p 121 7.31 of menke1984).  Use with
% subtraction method of calc ntz and nts (21-12-2007)
% A_all2=proportion of variance explained by all 3 terms

warning off
[rows cols]=size(tz);
tz2=NaN(rows,cols);
tn2=NaN(rows,cols);
sz2=NaN(rows,cols);
sn2=NaN(rows,cols);
Ntz2=NaN(rows,cols);
Nsz2=NaN(rows,cols);
A_all2=NaN(rows,cols);

for c=1:cols % 
        for r=1:rows % 
            x=1/(rp(r,c)-1);
%             x=[sz(r,c)-tz(r,c)]/(rp(r,c)-1); % added as test (18-12-2007)
            M=[-(rp(r,c)-1) 0 -rp(r,c);...
                1 rp(r,c) 0;...
                rp(r,c) rp(r,c) rp(r,c);...
                0 rp(r,c)-1 -1;...
                1 rp(r,c) 0;...
                1 1 1];
            M=x*M; % added coefficient x out the front (24-9-2007)
            d=[tz(r,c); tn(r,c); Ntz(r,c); sz(r,c); sn(r,c); Nsz(r,c)];
            Aw(r,c)=M(:,1)\d;    % (Af=Ah=0)
            Af(r,c)=M(:,2)\d;    % (Aw=Ah=0)
            Ah(r,c)=M(:,3)\d;    % (Aw=Af=0)
            
            % SVD to check closure:  p==2 
            if isnan(rp(r,c))==0; p=1;
                [U,S,V]=svd(M);
                A_all=V(:,1:p)*inv(S(1:p,1:p))*U(:,1:p)'*d; 
                A_allest=M*A_all;
                    tz2(r,c)=A_allest(1);
                    tn2(r,c)=A_allest(2);
                    Ntz2(r,c)=A_allest(3);
                    sz2(r,c)=A_allest(4);
                    sn2(r,c)=A_allest(5);
                    Nsz2(r,c)=A_allest(6);
                A_all2(r,c)=(A_allest'*A_allest)/(d'*d);
            end
        end
end

