Skip to content
Snippets Groups Projects
Select Git revision
  • 84c33e44c878aceab74b28b8ad810eef1d326a50
  • master default protected
2 results

potential_EBAs.m

Blame
  • Axel Schnitger's avatar
    Axel Schnitger authored
    m-files in the release directory are now executable from the release
    diretory. No need to change the path to the repo's root level.
    51ae6426
    History
    potential_EBAs.m 1.24 KiB
    function  v12sst = potential_EBAs(rho_dot,s_ea,s_eb,we)
    % POTENTIAL_EBAS.m calculates the potential along the GRACE simulated orbit. 
    
    % rho_dot is the range rate
    % s_ea and s_eb are the state vectors, position and velocity in itrf
    % we is the mean Earth rotation rate
    %
    % Written by Neda Darbeheshti, AEI Hannover, 2018-12.
    
    global e12 en er
    
    xdot12=s_eb(:,4:6)-s_ea(:,4:6);
    sumxdot=s_eb(:,4:6)+s_ea(:,4:6);
    
    %unit vectors
    
    e12=s_eb(:,1:3)-s_ea(:,1:3);
    e12 = bsxfun(@rdivide,e12,sqrt(sum(e12.^2,2)));
    
    en = [s_ea(:,2).*s_eb(:,3)-s_ea(:,3).*s_eb(:,2) s_ea(:,3).*s_eb(:,1)-s_ea(:,1).*s_eb(:,3) s_ea(:,1).*s_eb(:,2)-s_ea(:,2).*s_eb(:,1)];
    en = bsxfun(@rdivide,en,sqrt(sum(en.^2,2)));
    
    er = [e12(:,2).*en(:,3)-e12(:,3).*en(:,2) e12(:,3).*en(:,1)-e12(:,1).*en(:,3) e12(:,1).*en(:,2)-e12(:,2).*en(:,1)];
    er = bsxfun(@rdivide,er,sqrt(sum(er.^2,2)));
    
    e3=[0 0 1];
    
    % calculate potential along each vector
    for i=1:length(rho_dot)
        
    v12sst_rho(i)=.5*rho_dot(i)*e12(i,:)*sumxdot(i,:)';
    v12sst_n(i)=.5*en(i,:)*xdot12(i,:)'*en(i,:)*sumxdot(i,:)';
    v12sst_r(i)=.5*er(i,:)*xdot12(i,:)'*er(i,:)*sumxdot(i,:)';
    
    v12sst_rota(i)=-we^2/2*(norm(cross(e3,s_eb(i,1:3)))^2-norm(cross(e3,s_ea(i,1:3)))^2);
    
    end
    
    % sum all potentials 
    v12sst=v12sst_n+v12sst_r+v12sst_rho+v12sst_rota;
    end