transformation_module_test Subroutine

public subroutine transformation_module_test()

Uses

  • proc~~transformation_module_test~~UsesGraph proc~transformation_module_test transformation_module::transformation_module_test module~jpl_ephemeris_module jpl_ephemeris_module proc~transformation_module_test->module~jpl_ephemeris_module iso_fortran_env iso_fortran_env module~jpl_ephemeris_module->iso_fortran_env module~ephemeris_module ephemeris_module module~jpl_ephemeris_module->module~ephemeris_module module~kind_module kind_module module~jpl_ephemeris_module->module~kind_module module~ephemeris_module->module~kind_module module~celestial_body_module celestial_body_module module~ephemeris_module->module~celestial_body_module module~kind_module->iso_fortran_env module~celestial_body_module->module~kind_module module~base_class_module base_class_module module~celestial_body_module->module~base_class_module module~numbers_module numbers_module module~celestial_body_module->module~numbers_module module~numbers_module->module~kind_module

Transformation units test

Arguments

None

Calls

proc~~transformation_module_test~~CallsGraph proc~transformation_module_test transformation_module::transformation_module_test proc~get_c_cdot_ecliptic transformation_module::ecliptic_frame%get_c_cdot_ecliptic proc~transformation_module_test->proc~get_c_cdot_ecliptic proc~initialize_ephemeris jpl_ephemeris_module::jpl_ephemeris%initialize_ephemeris proc~transformation_module_test->proc~initialize_ephemeris proc~transform transformation_module::reference_frame%transform proc~transformation_module_test->proc~transform proc~equatorial_to_mean_ecliptic_rotmat obliquity_module::equatorial_to_mean_ecliptic_rotmat proc~get_c_cdot_ecliptic->proc~equatorial_to_mean_ecliptic_rotmat proc~mean_ecliptic_to_equatorial_rotmat obliquity_module::mean_ecliptic_to_equatorial_rotmat proc~get_c_cdot_ecliptic->proc~mean_ecliptic_to_equatorial_rotmat get_c_cdot get_c_cdot proc~transform->get_c_cdot proc~rvcto_rvcfrom_icrf transformation_module::rvcto_rvcfrom_icrf proc~transform->proc~rvcto_rvcfrom_icrf proc~equatorial_to_mean_ecliptic_rotmat->proc~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 get_rv get_rv proc~rvcto_rvcfrom_icrf->get_rv proc~from_primary_to_center transformation_module::two_body_rotating_frame%from_primary_to_center proc~rvcto_rvcfrom_icrf->proc~from_primary_to_center proc~from_primary_to_center->get_rv

Source Code

    subroutine transformation_module_test()

    use jpl_ephemeris_module, only: jpl_ephemeris

    implicit none

    real(wp),dimension(6),parameter :: initial_state = [10000.0_wp,&
                                                        0.0_wp,&
                                                        0.0_wp,&
                                                        1.0_wp,&
                                                        2.0_wp,&
                                                        3.0_wp] !! km, km/s
    real(wp),parameter :: et = zero           !! ephemeris time [sec]
    real(wp),parameter :: scale = 384400.0_wp !! scale factor [km]
    character(len=*),parameter :: ephemeris_file_421 = './eph/JPLEPH.421' !! JPL DE421 ephemeris file

    type(icrf_frame) :: from
    type(two_body_rotating_pulsating_frame) :: to
    type(jpl_ephemeris) :: eph421
    logical :: status_ok
    real(wp),dimension(6) :: rv_out
    type(ecliptic_frame) :: ecliptic_f
    real(wp),dimension(3,3) :: c
    integer :: i !! counter

    write(*,*) ''
    write(*,*) '---------------'
    write(*,*) ' transformation_module_test'
    write(*,*) '---------------'
    write(*,*) ''

    !open the ephemeris:
    call eph421%initialize(filename=ephemeris_file_421,status_ok=status_ok)
    if (.not. status_ok) error stop 'Error initializing DE421 ephemeris.'

    !initialize frames:
    !  [NOTE: need to make constructors for the frames]
    ! note: default is Earth-Moon-barycenter:
    from = icrf_frame(et=et)
    to   = two_body_rotating_pulsating_frame(et=et,scale=scale)

    call from%transform(initial_state,to,et,eph421,rv_out,status_ok)
    if (.not. status_ok) error stop 'Error in state transformation.'

    !results:
    write(*,*) ''
    write(*,'(A/,*(E30.16/))') 'initial state (J2000-Earth):', initial_state
    write(*,'(A/,*(E30.16/))') 'final state (Earth-Moon rotating, centered at barycenter, scale=384400):', rv_out
    write(*,*) ''

    ! ecliptic frame:
    ecliptic_f = ecliptic_frame(b=body_earth)
    call ecliptic_f%get_c_cdot(eph421,to_icrf=.true.,c=c,status_ok=status_ok)
    write(*,*) ''
    write(*,*) 'ecliptic to j2000:'
    do i=1,3
        write(*,*) c(i,:)
    end do

    ! from SPICE:
    ! rot = [1.0000000000000000E+00, 0.0000000000000000E+00, 0.0000000000000000E+00,
    !        0.0000000000000000E+00, 9.1748206206918181E-01, -3.9777715593191371E-01,
    !        0.0000000000000000E+00, 3.9777715593191371E-01, 9.1748206206918181E-01]
    !
    ! from FAT:
    !        1.0000000000000000        0.0000000000000000        0.0000000000000000
    !        0.0000000000000000       0.91748206206918181      -0.39777715593191371
    !        0.0000000000000000       0.39777715593191371       0.91748206206918181

    !close the ephemeris:
    call eph421%close()

    end subroutine transformation_module_test