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