Star proper motion: update star catalog data for space motion.
Status: support routine.
The starting and ending TDB epochs EP1A+EP1B and EP2A+EP2B are Julian Dates, apportioned in any convenient way between the two parts (A and B). For example, JD(TDB)=2450123.7 could be expressed in any of these ways, among others:
EPnA EPnB
2450123.7D0 0D0 (JD method)
2451545D0 -1421.3D0 (J2000 method)
2400000.5D0 50123.2D0 (MJD method)
2450123.5D0 0.2D0 (date & time method)
The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
In accordance with normal star-catalog conventions, the object's right ascension and declination are freed from the effects of secular aberration. The frame, which is aligned to the catalog equator and equinox, is Lorentzian and centered on the SSB.
The proper motions are the rate of change of the right ascension and declination at the catalog epoch and are in radians per TDB Julian year.
The parallax and radial velocity are in the same frame.
Care is needed with units. The star coordinates are in radians and the proper motions in radians per Julian year, but the parallax is in arcseconds.
The RA proper motion is in terms of coordinate angle, not true angle. If the catalog uses arcseconds for both RA and Dec proper motions, the RA proper motion will need to be divided by cos(Dec) before use.
Straight-line motion at constant speed, in the inertial frame, is assumed.
An extremely small (or zero or negative) parallax is interpreted to mean that the object is on the "celestial sphere", the radius of which is an arbitrary (large) value (see the STARPV routine for the value used). When the distance is overridden in this way, the status, initially zero, has 1 added to it.
If the space velocity is a significant fraction of c (see the constant VMAX in the routine STARPV), it is arbitrarily set to zero. When this action occurs, 2 is added to the status.
The relativistic adjustment carried out in the STARPV routine involves an iterative calculation. If the process fails to converge within a set number of iterations, 4 is added to the status.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=wp), | intent(in) | :: | ra1 | right ascension (radians), before |
||
| real(kind=wp), | intent(in) | :: | dec1 | declination (radians), before |
||
| real(kind=wp), | intent(in) | :: | pmr1 | RA proper motion (radians/year), before |
||
| real(kind=wp), | intent(in) | :: | pmd1 | Dec proper motion (radians/year), before |
||
| real(kind=wp), | intent(in) | :: | px1 | parallax (arcseconds), before |
||
| real(kind=wp), | intent(in) | :: | rv1 | radial velocity (km/s, +ve = receding), before |
||
| real(kind=wp), | intent(in) | :: | ep1a | "before" epoch, part A (Note 1) |
||
| real(kind=wp), | intent(in) | :: | ep1b | "before" epoch, part B (Note 1) |
||
| real(kind=wp), | intent(in) | :: | ep2a | "after" epoch, part A (Note 1) |
||
| real(kind=wp), | intent(in) | :: | ep2b | "after" epoch, part B (Note 1) |
||
| real(kind=wp), | intent(out) | :: | ra2 | right ascension (radians), after |
||
| real(kind=wp), | intent(out) | :: | dec2 | declination (radians), after |
||
| real(kind=wp), | intent(out) | :: | pmr2 | RA proper motion (radians/year), after |
||
| real(kind=wp), | intent(out) | :: | pmd2 | Dec proper motion (radians/year), after |
||
| real(kind=wp), | intent(out) | :: | px2 | parallax (arcseconds), after |
||
| real(kind=wp), | intent(out) | :: | rv2 | radial velocity (km/s, +ve = receding), after |
||
| integer, | intent(out) | :: | j |
|
subroutine STARPM ( ra1, dec1, pmr1, pmd1, px1, rv1, &
ep1a, ep1b, ep2a, ep2b, &
ra2, dec2, pmr2, pmd2, px2, rv2, j )
implicit none
real(wp),intent(in) :: ra1 !! right ascension (radians), before
real(wp),intent(in) :: dec1 !! declination (radians), before
real(wp),intent(in) :: pmr1 !! RA proper motion (radians/year), before
real(wp),intent(in) :: pmd1 !! Dec proper motion (radians/year), before
real(wp),intent(in) :: px1 !! parallax (arcseconds), before
real(wp),intent(in) :: rv1 !! radial velocity (km/s, +ve = receding), before
real(wp),intent(in) :: ep1a !! "before" epoch, part A (Note 1)
real(wp),intent(in) :: ep1b !! "before" epoch, part B (Note 1)
real(wp),intent(in) :: ep2a !! "after" epoch, part A (Note 1)
real(wp),intent(in) :: ep2b !! "after" epoch, part B (Note 1)
real(wp),intent(out) :: ra2 !! right ascension (radians), after
real(wp),intent(out) :: dec2 !! declination (radians), after
real(wp),intent(out) :: pmr2 !! RA proper motion (radians/year), after
real(wp),intent(out) :: pmd2 !! Dec proper motion (radians/year), after
real(wp),intent(out) :: px2 !! parallax (arcseconds), after
real(wp),intent(out) :: rv2 !! radial velocity (km/s, +ve = receding), after
integer,intent(out) :: j !! status:
!! * -1 = system error (should not occur)
!! * 0 = no warnings or errors
!! * 1 = distance overridden (Note 6)
!! * 2 = excessive velocity (Note 7)
!! * 4 = solution didn't converge (Note 8)
!! * else = binary logical OR of the above warnings
! Astronomical unit (m, IAU 2012)
real(wp),parameter :: aum = 149597870.7d3
! Speed of light (au per day)
real(wp),parameter :: c = d2s*cmps/aum
real(wp) :: pv1(3,2), r, tl1, dt, pv(3,2), r2, rdv, v2, &
c2mv2, tl2, pv2(3,2)
integer :: j1, j2
! RA,Dec etc. at the "before" epoch to space motion pv-vector.
call STARPV ( ra1, dec1, pmr1, pmd1, px1, rv1, pv1, j1 )
! Light time when observed (days).
call PM ( pv1, r )
tl1 = r / c
! Time interval, "before" to "after" (days).
dt = ( ep2a-ep1a ) + ( ep2b-ep1b )
! Move star along track from the "before" observed position to the
! "after" geometric position.
call PVU ( dt+tl1, pv1, pv )
! From this geometric position, deduce the observed light time (days)
! at the "after" epoch (with theoretically unneccessary error check).
call PDP ( pv(1,1), pv(1,1), r2 )
call PDP ( pv(1,1), pv(1,2), rdv )
call PDP ( pv(1,2), pv(1,2), v2 )
c2mv2 = c*c - v2
if ( c2mv2 <= 0.0_wp ) then
j = -1
return
end if
tl2 = ( - rdv + sqrt(rdv*rdv + c2mv2*r2) ) / c2mv2
! Move the position along track from the observed place at the
! "before" epoch to the observed place at the "after" epoch.
call PVU ( dt + ( tl1-tl2 ), pv1, pv2 )
! Space motion pv-vector to RA,Dec etc. at the "after" epoch.
call PVSTAR ( pv2, ra2, dec2, pmr2, pmd2, px2, rv2, j2 )
! Return the status.
if ( j2 /= 0 ) j1 = -1
j = j1
end subroutine STARPM