Re: SDP4 failure example

William Penner (bpenner@ix.netcom.com)
Sun, 30 Mar 1997 11:27:55 -0800

Bruno Tilgner wrote:
> 
> Well, it seems to depend on the implementation of the algorithm. The
following
> is a snippet from TRAKSAT 2.8 with some unnecessary numbers removed:
> 
> 
>     Date      Time (UTC)       Lat    Long    Alt   
>                HR:MN:SEC        Deg    Deg     km     
>  --------------------------------------------------                      
      
>  Tue  10Jan95  21:40:00.0      3.29    20.77  36510  
>  Tue  10Jan95  21:40:30.0      3.29    20.71  36493  
>  Tue  10Jan95  21:41:00.0      3.29    20.65  36476  
>  Tue  10Jan95  21:41:30.0      3.29    20.59  36459  
>  Tue  10Jan95  21:42:00.0      3.29    20.53  36441  
>  Tue  10Jan95  21:42:30.0      3.29    20.47  36424  
>  Tue  10Jan95  21:43:00.0      3.28    20.41  36406  
>  Tue  10Jan95  21:43:30.0      3.28    20.35  36388  
>  Tue  10Jan95  21:44:00.0      3.28    20.29  36370  
>  Tue  10Jan95  21:44:30.0      3.28    20.23  36352  
>  Tue  10Jan95  21:45:00.0      3.28    20.17  36334  
>  Tue  10Jan95  21:45:30.0      3.28    20.11  36315  <--
>  Tue  10Jan95  21:46:00.0      3.28    20.64  36297  <--
>  Tue  10Jan95  21:46:30.0      3.28    20.58  36278  
>  Tue  10Jan95  21:47:00.0      3.28    20.52  36259  
>  Tue  10Jan95  21:47:30.0      3.28    20.46  36240  
>  Tue  10Jan95  21:48:00.0      3.28    20.40  36221  
>  Tue  10Jan95  21:48:30.0      3.28    20.34  36202  
>  Tue  10Jan95  21:49:00.0      3.28    20.28  36182  
>  
> As can be seen, the (EAST) longitude of the subsatellite point is
decreasing
> until 21:45:30 by 0.06 degrees per 30 seconds. Therefore, at 21:46:00 it
> should be 20.05 deg but it is 20.64 deg, i.e. there is a jump of 0.59 deg
and
> in the opposite direction. This is clearly wrong.
> 
> 
> By contrast, my own program does not exhibit this jump (I am using WEST
> longitudes, hence the minus sign):
> 
> 
> Date    Time (UTC)   Geocentric    Subsatellite Point
>                      RA     Dec      Lat   West Long
> ----------------------------------------------------------------------
> 10Jan95 21:40:00    +06:23 +03:16   +3.26    -20.77
> 10Jan95 21:40:30    +06:23 +03:16   +3.26    -20.71
> 10Jan95 21:41:00    +06:23 +03:16   +3.26    -20.65
> 10Jan95 21:41:30    +06:24 +03:16   +3.26    -20.59
> 10Jan95 21:42:00    +06:24 +03:15   +3.26    -20.53
> 10Jan95 21:42:30    +06:24 +03:15   +3.26    -20.47
> 10Jan95 21:43:00    +06:24 +03:15   +3.26    -20.41
> 10Jan95 21:43:30    +06:25 +03:15   +3.26    -20.35
> 10Jan95 21:44:00    +06:25 +03:15   +3.26    -20.29
> 10Jan95 21:44:30    +06:25 +03:15   +3.26    -20.23
> 10Jan95 21:45:00    +06:25 +03:15   +3.26    -20.17
> 10Jan95 21:45:30    +06:26 +03:15   +3.26    -20.11 <--
> 10Jan95 21:46:00    +06:26 +03:15   +3.26    -20.05 <--
> 10Jan95 21:46:30    +06:26 +03:15   +3.26    -19.99
> 10Jan95 21:47:00    +06:26 +03:15   +3.26    -19.93
> 10Jan95 21:47:30    +06:27 +03:15   +3.26    -19.87
> 10Jan95 21:48:00    +06:27 +03:15   +3.26    -19.81
> 10Jan95 21:48:30    +06:27 +03:15   +3.26    -19.75
> 10Jan95 21:49:00    +06:27 +03:15   +3.26    -19.69
> 
> 
> My implementation of the SDP4 code is a translation to Microsoft
Professional
> Development System 7.1 (which is a close relative of MS QuickBasic 4.5)
from
> T.S. Kelso's PASCAL code which in itself is a translation from FORTRAN.
The
> deep space routines of the PASCAL code reflect some FORTRAN features, in
> particular calculated return addresses from GOSUB's, and these constructs
> can easily be misinterpreted or overlooked when converting the FORTRAN
code
> to other languages.
> 
> Apart from the jump, TRAKSAT and my own program agree very well in the
> longitude of the subsatellite point, but there is a small discrepancy in
the
> latitude for which I have no explanation. I take the Earth's flattening
into
> account, although not with very high precision, and perhaps TRAKSAT does
not. 
> 
> At any rate, the comparision of results from different programs which
have
> undergone several coordinate transformations is problematic. More direct
> comparisons could be made with the cartesian coordinates which the
various
> NORAD models produce as their native output.
> 
> Finally, a consolation for the users of satellite prediction software:
The deep
> space model is only used for satellites with a period of more than 225
minutes
> i.e. less than 6.4 revolutions per day. Most observers will be concerned
with
> faster moving objects, hence the flaws in the SDP4 code are of no
practical
> consequences for them, but program developers should of course pay
attention
> to this problem. 
> 
> 
> Bruno Tilgner
> Paris, France
> 

Well, this is very interesting.  I have also converted the SDP4 FORTRAN
routines to C++ and have similar results the TRAKSAT 2.8 numbers (except
the Subsatellite Latitude which is closer to Bruno Tilgner's results). 
Here are the numbers that I generated in my implementation:

Satellite: 92021N 23243/92-021N  

Time        Lat      Long       Alt
 UTC         N         W        km

>>>>> Date: 01/10/1995 -- Orbit: 973 <<<<<

21:40:00    3.265   -20.771    36510
21:40:30    3.265   -20.710    36493
21:41:00    3.264   -20.649    36476
21:41:30    3.264   -20.589    36458
21:42:00    3.264   -20.528    36441
21:42:30    3.263   -20.467    36423
21:43:00    3.263   -20.407    36406
21:43:30    3.262   -20.346    36388
21:44:00    3.262   -20.286    36370
21:44:30    3.262   -20.225    36352
21:45:00    3.261   -20.165    36334
21:45:30    3.261   -20.105    36315
21:46:00    3.260   -20.044    36297
21:46:30    3.260   -19.984 <--36278--- The infamous
21:47:00    3.255   -20.515 <--36259--- SDP4 step
21:47:30    3.255   -20.455    36240
21:48:00    3.254   -20.395    36221
21:48:30    3.254   -20.335    36202
21:49:00    3.253   -20.275    36182
21:49:30    3.253   -20.216    36163

Since various implementations of the same routine yield similar results,
and since it is unlikely that we all mis-understand FORTRAN in the same
way, I would start looking at differences between the PASCAL and FORTRAN
implementations.

Bill Penner