1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
|
/******************************************************************************
* Johns Hopkins University
* Center For Astrophysical Sciences
* FUSE
******************************************************************************
*
* Synopsis: void state_geod(double pos,
* double mjdate, double lon, double lat)
*
* Description: Computes the geocentric longitude and latitude of FUSE
* from the given state 6-vector.
*
* Arguments: double pos (km) 3-vector xyz position of satellite
* double mjdate Modified Julian date of observation
*
* Returns: double lon, lat (deg) Geocentric longitude and latitude
*
* Calls: slaPreces preces.f
*
* History: 03/04/98 E. Murphy Begin work.
* 03/04/98 E. Murphy Initial version working
* 03/15/99 E. murphy Changed function definition to void
* (was double)
* 07/07/99 E. Murphy Changed call to use pos instead of x,y,z
* 05/31/00 peb Implemented cfortran.h calls for slalib
* functions.
* 12/18/03 bjg Change calfusettag.h to calfuse.h
* 06/17/04 bjg 1.4 Corrected cfortran call to sla_preces
*
* Ake, T. 1998 in The Scientific Impact of the Goddard
* High Resolution Spectrograph, ed. J. C. Brandt et al.,
* ASP Conference Series, in preparation.
*****************************************************************************/
#include <stdio.h>
#include <math.h>
#ifdef CFORTRAN
#include "cfortran.h"
PROTOCCALLSFSUB5(SLA_PRECES, sla_preces, STRING, DOUBLE, DOUBLE, PDOUBLE, \
PDOUBLE)
#define slaPreces(SYSTEM, EP0, EP1, RA, DC) \
CCALLSFSUB5(SLA_PRECES, sla_preces, STRING, DOUBLE, DOUBLE, PDOUBLE, \
PDOUBLE, SYSTEM, EP0, EP1, RA, DC)
#else
#include "slalib.h"
#include "slamac.h"
#endif
#include "calfuse.h"
void state_geod(double pos[3], double mjdate,
double *lon, double *lat, double *gmst)
{
double ra, dec, epoch2, c1, r, intmjd, frcmjd;
char fk5_st[10]="FK5";
/* pos[0]=x
pos[1]=y
pos[2]=z
*/
r=sqrt(pos[0]*pos[0]+pos[1]*pos[1]+pos[2]*pos[2]);
ra=atan2(pos[1],pos[0]);
dec=asin(pos[2]/r);
/* We must precess the J2000 RA and DEC to date of observation. This
* is a very small correction, and could probably be ignored.
*/
epoch2=2000.0-((51544.0-mjdate)/365.25);
#ifdef CFORTRAN
slaPreces(fk5_st, 2000.0, epoch2, ra, dec);
#else
slaPreces(fk5_st, 2000.0, epoch2, &ra, &dec);
#endif
/* Calculate the Greenwich Mean Sidereal Time in seconds
* Astronomical Almanac, 1998, p. B6.
*/
frcmjd=modf(mjdate, &intmjd);
c1=(intmjd-51544.5)/36525.0;
/* The following line gives the GMST at 0 hr UT */
*gmst=24110.54841 + 8640184.812866*c1+0.093104*c1*c1-6.2E-6*c1*c1*c1;
/* The following adds on the number of seconds since the UT above */
*gmst+=frcmjd*1.00273791*86400.0;
*gmst*=(2*M_PI)/86400.0;
while (*gmst < 0.0) *gmst+=2*M_PI;
while (*gmst > 2*M_PI) *gmst-=2*M_PI;
*lon=(ra-*gmst)/RADIAN;
while (*lon > 180.0) *lon-=360.0;
while (*lon < -180.0) *lon+=360.0;
*lat=dec/RADIAN;
*gmst*=24.0/(2*M_PI);
}
|
