cartesian_to_equinoctial Subroutine

public subroutine cartesian_to_equinoctial(mu, rv, evec)

Convert Cartesian coordinates to modified equinoctial elements (posigrade formulation).

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: mu

central body gravitational parameter ()
real(kind=wp), intent(in), dimension(6) :: rv

Cartesian state vector
real(kind=wp), intent(out), dimension(6) :: evec

Modified equinoctial element vector

Calls

proc~~cartesian_to_equinoctial~~CallsGraph proc~cartesian_to_equinoctial modified_equinoctial_module::cartesian_to_equinoctial proc~cross vector_module::cross proc~cartesian_to_equinoctial->proc~cross proc~unit vector_module::unit proc~cartesian_to_equinoctial->proc~unit

Called by

proc~~cartesian_to_equinoctial~~CalledByGraph proc~cartesian_to_equinoctial modified_equinoctial_module::cartesian_to_equinoctial proc~modified_equinoctial_test modified_equinoctial_module::modified_equinoctial_test proc~modified_equinoctial_test->proc~cartesian_to_equinoctial

Source Code

    subroutine cartesian_to_equinoctial(mu,rv,evec)

    implicit none

    real(wp),intent(in)               :: mu    !! central body gravitational parameter (\(km^3/s^2\))
    real(wp),dimension(6),intent(in)  :: rv    !! Cartesian state vector
    real(wp),dimension(6),intent(out) :: evec  !! Modified equinoctial element vector

    real(wp),dimension(3) :: r,v,hvec,hhat,ecc,fhat,ghat,rhat,vhat
    real(wp) :: hmag,rmag,p,f,g,h,k,L,kk,hh,s2,tkh,rdv

    r = rv(1:3)
    v = rv(4:6)

    rdv      = dot_product(r,v)
    rhat     = unit(r)
    rmag     = norm2(r)
    hvec     = cross(r,v)
    hmag     = norm2(hvec)
    hhat     = unit(hvec)
    vhat     = (rmag*v - rdv*rhat)/hmag
    p        = hmag*hmag / mu
    k        = hhat(1)/(one + hhat(3))
    h        = -hhat(2)/(one + hhat(3))
    kk       = k*k
    hh       = h*h
    s2       = one+hh+kk
    tkh      = two*k*h
    ecc      = cross(v,hvec)/mu - rhat
    fhat(1)  = one-kk+hh
    fhat(2)  = tkh
    fhat(3)  = -two*k
    ghat(1)  = tkh
    ghat(2)  = one+kk-hh
    ghat(3)  = two*h
    fhat     = fhat/s2
    ghat     = ghat/s2
    f        = dot_product(ecc,fhat)
    g        = dot_product(ecc,ghat)
    L        = atan2(rhat(2)-vhat(1),rhat(1)+vhat(2))

    evec = [p,f,g,h,k,L]

    end subroutine cartesian_to_equinoctial