Post Reply 
(41C) Tracking Satellites in Elliptical Orbits
11-23-2021, 09:08 PM (This post was last modified: 11-23-2021 09:13 PM by twdeckard.)
Post: #10
RE: (41C) Tracking Satellites in Elliptical Orbits
Greetings Angel,

thank you for your gracious follow up. I never made much progress on the project. All my energy was squandered on the Arduino driven antenna rotor drive and alas that is quite dusty as well.

Here is an early test:

https://www.youtube.com/watch?v=VG5L5oP2...j_&index=4

Here is the unit running in sky survey mode:

https://www.youtube.com/watch?v=3rMM7kyH..._&index=14

The telescope provided quite some amusement for my neighbor across the street when I teased him that based on the 300K noise I was picking up was throwing off my survey of the cosmic background radiation and probably a sign of some alien presence -- an elaborate aluminum foil ET appeared in his window a few days later.

Here is an actual attempt at a contact over an amateur radio satellite. This one was funny. There are several attempts that worked.

https://www.youtube.com/watch?v=O5nFwU_w...j_&index=3

You are correct the INPUTSTA is flawed. Early in the projects I generally "rough" in the application to see if it will even fit in the mighty HP41 calculator before I really dive in.

Here was the last variant of the code but again, I never carried it across the finish line. Thank you for your inquiry -- I really enjoy the forum and the community.

Code:
; SATTRK 
; DERIVED ENTIRELY FROM PLAN-13 BY J.R. Miller G3RUH
; REQUIRES 41CX WITH ADVANTAGE ROM 
;
; REGISTERS AND MNEMONICS BELOW:

00      4^TMP        ; 00-13 ALSO SCRATCH SPACE FOR ADV SOLVE FUNCTION
01    4^TMP
02    4^TMP
03    4^TMP        
04    4^SX        ; DAMN YOU ADVANTAGE ROM, VECTORS ON Rn+1 BOUNDARIES
05    4^SY
06    4^SZ
07    4^OX
08    4^OY
09    4^OZ
10     4^VOx
11    4^VOy
12    4^VOz
13    4^RX
14    4^RY
15     4^RZ
16    4ELAPST        ; (COMPUTED) ELAPSED TIME (DAYS) SINCE EPOCH [EPTIME] 
17    4CURORB        ; (COMPUTED) CURRENT ORBIT NUMBER
18     4OBJECT        ; (INPUT) SATELLITE OBJECT BEING TRACKED
19    1EPYEAR        ; (TLE) EPOCH YEAR                 1:19-20        
20    2EPTIME        ; (TLE) EPOCH TIME                    1:21-32
21    3INCLIN        ; (TLE) ORBITAL INCLINATION             2:09-16        
22    2R.A.A.N    ; (TLE) ORBITAL RIGHT ASCENSION OF ASCEND. NODE 2:27-33
23    4ECCENT        ; (TLE) ORBITAL ECCENTRICITY
24    2PERIGE        ; (TLE) ORBITAL PERIGEE             2:35-42
25    2MEANAN        ; (TLE) ORBITAL MEAN ANOMOLY             2:44-51
26    3MMDOT1        ; (TLE) ORBITAL ANGULAR MOMENTUM         1:34-43
27    3MMDOT2        ; (TLE) ORBITAL ANGULAR MOMENTUM RATE         1:45-52
28    4DECAYR        ; (TLE) ORBITAL DECAY RATE             1:54-61
29    1ORBITN        ; (TLE) ORBIT NUMBER                
30    6OBSLAT        ; (INPUT) OBSERVERS LATITUDE
31    6OBSLON        ; (INPUT) OBSERVERS LONGITUDE
32    4LAMBDA        ; (COMPUTED) LAMBDA
33    1SATRNG        ; (COMPUTED) SATELLITE RANGE
34    4DOPPLR        ; (COMPUTED) DOPPLER CORRECTION
35    6SATLAT        ; (COMPUTED) SATELLITE LATITUDE        
36    6SATLON        ; (COMPUTED) SATELLITE LONGITUDE
37    1AZIMUT        ; (COMPUTED) SATELLITE AZIMUTH
38    1ELEVAT        ; (COMPUTED) SATELLITE ELEVATION    
39    1OBSHGT        ; (INPUT) OBSERVERS HEIGHT
40     4MEANAT        ; (COMPUTED) MEAN ANOMOLY AT TIME T
41    4KD        ; (COMPUTED) LINEAR DRAG TERM             
42    4KDP        ; (COMPUTED) LINEAR DRAG TERM
;
;
; EXTENDED MEMORY FILES
;
; ISS        ALPHA        TWO-LINE-ELEMENTS
; SATREG    ALPHA        REGISTER MNEMONICS
;
; ENTRY POINTS
;
;SATTRK    SATELLITE TRACKING PROGRAM HP41
;DISTLE    DISPLAY TLE STORED IN EXTENDED MEMORY
;GETTLE    READ ORBITAL PARAMETERS FROM TLE IN EM
;KEPLER    COMPUTE LAMBDA 
;SATPOS    COMPUTE SATELLITE LATITUDE AND LONGITUDE
;AZELEV    COMPUTE AZIMUTH AND ELEVATION
;SATRNG    COMPUTE RANGE TO SATELLITE
;DOPPLR    COMPUTE DOPPLER SHIFT
;OUTPUT AZIMUTH AND ELEVATION
;EVRY20    ACTIVATE TIMER ROUTINE TO COMPUTE AND DISPLAY EVERY TWENTY SECONDS
;UPDATE COMPUTE SATPOS, AZIMUTH AND ELEVATION
;NEXT24    AZIMUTH AND ELEVATION OVER NEXT 24 HOURS
;NEXT12 AZIMUTH AND ELEVATION OVER NEXT 12 HOURS
;DISPLA    PRINT OUT REGISTERS WITH MNEMONICS
;DEORBI END PROGRAM

LBL "PLAN13"
000.04201 
CLRGX            ; CLEAR REGISTERS
"A"            ; CLEAR EXTENDED MEMORY FILES
XEQ 40
"B"
XEQ 40
"C" 
XEQ 40
"D" 
XEQ 40
"RAAN"
XEQ 40
"IN"
XEQ 40
XEQ INPSTA        ; PROMPT FOR OBJECT, OBS LAT/LON/HGT
GTO a            ; ** TEMPORARY SKIP TLE FOR DEBUG
XEQ GETTLE         ; PARSE TLE FILE IN EM
LBL a
XEQ **TEST
GTO b            ; TEMPORARY SKIP TLE DISPLAY FOR DEBUG
XEQ DISTLE        ; DISPLAY TWO LINE ELEMENTS
LBL b
XEQ *ECFENU        ; TRANSFORM MATRIX GEOCENTRIC TO LOCAL CARTESIAN    
XEQ *RADIUS        ; COMPUTE EARTH RADIUS WGS-84 AND SET UP OBS X,Y,Z CALCULATION
XEQ *OBSXYZ        ; COMPUTE OBSERVERS X,Y,Z COORDINATES        
XEQ *OBSVEL        ; COMPUTE OBSERVERS VELOCITY
;
; DEBUG VERSION ONLY PERFORMS CALCULATION         
; ONCE AND USES CURRENT SYSTEM TIME
;
XEQ *MA+ORB        ; COMPUTE CURRENT ORBIT NUMBER             
XEQ *SOLVEA        ; SOLVE KEPLERS EQUATION
XEQ *CELEST        ; [C] = [RAAN]*[IN]*[AP]        
XEQ *SATCEL        ; COMPUTE SATELLITE CELESTIAL COORDINATES
XEQ *SATGEO        ; COMPUTE SATELLITE GEOCENTRIC COORDINATES
XEQ *^RANGE        ; COMPUTE RANGE VECTOR
XEQ *AZIELE        ; COMPUTE AZIMUTH ELEVATION
STOP
;
; READ TWO LINE ELEMENTS FROM EM FILE NAMED IN [OBJECT] 
;
LBL "GETTLE"
0.000
CLA
ARCL 18
SEEKPTA
19
0.0017018
XEQ 20            ; <*GETEL>
XEQ 30            ; <*PARSE>
2000
ST+ 19
20
0.019031
XEQ 20            ; <*GETEL>
XEQ 30            ; <*PARSE>
21
1.007015
XEQ 20            ; <*GETEL>
XEQ 30            ; <*PARSE>
22
1.016024
XEQ 20            ; <*GETEL>
XEQ 30            ; <*PARSE>
23
1.025032
XEQ 20            ; <*GETEL>
>"."
-1
AROT            ; LEADING DECIMAL PRESUMED
RDN
XEQ 30            ; <*PARSE>
24
1.033041
XEQ 20            ; <*GETEL>
XEQ 30            ; <*PARSE>
25
1.042050
XEQ 20            ; <*GETEL>
XEQ 30            ; <*PARSE>
26
0.032042
XEQ 20            ; <*GETEL>
XEQ 30            ; <*PARSE>
27
0.043051
XEQ 20            ; <*GETEL>
XEQ 30            ; <*PARSE>
28
0.052060
XEQ 20            ; <*GETEL>
XEQ 30            ; <*PARSE>
RTN
LBL "*GETEL"
LBL 20
CLA
ARCL 18
SEEKPTA
GETREC
RTN
LBL "*PARSE"
LBL 30
FRC
1000
*
ENTER
INT
X<>Y
FRC
1000
*
X<>Y
-
2
+
42
XEQ "SC"        ; FROM PPC, SUBSTITUTE CHARACTER
ANUM            ; USED TO INSERT NON-NUMERIC TERMINATOR IN ALPHA 
STO IND T        ; TO PARSE NUMBERS FROM TWO LINE ELEMENTS
RTN
LBL "SC"
X<>Y
AROT
-1
AROT
RDN
ATOX
CLX
RDN
X<>Y
XTOA
X<>Y
CHS
AROT
CHS
X<>Y
RTN
;
; PROMPT FOR SATELLITE OBJECT AND OBSERVER COORDINATES
; DEBUG HARD CODED FOR INTERNATIONAL SPACE STATION AND EN34FV 
;
LBL "INPSTA"            
ISS
ASTO 18
44.885077         ; ARLINGTON CT
55.21            ; TEST LATITUDE IN PLAN13    
STO 30
-93.551546        ; ARLINGTON CT
.06            ; TEST LONGITUDE IN PLAN13
STO 31
800            ; ARLINGTON CT
79            ; TEST HEIGHT IN PLAN13
STO 39
RTN
LBL "SATPOS"
RTN
LBL "AZELEV"
RTN
LBL "SATRNG"
RTN
LBL "*SATGEO"
RTN
LBL "*^RANGE"
RTN
LBL "*AZIELE"
RTN
LBL "OUTPUT"
RTN
LBL "NEXT24"
RTN
LBL "NEXT12"
RTN
LBL 40            ; DELETE EM FILE NAMED IN ALPHA
SF 25
PURFL
RTN
LBL "DISTLE"        ; NOTE: CREATE A COMMON SUBROUTINE TO PRINT MNEUMONICS
018.02901
XEQ 36
RTN
LBL "DISALL"
000.03501
XEQ 36
RTN
LBL 36
STO 00
LBL 46
RCL 00
INT
XEQ 47
ISG 00
GTO 46
RTN
LBL 47
"SATREG"
SEEKPTA
GETREC
+":"
RCL IND X
ATOX
ASTO Z
X<>Y
FIX 0
XEQ IND Y
AVIEW
RTN
LBL 48
ARCL X
RTN
LBL 49
FIX 0
ARCL X
RTN
LBL 50
FIX 1
ARCL X
RTN
LBL 51
LBL 52
FIX 2
ARCL X
RTN
LBL 53
CLK24
FIX 4
ATIME24
RTN
LBL 54
FIX 4
HMS
ATIME24
RTN

LBL *RAAN        ; CREATE ROTATION ABOUT Z OF RAAN
"RAAN"
RCL 22            ; [RAAN]
XEQ *ROTZ        
RTN
LBL *IN0        ; CREATE ROTATION ABOUT X OF INCLINATION              
"IN"
RCL 21            ; [INCL]
XEQ *ROTX        
RTN
LBL *AP            ; CREATE ROTATION ABOUT Z OF ARGUMENT OF PERIGEE 
"AP"
RCL 24            ; [PERIGEE]
XEQ ROTZ
RTN
LBL *CELEST        ; ORBITAL PLANE TO CELESTIAL COORDINATE TRANSFORMATION
3.003
"C"
MATDIM
"B"
MATDIM
XEQ "*RAAN"
XEQ "*IN0"
"RAAN,IN,B"
M*M
"RAAN"
XEQ 40
"IN"
XEQ 40
XEQ "*AP"
"B,AP,C"
M*M            ; [C] = [RAAN]*[IN]*[AP]
"AP"
XEQ 40
"B"
XEQ 40
RTN    
;
; TRANSFORMATION FROM EAST NORTH UP TO EARTH CENTERED EARTH FIXED COORDINATES
;
; clockwise rotation over east-axis by angle (90-phi) to align up-axis with z-axis

; clockwise rotation over z-axis by angle (90+lambda) to align east axis with x-axis

; >>>>>USING THIS INSTEAD R3(PI)*R2(PHI-PI/2)*R3(LAMBDA-PI)

; http://gis.stackexchange.com/questions/82998/trasformation-from-ecef-to-enu
;
LBL "*ECFENU"        ;
3.003
ENU
MATDIM
"C"
MATDIM
-90
"A"
XEQ "*ROT3"
RCL 30
90
-
"B"
XEQ "*ROT2"
"A,B,C"
M*M
"A"
PURFL
"B" 
PURFL
RCL 31
180
-
"D"
XEQ "*ROT3"
"C,D,ENU"
M*M
"C"
PURFL
"D"
PURFL 
RTN
; ENU MATRIX DOES NOT LOOK CORRECT
;
; TEST RESULTS
; Ux=0.6127688125485028
; Uy=0.0006416902235684213
; Uz=0.7902619632764536

; Ex=-0.001047197341936584
; Ey=0.9999994516887132
; Ez=0

; Nx=-0.7902615299668996
; Ny=-0.0008275602273766885
; Nz=0.6127691485367431
;
;
; A*B=ENU (-LAT,-LON)
;
; 1:1= 0.79026=Uz   1:2= 0.00064=Uy  1:3= 0.61277=Ux
; 2:1= 0.61277=Nz   2:2=-0.00083=Ny  2:3=-0.79026=Nx
; 3:1= 0.00000=Ez   3:2= 1.00000=Ey  3:3=-0.00105=Ex
;
; 1:1= 0.79026=Uz   2:1= 0.61277=Nz  3:1= 0.00000=Ez
; 1:2= 0.00064=Uy   2:2=-0.00083=Ny  3:2= 1.00000=Ey
; 1:3= 0.61277=Ux   2:3=-0.79026=Nx  3:3=-0.00105=Ex
;

LBL *ROTX        ;         |  1     0     0     | 
STO 00            ;        |                    |        
3.003            ; Rx(T) =     |  0  COS(T) -SIN(T) |    
MATDIM            ;         |                    |
0            ;         |  0  SIN(T)  COS(T) |
MSIJ
1
MSC+    ; MSR+    ; OOOPS
0
MSC+    ; MSR+
MSC+    ; MSR+
MSC+     ; MSR+
RCL 00
ENTER
SIN
X<>Y
COS
MSC+    ; MSR+
X<>Y
MSC+     ; MSR+
0
MSC+     ; MSR+
RDN
CHS
MSC+    ; MSR+
X<>Y
MS
RTN            ; VERIFIED 9/17/16


LBL *ROTZ        ;        | COS(T) -SIN(T)  0  |
STO 00            ;        |                    |    
3.003            ; Rz(T) =       | SIN(T)  COS(T)  0  |
MATDIM            ;        |                    |    
0            ;        |   0       0     1  |
MSIJ
RCL 00
ENTER
SIN
X<>Y
COS
MSC+    ;    MSR+
X<>Y
MSC+    ;    MSR+
0
MSC+    ;    MSR+
RDN
CHS
MSC+    ;    MSR+
X<>Y
MSC+    ;    MSR+
0
MSC+    ;    MSR+
MSC+    ;    MSR+
MSC+    ;    MSR+
1
MS
RTN            ; VERIFIED 9/17/2016

LBL *ROTY        ;        | COS(T)   0   SIN(T) |
STO 00            ;        |                     |
3.003            ; Ry(T) =       |    0     1     0    |
MATDIM            ;         |                     |
0            ;               |-SIN(T)   0   COS(T) |
MSIJ
RCL 00
ENTER
SIN
X<>Y
COS
MSC+
0
MSC+
X<>Y
MSC+
X<>Z
MSC+
1
MSC+
RDN
MSC+
X<>Z
MSC+
RDN
MS
RTN

LBL *ROT1        ;         |  1     0     0     | 
STO 00            ;        |                    |        
3.003            ; R1(T) =     |  0  COS(T)  SIN(T) |    
MATDIM            ;         |                    |
0            ;         |  0 -SIN(T)  COS(T) |
MSIJ
1
MSC+    ;
0
MSC+    ; 
MSC+    ; 
MSC+     ; 
RCL 00
ENTER
SIN
CHS
X<>Y
COS
MSC+    ; 
X<>Y    
MSC+     ; 
0
MSC+     ; 
RDN
CHS
MSC+    ; 
X<>Y
MS
RTN    

LBL *ROT3        ;        | COS(T)  SIN(T)  0  |
STO 00            ;        |                    |    
3.003            ; R3(T) =       | -SIN(T) COS(T)  0  |
MATDIM            ;        |                    |    
0            ;        |   0       0     1  |
MSIJ
RCL 00
ENTER
SIN
CHS
X<>Y
COS
MSC+    ;    
X<>Y
MSC+    ;
0
MSC+    ;    
RDN
CHS
MSC+    ;    
X<>Y
MSC+    ;    
0
MSC+    ;    
MSC+    ;    
MSC+    ;    
1
MS
RTN
;
;
;            
LBL *ROT2        ;        | COS(T)   0   -SIN(T) |
STO 00            ;        |                      |
3.003            ; R2(T) =       |    0     1     0     |
MATDIM            ;         |                      |
0            ;               |SIN(T)   0   COS(T)   |
MSIJ
RCL 00
ENTER
SIN
;;  CHS
X<>Y
COS
MSC+
0
MSC+
X<>Z
MSC+
X<>Z
MSC+
1
MSC+
RDN
MSC+
X<>Z
CHS
MSC+
X<>Z
MSC+
RDN
MS
RTN

;
; RE = 6378.137
; RP = 6256.752
;
; D  = SQRT(RE*RE*CL*CL + RP*RP*SL*SL)
; Rx = RE*RE/D + HT: Rz = RP*RP/D + HT
;
LBL *RADIUS
6378.137        ; RADIUS OF EARTH AT EQUATOR 
X^2            ; X=XX
STO 07            ; [OX]
STO 08            ; [OY]
6256.752        ; RADIUS OF EARTH THRU POLES 
X^2            ; 
STO 09            ; [OZ]          
RCL 30            ; [OBSLAT]
COS
X^2
RCL 07
*
RCL 30
SIN
X^2
RCL 09
+
SQRT
ST/ 07
ST/ 08
ST/ 09
RCL 40            ; [OBSHGT]
1000
/
ST+ 07            ; [OX]
ST+ 08            ; [OY]
ST+ 09            ; [OZ]
RTN
LBL *OBSXYZ
"ENU"
3.001
MSIJ
MRR+            ; UX
ST* 07            ; OPOS,X
MRR+            ; UY
ST* 08            ; OPOS,Y
MRR+            ; UZ
ST* 09            ; OPOS,Z
RTN
LBL *OBSVEL
RCL 08
STO 10
RCL 07
STO 11
CLX
STO 12
PI                 ; ERROR!!! CALC IN DEGREES


365.2421874 


ST* 11 


86400 
ST* 11 
CHS 
ST* 10 
RTN 
; Calculate Satellite Position at DN,TN
; T  = (DN - DE) + (TN-TE) Elapsed T since epoch, days

LBL *MA+ORB
RCL 19                ; EPYEAR
1E6
/
1.01
+
DATE
DDAYS
RCL 20                ; EPTIME
INT

TIME
24
/
RCL 20
FRC
-
+
STO 16                ; ELAPST SINCE EPOCH, DAYS
; DRAG COEFF = -2*M2/MM/3: (Angular momentum rate)/(Ang mom)  s^-1
RCL 27                ; MMDOT2 
-2
*
RCL 26                ; MMDOT1
/
3
/
*
; LINEAR DRAG TERMS KD = 1+4*DT:KDP= 1-7*DT
; Mean anomaly at YR,TN M = MA + MM*T*(1-3*DT)
STO 40                ; MEANAT
STO 41                ; KD
STO 42                ; KDP

ST* 41                ; KD
-7
ST* 42                ; KDP
-3
ST* 40                ; MEANAT

ST+ 41                ; KD
ST+ 42                ; KDP
ST+ 40                ; MEANAT
RCL 16                ; ELAPST
ST* 40                ; MEANAT
RCL 25                ; MMDOT1
ST* 40                ; MEANAT
RCL 24                ; MEANAN
ST+ 40                ; MEANAT
RCL 40                ; MEANAT
RTN
STOP
360
/
INT
ENTER
360
*
ST- 40                ; MEANAT
X<>Y
RCL 28                ; ORBITN
+
STO 17                ; CURORB
RTN

 2100 EA = M:                 REM Initial solution
 2110 REPEAT
 2120   C = COS(EA): S = SIN(EA): DNOM=1-EC*C
 2130   D = (EA-EC*S-M)/DNOM:   REM Change to EA for better solution
 2140   EA = EA - D:            REM by this amount
 2150 UNTIL ABS(D) < 1E-5:      REM Until converged

LBL *SOLVEA
RCL 40
STO 01
STO 02
"*KEPLE"
STO 06
XEQ SOL                ; IMPORTANT - MY HP41CX HAS BOTH MATH/STAT AND 
                ; ADVANTAGE.  REMOVE MATH MODULE WHEN ENTERING
                ; INSTRUCTION.  SHOULD BE XROM XX

                ; DISREGARD USING MATH SOL
RTN
LBL "*KEPLE"
;
;
LBL *KEPLER
;
;
; M = EA - EC*SIN(EA)
ENTER
SIN
X<>Y
RCL 23 ; ECCENTRICITY
*
-
RCL 40 ; MEANAT 
-
RTN
;
; ADVANTAGE <SOLVE> FUNCTION CONSUMES 13 REGISTERS
; MATH      <SOLVE> FUNCTION CONSUMES 7  REGISTERS
; NOTE: CHECK WHEN VECTORS COMPUTED

;    RCL 41                ; [M]
;    STO 32                ; [LAMBDA]
;    
;    LBL 60
;    RCL 32                ; [LAMBDA]
;    RCL 32                ; [LAMBDA]
;    SIN
;    RCL 23                ; [ECCENT]
;    -
;    RCL 41
;    -
;    1
;    RCL 23                ; [ECCENT]
;    RCL 32                ; [LAMBDA]
;    COS
;    *
;    -
;    /
;    ST- 32                ; [LAMBDA]
;    ABS
;    1E-5
;    X>Y?
;    GTO 60
;    RTN
;
; Cambridge UK 52.2053? N, 0.1218? E
;
;
LBL **TEST
1990 
STO 19                ; EPYEAR
191.145409
STO 20                ; EPOTIM
59.9975
STO 21                ; INCLIN
146.4527
STO 22                ; RAAN
0.6986        
STO 23                ; ECCENT
231.0027
STO 24                ; PERIGE
43.2637 
STO 25                ; MEANAN
2.09695848
STO 26                ; MMDOT1
1E-8 
STO 27                ; MMDOT2
1585
STO 29                ; ORBITN
52.21
STO 30
0.06
STO 31
79
STO 39
RTN
Find all posts by this user
Quote this message in a reply
Post Reply 


Messages In This Thread
RE: (41C) Tracking Satellites in Elliptical Orbits - twdeckard - 11-23-2021 09:08 PM



User(s) browsing this thread: 2 Guest(s)