/* Diurnal parallax, AA page D3
*/
#include "kep.h"
/* Earth radii per au. */
static double DISFAC;
extern double aearth, au;
/* distance to center of Earth, in Earth radii */
extern double trho;
/* observer's geocentric latitude, degrees */
extern double tlat;
int diurpx( last, ra, dec, dist )
double last; /* local apparent sidereal time, radians */
double *ra; /* right ascension, radians */
double *dec; /* declination, radians */
double dist; /* Earth - object distance, au */
{
double cosdec, sindec, coslat, sinlat;
double p[3], dp[3], x, y, z, D;
/* Don't bother with this unless the equatorial horizontal parallax
* is at least 0.005"
*/
if( dist > 1758.8 )
return(-1);
DISFAC = au / (0.001 * aearth);
cosdec = cos(*dec);
sindec = sin(*dec);
/* Observer's astronomical latitude
*/
x = tlat * DTR;
coslat = cos(x);
sinlat = sin(x);
/* Convert to equatorial rectangular coordinates
* in which unit distance = earth radius
*/
D = dist * DISFAC;
p[0] = D*cosdec*cos(*ra);
p[1] = D*cosdec*sin(*ra);
p[2] = D*sindec;
dp[0] = -trho*coslat*cos(last);
dp[1] = -trho*coslat*sin(last);
dp[2] = -trho*sinlat;
x = p[0] + dp[0];
y = p[1] + dp[1];
z = p[2] + dp[2];
D = x*x + y*y + z*z;
D = sqrt(D); /* topocentric distance */
/* recompute ra and dec */
*ra = zatan2(x,y);
*dec = asin(z/D);
showcor( "diurnal parallax", p, dp );
return(0);
}
syntax highlighted by Code2HTML, v. 0.9.1