1 /*
2 *+
3 * Name:
4 * palOne2One
5
6 * Purpose:
7 * File containing simple PAL wrappers for SLA routines that are identical in SOFA
8
9 * Invocation:
10 * Matches SLA API
11
12 * Language:
13 * Starlink ANSI C
14
15 * Type of Module:
16 * Library routine
17
18 * Description:
19 * Some SOFA/ERFA routines are identical to their SLA counterparts. PAL provides
20 * direct counterparts to these although it is generally a better idea to
21 * use the SOFA/ERFA routine directly in new code.
22 *
23 * The PAL routines with direct equivalents in SOFA/ERFA are:
24 * - palCldj
25 * - palDbear
26 * - palDaf2r
27 * - palDav2m
28 * - palDcc2s
29 * - palDcs2c
30 * - palDd2tf
31 * - palDimxv
32 * - palDm2av
33 * - palDjcl
34 * - palDmxm
35 * - palDmxv
36 * - palDpav
37 * - palDr2af
38 * - palDr2tf
39 * - palDranrm
40 * - palDsep
41 * - palDsepv
42 * - palDtf2d
43 * - palDtf2r
44 * - palDvdv
45 * - palDvn
46 * - palDvxv
47 * - palEpb
48 * - palEpb2d
49 * - palEpj
50 * - palEpj2d
51 * - palEqeqx
52 * - palFk5hz
53 * - palGmst
54 * - palGmsta
55 * - palHfk5z
56 * - palRefcoq
57
58 * Authors:
59 * TIMJ: Tim Jenness (JAC, Hawaii)
60 * DSB: David S Berry (JAC, Hawaii)
61 * {enter_new_authors_here}
62
63 * Notes:
64 * - Do not call these functions from other PAL functions. Always use
65 * the SOFA/ERFA routines directly in new code.
66 * - These are implemented as real functions rather than C preprocessor
67 * macros so there may be a performance penalty in using the PAL
68 * version instead of the SOFA/ERFA version.
69 * - Routines that take MJDs have SOFA/ERFA equivalents that have an explicit
70 * MJD offset included.
71 * - palEqeqx, palGmst and palGmsta use the IAU 2006 precession model.
72
73 * History:
74 * 2012-02-10 (TIMJ):
75 * Initial version
76 * Adapted with permission from the Fortran SLALIB library.
77 * 2012-03-23 (TIMJ):
78 * Update prologue.
79 * 2012-05-09 (DSBJ):
80 * Move palDrange into a separate file.
81 * 2014-07-15 (TIMJ):
82 * SOFA now has palRefcoq equivalent.
83 * {enter_further_changes_here}
84
85 * Copyright:
86 * Copyright (C) 2014 Tim Jenness
87 * Copyright (C) 2012 Science and Technology Facilities Council.
88 * All Rights Reserved.
89
90 * Licence:
91 * This program is free software: you can redistribute it and/or
92 * modify it under the terms of the GNU Lesser General Public
93 * License as published by the Free Software Foundation, either
94 * version 3 of the License, or (at your option) any later
95 * version.
96 *
97 * This program is distributed in the hope that it will be useful,
98 * but WITHOUT ANY WARRANTY; without even the implied warranty of
99 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
100 * GNU Lesser General Public License for more details.
101 *
102 * You should have received a copy of the GNU Lesser General
103 * License along with this program. If not, see
104 * <http://www.gnu.org/licenses/>.
105
106 * Bugs:
107 * {note_any_bugs_here}
108 *-
109 */
110
111 #include "pal.h"
112 #include "palmac.h"
113 #include "pal1sofa.h"
114
115 /*
116 *+
117 * Name:
118 * palCldj
119
120 * Purpose:
121 * Gregorian Calendar to Modified Julian Date
122
123 * Language:
124 * Starlink ANSI C
125
126 * Type of Module:
127 * Library routine
128
129 * Invocation:
130 * palCldj( int iy, int im, int id, double *djm, int *j );
131
132 * Arguments:
133 * iy = int (Given)
134 * Year in Gregorian calendar
135 * im = int (Given)
136 * Month in Gregorian calendar
137 * id = int (Given)
138 * Day in Gregorian calendar
139 * djm = double * (Returned)
140 * Modified Julian Date (JD-2400000.5) for 0 hrs
141 * j = int * (Returned)
142 * status: 0 = OK, 1 = bad year (MJD not computed),
143 * 2 = bad month (MJD not computed), 3 = bad day (MJD computed).
144
145 * Description:
146 * Gregorian calendar to Modified Julian Date.
147
148 * Notes:
149 * - Uses eraCal2jd(). See SOFA/ERFA documentation for details.
150
151 *-
152 */
153
palCldj(int iy,int im,int id,double * djm,int * j)154 void palCldj ( int iy, int im, int id, double *djm, int *j ) {
155 double djm0;
156 *j = eraCal2jd( iy, im, id, &djm0, djm );
157 }
158
159 /*
160 *+
161 * Name:
162 * palDbear
163
164 * Purpose:
165 * Bearing (position angle) of one point on a sphere relative to another
166
167 * Language:
168 * Starlink ANSI C
169
170 * Type of Module:
171 * Library routine
172
173 * Invocation:
174 * pa = palDbear( double a1, double b1, double a2, double b2 );
175
176 * Arguments:
177 * a1 = double (Given)
178 * Longitude of point 1 (e.g. RA) in radians.
179 * b1 = double (Given)
180 * Latitude of point 1 (e.g. Dec) in radians.
181 * a2 = double (Given)
182 * Longitude of point 2 in radians.
183 * b2 = double (Given)
184 * Latitude of point 2 in radians.
185
186 * Returned Value:
187 * The result is the bearing (position angle), in radians, of point
188 * A2,B2 as seen from point A1,B1. It is in the range +/- pi. If
189 * A2,B2 is due east of A1,B1 the bearing is +pi/2. Zero is returned
190 * if the two points are coincident.
191
192 * Description:
193 * Bearing (position angle) of one point in a sphere relative to another.
194
195 * Notes:
196 * - Uses eraPas(). See SOFA/ERFA documentation for details.
197
198 *-
199 */
200
palDbear(double a1,double b1,double a2,double b2)201 double palDbear ( double a1, double b1, double a2, double b2 ) {
202 return eraPas( a1, b1, a2, b2 );
203 }
204
205 /*
206 *+
207 * Name:
208 * palDaf2r
209
210 * Purpose:
211 * Convert degrees, arcminutes, arcseconds to radians
212
213 * Language:
214 * Starlink ANSI C
215
216 * Type of Module:
217 * Library routine
218
219 * Invocation:
220 * palDaf2r( int ideg, int iamin, double asec, double *rad, int *j );
221
222 * Arguments:
223 * ideg = int (Given)
224 * Degrees.
225 * iamin = int (Given)
226 * Arcminutes.
227 * iasec = double (Given)
228 * Arcseconds.
229 * rad = double * (Returned)
230 * Angle in radians.
231 * j = int * (Returned)
232 * Status: 0 = OK, 1 = "ideg" out of range 0-359,
233 * 2 = "iamin" outside of range 0-59,
234 * 2 = "asec" outside range 0-59.99999
235
236 * Description:
237 * Convert degrees, arcminutes, arcseconds to radians.
238
239 * Notes:
240 * - Uses eraAf2a(). See SOFA/ERFA documentation for details.
241
242 *-
243 */
244
245 /* Arguments differ slightly. Assumes that the sign is always positive
246 and dealt with externally. */
palDaf2r(int ideg,int iamin,double asec,double * rad,int * j)247 void palDaf2r ( int ideg, int iamin, double asec, double *rad, int *j ) {
248 *j = eraAf2a( ' ', ideg, iamin, asec, rad );
249 }
250
251 /*
252 *+
253 * Name:
254 * palDav2m
255
256 * Purpose:
257 * Form the rotation matrix corresponding to a given axial vector.
258
259 * Language:
260 * Starlink ANSI C
261
262 * Type of Module:
263 * Library routine
264
265 * Invocation:
266 * palDav2m( double axvec[3], double rmat[3][3] );
267
268 * Arguments:
269 * axvec = double [3] (Given)
270 * Axial vector (radians)
271 * rmat = double [3][3] (Returned)
272 * Rotation matrix.
273
274 * Description:
275 * A rotation matrix describes a rotation about some arbitrary axis,
276 * called the Euler axis. The "axial vector" supplied to this routine
277 * has the same direction as the Euler axis, and its magnitude is the
278 * amount of rotation in radians.
279
280 * Notes:
281 * - Uses eraRv2m(). See SOFA/ERFA documentation for details.
282
283 *-
284 */
285
palDav2m(double axvec[3],double rmat[3][3])286 void palDav2m ( double axvec[3], double rmat[3][3] ) {
287 eraRv2m( axvec, rmat );
288 }
289
290 /*
291 *+
292 * Name:
293 * palDcc2s
294
295 * Purpose:
296 * Cartesian to spherical coordinates
297
298 * Language:
299 * Starlink ANSI C
300
301 * Type of Module:
302 * Library routine
303
304 * Invocation:
305 * palDcc2s( double v[3], double *a, double *b );
306
307 * Arguments:
308 * v = double [3] (Given)
309 * x, y, z vector.
310 * a = double * (Returned)
311 * Spherical coordinate (radians)
312 * b = double * (Returned)
313 * Spherical coordinate (radians)
314
315 * Description:
316 * The spherical coordinates are longitude (+ve anticlockwise looking
317 * from the +ve latitude pole) and latitude. The Cartesian coordinates
318 * are right handed, with the x axis at zero longitude and latitude, and
319 * the z axis at the +ve latitude pole.
320
321 * Notes:
322 * - Uses eraC2s(). See SOFA/ERFA documentation for details.
323
324 *-
325 */
326
palDcc2s(double v[3],double * a,double * b)327 void palDcc2s ( double v[3], double *a, double *b ) {
328 eraC2s( v, a, b );
329 }
330
331 /*
332 *+
333 * Name:
334 * palDcs2c
335
336 * Purpose:
337 * Spherical coordinates to direction cosines
338
339 * Language:
340 * Starlink ANSI C
341
342 * Type of Module:
343 * Library routine
344
345 * Invocation:
346 * palDcs2c( double a, double b, double v[3] );
347
348 * Arguments:
349 * a = double (Given)
350 * Spherical coordinate in radians (ra, long etc).
351 * b = double (Given)
352 * Spherical coordinate in radians (dec, lat etc).
353 * v = double [3] (Returned)
354 * x, y, z vector
355
356 * Description:
357 * The spherical coordinates are longitude (+ve anticlockwise looking
358 * from the +ve latitude pole) and latitude. The Cartesian coordinates
359 * are right handed, with the x axis at zero longitude and latitude, and
360 * the z axis at the +ve latitude pole.
361
362 * Notes:
363 * - Uses eraS2c(). See SOFA/ERFA documentation for details.
364
365 *-
366 */
367
palDcs2c(double a,double b,double v[3])368 void palDcs2c ( double a, double b, double v[3] ) {
369 eraS2c( a, b, v );
370 }
371
372 /*
373 *+
374 * Name:
375 * palDd2tf
376
377 * Purpose:
378 * Convert an interval in days into hours, minutes, seconds
379
380 * Language:
381 * Starlink ANSI C
382
383 * Type of Module:
384 * Library routine
385
386 * Invocation:
387 * palDd2tf( int ndp, double days, char *sign, int ihmsf[4] );
388
389 * Arguments:
390 * ndp = int (Given)
391 * Number of decimal places of seconds
392 * days = double (Given)
393 * Interval in days
394 * sign = char * (Returned)
395 * '+' or '-' (single character, not string)
396 * ihmsf = int [4] (Returned)
397 * Hours, minutes, seconds, fraction
398
399 * Description:
400 * Convert and interval in days into hours, minutes, seconds.
401
402 * Notes:
403 * - Uses eraD2tf(). See SOFA/ERFA documentation for details.
404
405 *-
406 */
407
palDd2tf(int ndp,double days,char * sign,int ihmsf[4])408 void palDd2tf ( int ndp, double days, char *sign, int ihmsf[4] ) {
409 eraD2tf( ndp, days, sign, ihmsf );
410 }
411
412 /*
413 *+
414 * Name:
415 * palDimxv
416
417 * Purpose:
418 * Perform the 3-D backward unitary transformation
419
420 * Language:
421 * Starlink ANSI C
422
423 * Type of Module:
424 * Library routine
425
426 * Invocation:
427 * palDimxv( double dm[3][3], double va[3], double vb[3] );
428
429 * Arguments:
430 * dm = double [3][3] (Given)
431 * Matrix
432 * va = double [3] (Given)
433 * vector
434 * vb = double [3] (Returned)
435 * Result vector
436
437 * Description:
438 * Perform the 3-D backward unitary transformation.
439
440 * Notes:
441 * - Uses eraTrxp(). See SOFA/ERFA documentation for details.
442
443 *-
444 */
445
palDimxv(double dm[3][3],double va[3],double vb[3])446 void palDimxv ( double dm[3][3], double va[3], double vb[3] ) {
447 eraTrxp( dm, va, vb );
448 }
449
450 /*
451 *+
452 * Name:
453 * palDm2av
454
455 * Purpose:
456 * From a rotation matrix, determine the corresponding axial vector
457
458 * Language:
459 * Starlink ANSI C
460
461 * Type of Module:
462 * Library routine
463
464 * Invocation:
465 * palDm2av( double rmat[3][3], double axvec[3] );
466
467 * Arguments:
468 * rmat = double [3][3] (Given)
469 * Rotation matrix
470 * axvec = double [3] (Returned)
471 * Axial vector (radians)
472
473 * Description:
474 * A rotation matrix describes a rotation about some arbitrary axis,
475 * called the Euler axis. The "axial vector" returned by this routine
476 * has the same direction as the Euler axis, and its magnitude is the
477 * amount of rotation in radians. (The magnitude and direction can be
478 * separated by means of the routine palDvn.)
479
480 * Notes:
481 * - Uses eraRm2v(). See SOFA/ERFA documentation for details.
482
483 *-
484 */
485
palDm2av(double rmat[3][3],double axvec[3])486 void palDm2av ( double rmat[3][3], double axvec[3] ) {
487 eraRm2v( rmat, axvec );
488 }
489
490 /*
491 *+
492 * Name:
493 * palDjcl
494
495 * Purpose:
496 * Modified Julian Date to Gregorian year, month, day and fraction of day
497
498 * Language:
499 * Starlink ANSI C
500
501 * Type of Module:
502 * Library routine
503
504 * Invocation:
505 * palDjcl( double djm, int *iy, int *im, int *id, double *fd, int *j );
506
507 * Arguments:
508 * djm = double (Given)
509 * modified Julian Date (JD-2400000.5)
510 * iy = int * (Returned)
511 * year
512 * im = int * (Returned)
513 * month
514 * id = int * (Returned)
515 * day
516 * fd = double * (Returned)
517 * Fraction of day.
518
519 * Description:
520 * Modified Julian Date to Gregorian year, month, day and fraction of day.
521
522 * Notes:
523 * - Uses eraJd2cal(). See SOFA/ERFA documentation for details.
524
525 *-
526 */
527
528 /* Requires additional SLA MJD reference date */
palDjcl(double djm,int * iy,int * im,int * id,double * fd,int * j)529 void palDjcl ( double djm, int *iy, int *im, int *id, double *fd, int *j ) {
530 *j = eraJd2cal( PAL__MJD0, djm, iy, im, id, fd );
531 }
532
533 /*
534 *+
535 * Name:
536 * palDmxm
537
538 * Purpose:
539 * Product of two 3x3 matrices
540
541 * Language:
542 * Starlink ANSI C
543
544 * Type of Module:
545 * Library routine
546
547 * Invocation:
548 * palDmxm( double a[3][3], double b[3][3], double c[3][3] );
549
550 * Arguments:
551 * a = double [3][3] (Given)
552 * Matrix
553 * b = double [3][3] (Given)
554 * Matrix
555 * c = double [3][3] (Returned)
556 * Matrix result
557
558 * Description:
559 * Product of two 3x3 matrices.
560
561 * Notes:
562 * - Uses eraRxr(). See SOFA/ERFA documentation for details.
563
564 *-
565 */
566
palDmxm(double a[3][3],double b[3][3],double c[3][3])567 void palDmxm ( double a[3][3], double b[3][3], double c[3][3] ) {
568 eraRxr( a, b, c );
569 }
570
571 /*
572 *+
573 * Name:
574 * palDmxv
575
576 * Purpose:
577 * Performs the 3-D forward unitary transformation
578
579 * Language:
580 * Starlink ANSI C
581
582 * Type of Module:
583 * Library routine
584
585 * Invocation:
586 * palDmxv( double dm[3][3], double va[3], double vb[3] );
587
588 * Arguments:
589 * dm = double [3][3] (Given)
590 * matrix
591 * va = double [3] (Given)
592 * vector
593 * dp = double [3] (Returned)
594 * result vector
595
596 * Description:
597 * Performs the 3-D forward unitary transformation.
598
599 * Notes:
600 * - Uses eraRxp(). See SOFA/ERFA documentation for details.
601
602 *-
603 */
604
palDmxv(double dm[3][3],double va[3],double vb[3])605 void palDmxv ( double dm[3][3], double va[3], double vb[3] ) {
606 eraRxp( dm, va, vb );
607 }
608
609 /*
610 *+
611 * Name:
612 * palDpav
613
614 * Purpose:
615 * Position angle of one celestial direction with respect to another
616
617 * Language:
618 * Starlink ANSI C
619
620 * Type of Module:
621 * Library routine
622
623 * Invocation:
624 * pa = palDpav( double v1[3], double v2[3] );
625
626 * Arguments:
627 * v1 = double [3] (Given)
628 * direction cosines of one point.
629 * v2 = double [3] (Given)
630 * direction cosines of the other point.
631
632 * Returned Value:
633 * The result is the bearing (position angle), in radians, of point
634 * V2 with respect to point V1. It is in the range +/- pi. The
635 * sense is such that if V2 is a small distance east of V1, the
636 * bearing is about +pi/2. Zero is returned if the two points
637 * are coincident.
638
639 * Description:
640 * Position angle of one celestial direction with respect to another.
641
642 * Notes:
643 * - The coordinate frames correspond to RA,Dec, Long,Lat etc.
644 * - Uses eraPap(). See SOFA/ERFA documentation for details.
645
646 *-
647 */
648
palDpav(double v1[3],double v2[3])649 double palDpav ( double v1[3], double v2[3] ) {
650 return eraPap( v1, v2 );
651 }
652
653 /*
654 *+
655 * Name:
656 * palDr2af
657
658 * Purpose:
659 * Convert an angle in radians to degrees, arcminutes, arcseconds
660
661 * Language:
662 * Starlink ANSI C
663
664 * Type of Module:
665 * Library routine
666
667 * Invocation:
668 * palDr2af( int ndp, double angle, char *sign, int idmsf[4] );
669
670 * Arguments:
671 * ndp = int (Given)
672 * number of decimal places of arcseconds
673 * angle = double (Given)
674 * angle in radians
675 * sign = char * (Returned)
676 * '+' or '-' (single character)
677 * idmsf = int [4] (Returned)
678 * Degrees, arcminutes, arcseconds, fraction
679
680 * Description:
681 * Convert an angle in radians to degrees, arcminutes, arcseconds.
682
683 * Notes:
684 * - Uses eraA2af(). See SOFA/ERFA documentation for details.
685
686 *-
687 */
688
palDr2af(int ndp,double angle,char * sign,int idmsf[4])689 void palDr2af ( int ndp, double angle, char *sign, int idmsf[4] ) {
690 eraA2af( ndp, angle, sign, idmsf );
691 }
692
693 /*
694 *+
695 * Name:
696 * palDr2tf
697
698 * Purpose:
699 * Convert an angle in radians to hours, minutes, seconds
700
701 * Language:
702 * Starlink ANSI C
703
704 * Type of Module:
705 * Library routine
706
707 * Invocation:
708 * palDr2tf ( int ndp, double angle, char *sign, int ihmsf[4] );
709
710 * Arguments:
711 * ndp = int (Given)
712 * number of decimal places of arcseconds
713 * angle = double (Given)
714 * angle in radians
715 * sign = char * (Returned)
716 * '+' or '-' (single character)
717 * idmsf = int [4] (Returned)
718 * Hours, minutes, seconds, fraction
719
720 * Description:
721 * Convert an angle in radians to hours, minutes, seconds.
722
723 * Notes:
724 * - Uses eraA2tf(). See SOFA/ERFA documentation for details.
725
726 *-
727 */
728
palDr2tf(int ndp,double angle,char * sign,int ihmsf[4])729 void palDr2tf( int ndp, double angle, char *sign, int ihmsf[4] ) {
730 eraA2tf( ndp, angle, sign, ihmsf );
731 }
732
733 /*
734 *+
735 * Name:
736 * palDranrm
737
738 * Purpose:
739 * Normalize angle into range 0-2 pi
740
741 * Language:
742 * Starlink ANSI C
743
744 * Type of Module:
745 * Library routine
746
747 * Invocation:
748 * norm = palDranrm( double angle );
749
750 * Arguments:
751 * angle = double (Given)
752 * angle in radians
753
754 * Returned Value:
755 * Angle expressed in the range 0-2 pi
756
757 * Description:
758 * Normalize angle into range 0-2 pi.
759
760 * Notes:
761 * - Uses eraAnp(). See SOFA/ERFA documentation for details.
762
763 *-
764 */
765
palDranrm(double angle)766 double palDranrm ( double angle ) {
767 return eraAnp( angle );
768 }
769
770 /*
771 *+
772 * Name:
773 * palDsep
774
775 * Purpose:
776 * Angle between two points on a sphere
777
778 * Language:
779 * Starlink ANSI C
780
781 * Type of Module:
782 * Library routine
783
784 * Invocation:
785 * ang = palDsep( double a1, double b1, double a2, double b2 );
786
787 * Arguments:
788 * a1 = double (Given)
789 * Spherical coordinate of one point (radians)
790 * b1 = double (Given)
791 * Spherical coordinate of one point (radians)
792 * a2 = double (Given)
793 * Spherical coordinate of other point (radians)
794 * b2 = double (Given)
795 * Spherical coordinate of other point (radians)
796
797 * Returned Value:
798 * Angle, in radians, between the two points. Always positive.
799
800 * Description:
801 * Angle between two points on a sphere.
802
803 * Notes:
804 * - The spherical coordinates are [RA,Dec], [Long,Lat] etc, in radians.
805 * - Uses eraSeps(). See SOFA/ERFA documentation for details.
806
807 *-
808 */
809
palDsep(double a1,double b1,double a2,double b2)810 double palDsep ( double a1, double b1, double a2, double b2 ) {
811 return eraSeps( a1, b1, a2, b2 );
812 }
813
814 /*
815 *+
816 * Name:
817 * palDsepv
818
819 * Purpose:
820 * Angle between two vectors
821
822 * Language:
823 * Starlink ANSI C
824
825 * Type of Module:
826 * Library routine
827
828 * Invocation:
829 * ang = palDsepv( double v1[3], double v2[3] );
830
831 * Arguments:
832 * v1 = double [3] (Given)
833 * First vector
834 * v2 = double [3] (Given)
835 * Second vector
836
837 * Returned Value:
838 * Angle, in radians, between the two points. Always positive.
839
840 * Description:
841 * Angle between two vectors.
842
843 * Notes:
844 * - Uses eraSepp(). See SOFA/ERFA documentation for details.
845
846 *-
847 */
848
palDsepv(double v1[3],double v2[3])849 double palDsepv ( double v1[3], double v2[3] ) {
850 return eraSepp( v1, v2 );
851 }
852
853 /*
854 *+
855 * Name:
856 * palDtf2d
857
858 * Purpose:
859 * Convert hours, minutes, seconds to days
860
861 * Language:
862 * Starlink ANSI C
863
864 * Type of Module:
865 * Library routine
866
867 * Invocation:
868 * palDtf2d( int ihour, int imin, double sec, double *days, int *j );
869
870 * Arguments:
871 * ihour = int (Given)
872 * Hours
873 * imin = int (Given)
874 * Minutes
875 * sec = double (Given)
876 * Seconds
877 * days = double * (Returned)
878 * Interval in days
879 * j = int * (Returned)
880 * status: 0 = ok, 1 = ihour outside range 0-23,
881 * 2 = imin outside range 0-59, 3 = sec outside range 0-59.999...
882
883 * Description:
884 * Convert hours, minutes, seconds to days.
885
886 * Notes:
887 * - Uses eraTf2d(). See SOFA/ERFA documentation for details.
888
889 *-
890 */
891
892 /* Assumes that the sign is always positive and is dealt with externally */
palDtf2d(int ihour,int imin,double sec,double * days,int * j)893 void palDtf2d ( int ihour, int imin, double sec, double *days, int *j ) {
894 *j = eraTf2d( ' ', ihour, imin, sec, days );
895 }
896
897 /*
898 *+
899 * Name:
900 * palDtf2r
901
902 * Purpose:
903 * Convert hours, minutes, seconds to radians
904
905 * Language:
906 * Starlink ANSI C
907
908 * Type of Module:
909 * Library routine
910
911 * Invocation:
912 * palDtf2r( int ihour, int imin, double sec, double *rad, int *j );
913
914 * Arguments:
915 * ihour = int (Given)
916 * Hours
917 * imin = int (Given)
918 * Minutes
919 * sec = double (Given)
920 * Seconds
921 * days = double * (Returned)
922 * Angle in radians
923 * j = int * (Returned)
924 * status: 0 = ok, 1 = ihour outside range 0-23,
925 * 2 = imin outside range 0-59, 3 = sec outside range 0-59.999...
926
927 * Description:
928 * Convert hours, minutes, seconds to radians.
929
930 * Notes:
931 * - Uses eraTf2a(). See SOFA/ERFA documentation for details.
932
933 *-
934 */
935
936 /* Assumes that the sign is dealt with outside this routine */
palDtf2r(int ihour,int imin,double sec,double * rad,int * j)937 void palDtf2r ( int ihour, int imin, double sec, double *rad, int *j ) {
938 *j = eraTf2a( ' ', ihour, imin, sec, rad );
939 }
940
941 /*
942 *+
943 * Name:
944 * palDvdv
945
946 * Purpose:
947 * Scalar product of two 3-vectors
948
949 * Language:
950 * Starlink ANSI C
951
952 * Type of Module:
953 * Library routine
954
955 * Invocation:
956 * prod = palDvdv ( double va[3], double vb[3] );
957
958 * Arguments:
959 * va = double [3] (Given)
960 * First vector
961 * vb = double [3] (Given)
962 * Second vector
963
964 * Returned Value:
965 * Scalar product va.vb
966
967 * Description:
968 * Scalar product of two 3-vectors.
969
970 * Notes:
971 * - Uses eraPdp(). See SOFA/ERFA documentation for details.
972
973 *-
974 */
975
palDvdv(double va[3],double vb[3])976 double palDvdv ( double va[3], double vb[3] ) {
977 return eraPdp( va, vb );
978 }
979
980 /*
981 *+
982 * Name:
983 * palDvn
984
985 * Purpose:
986 * Normalizes a 3-vector also giving the modulus
987
988 * Language:
989 * Starlink ANSI C
990
991 * Type of Module:
992 * Library routine
993
994 * Invocation:
995 * palDvn( double v[3], double uv[3], double *vm );
996
997 * Arguments:
998 * v = double [3] (Given)
999 * vector
1000 * uv = double [3] (Returned)
1001 * unit vector in direction of "v"
1002 * vm = double * (Returned)
1003 * modulus of "v"
1004
1005 * Description:
1006 * Normalizes a 3-vector also giving the modulus.
1007
1008 * Notes:
1009 * - Uses eraPn(). See SOFA/ERFA documentation for details.
1010
1011 *-
1012 */
1013
1014 /* Note that the arguments are flipped */
palDvn(double v[3],double uv[3],double * vm)1015 void palDvn ( double v[3], double uv[3], double *vm ) {
1016 eraPn( v, vm, uv );
1017 }
1018
1019 /*
1020 *+
1021 * Name:
1022 * palDvxv
1023
1024 * Purpose:
1025 * Vector product of two 3-vectors
1026
1027 * Language:
1028 * Starlink ANSI C
1029
1030 * Type of Module:
1031 * Library routine
1032
1033 * Invocation:
1034 * palDvxv( double va[3], double vb[3], double vc[3] );
1035
1036 * Arguments:
1037 * va = double [3] (Given)
1038 * First vector
1039 * vb = double [3] (Given)
1040 * Second vector
1041 * vc = double [3] (Returned)
1042 * Result vector
1043
1044 * Description:
1045 * Vector product of two 3-vectors.
1046
1047 * Notes:
1048 * - Uses eraPxp(). See SOFA/ERFA documentation for details.
1049
1050 *-
1051 */
1052
palDvxv(double va[3],double vb[3],double vc[3])1053 void palDvxv ( double va[3], double vb[3], double vc[3] ) {
1054 eraPxp( va, vb, vc );
1055 }
1056
1057 /*
1058 *+
1059 * Name:
1060 * palEpb
1061
1062 * Purpose:
1063 * Conversion of modified Julian Data to Besselian Epoch
1064
1065 * Language:
1066 * Starlink ANSI C
1067
1068 * Type of Module:
1069 * Library routine
1070
1071 * Invocation:
1072 * epb = palEpb ( double date );
1073
1074 * Arguments:
1075 * date = double (Given)
1076 * Modified Julian Date (JD - 2400000.5)
1077
1078 * Returned Value:
1079 * Besselian epoch.
1080
1081 * Description:
1082 * Conversion of modified Julian Data to Besselian Epoch.
1083
1084 * Notes:
1085 * - Uses eraEpb(). See SOFA/ERFA documentation for details.
1086
1087 *-
1088 */
1089
1090 /* Requires additional SLA MJD reference date */
palEpb(double date)1091 double palEpb ( double date ) {
1092 return eraEpb( PAL__MJD0, date );
1093 }
1094
1095 /*
1096 *+
1097 * Name:
1098 * palEpb2d
1099
1100 * Purpose:
1101 * Conversion of Besselian Epoch to Modified Julian Date
1102
1103 * Language:
1104 * Starlink ANSI C
1105
1106 * Type of Module:
1107 * Library routine
1108
1109 * Invocation:
1110 * mjd = palEpb2d ( double epb );
1111
1112 * Arguments:
1113 * epb = double (Given)
1114 * Besselian Epoch
1115
1116 * Returned Value:
1117 * Modified Julian Date (JD - 2400000.5)
1118
1119 * Description:
1120 * Conversion of Besselian Epoch to Modified Julian Date.
1121
1122 * Notes:
1123 * - Uses eraEpb2jd(). See SOFA/ERFA documentation for details.
1124
1125 *-
1126 */
1127
1128
palEpb2d(double epb)1129 double palEpb2d ( double epb ) {
1130 double djm0, djm;
1131 eraEpb2jd( epb, &djm0, &djm );
1132 return djm;
1133 }
1134
1135 /*
1136 *+
1137 * Name:
1138 * palEpj
1139
1140 * Purpose:
1141 * Conversion of Modified Julian Date to Julian Epoch
1142
1143 * Language:
1144 * Starlink ANSI C
1145
1146 * Type of Module:
1147 * Library routine
1148
1149 * Invocation:
1150 * epj = palEpj ( double date );
1151
1152 * Arguments:
1153 * date = double (Given)
1154 * Modified Julian Date (JD - 2400000.5)
1155
1156 * Returned Value:
1157 * The Julian Epoch.
1158
1159 * Description:
1160 * Conversion of Modified Julian Date to Julian Epoch.
1161
1162 * Notes:
1163 * - Uses eraEpj(). See SOFA/ERFA documentation for details.
1164
1165 *-
1166 */
1167
1168 /* Requires additional SLA MJD reference date */
palEpj(double date)1169 double palEpj ( double date ) {
1170 return eraEpj( PAL__MJD0, date );
1171 }
1172
1173 /*
1174 *+
1175 * Name:
1176 * palEpj2d
1177
1178 * Purpose:
1179 * Conversion of Julian Epoch to Modified Julian Date
1180
1181 * Language:
1182 * Starlink ANSI C
1183
1184 * Type of Module:
1185 * Library routine
1186
1187 * Invocation:
1188 * mjd = palEpj2d ( double epj );
1189
1190 * Arguments:
1191 * epj = double (Given)
1192 * Julian Epoch.
1193
1194 * Returned Value:
1195 * Modified Julian Date (JD - 2400000.5)
1196
1197 * Description:
1198 * Conversion of Julian Epoch to Modified Julian Date.
1199
1200 * Notes:
1201 * - Uses eraEpj2d(). See SOFA/ERFA documentation for details.
1202
1203 *-
1204 */
palEpj2d(double epj)1205 double palEpj2d ( double epj ) {
1206 double djm0, djm;
1207 eraEpj2jd( epj, &djm0, &djm );
1208 return djm;
1209 }
1210
1211 /*
1212 *+
1213 * Name:
1214 * palEqeqx
1215
1216 * Purpose:
1217 * Equation of the equinoxes (IAU 2000/2006)
1218
1219 * Language:
1220 * Starlink ANSI C
1221
1222 * Type of Module:
1223 * Library routine
1224
1225 * Invocation:
1226 * palEqeqx( double date );
1227
1228 * Arguments:
1229 * date = double (Given)
1230 * TT as Modified Julian Date (JD-400000.5)
1231
1232 * Description:
1233 * Equation of the equinoxes (IAU 2000/2006).
1234
1235 * Notes:
1236 * - Uses eraEe06a(). See SOFA/ERFA documentation for details.
1237
1238 *-
1239 */
1240
1241 /* Requires additional SLA MJD reference date */
palEqeqx(double date)1242 double palEqeqx ( double date ) {
1243 return eraEe06a( PAL__MJD0, date );
1244 }
1245
1246 /*
1247 *+
1248 * Name:
1249 * palFk5hz
1250
1251 * Purpose:
1252 * Transform an FK5 (J2000) star position into the frame of the
1253 * Hipparcos catalogue.
1254
1255 * Language:
1256 * Starlink ANSI C
1257
1258 * Type of Module:
1259 * Library routine
1260
1261 * Invocation:
1262 * palFk5hz ( double r5, double d5, double epoch,
1263 * double *rh, double *dh );
1264
1265 * Arguments:
1266 * r5 = double (Given)
1267 * FK5 RA (radians), equinox J2000, epoch "epoch"
1268 * d5 = double (Given)
1269 * FK5 dec (radians), equinox J2000, epoch "epoch"
1270 * epoch = double (Given)
1271 * Julian epoch
1272 * rh = double * (Returned)
1273 * RA (radians)
1274 * dh = double * (Returned)
1275 * Dec (radians)
1276
1277 * Description:
1278 * Transform an FK5 (J2000) star position into the frame of the
1279 * Hipparcos catalogue.
1280
1281 * Notes:
1282 * - Assumes zero Hipparcos proper motion.
1283 * - Uses eraEpj2jd() and eraFk5hz.
1284 * See SOFA/ERFA documentation for details.
1285
1286 *-
1287 */
1288
palFk5hz(double r5,double d5,double epoch,double * rh,double * dh)1289 void palFk5hz ( double r5, double d5, double epoch,
1290 double *rh, double *dh ) {
1291 /* Need to convert epoch to Julian date first */
1292 double date1, date2;
1293 eraEpj2jd( epoch, &date1, &date2 );
1294 eraFk5hz( r5, d5, date1, date2, rh, dh );
1295 }
1296
1297 /*
1298 *+
1299 * Name:
1300 * palGmst
1301
1302 * Purpose:
1303 * Greenwich mean sidereal time (consistent with IAU 2006 precession).
1304
1305 * Language:
1306 * Starlink ANSI C
1307
1308 * Type of Module:
1309 * Library routine
1310
1311 * Invocation:
1312 * mst = palGmst ( double ut1 );
1313
1314 * Arguments:
1315 * ut1 = double (Given)
1316 * Universal time (UT1) expressed as modified Julian Date (JD-2400000.5)
1317
1318 * Returned Value:
1319 * Greenwich mean sidereal time
1320
1321 * Description:
1322 * Greenwich mean sidereal time (consistent with IAU 2006 precession).
1323
1324 * Notes:
1325 * - Uses eraGmst06(). See SOFA/ERFA documentation for details.
1326
1327 *-
1328 */
1329
1330 /* Note that SOFA/ERFA has more accurate time arguments
1331 and we use the 2006 precession model */
palGmst(double ut1)1332 double palGmst ( double ut1 ) {
1333 return eraGmst06( PAL__MJD0, ut1, PAL__MJD0, ut1 );
1334 }
1335
1336 /*
1337 *+
1338 * Name:
1339 * palGmsta
1340
1341 * Purpose:
1342 * Greenwich mean sidereal time (consistent with IAU 2006 precession).
1343
1344 * Language:
1345 * Starlink ANSI C
1346
1347 * Type of Module:
1348 * Library routine
1349
1350 * Invocation:
1351 * mst = palGmsta ( double date, double ut1 );
1352
1353 * Arguments:
1354 * date = double (Given)
1355 * UT1 date (MJD: integer part of JD-2400000.5)
1356 * ut1 = double (Given)
1357 * UT1 time (fraction of a day)
1358
1359 * Returned Value:
1360 * Greenwich mean sidereal time (in range 0 to 2 pi)
1361
1362 * Description:
1363 * Greenwich mean sidereal time (consistent with IAU 2006 precession).
1364
1365 * Notes:
1366 * - For best accuracy use eraGmst06() directly.
1367 * - Uses eraGmst06(). See SOFA/ERFA documentation for details.
1368
1369 *-
1370 */
1371
1372 /* Slightly better but still not as accurate as SOFA/ERFA */
1373
palGmsta(double date,double ut)1374 double palGmsta( double date, double ut ) {
1375 date += PAL__MJD0;
1376 return eraGmst06( date, ut, date, ut );
1377 }
1378
1379 /*
1380 *+
1381 * Name:
1382 * palHfk5z
1383
1384 * Purpose:
1385 * Hipparcos star position to FK5 J2000
1386
1387 * Language:
1388 * Starlink ANSI C
1389
1390 * Type of Module:
1391 * Library routine
1392
1393 * Invocation:
1394 * palHfk5z( double rh, double dh, double epoch,
1395 * double *r5, double *d5, double *dr5, double *dd5 );
1396
1397 * Arguments:
1398 * rh = double (Given)
1399 * Hipparcos RA (radians)
1400 * dh = double (Given)
1401 * Hipparcos Dec (radians)
1402 * epoch = double (Given)
1403 * Julian epoch (TDB)
1404 * r5 = double * (Returned)
1405 * RA (radians, FK5, equinox J2000, epoch "epoch")
1406 * d5 = double * (Returned)
1407 * Dec (radians, FK5, equinox J2000, epoch "epoch")
1408
1409 * Description:
1410 * Transform a Hipparcos star position into FK5 J2000, assuming
1411 * zero Hipparcos proper motion.
1412
1413 * Notes:
1414 * - Uses eraEpj2jd and eraHfk5z(). See SOFA/ERFA documentation for details.
1415
1416 *-
1417 */
1418
palHfk5z(double rh,double dh,double epoch,double * r5,double * d5,double * dr5,double * dd5)1419 void palHfk5z ( double rh, double dh, double epoch,
1420 double *r5, double *d5, double *dr5, double *dd5 ) {
1421 /* Need to convert epoch to Julian date first */
1422 double date1, date2;
1423 eraEpj2jd( epoch, &date1, &date2 );
1424 eraHfk5z( rh, dh, date1, date2, r5, d5, dr5, dd5 );
1425 }
1426
1427 /*
1428 *+
1429 * Name:
1430 * palRefcoq
1431
1432 * Purpose:
1433 * Determine the constants A and B in the atmospheric refraction model
1434
1435 * Language:
1436 * Starlink ANSI C
1437
1438 * Type of Module:
1439 * Library routine
1440
1441 * Invocation:
1442 * palRefcoq( double tdk, double pmb, double rh, double wl,
1443 * double *refa, double *refb );
1444
1445 * Arguments:
1446 * tdk = double (Given)
1447 * Ambient temperature at the observer (K)
1448 * pmb = double (Given)
1449 * Pressure at the observer (millibar)
1450 * rh = double (Given)
1451 * Relative humidity at the observer (range 0-1)
1452 * wl = double (Given)
1453 * Effective wavelength of the source (micrometre).
1454 * Radio refraction is chosen by specifying wl > 100 micrometres.
1455 * refa = double * (Returned)
1456 * tan Z coefficient (radian)
1457 * refb = double * (Returned)
1458 * tan**3 Z coefficient (radian)
1459
1460 * Description:
1461 * Determine the constants A and B in the atmospheric refraction
1462 * model dZ = A tan Z + B tan**3 Z. This is a fast alternative
1463 * to the palRefco routine.
1464 *
1465 * Z is the "observed" zenith distance (i.e. affected by refraction)
1466 * and dZ is what to add to Z to give the "topocentric" (i.e. in vacuo)
1467 * zenith distance.
1468
1469 * Notes:
1470 * - Uses eraRefco(). See SOFA/ERFA documentation for details.
1471 * - Note that the SOFA/ERFA routine uses different order of
1472 * of arguments and uses deg C rather than K.
1473
1474 *-
1475 */
1476
palRefcoq(double tdk,double pmb,double rh,double wl,double * refa,double * refb)1477 void palRefcoq ( double tdk, double pmb, double rh, double wl,
1478 double *refa, double *refb ) {
1479 /* Note that SLA (and therefore PAL) uses units of kelvin
1480 but SOFA/ERFA uses deg C */
1481 eraRefco( pmb, tdk - 273.15, rh, wl, refa, refb );
1482 }
1483