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
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
|
/*
*+
* Name:
* palDtps2c
* Purpose:
* Determine RA,Dec of tangent point from coordinates
* Language:
* Starlink ANSI C
* Type of Module:
* Library routine
* Invocation:
* palDtps2c( double xi, double eta, double ra, double dec,
* double * raz1, double decz1,
* double * raz2, double decz2, int *n);
* Arguments:
* xi = double (Given)
* First rectangular coordinate on tangent plane (radians)
* eta = double (Given)
* Second rectangular coordinate on tangent plane (radians)
* ra = double (Given)
* RA spherical coordinate of star (radians)
* dec = double (Given)
* Dec spherical coordinate of star (radians)
* raz1 = double * (Returned)
* RA spherical coordinate of tangent point, solution 1 (radians)
* decz1 = double * (Returned)
* Dec spherical coordinate of tangent point, solution 1 (radians)
* raz2 = double * (Returned)
* RA spherical coordinate of tangent point, solution 2 (radians)
* decz2 = double * (Returned)
* Dec spherical coordinate of tangent point, solution 2 (radians)
* n = int * (Returned)
* number of solutions: 0 = no solutions returned (note 2)
* 1 = only the first solution is useful (note 3)
* 2 = both solutions are useful (note 3)
* Description:
* From the tangent plane coordinates of a star of known RA,Dec,
* determine the RA,Dec of the tangent point.
* Authors:
* PTW: Pat Wallace (STFC)
* TIMJ: Tim Jenness (JAC, Hawaii)
* {enter_new_authors_here}
* Notes:
* - The RAZ1 and RAZ2 values are returned in the range 0-2pi.
* - Cases where there is no solution can only arise near the poles.
* For example, it is clearly impossible for a star at the pole
* itself to have a non-zero XI value, and hence it is
* meaningless to ask where the tangent point would have to be
* to bring about this combination of XI and DEC.
* - Also near the poles, cases can arise where there are two useful
* solutions. The argument N indicates whether the second of the
* two solutions returned is useful. N=1 indicates only one useful
* solution, the usual case; under these circumstances, the second
* solution corresponds to the "over-the-pole" case, and this is
* reflected in the values of RAZ2 and DECZ2 which are returned.
* - The DECZ1 and DECZ2 values are returned in the range +/-pi, but
* in the usual, non-pole-crossing, case, the range is +/-pi/2.
* - This routine is the spherical equivalent of the routine sla_DTPV2C.
* History:
* 2012-02-08 (TIMJ):
* Initial version with documentation taken from Fortran SLA
* 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 Lesser 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 Lesser General Public License for more details.
*
* You should have received a copy of the GNU Lesser General
* License along with this program. If not, see
* <http://www.gnu.org/licenses/>.
* Bugs:
* {note_any_bugs_here}
*-
*/
#include "pal.h"
#include "pal1sofa.h"
#include <math.h>
void
palDtps2c( double xi, double eta, double ra, double dec,
double * raz1, double * decz1,
double * raz2, double * decz2, int *n) {
double x2;
double y2;
double sd;
double cd;
double sdf;
double r2;
x2 = xi * xi;
y2 = eta * eta;
sd = sin(dec);
cd = cos(dec);
sdf = sd * sqrt(x2 + 1. + y2);
r2 = cd * cd * (y2 + 1.) - sd * sd * x2;
if (r2 >= 0.) {
double r;
double s;
double c;
r = sqrt(r2);
s = sdf - eta * r;
c = sdf * eta + r;
if (xi == 0. && r == 0.) {
r = 1.;
}
*raz1 = eraAnp(ra - atan2(xi, r));
*decz1 = atan2(s, c);
r = -r;
s = sdf - eta * r;
c = sdf * eta + r;
*raz2 = eraAnp(ra - atan2(xi, r));
*decz2 = atan2(s, c);
if (fabs(sdf) < 1.) {
*n = 1;
} else {
*n = 2;
}
} else {
*n = 0;
}
return;
}
|
