# hp41programs

Coordinates Transformation of Coordinates and Precession for the HP-41

Overview:

1°) 3 Subroutines

a) Rectangular-Spherical conversion
b) Spherical-Rectangular conversion
c) A very useful subroutine: "EE"

2°) Equatorial & Ecliptic Coordinates

a) Equatorial >>> Ecliptic
b) Ecliptic >>> Equatorial

3°) Equatorial & Azimuthal Coordinates

a) Equatorial >>> Azimuthal
b) Azimuthal >>> Equatorial

4°) Galactic Coordinates

a) Equatorial >>> Galactic
b) Galactic >>> Equatorial

5°) Precession

a) Equatorial coordinates
b) Ecliptic coordinates, program#1
c) Ecliptic coordinates, program#2

1°)  3 Subroutines

a) Rectangular-Spherical conversion:

x = r cos b cos l
y = r cos b sin l
z = r sin b

where    x , y , z = rectangular coordinates,    r  = ( x2 + y2 + z2 )1/2  ,   l  = longitude ,  b = latitude

-However, the results can be obtained more easily by the R-P and P-R functions:
T-register is saved and no data register is used!

 01  LBL "R-S"  02  R-P  03  X<>Y  04  RDN  05  R-P  06  R^  07  X<>Y  08  END

( 16 bytes / SIZE 000 )

 STACK INPUTS OUTPUTS T T T Z z b Y y l X x r L / (x2+y2)1/2

Example:     x = 3 ; y = 4 ; z = -7   find the spherical coordinates ( in DEG mode )

-7  ENTER^
4  ENTER^
3  XEQ "R-S"    r  =  8.602325267
RDN     = 53.13010235°
RDN     b = -54.46232221°

b) Spherical-Rectangular conversion:

 01  LBL "S-R"  02  X<>Y  03  RDN  04  P-R  05  R^  06  X<>Y  07  P-R  08  END

( 16 bytes / SIZE 000 )

 STACK INPUTS OUTPUTS T T T Z b z Y l y X r x L / (x2+y2)1/2

Example:     r = 10 ; l = 124° ; b = 37°   find    x ; y ; z   ( in DEG mode )

37   ENTER^
124  ENTER^
10  XEQ "S-R"   x = -4.465913097
RDN   y =  6.620988446
RDN   z =  6.018150232

- "R-S" and "S-R" work in all angular modes

c) A very useful subroutine: "EE"

-Many transformations use the same type of formulae which appear in the equatorial-ecliptic conversion, namely:

sin  b = cos e  sin d - sin e  cos d  sin a
cos cos l = cos d  cos a
cos b  sin  l = sin e  sin d  +  cos e  cos d  sin a

- But once again, P-R and R-P lead to a shorter and faster algorithm.

 01  LBL "EE"  02  1  03  XEQ "S-R"  04  RDN  05  R-P  06  X<> Z  07  ST- Y  08  X<> Z  09  P-R  10  R^  11  XEQ "R-S"  12  RDN  13  END

( 31 bytes / SIZE 000 )

-If "EE" is useful for you but "R-S" and "S-R" are not, here is another version of this program that doesn't need any subroutine:

 01  LBL "EE"  02  1  03  X<>Y  04  RDN  05  P-R  06  R^  07  X<>Y  08  P-R  09  RDN  10  R-P  11  X<> Z  12  ST- Y  13  X<> Z  14  P-R  15  R^  16  R-P  17  X<>Y  18  RDN  19  R-P  20  X<> T  21  END

( 32 bytes / SIZE 000 )

 STACK INPUTS OUTPUTS Z e e Y decl b X RA l

where  e = obliquity of the ecliptic

Example:     if right-ascension = RA = 116.328942 , declination = decl = 28.026183 and  e = 23.4392911

23.4392911  ENTER^
28.026183    ENTER^
116.328942    XEQ "EE"  >>>> l  =  113.215630°   RDN   b  =  6.684170°

-Like "R-S" and "S-R" , "EE"  works in all angular modes.
-But here are more easy-to-use programs:

2°) Equatorial & Ecliptic Coordinates

a) Equatorial >>> Ecliptic

Data Registers:    R00 = number of millenia since 2000/01/01  0h        R01 & R02: temp

Flag:   F01   Set flag F01 if you are using apparent coordinates
Clear flag F01 ---------------   mean     -----------

Subroutines:  "J0"      ( cf "Julian & Gregorian Calendar for the HP-41" )
"OBL"  ( cf "Nutation - Obliquity of the Ecliptic - Sidereal time for the HP-41" )
"NUT"  ( idem )  if SF 01  ;  "EE"

 01  LBL "EQ-ECL" 02  STO 01 03  RDN 04  STO 02 05  X<>Y 06  XEQ "J0" 07  365250 08  / 09  STO 00 10  XEQ "OBL" 11  LASTX 12  RCL 02            13  HR 14  RCL 01 15  HR 16  15 17  * 18  XEQ "EE"        19  X<>Y 20  HMS 21  X<>Y              22  360 23  MOD 24  HMS 25  END

( 54 bytes SIZE 003 )

 STACK INPUTS OUTPUTS Z YYYY.MNDDdd e or em Y decl ( ° . ' " ) b ( ° . ' " ) X RA (hh.mnss) l ( ° . ' " )

where   l  = ecliptic longitude , b  = ecliptic latitude

Example1:     On  2134/04/04 at 0h  if  Right-Ascension = 12h34m56s  and  Declination =  25°12'49"  ( mean coordinates )

CF 01
2134.0404  ENTER^
25.1249  ENTER^
12.3456  XEQ "EQ-ECL"  >>>>  l  =  177°13'44"69    RDN    b  =  26°27'18"71

Example2:     On  2134/04/04 at 0h  if  RA = 12h34m56s  and  decl =  25°12'49"  ( apparent coordinates )

SF 01
2134.0404  ENTER^
25.1249  ENTER^
12.3456  XEQ "EQ-ECL"  >>>>  l  =  177°13'41"39    RDN    b  =  26°27'18"39

b) Ecliptic >>> Equatorial

Data Registers:    R00 = number of millenia since 2000/01/01  0h        R01 & R02: temp

Flag:   F01   Set flag F01 if you are using apparent coordinates
Clear flag F01 ---------------   mean     -----------

Subroutines:  "J0"      ( cf "Julian & Gregorian Calendar for the HP-41" )
"OBL"  ( cf "Nutation - Obliquity of the Ecliptic - Sidereal time for the HP-41" )
"NUT"  ( idem )  if SF 01 ;  "EE"

 01  LBL "ECL-EQ" 02  STO 01 03  RDN 04  STO 02 05  X<>Y 06  XEQ "J0" 07  365250 08  / 09  STO 00 10  XEQ "OBL"     11  LASTX 12  CHS 13  RCL 02 14  HR 15  RCL 01 16  HR 17  XEQ "EE"        18  X<>Y 19  HMS 20  X<>Y              21  15 22  / 23  24 24  MOD 25  HMS               26  END

( 54 bytes / SIZE 003 )

 STACK INPUTS OUTPUTS Z YYYY.MNDDdd e or em Y decl ( ° . ' " ) b ( ° . ' " ) X RA (hh.mnss) l ( ° . ' " )

where   l  = ecliptic longitude , b  = ecliptic latitude

Example1:     On  2134/04/04 at 0h  if  l  =  177°13'44"69  and b  =  26°27'18"71  ( mean coordinates )

CF 01
2134.0404      ENTER^
26.271871  ENTER^
177.134469  XEQ "ECL-EQ"  >>>>  RA = 12h34m56s00    RDN    decl =  25°12'49"00

Example2:     On  2134/04/04 at 0h  if  l  =  177°13'41"39  and b  =  26°27'18"39  ( apparent coordinates )

SF 01
2134.0404  ENTER^
26.271839  ENTER^
177.134139  XEQ "ECL-EQ"  >>>>  RA = 12h34m56s00    RDN    decl =  25°12'49"00

Note:    If you know the mean equatorial coordinates and you want to find the apparent equatorial coordinates,

a)  Execute "EQ-ECL" with CF 01 to calculate the mean ecliptic coordinates
b)  Add the nutation in longitude to the mean ecliptic longitude
c)  Execute "ECL-EQ" with SF 01

3°) Equatorial & Azimuthal Coordinates

a) Equatorial >>> Azimuthal

Data Registers:    R00 = number of millenia since 2000/01/01  0h        R01 & R05: temp   R03 = local sidereal time, R04 = hour angle

Flag:   F01   Set flag F01 if you are using apparent coordinates
Clear flag F01 ---------------   mean     -----------

Subroutines:  "J0"      ( cf "Julian & Gregorian Calendar for the HP-41" )
"OBL"  ( cf "Nutation - Obliquity of the Ecliptic - Sidereal time for the HP-41" )
"ST"      ( idem )
"NUT"  ( idem )  if SF 01 ;  "EE"

-Line 07 = the East longitude of the US Naval Observatory.  Change this line according to your location or store your longitude in a data register.
-Line 20 = the latitude of the US Naval Observatory.  Change this line according to your location or store your latitude in a data register.

 01  LBL "EQ-AZ" 02  STO 01 03  RDN 04  STO 02 05  RDN 06  XEQ "ST" 07  -77.0356 08  HR 09  15 10  / 11  HMS 12  HMS+ 13  STO 03           14  RCL 01 15  HMS- 16  STO 04 17  HR 18  15 19  * 20  38.55172 21  HR 22  RCL 02           23  HR 24  90 25  ST- Z 26  R^ 27  - 28  XEQ "EE" 29  CHS 30  90 31  + 32  X<>Y             33  HMS 34  X<>Y 35  360  36  MOD 37  180  38  X

( 84 bytes / SIZE 005 )

 STACK INPUTS OUTPUTS T YYYY.MNDD / Z HH.MNSS (UT) / Y decl ( ° . ' " ) altitude ( ° . ' " ) X RA (hh.mnss) azimuth ( ° . ' " )

-Azimuths are measured clockwise ( westwards ) from South
-If you reckon the azimuths clockwise from North, as the navigators prefer, replace the + ( line 31 ) with a  -

Example1:    On  2005/12/12  at  20h51m29s (UT)  we have  R.A. = 7h41m16s  &  decl = 60°21'37"   ( mean coordinates ),
compute the azimuthal coordinates at the US Naval Observatory.

CF 01
2005.1212  ENTER^
20.5129  ENTER^
60.2137  ENTER^
7.4116  XEQ "EQ-AZ"  >>>>   Azimuth = -169°03'28"33   RDN   Altitude = 10°55'57"90

Example2:    On  2005/12/12  at  20h51m29s (UT)  we have  R.A. = 7h41m16s  &  decl = 60°21'37"   ( apparent coordinates )

SF 01
2005.1212  ENTER^
20.5129  ENTER^
60.2137  ENTER^
7.4116  XEQ "EQ-AZ"  >>>>   Azimuth = -169°03'29"87   RDN   Altitude = 10°55'57"43

-If you want to take the refraction into account, cf for instance "Atmospheric Refraction for the HP-41"

b) Azimuthal >>> Equatorial

Data Registers:    R00 = number of millenia since 2000/01/01  0h        R01 & R05: temp   R03 = local sidereal time

Flag:   F01   Set flag F01 if you are using apparent coordinates
Clear flag F01 ---------------   mean     -----------

Subroutines:  "J0"      ( cf "Julian & Gregorian Calendar for the HP-41" )
"OBL"  ( cf "Nutation - Obliquity of the Ecliptic - Sidereal time for the HP-41" )
"ST"      ( idem )
"NUT"  ( idem )  if SF 01 ;  "EE"

-Line 07 = the East longitude of the USNO.  Change this line according to your location or store your longitude in a data register.
-Line 14 = the latitude of the USNO.  Change this line according to your location or store your latitude in a data register.
-Add   CHS  after line 20 if you measure the azimuths from North.

 01  LBL "AZ-EQ" 02  STO 01 03  RDN 04  STO 02 05  RDN 06  XEQ "ST" 07  -77.0356 08  HR 09  15 10  / 11  HMS 12  HMS+ 13  STO 03 14  38.55172        15  CHS 16  HR 17  RCL 02 18  HR 19  90 20  ST+ Z 21  RCL 01 22  HR 23  - 24  XEQ "EE" 25  90 26  - 27  X<>Y 28  HMS 29  X<>Y 30  15 31  / 32  RCL 03 33  HR 34  + 35  24 36  MOD              37  HMS 38  END

( 74 bytes / SIZE 004 )

 STACK INPUTS OUTPUTS T YYYY.MNDD / Z HH.MNSS (UT) / Y altitude ( ° . ' " ) decl ( ° . ' " ) X azimuth ( ° . ' " ) RA (hh.mnss)

Example1:    On  2005/12/12  at  20h51m29s (UT)  we have  Azimuth = 190°56'31"67   &   Altitude = 10°55'57"90   ( mean coordinates ),
compute the equatorial coordinates at the US Naval Observatory.

CF 01
2005.1212      ENTER^
20.5129      ENTER^
10.555790  ENTER^
190.563167  XEQ "AZ-EQ"  >>>>   R.A. = 7h41m16s00   RDN   decl = 60°21'37"00

Example2:    On  2005/12/12  at  20h51m29s (UT)  we have  Azimuth = 190°56'30"13   &  Altitude = 10°55'57"43   ( apparent coordinates )

SF 01
2005.1212      ENTER^
20.5129      ENTER^
10.555743  ENTER^
190.563013  XEQ "AZ-EQ"  >>>>   R.A. = 7h41m16s00   RDN   decl = 60°21'37"00

4°) Galactic Coordinates

a) Equatorial >>> Galactic

Data Registers:  /
Flags:  /
Subroutine:   "EE"

 01  LBL "EQ-GAL" 02  DEG 03  HR 04  15 05  * 06  77.140702 07  + 08  X<>Y 09  HR 10  62.871732 11  X<> Z 12  XEQ "EE" 13  X<>Y 14  HMS 15  X<>Y 16  32.93192 17  + 18  360 19  MOD 20  HMS 21  END

( 62 bytes / SIZE 000 )

 STACK INPUTS OUTPUTS Y decl ( ° . ' " ) b ( ° . ' " ) X RA (hh.mnss) l ( ° . ' " )

where   l  = galactic longitude , b  = galactic latitude

Example:    The coordinates of Procyon are  Right Ascension = 7h39m18s1  ,  Declination = 5°13'30"  for the equinox of  J2000.0

5.1330    ENTER^
7.39181  XEQ "EQ-GAL"  >>>>  l  = 213°42'08"19   RDN   b  =  13°01'10"16

-In this program, right ascensions and declinations must be referred to J2000.0  ( Julian day = 2451545 )
-If  R.A. & decl are referred to B1950.0 ( Julian day = 2433282.4235 ),  replace lines 06-10-16  with  77.75  ;  62.6  ;  33  respectively.

b) Galactic >>> Equatorial

Data Registers:  /
Flags:  /
Subroutine:   "EE"

 01  LBL "GAL-EQ" 02  DEG 03  HR 04  32.93192 05  - 06  X<>Y 07  HR 08  62.871732 09  CHS 10  X<> Z 11  XEQ "EE" 12  X<>Y 13  HMS 14  X<>Y 15  77.140702 16  - 17  15 18  / 19  24 20  MOD 21  HMS 22  END

( 62 bytes / SIZE 000 )

 STACK INPUTS OUTPUTS Y b ( ° . ' " ) decl ( ° . ' " ) X l  ( ° . ' " ) RA ( hh.mnss )

where   l  = galactic longitude , b  = galactic latitude

Example:    l  = 213°42'08"19   &   b  =  13°01'10"16   for the equinox of  J2000.0

13.011016    ENTER^
213.420819    XEQ "GAL-EQ"  >>>>  RA  = = 7h39m18s10   RDN   Declination = 5°13'30"00

-Right ascensions and declinations are referred to J2000.0  ( Julian day = 2451545 )
-If  R.A. & decl are referred to B1950.0 ( Julian day = 2433282.4235 ),  replace lines 04-08-15  with  33 ;  62.6  ;  77.75  respectively.

5°) Precession

a) Equatorial coordinates

Data Registers:    R00 thru R03: temp
Flag:   F00
Subroutines:  "J0"      ( cf "Julian & Gregorian Calendar for the HP-41" ) , "EE"

-Line 95 is a three-byte GTO 00

 01  LBL "PREQ"   02  DEG   03  HR   04  STO 01   05  RDN   06  HR   07  STO 02   08  15   09  ST* 01   10  R^   11  STO 00   12  R^   13  SF 00   14  LBL 00   15  X<> 00   16  XEQ "J0"   17  .5   18  -   19  365250   20  /   21  17 22  RCL Y   23  9   24  *   25  +   26  *   27  5005   28  -   29  *   30  8301   31  -   32  *   33  6405787         34  -   35  *   36  CHS   37  736   38  +   39  STO 03   40  CLX   41  8   42  * 43  79   44  +   45  *   46  5075   47  -   48  *   49  30354   50  -   51  *   52  6405770         53  -   54  *   55  736   56  +   57   E6   58  ST/ 03   59  /   60  FC? 00   61  X<> 03   62  90   63  - 64  ST+ 01   65  CLX   66  5   67  +   68  *   69  4   70  *   71  11617   72  +   73  *   74  11930   75  +   76  *   77   E6   78  /   79  5.5672            80  -   81  *   82  FC? 00   83  CHS   84  RCL 02 85  RCL 01   86  XEQ "EE"   87  90   88  +   89  RCL 03   90  -   91  STO 01   92  X<>Y   93  STO 02   94  FS?C 00   95  GTO 00          96  HMS   97  X<>Y   98  15   99  / 100  24 101  MOD 102  HMS 103  END

( 187 bytes / SIZE 004 )

 STACK INPUTS OUTPUTS T YYYY.MNDDdd1 / Z YYYY.MNDDdd2 / Y decl1 ( ° . ' " ) decl2 ( ° . ' " ) X RA1 (hh.mnss) RA2 (hh.mnss)

Example:   The mean coordinates of Sirius on 1600/04/04 0h are  AR = 6h27m17s88 ;  Decl = -16°21'56"34
Calculate the equatorial coordinates on 2134/12/12 0h

1600.0404    ENTER^
2134.1212    ENTER^
-16.215634   ENTER^
6.271788   XEQ "PREQ"  >>>>  RA = 6h51m10s79   RDN   Decl =  -16°52'22"05

b) Ecliptic coordinates, program#1

-One could use similar expressions for the ecliptic coordinates ( see the next program below ).
-But the following program converts ecliptic coordinates into equatorial coordinates,
then "PREQ" is called and finally, "EQ-ECL" produces the required results.
-This method has at least one advantage: it saves bytes!
-However, it's also slower than direct formulae.

Data Registers:    R00 thru R04: temp
Flags:   F00 & F01
Subroutines:  "J0"  ( cf "Julian & Gregorian Calendar for the HP-41" ) "ECL-EQ"  "PREQ"  "EQ-ECL"  "EE"

 01  LBL "PREC"  02  CF 01  03  R^  04  STO 03  05  X<> T  06  STO 04  07  RDN  08  XEQ "ECL-EQ"  09  RCL 03  10  RCL 04  11  R^  12  R^  13  XEQ "PREQ"  14  X<>Y  15  RCL 04  16  X<> Z  17  XEQ "EQ-ECL"  18  END

( 49 bytes / SIZE 005 )

 STACK INPUTS OUTPUTS T YYYY.MNDDdd1 / Z YYYY.MNDDdd2 / Y b1 ( ° . ' " ) b2 ( ° . ' " ) X l1  ( ° . ' " ) l2  ( ° . ' " )

where   l  = ecliptic longitudes , b  = ecliptic latitudes

Example:   The mean ecliptic coordinates of Sirius on 1600/04/04 0h are  l1 = 98°30'58"32  ;  b1 =  -39°39'17"79
Calculate the ecliptic coordinates on 2134/12/12 0h

1600.0404    ENTER^
2134.1212    ENTER^
-39.391779   ENTER^
98.305832   XEQ "PREC"  >>>> l2=  105°57'44"15   RDN b1=  -39°35'19"15

c) Ecliptic coordinates, program#2

-Now, we use the specific formulae for ecliptic coordinates.

Data Registers:    R00 thru R04: temp
Flag:   F00
Subroutines:  "J0"      ( cf "Julian & Gregorian Calendar for the HP-41" )  "EE"

-Line 82 is a three-byte GTO 00

 01  LBL "PREC2" 02  DEG 03  HR 04  STO 01 05  X<>Y 06  HR 07  STO 02 08  R^ 09  STO 00 10  R^ 11  SF 00 12  LBL 00 13  X<> 00 14  XEQ "J0" 15  .5 16  - 17  365250 18  / 19  133 20  CHS 21  RCL Y 22  ENTER^ 23  + 24  + 25  * 26  149 27  - 28  * 29  4389 30  + 31  * 32  2410993 33  - 34  * 35  174874109     36  + 37   E6 38  / 39  STO 03 40  ST- 01 41  RDN 42  66 43  - 44  * 45  22 46  + 47  * 48  30707 49  + 50  * 51  13968878      52  + 53  * 54  CHS 55   E6 56  / 57  STO 04 58  FS? 00 59  ST+ 01 60  CLX 61  35 62  * 63  930 64  + 65  * 66  130553          67  - 68  * 69   E6 70  / 71  FC? 00 72  CHS 73  RCL02 74  RCL 01 75  XEQ "EE" 76  RCL 03 77  + 78  STO 01 79  X<>Y 80  STO 02 81  FS?C 00 82  GTO 00          83  HMS 84  X<>Y 85  RCL 04 86  - 87  360 88  MOD 89  HMS 90  END

( 169 bytes / SIZE 005 )

 STACK INPUTS OUTPUTS T YYYY.MNDDdd1 / Z YYYY.MNDDdd2 / Y b1 ( ° . ' " ) b2 ( ° . ' " ) X l1  ( ° . ' " ) l2  ( ° . ' " )

where   l  = ecliptic longitudes , b  = ecliptic latitudes

Example:   The mean ecliptic coordinates of Sirius on 1600/04/04 0h are  l1 = 98°30'58"32  ;  b1 =  -39°39'17"79
Calculate the ecliptic coordinates on 2134/12/12 0h

1600.0404    ENTER^
2134.1212    ENTER^
-39.391779   ENTER^
98.305832   XEQ "PREC2"  >>>> l2=  105°57'44"15   RDN b1=  -39°35'19"15

-The polynomial coefficients of the precession angles ( expressed in degrees ) are rounded to 6 decimals, with T expressed in millenia from J2000.
-So, the accuracy is of the order of 0.01 arcsecond over the time span [ 1000 ; 3000 ]

References:

[1] Astronomical Algorithms - Jean Meeus - Willmann-Bell   ISBN 0-943396-35-2
[2] Spherical Astronomy - Robin M.Green - Cambridge University Press    ISBN  0-521-31779-7
[3] Introduction aux Ephemerides Astronomiques - EDP Sciences     ISBN  2-86883-298-9  ( in French )
[4] United States Naval Observatory, Circular n° 179 http://aa.usno.navy.mil/publications/docs/circular_179.html
[5] Report of the IAU division I working group on precession and the ecliptic.