source: trunk/FACT++/pal/palDtps2c.c@ 20022

Last change on this file since 20022 was 18347, checked in by tbretz, 9 years ago
File size: 4.6 KB
Line 
1/*
2*+
3* Name:
4* palDtps2c
5
6* Purpose:
7* Determine RA,Dec of tangent point from coordinates
8
9* Language:
10* Starlink ANSI C
11
12* Type of Module:
13* Library routine
14
15* Invocation:
16* palDtps2c( double xi, double eta, double ra, double dec,
17* double * raz1, double decz1,
18* double * raz2, double decz2, int *n);
19
20* Arguments:
21* xi = double (Given)
22* First rectangular coordinate on tangent plane (radians)
23* eta = double (Given)
24* Second rectangular coordinate on tangent plane (radians)
25* ra = double (Given)
26* RA spherical coordinate of star (radians)
27* dec = double (Given)
28* Dec spherical coordinate of star (radians)
29* raz1 = double * (Returned)
30* RA spherical coordinate of tangent point, solution 1 (radians)
31* decz1 = double * (Returned)
32* Dec spherical coordinate of tangent point, solution 1 (radians)
33* raz2 = double * (Returned)
34* RA spherical coordinate of tangent point, solution 2 (radians)
35* decz2 = double * (Returned)
36* Dec spherical coordinate of tangent point, solution 2 (radians)
37* n = int * (Returned)
38* number of solutions: 0 = no solutions returned (note 2)
39* 1 = only the first solution is useful (note 3)
40* 2 = both solutions are useful (note 3)
41
42
43* Description:
44* From the tangent plane coordinates of a star of known RA,Dec,
45* determine the RA,Dec of the tangent point.
46
47* Authors:
48* PTW: Pat Wallace (STFC)
49* TIMJ: Tim Jenness (JAC, Hawaii)
50* {enter_new_authors_here}
51
52* Notes:
53* - The RAZ1 and RAZ2 values are returned in the range 0-2pi.
54* - Cases where there is no solution can only arise near the poles.
55* For example, it is clearly impossible for a star at the pole
56* itself to have a non-zero XI value, and hence it is
57* meaningless to ask where the tangent point would have to be
58* to bring about this combination of XI and DEC.
59* - Also near the poles, cases can arise where there are two useful
60* solutions. The argument N indicates whether the second of the
61* two solutions returned is useful. N=1 indicates only one useful
62* solution, the usual case; under these circumstances, the second
63* solution corresponds to the "over-the-pole" case, and this is
64* reflected in the values of RAZ2 and DECZ2 which are returned.
65* - The DECZ1 and DECZ2 values are returned in the range +/-pi, but
66* in the usual, non-pole-crossing, case, the range is +/-pi/2.
67* - This routine is the spherical equivalent of the routine sla_DTPV2C.
68
69* History:
70* 2012-02-08 (TIMJ):
71* Initial version with documentation taken from Fortran SLA
72* Adapted with permission from the Fortran SLALIB library.
73* {enter_further_changes_here}
74
75* Copyright:
76* Copyright (C) 1995 Rutherford Appleton Laboratory
77* Copyright (C) 2012 Science and Technology Facilities Council.
78* All Rights Reserved.
79
80* Licence:
81* This program is free software: you can redistribute it and/or
82* modify it under the terms of the GNU Lesser General Public
83* License as published by the Free Software Foundation, either
84* version 3 of the License, or (at your option) any later
85* version.
86*
87* This program is distributed in the hope that it will be useful,
88* but WITHOUT ANY WARRANTY; without even the implied warranty of
89* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
90* GNU Lesser General Public License for more details.
91*
92* You should have received a copy of the GNU Lesser General
93* License along with this program. If not, see
94* <http://www.gnu.org/licenses/>.
95
96* Bugs:
97* {note_any_bugs_here}
98*-
99*/
100
101#include "pal.h"
102#include "pal1sofa.h"
103
104#include <math.h>
105
106void
107palDtps2c( double xi, double eta, double ra, double dec,
108 double * raz1, double * decz1,
109 double * raz2, double * decz2, int *n) {
110
111 double x2;
112 double y2;
113 double sd;
114 double cd;
115 double sdf;
116 double r2;
117
118 x2 = xi * xi;
119 y2 = eta * eta;
120 sd = sin(dec);
121 cd = cos(dec);
122 sdf = sd * sqrt(x2 + 1. + y2);
123 r2 = cd * cd * (y2 + 1.) - sd * sd * x2;
124 if (r2 >= 0.) {
125 double r;
126 double s;
127 double c;
128
129 r = sqrt(r2);
130 s = sdf - eta * r;
131 c = sdf * eta + r;
132 if (xi == 0. && r == 0.) {
133 r = 1.;
134 }
135 *raz1 = eraAnp(ra - atan2(xi, r));
136 *decz1 = atan2(s, c);
137 r = -r;
138 s = sdf - eta * r;
139 c = sdf * eta + r;
140 *raz2 = eraAnp(ra - atan2(xi, r));
141 *decz2 = atan2(s, c);
142 if (fabs(sdf) < 1.) {
143 *n = 1;
144 } else {
145 *n = 2;
146 }
147 } else {
148 *n = 0;
149 }
150 return;
151}
Note: See TracBrowser for help on using the repository browser.