mean_ecliptic_to_equatorial_rotmat Function

public pure function mean_ecliptic_to_equatorial_rotmat(obliquity_func) result(rot)

Rotation matrix from Mean Ecliptic to J2000.

Reference

  • https://en.wikipedia.org/wiki/Ecliptic_coordinate_system

Arguments

Type IntentOptional Attributes Name
procedure(mean_obliquity_func), optional :: obliquity_func

optional function to compute the mean obliquity. If not present, then mean_obliquity_of_ecliptic_iau1980 is used.

Return Value real(kind=wp), dimension(3,3)

rotation matrix


Calls

proc~~mean_ecliptic_to_equatorial_rotmat~~CallsGraph proc~mean_ecliptic_to_equatorial_rotmat obliquity_module::mean_ecliptic_to_equatorial_rotmat proc~mean_obliquity_of_ecliptic_iau1980 obliquity_module::mean_obliquity_of_ecliptic_iau1980 proc~mean_ecliptic_to_equatorial_rotmat->proc~mean_obliquity_of_ecliptic_iau1980

Called by

proc~~mean_ecliptic_to_equatorial_rotmat~~CalledByGraph proc~mean_ecliptic_to_equatorial_rotmat obliquity_module::mean_ecliptic_to_equatorial_rotmat proc~equatorial_to_mean_ecliptic_rotmat obliquity_module::equatorial_to_mean_ecliptic_rotmat proc~equatorial_to_mean_ecliptic_rotmat->proc~mean_ecliptic_to_equatorial_rotmat proc~get_c_cdot_ecliptic transformation_module::ecliptic_frame%get_c_cdot_ecliptic proc~get_c_cdot_ecliptic->proc~mean_ecliptic_to_equatorial_rotmat proc~get_c_cdot_ecliptic->proc~equatorial_to_mean_ecliptic_rotmat proc~transformation_module_test transformation_module::transformation_module_test proc~transformation_module_test->proc~get_c_cdot_ecliptic

Source Code

    pure function mean_ecliptic_to_equatorial_rotmat(obliquity_func) result(rot)

    implicit none

    real(wp),dimension(3,3) :: rot  !! rotation matrix
    procedure(mean_obliquity_func),optional :: obliquity_func !! optional function to compute
                                                              !! the mean obliquity. If not
                                                              !! present, then
                                                              !! [[mean_obliquity_of_ecliptic_iau1980]]
                                                              !! is used.

    real(wp) :: e !! mean obliquity at J2000 (rad)
    real(wp) :: s,c

    if (present(obliquity_func)) then
        e = obliquity_func(0.0_wp)
    else
        e = mean_obliquity_of_ecliptic_iau1980(0.0_wp)
    end if
    e = e * deg2rad
    s = sin(e)
    c = cos(e)

    rot(:,1) = [1.0_wp, 0.0_wp, 0.0_wp]
    rot(:,2) = [0.0_wp, c, s]
    rot(:,3) = [0.0_wp, -s, c]

    end function mean_ecliptic_to_equatorial_rotmat