-
Notifications
You must be signed in to change notification settings - Fork 12
/
palPreces.c
126 lines (104 loc) · 3.55 KB
/
palPreces.c
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
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
/*
*+
* Name:
* palPreces
* Purpose:
* Precession - either FK4 or FK5 as required.
* Language:
* Starlink ANSI C
* Type of Module:
* Library routine
* Invocation:
* void palPreces ( const char sys[3], double ep0, double ep1,
* double *ra, double *dc );
* Arguments:
* sys = const char [3] (Given)
* Precession to be applied: FK4 or FK5. Case insensitive.
* ep0 = double (Given)
* Starting epoch.
* ep1 = double (Given)
* Ending epoch
* ra = double * (Given & Returned)
* On input the RA mean equator & equinox at epoch ep0. On exit
* the RA mean equator & equinox of epoch ep1.
* dec = double * (Given & Returned)
* On input the dec mean equator & equinox at epoch ep0. On exit
* the dec mean equator & equinox of epoch ep1.
* Description:
* Precess coordinates using the appropriate system and epochs.
* Authors:
* PTW: Patrick T. Wallace
* TIMJ: Tim Jenness (JAC, Hawaii)
* {enter_new_authors_here}
* Notes:
* - Uses palPrec for FK5 data and palPrebn for FK4 data.
* - The epochs are Besselian if SYSTEM='FK4' and Julian if 'FK5'.
* For example, to precess coordinates in the old system from
* equinox 1900.0 to 1950.0 the call would be:
* palPreces( "FK4", 1900.0, 1950.0, &ra, &dc );
* - This routine will NOT correctly convert between the old and
* the new systems - for example conversion from B1950 to J2000.
* For these purposes see palFk425, palFk524, palFk45z and
* palFk54z.
* - If an invalid SYSTEM is supplied, values of -99D0,-99D0 will
* be returned for both RA and DC.
* History:
* 2012-03-02 (TIMJ):
* Initial version
* Adapted with permission from the Fortran SLALIB library.
* {enter_further_changes_here}
* Copyright:
* Copyright (C) 1995 Rutherford Appleton Laboratory
* Copyright (C) 2012 Science and Technology Facilities Council.
* All Rights Reserved.
* Licence:
* This program is free software; you can redistribute it and/or
* modify it under the terms of the GNU General Public License as
* published by the Free Software Foundation; either version 3 of
* the License, or (at your option) any later version.
*
* This program is distributed in the hope that it will be
* useful, but WITHOUT ANY WARRANTY; without even the implied
* warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
* PURPOSE. See the GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program; if not, write to the Free Software
* Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston,
* MA 02110-1301, USA.
* Bugs:
* {note_any_bugs_here}
*-
*/
#if HAVE_CONFIG_H
# include <config.h>
#endif
#include "pal.h"
#include "pal1sofa.h"
#include <string.h>
#if HAVE__STRNICMP || defined(_WIN32) || defined(_WIN64)
# define strncasecmp _strnicmp
#endif
void palPreces ( const char sys[3], double ep0, double ep1,
double *ra, double *dc ) {
double pm[3][3];
double v1[3];
double v2[3];
/* Generate appropriate precession matrix */
if ( strncasecmp( "FK4", sys, 3 ) == 0 ) {
palPrebn( ep0, ep1, pm );
} else if (strncasecmp( "FK5", sys, 3 ) == 0 ) {
palPrec( ep0, ep1, pm );
} else {
*ra = -99.0;
*dc = -99.0;
return;
}
/* Convert RA,Dec to x,y,z */
eraS2c( *ra, *dc, v1 );
/* Precess */
eraRxp( pm, v1, v2 );
/* Back to RA,Dec */
eraC2s( v2, ra, dc );
*ra = eraAnp( *ra );
}