Refraction2

# Astronomical Refraction(2) for the HP-41

Overview

1°)  Program#1-One Polynomial of degree 7
2°)  Program#2-One Polynomial of degree 13
3°)  Program#3-Three Polynomials of degree 6
4°)  Program#4-With XM-Data File
5°)  Auer & Standish Model of the Atmosphere
6°)  Faster Programs

a) REFR
c) PH0-H & PH-H0

-The programs listed in paragraphs 1 to 4 employs the data given in https://ccmc.gsfc.nasa.gov/modelweb/models/msis_vitmo.php
-After fitting these data to functions of the type:

n - 1 = ( n0 - 1 ) exp [ p(x)]     where  n is the refractive index of the atmosphere, p(x) a polynomial and x = the altitude

the refraction is calculated by a numerical integration, as described by Auer & Standish ( cf reference [2] for the details )

R = § ( r dn/dr ) / ( n + r dn/dr ) . d(psi)

-So, you can use these programs for h0 < 0

-Gauss-Legendre 3-point formula is used to evaluate this integral.
-The interval of integration is divided in 3 parts, corresponding to the altitudes [0,11]  [11,28] & [28,87]  ( in kilometers )

-In paragraph 5, the atmospheric model given by Auer & Standish is used - with small modifications.

-The value of ( n0 - 1 ) is computed by Owens' simplified formula
-The azimuth is also taken into account to calculate the radius of curvature Rho: Fermat's formula yields:

Rho = 1 / [ ( Cos2 Az ) / r1 + ( Sin2 Az ) / r2 ]

where  r1 = a ( 1 - e2 ) ( 1 - e2 Sin2 b ) -3/2    and    r2 = a  ( 1 - e2 Sin2 b ) -1/2

(  a = semi-major axis of the Earth , e = eccentricity ,  b = latitude of the observer )

-The refractivity at an altitude 0 ( or at the altitude of the observer ) is given by Owens' formula:

( n0 - 1 ) 108 = [ 2371.34 + 683939.7 ( 130 - 1/L2 ) -1 + 4547.3 ( 38.9 - 1/L2 ) -1 ] DS
+ [ 6487.31 + 58.058 / L2 - 0.71150 / L4 + 0.08851 / L6 ] DW

where    DS = ( PS / T ) [ 1 + PS ( 57.90 10 -8 - 9.325 10 -4 / T + 0.25844 / T2 ) ]
and      DW = ( F / T ) [ 1 + F ( 1 + 0.00037 F ) ( -0.00237321 + 2.23366 / T - 710.792 / T2 + 77514.1 / T3 ) ]

with   PS = pressure of dry air in mbar                     P = PS + F                              L = wavelength in µm
F  = pressure of water vapor in mbar             T = temperature in Kelvin

-In paragraph 6, simplified formulas are used to get faster - but less accurate - results.

1°)  Program#1-One Polynomial of degree 7

-Here, we assume that p(x) is a polynomial of degree 7 for x between 0 and 87 km

p(x) = c1 x + c2 x2 + ............... + c7 x7

-So,  n - 1 = ( n0 - 1 ) exp [ c1 x + c2 x2 + ............... + c7 x7 ]     ( approximately )   with:

c1 = t / 1250 - 0.109142               c4 = 8.89464 E-6               c7 = -3.63388 E-12
c2 = -1/97162 - t x 9 E-6              c5 = -1.53611 E-7
c3 = -204894 E-9                         c6 = 1.21088 E-9

-You have to initialize R01 thru R06:

t = temperature
P = total pressure
F = pressure of water-vapor

Data Registers:              R00 = dn/dx              ( Registers R01 thru R06 are to be initialized before executing "REF7" )

•  R01 = t  ( °C )           •  R04 = wave-length ( µm )
•  R02 = P ( mbar )       •  R05 = Latitude ( ° ' " )
•  R03 = F ( mbar )       •  R06 = Altitude of the observer ( in meters ) < 11000

R07 = c1  ....................  R13 = c7      R14 to R40: temp

Flag:  F01                      CF 01 = t , P , F  are measured for an altitude = 0
SF 01 = t , P , F  are measured at the observer's station
Subroutines: /

-Lines 08 & 218 are three-byte GTOs

 01 LBL "REF7"  02 DEG  03 HR  04 STO 22  05 COS  06 STO 20   07 X=0?  08 GTO 12  09 STO 32  10 X<>Y  11 HR  12 1  13 P-R  14 X^2  15 0.99330562      16 /  17 X<>Y  18 X^2  19 RCL 05   20 HR  21 SIN  22 X^2  23 669438 E-8  24 *  25 1  26 X<>Y  27 -  28 SQRT  29 ST* Z  30 ST* Z  31 ST* Z  32 *  33 +  34 6378.137  35 X<>Y  36 /  37 STO 23  38 RCL 06  39  E3  40 /  41 STO 20  42 +  43 STO 36  44 RCL 01  45 273.15  46 +  47 1/X  48 ENTER  49 ENTER  50 STO 00  51 77514.1  52 * 53 710.792  54 -  55 *  56 2.23366  57 +  58 *  59 237321 E-8  60 -  61 RCL 03   62 ST* Y  63 37 E-5  64 *  65 1  66 +  67 *  68 1  69 +  70 RCL 03            71 *  72 *  73 STO 07  74 RCL 04  75 X^2  76 1/X  77 ENTER  78 ENTER  79 STO 08  80 .08851  81 *  82 .7115  83 -  84 *  85 58.058  86 +  87 *  88 6487.31  89 +  90 ST* 07  91 RCL 00  92 .25844  93 *  94 9325 E-7  95 -  96 RCL 00  97 *  98 579 E-9  99 + 100 RCL 02 101 RCL 03 102 - 103 ST* 00 104 * 105 1 106 + 107 RCL 00 108 * 109 4547.3 110 38.9 111 RCL 08 112 - 113 / 114 130 115 RCL 08          116 - 117 683939.7 118 X<>Y 119 / 120 + 121 2371.34 122 + 123 * 124 RCL 07  125 + 126  E8 127 / 128 STO 39 129 1.6 130 STO 31 131 FRC 132 SQRT 133 STO 30 134 RCL 01 135 1250 136 / 137 .109142 138 - 139 STO 07 140 97162 141 1/X 142 CHS 143 RCL 01 144 9 E-6 145 * 146 - 147 STO 08 148 ST+ X 149 STO 14 150 204894 E-9 151 CHS 152 STO 09 153 3 154 * 155 STO 15 156 8.89464 E-6 157 STO 10 158 4 159 * 160 STO 16 161 1.53611 E-7 162 CHS 163 STO 11 164 5 165 * 166 STO 17 167 1.21088 E-9 168 STO 12 169 6 170 * 171 STO 18  172 3.63388 E-12 173 CHS 174 STO 13  175 7 176 * 177 STO 19 178 CLX 179 X<> 20 180 XEQ 04 181 LASTX 182 STO 00 183 FS? 01 184 RCL 39 185 RCL 36 186 ST* Y 187 + 188 RCL 32 189 * 190 STO 35 191 FC? 01 192 GTO 00 193 RCL 00 194 RCL 39 195 / 196 ST/ 39 197 LBL 00 198 3 E-6 199 STO 38 200 11 201 XEQ 02 202 STO 21 203 6 204 STO 40 205 XEQ 13 206 28 207 XEQ 02 208 X<> 21 209 STO 22 210 4 211 STO 40 212 XEQ 13 213 87 214 XEQ 02 215 X<> 21 216 STO 22  217 XEQ 13 218 GTO 12 219 LBL 13 220 RCL 22         221 RCL 21 222 STO 24 223 - 224 RCL 40 225 STO 28 226 / 227 STO 27  228 2 229 / 230 ST+ 24 231 RCL 30 232 * 233 STO 25 234 CLX 235 STO 26 236 LBL 14 237 RCL 24 238 RCL 25 239 - 240 XEQ 01 241 ST+ 26 242 RCL 24 243 XEQ 01 244 RCL 31 245 * 246 ST+ 26 247 RCL 24 248 RCL 25 249 + 250 XEQ 01 251 ST+ 26 252 RCL 27 253 ST+ 24 254 DSE 28 255 GTO 14 256 RCL 26 257 * 258 ST+ 20 259 RTN 260 LBL 02 261 STO 34 262 XEQ 04 263 X<>Y 264 RCL 23 265 + 266 * 267 RCL 35 268 X<>Y 269 / 270 ACOS 271 RTN 272 LBL 01 273 COS 274 STO 29          275 LBL 03 276 RCL 34  277 XEQ 04 278 STO 37 279 X<>Y 280 RCL 23 281 + 282 STO 33 283 * 284 RCL 35 285 RCL 29 286 / 287 - 288 RCL 00 289 RCL 33 290 * 291 RCL 37 292 + 293 / 294 ST- 34 295 ABS 296 RCL 38 297 X

( 689 bytes / SIZE 041 )

 STACK INPUTS OUTPUTS Y Az ( ° ' " ) R ( " ) X h0 ( ° ' " ) R ( ° ' " )

Where h0 = apparent height

Example:       t = 10°C  P = 1010 mbar  F = 6 mbar  Lambda = 0.577 µm  Latitude = 33°21'22"  Altitude = 1706 m

10    STO 01        0.577    STO 04
1010  STO 02      33.2122  STO 05
6     STO 03        1706     STO 06

1°)    CF 01    Az = 12°41'   h0 = 1°23'45"

12.41    ENTER^
1.2345   XEQ "REF7"  >>>>  R = 0°18'08"747                    ---Execution time = 7m24s---
RDN   R = 1088"747

2°)    SF 01    Az = 12°41'   h0 = 1°23'45"

12.41    ENTER^
1.2345       R/S     >>>>  R = 0°21'45"293
RDN   R = 1305"293

Notes:

-With  Az = 84° ,  we get  1090"366  &  1307"268  in the examples above.

-You can of course use your own constants:
-Change lines 134 to 174.

2°)  Program#2-One Polynomial of degree 13

-Here, we assume that p(x) is a polynomial of degree 13 for x between 0 and 87 km

p(x) = c1 x + c2 x2 + ............... + c13 x13

-So,  n - 1 = ( n0 - 1 ) exp [ c1 x + c2 x2 + ............... + c13 x13 ]     ( approximately )   with:

c1 = t / 1250 - 0.109671                  c4 = -1.553002 E-4           c7 = 1.012373 E-8            c10 = 2529916 E-20            c13 = -41167039 E-28
c2 = -0.0026952 - t x 95 E-7           c5 = 1.137826 E-5            c8 = -1054348 E-16           c11 = -31539538 E-23
c3 =  958131 E-9                             c6 = -4.532222 E-7           c9 = -3737867 E-19           c12 = 1805402 E-24

( lines 136 to 176 )

Data Registers:              R00 = dn/dx              ( Registers R01 thru R06 are to be initialized before executing "REF13" )

•  R01 = t  ( °C )           •  R04 = wave-length ( µm )
•  R02 = P ( mbar )       •  R05 = Latitude ( ° ' " )
•  R03 = F ( mbar )       •  R06 = Altitude of the observer ( in meters ) < 11000

R07 = c1  ....................  R19 = c13      R20 to R52: temp

Flags:  F00 & F01         CF 00 = The constants in R07 to R19 are those listed below
SF 00 = Store your own constants in R07 to R19 before executing "REF13"

CF 01 = t , P , F  are measured for an altitude = 0
SF 01 = t , P , F  are measured at the observer's station
Subroutines: /

-Lines 08-135-233 are three-byte GTOs

 01 LBL "REF13"  02 DEG  03 HR  04 STO 33   05 COS  06 STO 51  07 X=0?  08 GTO 12  09 STO 43   10 X<>Y  11 HR  12 1  13 P-R  14 X^2  15 .99330562  16 /  17 X<>Y  18 X^2  19 RCL 05  20 HR  21 SIN  22 X^2  23 .00669438  24 *  25 1  26 X<>Y  27 -  28 SQRT  29 ST* Z  30 ST* Z  31 ST* Z  32 *  33 +  34 6378.137  35 X<>Y  36 /  37 STO 34  38 RCL 06  39  E3  40 /  41 STO 51  42 +  43 STO 47  44 RCL 01  45 273.15  46 +  47 1/X  48 ENTER  49 ENTER  50 STO 00              51 77514.1  52 *  53 710.792  54 -  55 *  56 2.23366  57 +  58 *  59 .00237321  60 - 61 RCL 03   62 ST* Y  63 37 E-5  64 *  65 1  66 +  67 *  68 1  69 +  70 RCL 03  71 *  72 *  73 STO 21  74 RCL 04  75 X^2  76 1/X  77 ENTER  78 ENTER  79 STO 22  80 .08851  81 *  82 .7115  83 -  84 *  85 58.058  86 +  87 *  88 6487.31  89 +  90 ST* 21  91 RCL 00  92 .25844  93 *  94 9325 E-7  95 -  96 RCL 00  97 *  98 579 E-9  99 + 100 RCL 02 101 RCL 03 102 - 103 ST* 00 104 * 105 1 106 + 107 RCL 00 108 * 109 4547.3 110 38.9 111 RCL 22            112 - 113 / 114 130 115 RCL 22 116 - 117 683939.7 118 X<>Y 119 / 120 + 121 2371.34 122 + 123 * 124 RCL 21  125 + 126  E8 127 / 128 STO 50 129 1.6 130 STO 42  131 FRC 132 SQRT 133 STO 41 134 FS? 00 135 GTO 00 136 RCL 01 137 1250 138 / 139 .109671 140 - 141 STO 07 142 .0026952 143 CHS 144 RCL 01 145 95 E-7 146 * 147 - 148 STO 08 149 958131 E-9 150 STO 09 151 1.553002 E-4 152 CHS 153 STO 10 154 1.137826 E-5 155 STO 11 156 4.532222 E-7 157 CHS 158 STO 12 159 1.012373 E-8 160 STO 13 161 1054348 E-16 162 CHS 163 STO 14 164 3737867 E-19 165 CHS 166 STO 15 167 2529916 E-20 168 STO 16 169 31539538 E-23 170 CHS 171 STO 17 172 1805402 E-24 173 STO 18 174 41167039 E-28 175 CHS 176 STO 19 177 LBL 00 178 31.019 179 STO 48 180 19 181 13 182 LBL 07 183 RCL IND Y 184 RCL Y 185 * 186 STO IND 48 187 CLX 188 SIGN 189 ST- Z 190 - 191 DSE 48 192 GTO 07 193 CLX 194 X<> 51 195 XEQ 04 196 LASTX 197 STO 00     198 FS? 01 199 RCL 50 200 RCL 47  201 ST* Y 202 + 203 RCL 43 204 * 205 STO 46 206 FC? 01 207 GTO 00 208 RCL 00 209 RCL 50 210 / 211 ST/ 50 212 LBL 00 213 3 E-6 214 STO 49 215 11 216 XEQ 02 217 STO 32 218 6 219 STO 52 220 XEQ 13 221 28 222 XEQ 02 223 X<> 32 224 STO 33 225 4 226 STO 52 227 XEQ 13 228 87 229 XEQ 02 230 X<> 32 231 STO 33 232 XEQ 13 233 GTO 12           234 LBL 13 235 RCL 33 236 RCL 32 237 STO 35 238 - 239 RCL 52 240 STO 39 241 / 242 STO 38 243 2 244 / 245 ST+ 35 246 RCL 41 247 * 248 STO 36 249 CLX 250 STO 37 251 LBL 14 252 RCL 35 253 RCL 36 254 - 255 XEQ 01 256 ST+ 37 257 RCL 35 258 XEQ 01 259 RCL 42  260 * 261 ST+ 37 262 RCL 35  263 RCL 36 264 + 265 XEQ 01 266 ST+ 37 267 RCL 38 268 ST+ 35 269 DSE 39 270 GTO 14 271 RCL 37 272 * 273 ST+ 51 274 RTN 275 LBL 02 276 STO 45 277 XEQ 04 278 X<>Y 279 RCL 34 280 + 281 * 282 RCL 46 283 X<>Y 284 / 285 ACOS 286 RTN 287 LBL 01 288 COS 289 STO 40 290 LBL 03 291 RCL 45 292 XEQ 04 293 STO 48           294 X<>Y 295 RCL 34 296 + 297 STO 44 298 * 299 RCL 46 300 RCL 40 301 / 302 - 303 RCL 00 304 RCL 44 305 * 306 RCL 48 307 + 308 / 309 ST- 45 310 ABS 311 RCL 49 312 X

( 843 bytes / SIZE 053 )

 STACK INPUTS OUTPUTS Y Az ( ° ' " ) R ( " ) X h0 ( ° ' " ) R ( ° ' " )

Example:       t = 10°C  P = 1010 mbar  F = 6 mbar  Lambda = 0.577 µm  Latitude = 33°21'22"  Altitude = 1706 m

10    STO 01        0.577    STO 04
1010  STO 02      33.2122  STO 05
6     STO 03        1706     STO 06

1°)     CF 00  CF 01     Az = 12°41'   h0 = 1°23'45"

12.41    ENTER^
1.2345   XEQ "REF13"  >>>>  R = 0°17'56"321                    ---Execution time = 10m39s---
RDN   R = 1076"321

2°)     CF 00  SF 01     Az = 12°41'   h0 = 1°23'45"

12.41    ENTER^
1.2345      R/S       >>>>  R = 0°21'35"834
RDN   R = 1295"834

Note:

-With  Az = 84° ,  we get  1077"921  &  1297"595  in the examples above.

3°)  Program#3-Three Polynomials of degree 6

-Here, we assume that the logarithm of the densities d(x) is a polynomial of degree 6 for 0 < x < 11 km
-Another polynomial of degree 6 for 11 < x < 28 km
-And a 3rd polynomial of degree 6 for 28 < x < 87 km

d1(x) = exp [ c0 + c1 x + c2 x2 + ............... + c6 x6 ]     ( 0 < x < 11 )
d2(x) = exp [ c'0 + c'1 x + c'2 x2 + ............... + c'6 x6 ]     ( 11 < x < 28 )
d3(x) = exp [ c"0 + c"1 x + c"2 x2 + ............... + c"6 x6 ]     ( 28 < x < 87 )     with:

c0 = -6.704085          c3 = -5.19398 E-4          c6 = 2.811385 E-8            ( lines 136 to 153 )
c1 = -0.111511          c4 = 2.309197 E-5
c2 = 3.835206 E-3     c5 = -9.619965 E-7

c'0 = -6.05731            c'3 = -0.00341631           c'6 = 2.638265 E-8            ( lines 198 to 215 )
c'1 = -0.3472742        c'4 = 1.461133 E-4
c'2 = 0.03978828       c'5 = -3.109464 E-6

c"0 = -27.374327        c"3 = 0.003395521           c"6 = -9.3999402 E-10        ( lines 231 to 248 )
c"1 = 2.5006517          c"4 = -4.655209 E-5
c"2 = -0.1331043        c"5 = 3.284335 E-7

-You can use the constants listed below ( in this case, clear flag F00 )
-Or your own coefficients ( store these 21 numbers in a DATA file and set flag F00 )

Data Registers:              R00 = dn/dx              ( Registers R01 thru R06 are to be initialized before executing "REF6+" )

•  R01 = t  ( °C )           •  R04 = wave-length ( µm )
•  R02 = P ( mbar )       •  R05 = Latitude ( ° ' " )
•  R03 = F ( mbar )       •  R06 = Altitude of the observer ( in meters ) < 11000

R07 = c0  ....................  R13 = c6      R14 to R40: temp

Flags:  F00 & F01         CF 00 = The constants in R07 to R13 are those listed below
SF 00 = Store the 3x7 coefficients in a DATA file before executing "REF6+"

CF 01 = t , P , F  are measured for an altitude = 0
SF 01 = t , P , F  are measured at the observer's station
Subroutines: /

-Lines 08 & 260 are three-byte GTOs

 01 LBL "REF6+"  02 DEG  03 HR  04 STO 22  05 COS  06 STO 19  07 X=0?  08 GTO 12  09 STO 32   10 X<>Y  11 HR  12 1  13 P-R  14 X^2  15 0.99330562  16 /  17 X<>Y  18 X^2  19 RCL 05   20 HR  21 SIN  22 X^2  23 669438 E-8         24 *  25 1  26 X<>Y  27 -  28 SQRT  29 ST* Z  30 ST* Z  31 ST* Z  32 *  33 +  34 6378.137  35 X<>Y  36 /  37 STO 23  38 RCL 06  39  E3  40 /  41 STO 19  42 +  43 STO 36  44 RCL 01  45 273.15  46 +  47 1/X  48 ENTER  49 ENTER  50 STO 00  51 77514.1  52 *  53 710.792  54 -  55 *  56 2.23366  57 +  58 *  59 237321 E-8  60 - 61 RCL 03  62 ST* Y  63 37 E-5  64 *  65 1  66 +  67 *  68 1  69 +  70 RCL 03   71 *  72 *  73 STO 07   74 RCL 04  75 X^2  76 1/X  77 ENTER  78 ENTER  79 STO 08  80 .08851  81 *  82 .7115  83 -  84 *  85 58.058  86 +  87 *  88 6487.31  89 +  90 ST* 07  91 RCL 00               92 .25844  93 *  94 9325 E-7  95 -  96 RCL 00  97 *  98 579 E-9  99 + 100 RCL 02 101 RCL 03 102 - 103 ST* 00 104 * 105 1 106 + 107 RCL 00 108 * 109 4547.3 110 38.9 111 RCL 08 112 - 113 / 114 130 115 RCL 08 116 - 117 683939.7 118 X<>Y 119 / 120 + 121 2371.34 122 + 123 * 124 RCL 07 125 + 126  E8 127 / 128 STO 39 129 1.6 130 STO 31  131 FRC 132 SQRT 133 STO 30  134 FS? 00 135 GTO 00 136 6.704085 137 CHS 138 STO 07 139 .111511 140 CHS 141 STO 08 142 3.835206 E-3 143 STO 09 144 5.19398 E-4 145 CHS 146 STO 10 147 2.309197 E-5 148 STO 11 149 9.619965 E-7 150 CHS 151 STO 12 152 2.811385 E-8    153 STO 13 154 GTO 06 155 LBL 00 156 CLX 157 SEEKPT 158 7.013 159 STO 42 160 GETRX 161 LBL 06 162 RCL 39 163 RCL 07 164 E^X 165 / 166 STO 40 167 XEQ 07 168 CLX 169 X<> 19 170 XEQ 04 171 LASTX 172 STO 00 173 FS? 01 174 RCL 39 175 RCL 36 176 ST* Y 177 + 178 RCL 32 179 * 180 STO 35 181 FC? 01 182 GTO 00 183 RCL 00 184 RCL 39 185 / 186 ST/ 40 187 LBL 00 188 3 E-6 189 STO 38  190 11 191 XEQ 02 192 STO 21  193 6 194 STO 20 195 XEQ 13 196 FS? 00 197 GTO 00 198 6.05731 199 CHS 200 STO 07 201 .3472742 202 CHS 203 STO 08 204 .03978828 205 STO 09 206 .00341631 207 CHS 208 STO 10 209 1.461133 E-4 210 STO 11 211 3.109464 E-6    212 CHS 213 STO 12 214 2.638265 E-8 215 STO 13 216 GTO 06 217 LBL 00 218 RCL 42 219 GETRX 220 LBL 06 221 XEQ 07 222 28 223 XEQ 02 224 X<> 21 225 STO 22 226 4 227 STO 20 228 XEQ 13 229 FS? 00 230 GTO 00 231 27.374327 232 CHS 233 STO 07 234 2.5006517 235 STO 08 236 .1331043 237 CHS 238 STO 09 239 .003395521 240 STO 10 241 4.655209 E-5 242 CHS 243 STO 11 244 3.284335 E-7 245 STO 12 246 9.3999402 E-10 247 CHS 248 STO 13 249 GTO 06 250 LBL 00 251 RCL 42  252 GETRX 253 LBL 06 254 XEQ 07 255 87 256 XEQ 02 257 X<> 21 258 STO 22 259 XEQ 13 260 GTO 12 261 LBL 07 262 18.013 263 STO 41 264 13 265 6 266 LBL 08 267 RCL IND Y 268 RCL Y 269 * 270 STO IND 41     271 CLX 272 SIGN 273 ST- Z 274 - 275 DSE 41 276 GTO 08 277 RTN 278 LBL 13 279 RCL 22 280 RCL 21 281 STO 24 282 - 283 RCL 20 284 STO 28 285 / 286 STO 27 287 2 288 / 289 ST+ 24 290 RCL 30 291 * 292 STO 25 293 CLX 294 STO 26 295 LBL 14 296 RCL 24 297 RCL 25 298 - 299 XEQ 01 300 ST+ 26 301 RCL 24 302 XEQ 01 303 RCL 31 304 * 305 ST+ 26 306 RCL 24 307 RCL 25  308 + 309 XEQ 01 310 ST+ 26 311 RCL 27 312 ST+ 24 313 DSE 28 314 GTO 14 315 RCL 26  316 * 317 ST+ 19 318 RTN 319 LBL 02 320 STO 34 321 XEQ 04 322 X<>Y 323 RCL 23            324 + 325 * 326 RCL 35 327 X<>Y 328 / 329 ACOS 330 RTN 331 LBL 01 332 COS 333 STO 29 334 LBL 03 335 RCL 34 336 XEQ 04 337 STO 37 338 X<>Y 339 RCL 23 340 + 341 STO 33 342 * 343 RCL 35 344 RCL 29 345 / 346 - 347 RCL 00 348 RCL 33 349 * 350 RCL 37 351 + 352 / 353 ST- 34 354 ABS 355 RCL 38 356 X

( 905 bytes / SIZE 043 )

 STACK INPUTS OUTPUTS Y Az ( ° ' " ) R ( " ) X h0 ( ° ' " ) R ( ° ' " )

Example:       t = 10°C  P = 1010 mbar  F = 6 mbar  Lambda = 0.577 µm  Latitude = 33°21'22"  Altitude = 1706 m

10    STO 01        0.577    STO 04
1010  STO 02      33.2122  STO 05
6     STO 03        1706     STO 06

1°)     CF 01    Az = 12°41'   h0 = 1°23'45"

12.41    ENTER^
1.2345   XEQ "REF6+"  >>>>  R = 0°17'40"473                    ---Execution time = 7m12s---
RDN   R = 1060"473

2°)     SF 01    Az = 12°41'   h0 = 1°23'45"

12.41    ENTER^
1.2345      R/S     >>>>  R = 0°21'21"145
RDN   R = 1281"145

Notes:

-With  Az = 84° ,  we get  1062"060  &  1283"096  in the examples above.
-Register R42 is unused if CF 00.

4°)  Program#4-With XM-Data File

-Instead of fitting Ln d to a polynomial, the following program uses directly the densities at altitude 0km, 1km, ..... , 89km
-Lagrange interpolation & differentiation formulas are employed to evaluate n & dn/dr  ( 5 point-formulas, cf reference [6] )

-One of the numerous atmospheric models given in reference [4] gives the following values for the density as a function of the altitude ( in kilometers ):

0.0 1.226E-03
1.0 1.101E-03
2.0 9.932E-04
3.0 8.987E-04
4.0 8.136E-04
5.0 7.354E-04
6.0 6.627E-04
7.0 5.944E-04
8.0 5.303E-04
9.0 4.701E-04
10.0 4.138E-04
11.0 3.618E-04
12.0 3.144E-04
13.0 2.718E-04
14.0 2.340E-04
15.0 2.008E-04
16.0 1.718E-04
17.0 1.467E-04
18.0 1.251E-04
19.0 1.065E-04
20.0 9.071E-05
21.0 7.723E-05
22.0 6.576E-05
23.0 5.601E-05
24.0 4.771E-05
25.0 4.066E-05
26.0 3.466E-05
27.0 2.957E-05
28.0 2.523E-05
29.0 2.155E-05
30.0 1.842E-05
31.0 1.576E-05
32.0 1.349E-05
33.0 1.156E-05
34.0 9.909E-06
35.0 8.503E-06
36.0 7.306E-06
37.0 6.288E-06
38.0 5.423E-06
39.0 4.687E-06
40.0 4.062E-06
41.0 3.530E-06
42.0 3.077E-06
43.0 2.690E-06
44.0 2.360E-06
45.0 2.078E-06
46.0 1.836E-06
47.0 1.626E-06
48.0 1.444E-06
49.0 1.284E-06
50.0 1.143E-06
51.0 1.018E-06
52.0 9.066E-07
53.0 8.070E-07
54.0 7.177E-07
55.0 6.375E-07
56.0 5.654E-07
57.0 5.007E-07
58.0 4.427E-07
59.0 3.909E-07
60.0 3.447E-07
61.0 3.035E-07
62.0 2.669E-07
63.0 2.344E-07
64.0 2.056E-07
65.0 1.800E-07
66.0 1.575E-07
67.0 1.375E-07
68.0 1.200E-07
69.0 1.045E-07
70.0 9.088E-08
71.0 7.892E-08
72.0 6.843E-08
73.0 5.914E-08
74.0 5.119E-08
75.0 4.426E-08
76.0 3.821E-08
77.0 3.291E-08
78.0 2.827E-08
79.0 2.420E-08
80.0 2.065E-08
81.0 1.754E-08
82.0 1.485E-08
83.0 1.251E-08
84.0 1.050E-08
85.0 8.779E-09
86.0 7.308E-09
87.0 6.061E-09
88.0 5.009E-09
89.0 4.125E-09

-Store the logarithms of these densities in a DATA file of 90 numbers before executing "REFXM"

>>> This file must be the current file.

Data Registers:              R00 = dn/dx              ( Registers R01 thru R06 are to be initialized before executing "REFXM" )

•  R01 = t  ( °C )           •  R04 = wave-length ( µm )
•  R02 = P ( mbar )       •  R05 = Latitude ( ° ' " )
•  R03 = F ( mbar )       •  R06 = Altitude of the observer ( in meters ) < 11000

R07 to R46: temp

Flag:  F01                      CF 01 = t , P , F  are measured for an altitude = 0
SF 01 = t , P , F  are measured at the observer's station
Subroutines: /

-Lines 08 & 199 are three-byte GTOs

 01 LBL "REFXM"  02 DEG  03 HR  04 STO 17  05 COS  06 STO 36   07 X=0?  08 GTO 12  09 STO 27  10 X<>Y  11 HR  12 1  13 P-R  14 X^2  15 .99330562  16 /  17 X<>Y  18 X^2  19 RCL 05  20 HR  21 SIN  22 X^2  23 .00669438  24 *  25 1  26 X<>Y  27 -  28 SQRT  29 ST* Z  30 ST* Z  31 ST* Z  32 *  33 +  34 6378.137  35 X<>Y  36 /  37 STO 18             38 RCL 06  39  E3  40 /  41 STO 36  42 +  43 STO 31  44 RCL 01  45 273.15  46 +  47 1/X  48 ENTER  49 ENTER  50 STO 00  51 77514.1  52 *  53 710.792  54 -  55 *  56 2.23366  57 +  58 *  59 .00237321  60 -  61 RCL 03  62 ST* Y  63 37 E-5  64 *  65 1 66 +  67 *  68 1  69 +  70 RCL 03   71 *  72 *  73 STO 07  74 RCL 04  75 X^2  76 1/X  77 ENTER  78 ENTER  79 STO 08  80 .08851  81 *  82 .7115  83 -  84 *  85 58.058  86 +  87 *  88 6487.31  89 +  90 ST* 07  91 RCL 00  92 .25844  93 *  94 9325 E-7  95 -  96 RCL 00  97 *  98 579 E-9  99 + 100 RCL 02           101 RCL 03  102 - 103 ST* 00 104 * 105 1 106 + 107 RCL 00 108 * 109 4547.3 110 38.9 111 RCL 08 112 - 113 / 114 130 115 RCL 08 116 - 117 683939.7 118 X<>Y 119 / 120 + 121 2371.34 122 + 123 * 124 RCL 07 125 + 126  E8 127 / 128 STO 34 129 0 130 SEEKPT 131 CLX 132 GETX 133 E^X 134 / 135 STO 46  136 1.6 137 STO 26 138 FRC 139 SQRT 140 STO 25 141 SIGN 142 STO 08 143 2 144 STO 09 145 + 146 STO 10 147 4 148 STO 11 149 * 150 STO 45 151 5 152 STO 12 153 6 154 STO 13 155 8 156 STO 14 157 37.041 158 STO 07 159 CLX 160 X<> 36 161 XEQ 04 162 LASTX 163 STO 00 164 FS? 01 165 RCL 34           166 RCL 31  167 ST* Y 168 + 169 RCL 27 170 * 171 STO 30 172 FC? 01 173 GTO 00 174 RCL 00 175 RCL 34 176 / 177 ST/ 46 178 LBL 00 179 3 E-6 180 STO 33 181 11 182 XEQ 02 183 STO 16 184 6 185 STO 35 186 XEQ 13 187 28 188 XEQ 02 189 X<> 16 190 STO 17 191 4 192 STO 35 193 XEQ 13 194 87 195 XEQ 02 196 X<> 16 197 STO 17 198 XEQ 13 199 GTO 12 200 LBL 13 201 RCL 17  202 RCL 16 203 STO 19 204 - 205 RCL 35 206 STO 23 207 / 208 STO 22 209 2 210 / 211 ST+ 19 212 RCL 25 213 * 214 STO 20 215 CLX 216 STO 21 217 LBL 14 218 RCL 19 219 RCL 20 220 - 221 XEQ 01 222 ST+ 21 223 RCL 19 224 XEQ 01 225 RCL 26 226 * 227 ST+ 21 228 RCL 19 229 RCL 20 230 + 231 XEQ 01 232 ST+ 21 233 RCL 22           234 ST+ 19 235 DSE 23 236 GTO 14 237 RCL 21 238 * 239 ST+ 36 240 RTN 241 LBL 02 242 STO 29 243 XEQ 04 244 RCL 15 245 RCL 18 246 + 247 * 248 RCL 30 249 X<>Y 250 / 251 ACOS 252 RTN 253 LBL 01 254 COS 255 STO 24 256 LBL 03 257 RCL 29 258 XEQ 04 259 STO 32 260 RCL 15 261 RCL 18 262 + 263 STO 28 264 * 265 RCL 30  266 RCL 24 267 / 268 - 269 RCL 00 270 RCL 28 271 * 272 RCL 32 273 + 274 / 275 ST- 29 276 ABS 277 RCL 33 278 XY 379 RCL 09 380 - 381 * 382 * 383 RCL 37 384 * 385 STO 42 386 RDN 387 X^2 388 RCL 11 389 - 390 STO 44 391 * 392 X<>Y 393 RCL 08  394 - 395 * 396 RCL 11 397 * 398 RCL 38 399 * 400 ST- 42 401 CLX 402 RCL 43 403 RCL 44 404 * 405 RCL 13 406 * 407 RCL 39 408 * 409 ST+ 42 410 CLX 411 RCL 44 412 * 413 X<>Y 414 RCL 08 415 + 416 * 417 RCL 11 418 * 419 RCL 40 420 * 421 ST- 42 422 CLX 423 RCL 09  424 + 425 * 426 RCL 43 427 * 428 RCL 41           429 * 430 RCL 42 431 + 432 RCL 45 433 ST/ 00 434 ST+ X 435 / 436 E^X 437 RCL 46 438 * 439 ST* 00 440 SIGN 441 LASTX 442 + 443 RTN 444 LBL 12 445 RCL 36 446  E3 447 * 448 STO Y 449 3600 450 / 451 HMS 452 END

( 746 bytes / SIZE 047 )

 STACK INPUTS OUTPUTS Y Az ( ° ' " ) R ( " ) X h0 ( ° ' " ) R ( ° ' " )

Example:       t = 10°C  P = 1010 mbar  F = 6 mbar  Lambda = 0.577 µm  Latitude = 33°21'22"  Altitude = 1706 m

10    STO 01        0.577    STO 04
1010  STO 02      33.2122  STO 05
6     STO 03        1706     STO 06

1°)    CF 01    Az = 12°41'   h0 = 1°23'45"

12.41    ENTER^
1.2345   XEQ "REFXM"  >>>>  R = 0°17'40"358                    ---Execution time = 15m21s---
RDN   R = 1060"358

2°)    SF 01    Az = 12°41'   h0 = 1°23'45"

12.41    ENTER^
1.2345       R/S       >>>>    R = 0°21'19"933
RDN    R = 1279"933

Notes:

-The number of subintervals ( lines 184 & 191 ) is perhaps not quite enough.
-Try another number if need be.

-If you want to use the densities at altitudes 0 0.5 1 1.5 ... km, create a data file of 177 numbers and store

Ln d(0) Ln d(0.5) Ln d(1)  Ln d(1.5) .................. Ln d(88)

replace the first instructions in LBL 04  ( lines 288 to 304 )  by

LBL 04          RCL 09        RCL 07        /                               after replacing   lines  432 to 435  by
ENTER          -                   GETRX       -
STO 15          X<0?            CLX           ST+ X                       RCL 45      ST/ 00     /
ST+ X            CLX             RCL 09      ENTER                      RCL 09     RCL 11
INT                SEEKPT       ST+ Y       STO Z    ........             /                 *

-Here are the densities for the same atmospheric model:

0.0 1.226E-03
0.5 1.161E-03
1.0 1.101E-03
1.5 1.045E-03
2.0 9.932E-04
2.5 9.446E-04
3.0 8.987E-04
3.5 8.552E-04
4.0 8.136E-04
4.5 7.738E-04
5.0 7.354E-04
5.5 6.985E-04
6.0 6.627E-04
6.5 6.280E-04
7.0 5.944E-04
7.5 5.619E-04
8.0 5.303E-04
8.5 4.997E-04
9.0 4.701E-04
9.5 4.414E-04
10.0 4.138E-04
10.5 3.873E-04
11.0 3.618E-04
11.5 3.375E-04
12.0 3.144E-04
12.5 2.925E-04
13.0 2.718E-04
13.5 2.523E-04
14.0 2.340E-04
14.5 2.168E-04
15.0 2.008E-04
15.5 1.858E-04
16.0 1.718E-04
16.5 1.588E-04
17.0 1.467E-04
17.5 1.355E-04
18.0 1.251E-04
18.5 1.154E-04
19.0 1.065E-04
19.5 9.830E-05
20.0 9.071E-05
20.5 8.370E-05
21.0 7.723E-05
21.5 7.126E-05
22.0 6.576E-05
22.5 6.069E-05
23.0 5.601E-05
23.5 5.169E-05
24.0 4.771E-05
24.5 4.404E-05
25.0 4.066E-05
25.5 3.754E-05
26.0 3.466E-05
26.5 3.201E-05
27.0 2.957E-05
27.5 2.731E-05
28.0 2.523E-05
28.5 2.332E-05
29.0 2.155E-05
29.5 1.992E-05
30.0 1.842E-05
30.5 1.703E-05
31.0 1.576E-05
31.5 1.458E-05
32.0 1.349E-05
32.5 1.248E-05
33.0 1.156E-05
33.5 1.070E-05
34.0 9.909E-06
34.5 9.178E-06
35.0 8.503E-06
35.5 7.881E-06
36.0 7.306E-06
36.5 6.777E-06
37.0 6.288E-06
37.5 5.838E-06
38.0 5.423E-06
38.5 5.040E-06
39.0 4.687E-06
39.5 4.362E-06
40.0 4.062E-06
40.5 3.785E-06
41.0 3.530E-06
41.5 3.294E-06
42.0 3.077E-06
42.5 2.876E-06
43.0 2.690E-06
43.5 2.519E-06
44.0 2.360E-06
44.5 2.214E-06
45.0 2.078E-06
45.5 1.952E-06
46.0 1.836E-06
46.5 1.727E-06
47.0 1.626E-06
47.5 1.532E-06
48.0 1.444E-06
48.5 1.361E-06
49.0 1.284E-06
49.5 1.211E-06
50.0 1.143E-06
50.5 1.079E-06
51.0 1.018E-06
51.5 9.607E-07
52.0 9.066E-07
52.5 8.554E-07
53.0 8.070E-07
53.5 7.611E-07
54.0 7.177E-07
54.5 6.765E-07
55.0 6.375E-07
55.5 6.005E-07
56.0 5.654E-07
56.5 5.321E-07
57.0 5.007E-07
57.5 4.709E-07
58.0 4.427E-07
58.5 4.161E-07
59.0 3.909E-07
59.5 3.672E-07
60.0 3.447E-07
60.5 3.235E-07
61.0 3.035E-07
61.5 2.847E-07
62.0 2.669E-07
62.5 2.502E-07
63.0 2.344E-07
63.5 2.196E-07
64.0 2.056E-07
64.5 1.924E-07
65.0 1.800E-07
65.5 1.684E-07
66.0 1.575E-07
66.5 1.472E-07
67.0 1.375E-07
67.5 1.285E-07
68.0 1.200E-07
68.5 1.120E-07
69.0 1.045E-07
69.5 9.747E-08
70.0 9.088E-08
70.5 8.470E-08
71.0 7.892E-08
71.5 7.350E-08
72.0 6.843E-08
72.5 6.368E-08
73.0 5.914E-08
73.5 5.502E-08
74.0 5.119E-08
74.5 4.761E-08
75.0 4.426E-08
75.5 4.114E-08
76.0 3.821E-08
76.5 3.548E-08
77.0 3.291E-08
77.5 3.052E-08
78.0 2.827E-08
78.5 2.617E-08
79.0 2.420E-08
79.5 2.236E-08
80.0 2.065E-08
80.5 1.904E-08
81.0 1.754E-08
81.5 1.615E-08
82.0 1.485E-08
82.5 1.364E-08
83.0 1.251E-08
83.5 1.147E-08
84.0 1.050E-08
84.5 9.607E-09
85.0 8.779E-09
85.5 8.014E-09
86.0 7.308E-09
86.5 6.659E-09
87.0 6.061E-09
87.5 5.512E-09
88.0 5.009E-09

5°)  Auer & Standish Model of the Atmosphere

-In this model, the refractivity of air is given by formulas of the form:

n - 1 = ( n0 - 1 ) [ 1 + C ( rT / ( rT + x ) - 1 ) ]5                  if   the altitude  x < 11016 m
n - 1 = ( n0 - 1 ) Exp [ C' ( rT / ( rT + x ) - 1 ) ]                  if   the altitude  x > 11016 m

-Cf reference [2] for the details

-Here, there are a few modifications to take the latitude & the azimuth into account.

Data Registers:              R00 = dn/dx              ( Registers R01 thru R06 are to be initialized before executing "AUSTD" )

•  R01 = t  ( °C )           •  R04 = wave-length ( µm )
•  R02 = P ( mbar )       •  R05 = Latitude ( ° ' " )
•  R03 = F ( mbar )       •  R06 = Altitude of the observer ( in meters ) < 11019

R07 to R38: temp

Flag:  F01                      CF 01 = t , P , F  are measured for an altitude = 0
SF 01 = t , P , F  are measured at the observer's station
Subroutines: /

-Lines 08 & 218 are three-byte GTOs

 01 LBL "AUSTD"  02 DEG  03 HR  04 STO 22  05 COS  06 STO 20   07 X=0?  08 GTO 12  09 STO 32   10 X<>Y  11 HR  12 1  13 P-R  14 X^2  15 .99330562  16 /  17 X<>Y  18 X^2  19 RCL 05  20 HR  21 SIN  22 X^2  23 669438 E-8  24 *  25 1  26 X<>Y  27 -  28 SQRT  29 ST* Z  30 ST* Z  31 ST* Z  32 *  33 +  34 6378.137  35 X<>Y  36 /  37 STO 23            38 RCL 06  39  E3  40 /  41 STO 20  42 STO 38  43 +  44 STO 36  45 RCL 01  46 273.15  47 +  48 STO 15  49 1/X  50 ENTER  51 ENTER  52 STO 00  53 77514.1  54 * 55 710.792  56 -  57 *  58 2.23366  59 +  60 *  61 237321 E-8  62 -  63 RCL 03   64 ST* Y  65 37 E-5  66 *  67 1  68 +  69 *  70 1  71 +  72 RCL 03  73 *  74 *  75 STO 07  76 RCL 04  77 X^2  78 1/X  79 ENTER  80 ENTER  81 STO 08  82 .08851  83 *  84 .7115  85 -  86 *  87 58.058  88 +  89 *  90 6487.31  91 +  92 ST* 07  93 RCL 00            94 .25844  95 *  96 9325 E-7  97 -  98 RCL 00  99 * 100 579 E-9 101 + 102 RCL 02 103 RCL 03 104 - 105 ST* 00 106 * 107 1 108 + 109 RCL 00 110 * 111 4547.3 112 38.9 113 RCL 08 114 - 115 / 116 130 117 RCL 08  118 - 119 683939.7 120 X<>Y 121 / 122 + 123 2371.34 124 + 125 * 126 RCL 07 127 + 128  E8 129 / 130 STO 19 131 1.6 132 STO 31 133 FRC 134 SQRT 135 STO 30 136 5 137 CHS 138 STO 07 139 287053 E-9 140 STO 09 141 RCL 15 142 * 143 RCL 05 144 HR 145 SIN 146 X^2 147 STO 00           148 228 149 * 150 51631 151 + 152 RCL 00 153 * 154 9780325 155 + 156  E9 157 / 158 RCL 23 159 X^2 160 * 161 X<>Y 162 / 163 STO 11 164 6 165 / 166 STO 10 167 SF 10 168 CLX 169 X<> 20 170 XEQ 04 171 RCL 14  172 FS? 01 173 RCL 19 174 RCL 36  175 ST* Y 176 + 177 RCL 32 178 * 179 STO 35 180 FC? 01 181 GTO 00 182 RCL 14 183 RCL 19 184 / 185 ST/ 19 186 LBL 00 187 3 E-6 188 STO 18 189 11.019 190 XEQ 02 191 STO 21 192 RCL 14 193 STO 15 194 RCL 12 195 CF 10 196 XEQ 04 197 RCL 14 198 RCL 08 199 / 200 ST* 15 201 6 202 STO 17           203 SF 10 204 XEQ 13 205 CF 10 206 28 207 XEQ 02 208 X<> 21 209 STO 22 210 4 211 STO 17 212 XEQ 13 213  E2 214 XEQ 02 215 X<> 21 216 STO 22 217 XEQ 13 218 GTO 12 219 LBL 13 220 RCL 22 221 RCL 21 222 STO 24 223 - 224 RCL 17  225 STO 28 226 / 227 STO 27 228 2 229 / 230 ST+ 24 231 RCL 30  232 * 233 STO 25 234 CLX 235 STO 26 236 LBL 14 237 RCL 24 238 RCL 25 239 - 240 XEQ 01 241 ST+ 26 242 RCL 24 243 XEQ 01 244 RCL 31 245 * 246 ST+ 26 247 RCL 24 248 RCL 25 249 + 250 XEQ 01 251 ST+ 26 252 RCL 27 253 ST+ 24 254 DSE 28 255 GTO 14 256 RCL 26           257 * 258 ST+ 20 259 RTN 260 LBL 02 261 STO 34 262 XEQ 04 263 STO 16 264 RCL 12 265 RCL 23 266 + 267 * 268 RCL 35 269 X<>Y 270 / 271 ACOS 272 RTN 273 LBL 01 274 COS 275 STO 29 276 LBL 03 277 RCL 34  278 XEQ 04 279 STO 37 280 RCL 12 281 RCL 23 282 + 283 STO 33  284 * 285 RCL 35 286 RCL 29 287 / 288 - 289 RCL 00 290 RCL 33 291 * 292 RCL 37 293 + 294 / 295 ST- 34 296 ABS 297 RCL 18 298 X 00 338 RCL 10  339 * 340 RCL 07 341 * 342 RCL 13 343 X^2 344 / 345 X<> 00 346 RTN 347 LBL 05 348 RCL 11 349 * 350 E^X 351 RCL 15 352 * 353 SIGN 354 LASTX 355 STO 08 356 + 357 STO 00 358 LASTX 359 CHS 360 RCL 11 361 * 362 RCL 13            363 X^2 364 / 365 X<> 00 366 RTN 367 LBL 12 368 RCL 20 369  E3 370 * 371 STO Y 372 3600 373 / 374 HMS 375 END

( 669 bytes / SIZE 039 )

 STACK INPUTS OUTPUTS Y Az ( ° ' " ) R ( " ) X h0 ( ° ' " ) R ( ° ' " )

Example:       t = 10°C  P = 1010 mbar  F = 6 mbar  Lambda = 0.577 µm  Latitude = 33°21'22"  Altitude = 1706 m

10    STO 01        0.577    STO 04
1010  STO 02      33.2122  STO 05
6     STO 03        1706     STO 06

1°)    CF 01    Az = 12°41'   h0 = 1°23'45"

12.41    ENTER^
1.2345   XEQ "AUSTD"  >>>>  R = 0°17'54"337                    ---Execution time = 5m28s---
RDN   R = 1074"337

2°)    SF 01    Az = 12°41'   h0 = 1°23'45"

12.41    ENTER^
1.2345       R/S       >>>>    R = 0°21'28"454
RDN   R = 1288"454

6°)  Faster Programs

a) REFR

-"REFR" approximate the results given by "REF13" for

t = 10°C  P = 1010 mbar  F = 0 mbar  Lambda = 0.59 µm  Latitude = 45°  Altitude = 0

-The error ( compared to "REF13" ) are:

err < 0"003 over [5°,90°]
err < 0"004 over [1°,5°]
err < 0"005 over [0°,1°]
err < 0"011 over [-1°,0°]

-Such an accuracy is of course illusory near the horizon...
-And the precision is much smaller for other atmospheric conditions.

Formula for h0 > 5°

R ~ 57"91214 / Tan h0 - 0"06675061 / Tan3 h0 + 1"97745 E-4 / Tan5 h0 - 6"652813 E-7 Tan7 h0  + 1"306196 E-9 / Tan9 h0

Formula for h0 < 5°

R ~  Exp ( A + B h0 + C h02 + D h03 + E h04 + F h05 + G h06 + H h07 + I h08 + J h09 + K h010 + L h011 + M h012 )    with

A = +7.631589            E = - 0.010024199        I = - 4.8712695 E-4      M = - 3.7223778 E-8
B = - 0.3890402          F = +0.007237414        J = +1.0159107 E-4
C = +0.03649829        G = - 0.0039216984     K = - 1.3748284 E-5
D = +0.006352585      H = +0.0016179943     L = +1.0796128 E-6

Data Registers:              R00 = h0 ( deg )        ( Registers R01 thru R06 are to be initialized before executing "REFR" )

•  R01 = t  ( °C )           •  R04 = wave-length ( µm )
•  R02 = P ( mbar )       •  R05 = Latitude ( ° ' " )
•  R03 = F ( mbar )       •  R06 = Altitude of the observer ( in meters )

Flags: /
Subroutines: /

 01 LBL "REFR"  02 DEG  03 HR  04 STO 00   05 5  06 X<>Y  07 X>Y?  08 GTO 01   09 1.0796128 E-6  10 RCL Y  11 3.7223778 E-8  12 *  13 -  14 *  15 1.3748284 E-5  16 -  17 *  18 1.0159107 E-4  19 +  20 *  21 4.8712695 E-4  22 -  23 *  24 .0016179943  25 +  26 * 27 .0039216984  28 -  29 *  30 .007237414  31 +  32 *  33 .010024199  34 -  35 *  36 .006352585  37 +  38 *  39 .03649829  40 +  41 *  42 .3890402  43 -  44 *  45 7.631589          46 +  47 E^X  48 GTO 02   49 LBL 01  50 TAN  51 1/X  52 X^2 53 6.652813 E-7  54 CHS  55 RCL Y  56 1.306196 E-9  57 *  58 +  59 *  60 1.97745 E-4  61 +  62 *  63 .06675061  64 -  65 *  66 57.91214  67 +  68 X<>Y  69 SQRT  70 *  71 LBL 02               72 RCL 02   73 *  74 3.56701  75 /  76 RCL 01  77 273.15  78 + 79 /  80 1  81 RCL 03  82 18 E4  83 /  84 6579  85 1/X  86 +  87 RCL 03   88 *  89 -  90 *  91 5981  92 RCL 04  93 X^2  94 /  95 982818             96 +  97  E6  98 /  99 * 100 RCL 05  101 HR 102 ST+ X 103 COS 104 RCL 00 105 49 106 * 107 197 108 + 109 RCL 00  110 * 111 500 112 + 113 / 114 CHS 115 1 116 + 117 * 118 RCL 06            119 11000 120 / 121 CHS 122 E^X 123 * 124 STO Y 125 3600 126 / 127 HMS 128 END

( 350 bytes / SIZE 007 )

 STACK INPUTS OUTPUTS Y / R ( " ) X h0 ( ° ' " ) R ( ° ' " )

Example:       t = 10°C  P = 1010 mbar  F = 6 mbar  Lambda = 0.577 µm  Latitude = 33°21'22"  Altitude = 1706 m

10    STO 01        0.577    STO 04
1010  STO 02      33.2122  STO 05
6     STO 03        1706     STO 06

1°)      h0 = 1°23'45"

1.2345   XEQ "REFR"  >>>>  R = 0°18'20"742                    ---Execution time = 10s---
RDN   R = 1100"742

2°)     h0 = 12°34'56"

12.3456    R/S       >>>>    R = 0°03'37"253                    ---Execution time = 7s---
RDN    R = 217"253

Note:

-The results rapidly become meaningless for  h0 < -1°

-"RADAU" approximates the table of refraction given in reference [5] - pages A42 & A43 - over [-1°,+90°]
for a temperature = 0°C , atmospheric pressure = 760 mmHg , water-vapor pressure = 6 mmHg , latitude = 45° , altitude = 0

Formula:

R ~  A / Tan ( h0 + B / ( h0 + C / ( h0 + D / ( h0 + E / ( h0 + F ) ) ) ) )

with   A = 1°/59.79268   B = 3.68278  C = 7.37814  D = 15.08593  E = 64.96944  F = 13.55049

h =  h0 - R

precision ~ 0"06

 01 LBL "RADAU"  02 DEG  03 HR  04 64.96944  05 RCL Y  06 13.55049  07 +  08 /  09 +  10 15.08593  11 X<>Y  12 /  13 +  14 7.37814  15 X<>Y  16 /  17 +  18 3.68278  19 X<>Y  20 /  21 +  22 TAN  23 1/X  24 59.79268  25 /  26 ST- Y  27 HMS  28 X<>Y  29 HMS  30 END

( 82 bytes / SIZE 000 )

 STACK INPUTS OUTPUTS Y / R ( ° ' " ) X h0 ( ° ' " ) h ( ° ' " )

Example:        h0 = -0°12'34"

0.1234  CHS   XEQ "RADAU"  >>>>    h = -0°52'22"71                    ---Execution time = 3s---
RDN     R =  0°39'48"71

Formula:

R ~  A / Tan ( h + B / ( h + C / ( h + D / ( h + E / ( h + F ) ) ) ) )

with   A = 1°/59.76866   B = 4.67605  C = 7.93897  D = 16.24011  E = 73.68457  F = 14.61994

h0 =  h + R

precision ~ 0"10

 01 LBL "RADAU"  02 DEG  03 HR  04 73.68457  05 RCL Y  06 14.61994  07 +  08 /  09 +  10 16.24011  11 X<>Y  12 /  13 +  14 7.93897  15 X<>Y  16 /  17 +  18 4.67605  19 X<>Y  20 /  21 +  22 TAN  23 1/X  24 59.76866  25 /  26 ST+ Y  27 HMS  28 X<>Y  29 HMS  30 END

( 82 bytes / SIZE 000 )

 STACK INPUTS OUTPUTS Y / R ( ° ' " ) X h ( ° ' " ) h0 ( ° ' " )

Example:        h = -1°

1  CHS   XEQ "RADAU"  >>>>    h0 = -0°18'31"14                    ---Execution time = 3s---
RDN     R =  0°41'28"86

c) PH0-H & PH-H0

-"PH0-H" uses the refraction tables of Pulkovo for standard atmospheric conditions:

t = 15°C ,  P= 1013.25 mbar ,  F = 0 , lambda = 0.59 µ , latitude = 45° , altitude = 0

-It takes the apparent height of a star and returns the true height.

-PH-H0" performs the reverse transformation.

-With "PH0-H", the precision is ~ 0"26 over [0°,20°]

-Better than 0"02 over [20°,23°]
-Better than 0"01 over [23°,30°]
-Better than 0"002 over [30°,90°]

-With "PH-H0" the precision is ~ 0"11 over [0°,20°]

Formulae:

1°)  "PH0-H"

R ~  A / Tan ( h0 + B / ( h0 + C / ( h0 + D / ( h0 + E / ( h0 + F ) ) ) ) )           if  h0 < 20°

with   A = 1°/63.05561   B = 3.81451  C = 6.04529  D = 8.42681  E = 23.82074  F = 7.40780

R ~ 57"085 / Tan h0 - 0"0666 / Tan3 h0     if    h0 > 20°

2°)  "PH-H0"

R ~  A / Tan ( h + B / ( h + C / ( h + D / ( h + E / ( h + F ) ) ) ) )          if  h < 20°

with   A = 1°/62.93951   B = 4.80017  C = 6.90263  D = 10.06891  E = 31.76812  F = 8.87360

R ~ 57"0684 / Tan h - 0"081674 / Tan3 h     if     h > 20°

Data Registers: /
Flags: /
Subroutines: /

 01 LBL "PH0-H"  02 DEG  03 HR  04 STO M  05 20  06 X<>Y            07 XY            31 /  32 +  33 6.04529  34 X<>Y  35 /  36 +  37 3.81451  38 X<>Y  39 /  40 + 41 TAN  42 1/X  43 63.05561  44 /  45 LBL 02  46 0  47 X<> M  48 X<>Y            49 -  50 HMS  51 RTN  52 LBL "PH-H0"  53 DEG  54 HR  55 STO M  56 20  57 X<>Y  58 XY  82 /  83 +  84 6.90263  85 X<>Y  86 /  87 +  88 4.80017  89 X<>Y            90 /  91 +  92 TAN  93 1/X  94 62.93951  95 /  96 LBL 02  97 0  98 X<> M  99 + 100 HMS 101 END

( 233 bytes / SIZE 000 )

 STACK INPUT OUTPUT X h0 or h ( ° ' " ) h or h0 ( ° ' " )

Example1:       h0 = 1°23'45"

1.2345   XEQ "PH0-H"  >>>>   h = 1°02'51"39       R/S      >>>>    h0 = 1°23'45"02               ---Execution time = 2 or 3s---

Example2:       h = 24°12'57"

24.1257  XEQ "PH-H0"   >>>>   h0 = 24°15'02"99       R/S      >>>>   h = 24°12'57"00            ---Execution time = 2 or 3s---

Note:

-Synthetic register M is cleared.

References:

[1]  Jean Kovalevsky et al. - "Introduction aux Ephemerides Astronomiques" - EDP Sciences - ISBN 2-86883-298-9  ( in French )
[2]  Lawrence H. Auer & E. Myles Standish - "Astronomical Refraction, Computational method for all zenith angles" 2000 AJ, 119 , 2472
[3]  Andrew T. Young - "Understanding Astronomical Refraction" - 2006, The Observatory, Vol. 126, N° 1191, pp. 82-115
[4]  https://ccmc.gsfc.nasa.gov/modelweb/models/msis_vitmo.php
[5]  "Connaissance des Temps 1977" - Gauthier-Villars - Paris
[6]   Abramowitz and Stegun - "Handbook of Mathematical Functions" -  Dover Publications -  ISBN  0-486-61272-4