function [tz,tn,Ntz,sz,sn,Nsz]=kh_get_pure_terms_reversed(Aw,Af,Ah,rp,tz1,sz1)
%
% Script to calculate M*a_est=d_est.  dest=
%
% All input terms are density levels*profiles (ie 40*m)
%
% Written by Kieran Helm (18-12-2007)
% NOTE PROBLEM - CAN'T JUST ADD 3 PURE TERMS - SOLUTION WILL AWAYSBE TOO
% LARGE:  USE SVD ... KH_GET_PURE_TERM_SVDCHECK (21-12-2007)

warning off
[rows cols]=size(Aw);
tz=NaN(rows,cols);
tn=NaN(rows,cols);
Ntz=NaN(rows,cols);
sz=NaN(rows,cols);
sn=NaN(rows,cols);
Nsz=NaN(rows,cols);


for c=1:cols % 
        for r=1:rows % 
             x=1/(rp(r,c)-1);
%            x=[sz1(r,c)-tz1(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)
            dest=M*[Aw(r,c); Af(r,c); Ah(r,c)];
                tz(r,c)=dest(1);
                tn(r,c)=dest(2);
                Ntz(r,c)=dest(3);
                sz(r,c)=dest(4);
                sn(r,c)=dest(5);
                Nsz(r,c)=dest(6);
        end
end

