| 1 | #include "slalib.h"
|
|---|
| 2 | #include "slamac.h"
|
|---|
| 3 | void slaPdq2h ( double p, double d, double q,
|
|---|
| 4 | double *h1, int *j1, double *h2, int *j2 )
|
|---|
| 5 | /*
|
|---|
| 6 | ** - - - - - - - - -
|
|---|
| 7 | ** s l a P d q 2 h
|
|---|
| 8 | ** - - - - - - - - -
|
|---|
| 9 | **
|
|---|
| 10 | ** Hour Angle corresponding to a given parallactic angle
|
|---|
| 11 | **
|
|---|
| 12 | ** (double precision)
|
|---|
| 13 | **
|
|---|
| 14 | ** Given:
|
|---|
| 15 | ** p double latitude
|
|---|
| 16 | ** d double declination
|
|---|
| 17 | ** q double parallactic angle
|
|---|
| 18 | **
|
|---|
| 19 | ** Returned:
|
|---|
| 20 | ** *h1 double hour angle: first solution if any
|
|---|
| 21 | ** *j1 int flag: 0 = solution 1 is valid
|
|---|
| 22 | ** *h2 double hour angle: first solution if any
|
|---|
| 23 | ** *j2 int flag: 0 = solution 2 is valid
|
|---|
| 24 | **
|
|---|
| 25 | ** Called: slaDrange
|
|---|
| 26 | **
|
|---|
| 27 | ** Defined in slamac.h: DPI, DPIBY2
|
|---|
| 28 | **
|
|---|
| 29 | ** Last revision: 24 November 1994
|
|---|
| 30 | **
|
|---|
| 31 | ** Copyright P.T.Wallace. All rights reserved.
|
|---|
| 32 | */
|
|---|
| 33 |
|
|---|
| 34 | #define TINY 1e-12 /* Zone of avoidance around critical angles */
|
|---|
| 35 |
|
|---|
| 36 | {
|
|---|
| 37 | double pn, qn, dn, sq, cq, sqsd, qt, qb, hpt, t;
|
|---|
| 38 |
|
|---|
| 39 | /* Preset status flags to OK */
|
|---|
| 40 | *j1 = 0;
|
|---|
| 41 | *j2 = 0;
|
|---|
| 42 |
|
|---|
| 43 | /* Adjust latitude, azimuth, parallactic angle to avoid critical values */
|
|---|
| 44 | pn = slaDrange ( p );
|
|---|
| 45 | if ( fabs ( fabs ( pn ) - DPIBY2 ) < TINY ) {
|
|---|
| 46 | pn -= dsign ( TINY, pn);
|
|---|
| 47 | } else if ( fabs ( pn ) < TINY ) {
|
|---|
| 48 | pn = TINY;
|
|---|
| 49 | }
|
|---|
| 50 | qn = slaDrange ( q );
|
|---|
| 51 | if ( fabs ( fabs ( qn ) - DPI ) < TINY ) {
|
|---|
| 52 | qn -= dsign ( TINY, qn );
|
|---|
| 53 | } else if ( fabs ( qn ) < TINY ) {
|
|---|
| 54 | qn = TINY;
|
|---|
| 55 | }
|
|---|
| 56 | dn = slaDrange ( d );
|
|---|
| 57 | if ( fabs ( fabs ( dn ) - fabs ( p ) ) < TINY ) {
|
|---|
| 58 | dn -= dsign ( TINY, dn );
|
|---|
| 59 | } else if ( fabs ( fabs ( dn ) - DPIBY2 ) < TINY ) {
|
|---|
| 60 | dn -= dsign ( TINY, dn );
|
|---|
| 61 | } else if ( fabs ( dn ) < TINY ) {
|
|---|
| 62 | dn = TINY;
|
|---|
| 63 | }
|
|---|
| 64 |
|
|---|
| 65 | /* Useful functions */
|
|---|
| 66 | sq = sin ( qn );
|
|---|
| 67 | cq = cos ( qn );
|
|---|
| 68 | sqsd = sq * sin ( dn );
|
|---|
| 69 |
|
|---|
| 70 | /* Quotient giving sin(h+t) */
|
|---|
| 71 | qt = sin ( pn ) * sq * cos ( dn );
|
|---|
| 72 | qb = cos ( pn ) * sqrt ( cq * cq + sqsd * sqsd );
|
|---|
| 73 |
|
|---|
| 74 | /* Any solutions? */
|
|---|
| 75 | if ( fabs ( qt ) <= qb ) {
|
|---|
| 76 |
|
|---|
| 77 | /* Yes: find h+t and t */
|
|---|
| 78 | hpt = asin ( qt / qb );
|
|---|
| 79 | t = atan2 ( sqsd, cq );
|
|---|
| 80 |
|
|---|
| 81 | /* The two solutions */
|
|---|
| 82 | *h1 = slaDrange ( hpt - t );
|
|---|
| 83 | *h2 = slaDrange ( - hpt - ( t + DPI ) );
|
|---|
| 84 |
|
|---|
| 85 | /* Reject if h and Q different signs */
|
|---|
| 86 | if ( *h1 * qn < 0.0 ) *j1 = - 1;
|
|---|
| 87 | if ( *h2 * qn < 0.0 ) *j2 = - 1;
|
|---|
| 88 | } else {
|
|---|
| 89 | *j1 = - 1;
|
|---|
| 90 | *j2 = - 1;
|
|---|
| 91 | }
|
|---|
| 92 | }
|
|---|
