diff options
Diffstat (limited to 'ast/pal/palDtps2c.c')
-rw-r--r-- | ast/pal/palDtps2c.c | 151 |
1 files changed, 151 insertions, 0 deletions
diff --git a/ast/pal/palDtps2c.c b/ast/pal/palDtps2c.c new file mode 100644 index 0000000..ecb3090 --- /dev/null +++ b/ast/pal/palDtps2c.c @@ -0,0 +1,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; +} |