home *** CD-ROM | disk | FTP | other *** search
/ rtsi.com / 2014.01.www.rtsi.com.tar / www.rtsi.com / OS9 / OSK / TELECOM / OSKBox.lzh / MAILBOX / SAT / orbit.c < prev    next >
Text File  |  1990-08-27  |  18KB  |  671 lines
  1.  
  2. /* Copyright (c) 1986,1987,1988,1989,1990 Robert W. Berger N3EMO
  3.    May be freely distributed, provided this notice remains intact. */
  4.  
  5. /* Change Log
  6.     3/7/1990    v3.7 Make Phase III style phase (0-255) the default.
  7.             Ignore case in satellite names.
  8.  
  9.     12/19/1989    v3.6 Use more direct calculations for dates.
  10.             Calculate a new sidereal time reference for each run.
  11.  
  12.     12/8/1988    v3.5 Allow multiple overlapping modes in "mode.dat".
  13.  
  14.     6/28/1988    v3.4 Cleaned up Eclipse code. Fixed leap year handling
  15.             for centesimal years. Added a heuristic to GetDay to
  16.             allow 2 or 4 digit year specifications.
  17.                   
  18.     1/25/1988    v3.2 Rewrote orbitr.c to improve modularity,
  19.             efficiency, and accuracy. Adopted geocentric
  20.             cartesian coordinates as the standard representation
  21.             for position and velocity. Added direct calculation
  22.              of range-rate for better doppler predections.
  23.  
  24.     12/1/1988    v3.1 Allow spaces in satellite names. Provide
  25.             single character aliases for 62 satellites 
  26.             (up from 26).
  27.  
  28.     4/7/87        v3.0 Added eclipses.
  29.  
  30.     4/1/87        v2.4 Added "Flip" option in site file for
  31.             0-180 elevation support.
  32.  
  33.     3/24/87        v2.3 Adapted for new kepler.dat format.
  34.             Allow beacon frequencies in mode.dat.
  35.             Use decay rate for drag compensation.
  36.  
  37.     5/10/86        v2.2 Added single character aliases for satellite
  38.             names.
  39.  
  40.     4/30/86        v2.1 Print blank line if satellite dips below
  41.             horizon and reappears during same orbit and day
  42.  
  43.     4/29/86        v2.0 Changed GetSatelliteParams() to use AMSAT's
  44.             "kepler.dat" file. Moved schedule to "mode.dat" file.
  45.  
  46.         4/22/86         v1.3  Inserted N8FJB's suggestions for variable naming
  47.                         which maintain 8 character uniqueness.
  48.                         Also removed "include" file orbitr.h, which had two
  49.                         definitions of external functions defined in orbit.c
  50.                 -K3MC 
  51.  
  52.         4/1/86          v1.2  Corrected a scanf conversion to %d for an int
  53.                         type.    -K3MC
  54.  
  55.         3/19/86         v1.1  Changed GetSatelliteParams to not pass NULL
  56.                         to sscanf.
  57.                                                                         */
  58.  
  59. #define DRAG 1
  60.  
  61. #include <stdio.h>
  62. #include <math.h>
  63. #include <ctype.h>
  64.  
  65. extern double Kepler();
  66. extern long GetDayNum();
  67.  
  68. #define LC(c) (isupper(c) ? tolower(c) : (c))
  69.  
  70. #ifndef PI
  71. #define PI 3.14159265
  72. #endif
  73.  
  74. #ifdef PI2
  75. #undef PI2
  76. #endif
  77.  
  78. #define PI2 (PI*2)
  79.  
  80. #define MinutesPerDay (24*60.0)
  81. #define SecondsPerDay (60*MinutesPerDay)
  82. #define HalfSecond (0.5/SecondsPerDay)
  83. #define EarthRadius 6378.16             /* Kilometers           */
  84. #define C 2.997925e5                    /* Kilometers/Second    */
  85. #define TropicalYear 365.24199        /* Mean solar days    */
  86. #define EarthEccentricity 0.016713
  87. #define DegreesPerRadian (180/PI)
  88. #define RadiansPerDegree (PI/180)
  89. #define ABS(x) ((x) < 0 ? (-(x)) : (x))
  90. #define SQR(x) ((x)*(x))
  91.  
  92. #define MaxModes 10
  93. typedef struct {
  94.                 int MinPhase,MaxPhase;
  95.                 char ModeStr[20];
  96.                }  ModeRec;
  97.  
  98. char VersionStr[] = "N3EMO Orbit Simulator  v3.7";
  99.  
  100.     /*  Keplerian Elements and misc. data for the satellite              */
  101.     double  EpochDay;                   /* time of epoch                 */
  102.     double EpochMeanAnomaly;            /* Mean Anomaly at epoch         */
  103.     long EpochOrbitNum;                 /* Integer orbit # of epoch      */
  104.     double EpochRAAN;                   /* RAAN at epoch                 */
  105.     double epochMeanMotion;             /* Revolutions/day               */
  106.     double OrbitalDecay;                /* Revolutions/day^2             */
  107.     double EpochArgPerigee;             /* argument of perigee at epoch  */
  108.     double Eccentricity;
  109.     double Inclination;
  110.     char SatName[100];
  111.     int ElementSet;
  112.     double BeaconFreq;                  /* Mhz, used for doppler calc    */
  113.     double MaxPhase;                    /* Phase units in 1 orbit        */
  114.     double perigeePhase;
  115.     int NumModes;
  116.     ModeRec Modes[MaxModes];
  117.     int PrintApogee;
  118.     int PrintEclipses;
  119.     int Flip;
  120.  
  121.     /* Simulation Parameters */
  122.  
  123.     double StartTime,EndTime, StepTime; /* In Days, 1 = New Year        */
  124.                                         /*      of reference year       */
  125.  
  126.     /* Site Parameters */
  127.     char SiteName[100];
  128.     double SiteLat,SiteLong,SiteAltitude,SiteMinElev;
  129.  
  130.  
  131. /* List the satellites in kepler.dat, and return the number found */
  132. ListSatellites()
  133. {
  134.     char str[100];
  135.     FILE *InFile;
  136.     char satchar;
  137.     int NumSatellites;
  138.  
  139.     printf("Available satellites:\n");
  140.  
  141.     if ((InFile = fopen("kepler.dat","r")) == 0)
  142.         {
  143.     printf("\"kepler.dat\" not found\n");
  144.     exit(-1);
  145.     }
  146.  
  147.     satchar = 'a';
  148.     NumSatellites = 0;
  149.     while (fgets(str,100,InFile))
  150.     if (strncmp(str,"Satellite: ",11) == 0)
  151.         {
  152.         printf("    %c) %s",satchar,&str[11]);
  153.         if (satchar == 'z')
  154.         satchar = 'A';
  155.                else if (satchar == 'Z')
  156.            satchar = '0';
  157.                 else satchar++;
  158.         NumSatellites++;
  159.         }
  160.  
  161.     fclose(InFile);
  162.  
  163.     return NumSatellites;
  164. }
  165.  
  166. /* Match and skip over a string in the input file. Exits on failure. */
  167.  
  168. MatchStr(InFile,FileName,Target)
  169. FILE *InFile;
  170. char *FileName,*Target;
  171. {
  172.     char str[100];
  173.  
  174.     fgets(str,strlen(Target)+1,InFile);
  175.     if (strcmp(Target,str))
  176.        {
  177.        printf("%s: found \"%s\" while expecting \"%s\n\"",FileName,str,Target);
  178.        exit(-1);
  179.        }
  180. }
  181.  
  182. LetterNum(c)
  183. char c;
  184. {
  185.     if (c >= 'a' && c <= 'z')
  186.     return c - 'a' + 1;
  187.       else if (c >= 'A' && c <= 'Z')
  188.        return c - 'A'+ 27;
  189.     else if (c >= '0' && c <= '9')
  190.       return c - '0' + 53;
  191. }
  192.       
  193. /* Case insensitive strncmp */
  194. cstrncmp(str1,str2,l)
  195. char *str1,*str2;
  196. {
  197.     int i;
  198.  
  199.     for (i = 0; i < l; i++)
  200.     if (LC(str1[i]) != LC(str2[i]))
  201.         return 1;
  202.  
  203.     return 0;
  204. }
  205.  
  206.  
  207. cstrcmp(str1,str2)
  208. char *str1,*str2;
  209. {
  210.     int i,l;
  211.  
  212.     l = strlen(str1);
  213.     if (strlen(str2) != l)
  214.     return 1;
  215.  
  216.     for (i = 0; i < l; i++)
  217.     if (LC(str1[i]) != LC(str2[i]))
  218.         return 1;
  219.  
  220.     return 0;
  221. }
  222.  
  223.  
  224. GetSatelliteParams()
  225. {
  226.     FILE *InFile;
  227.     char str[100];
  228.     int EpochYear;
  229.     double EpochHour,EpochMinute,EpochSecond;
  230.     int found;
  231.     int i,NumSatellites;
  232.     char satchar;
  233.  
  234.     NumSatellites = ListSatellites();
  235.  
  236.     found = 0;
  237.  
  238.     while (!found)
  239.     {
  240.     printf("Letter or satellite name :");
  241.     gets(SatName);
  242.  
  243.     if ((InFile = fopen("kepler.dat","r")) == 0)
  244.         {
  245.         printf("kepler.dat not found\n");
  246.         exit(-1);
  247.         }
  248.  
  249.     if (strlen(SatName) == 1)
  250.         {            /* use single character label */
  251.         satchar = SatName[0];
  252.         if (LetterNum(satchar) > NumSatellites)
  253.             {
  254.             printf("'%c' is out of range\n",satchar);
  255.         fclose(InFile);
  256.         continue;
  257.         }
  258.  
  259.         for (i = 1; i <= LetterNum(satchar); i++)
  260.         {
  261.         do  /* find line beginning with "Satellite: " */
  262.             fgets(str,100,InFile);
  263.         while (strncmp(str,"Satellite: ",11) != 0);
  264.         }
  265.         found = 1;
  266.         strncpy(SatName,&str[11],strlen(str)-12);
  267.         }
  268.         
  269.      else 
  270.          {
  271.          while (!found)  /* use satellite name */
  272.                 {
  273.             if (! fgets(str,100,InFile))
  274.                 break;    /* EOF */
  275.  
  276.             if (strncmp(str,"Satellite: ",11) == 0)
  277.            if (cstrncmp(SatName,&str[11],strlen(SatName)) == 0)
  278.             found = 1;
  279.             }
  280.  
  281.         if (!found)
  282.         {
  283.         printf("Satellite %s not found\n",SatName);
  284.         fclose(InFile);
  285.         }
  286.         }
  287.     }
  288.  
  289.     BeaconFreq = 146.0;  /* Default value */
  290.  
  291.     fgets(str,100,InFile);    /* Skip line */
  292.  
  293.     MatchStr(InFile,"kepler.dat","Epoch time:");
  294.     fgets(str,100,InFile);
  295.     sscanf(str,"%lf",&EpochDay);
  296.  
  297.     EpochYear = EpochDay / 1000.0;
  298.     EpochDay -= EpochYear*1000.0;
  299.     EpochDay += GetDayNum(EpochYear,1,0);
  300.     fgets(str,100,InFile);
  301.  
  302.     if (sscanf(str,"Element set: %ld",&ElementSet) == 0)
  303.        {   /* Old style kepler.dat */
  304.        MatchStr(InFile,"kepler.dat","Element set:");
  305.        fgets(str,100,InFile);
  306.        sscanf(str,"%d",&ElementSet);
  307.        }
  308.  
  309.     MatchStr(InFile,"kepler.dat","Inclination:");
  310.     fgets(str,100,InFile);
  311.     sscanf(str,"%lf",&Inclination);
  312.     Inclination *= RadiansPerDegree;
  313.  
  314.     MatchStr(InFile,"kepler.dat","RA of node:");
  315.     fgets(str,100,InFile);
  316.     sscanf(str,"%lf",&EpochRAAN);
  317.     EpochRAAN *= RadiansPerDegree;
  318.  
  319.     MatchStr(InFile,"kepler.dat","Eccentricity:");
  320.     fgets(str,100,InFile);
  321.     sscanf(str,"%lf",&Eccentricity);
  322.  
  323.     MatchStr(InFile,"kepler.dat","Arg of perigee:");
  324.     fgets(str,100,InFile);
  325.     sscanf(str,"%lf",&EpochArgPerigee);
  326.     EpochArgPerigee *= RadiansPerDegree;
  327.  
  328.     MatchStr(InFile,"kepler.dat","Mean anomaly:");
  329.     fgets(str,100,InFile);
  330.     sscanf(str,"%lf",&EpochMeanAnomaly);
  331.     EpochMeanAnomaly *= RadiansPerDegree;
  332.  
  333.     MatchStr(InFile,"kepler.dat","Mean motion:");
  334.     fgets(str,100,InFile);
  335.     sscanf(str,"%lf",&epochMeanMotion);
  336.  
  337.     MatchStr(InFile,"kepler.dat","Decay rate:");
  338.     fgets(str,100,InFile);
  339.     sscanf(str,"%lf",&OrbitalDecay);
  340.  
  341.     MatchStr(InFile,"kepler.dat","Epoch rev:");
  342.     fgets(str,100,InFile);
  343.     sscanf(str,"%ld",&EpochOrbitNum);
  344.  
  345.     while (1)
  346.         {
  347.         if (! fgets(str,100,InFile))
  348.         break;    /* EOF */
  349.         if (strlen(str) <= 2)
  350.             break;  /* Blank line */
  351.         sscanf(str,"Beacon: %lf",&BeaconFreq);
  352.         }
  353.  
  354.     PrintApogee = (Eccentricity >= 0.3);
  355.  
  356.     perigeePhase = 0; MaxPhase = 256; /* Default values */
  357.     NumModes = 0;
  358.  
  359.     if ((InFile = fopen("mode.dat","r")) == 0)
  360.     return;
  361.  
  362.     found = 0;
  363.     while (!found)
  364.         {
  365.     if (! fgets(str,100,InFile))
  366.         break;    /* EOF */
  367.     if (sscanf(str,"Satellite: %s",str) == 1
  368.         && cstrcmp(SatName,str) == 0)
  369.         found = 1;
  370.     }
  371.     
  372.     if (found)
  373.     {
  374.     while (1)
  375.         {
  376.         if (! fgets(str,100,InFile))
  377.         break;    /* EOF */
  378.         if (strlen(str) <= 2)
  379.             break;  /* Blank line */
  380.         sscanf(str,"Beacon: %lf",&BeaconFreq);
  381.         sscanf(str,"Perigee phase: %lf",&perigeePhase);
  382.         sscanf(str,"Max phase: %lf",&MaxPhase);
  383.  
  384.         if (sscanf(str,"Mode: %20s from %d to %d",Modes[NumModes].ModeStr,
  385.         &Modes[NumModes].MinPhase,&Modes[NumModes].MaxPhase) == 3
  386.           && NumModes < MaxModes)
  387.           NumModes++;
  388.         }
  389.     fclose(InFile);
  390.     }
  391. }
  392.  
  393.  
  394. GetSiteParams()
  395. {
  396.     FILE *InFile;
  397.     char name[100],str[100];
  398.  
  399.     printf("Site name :");
  400.     gets(name);
  401.     strcat(name,".sit");
  402.  
  403.     if ((InFile = fopen(name,"r")) == 0)
  404.         {
  405.         printf("%s not found\n",name);
  406.         exit(-1);
  407.         }
  408.  
  409.     fgets(SiteName,100,InFile);
  410.  
  411.     fgets(str,100,InFile);
  412.     sscanf(str,"%lf",&SiteLat);
  413.     SiteLat *= RadiansPerDegree;
  414.  
  415.     fgets(str,100,InFile);
  416.     sscanf(str,"%lf",&SiteLong);
  417.     SiteLong *= RadiansPerDegree;
  418.  
  419.     fgets(str,100,InFile);
  420.     sscanf(str,"%lf",&SiteAltitude);
  421.     SiteAltitude /= 1000;   /* convert to km */
  422.  
  423.     fgets(str,100,InFile);
  424.     sscanf(str,"%lf",&SiteMinElev);
  425.     SiteMinElev *= RadiansPerDegree;
  426.  
  427.     Flip = PrintEclipses = 0;
  428.     while (fgets(str,100,InFile))
  429.     {
  430.     if (strncmp(str,"Flip",4) == 0)
  431.         Flip = 1;
  432.       else if (strncmp(str,"Eclipse",7) == 0)
  433.         PrintEclipses = 1;
  434.        else printf("\"%s\" unknown option: %s",name,str);
  435.     }
  436. }
  437.  
  438. GetSimulationParams()
  439. {
  440.     double hour,duration;
  441.     int Month,Day,Year;
  442.  
  443.     printf("Start date (UTC) (Month/Day/Year) :");
  444.     scanf("%d/%d/%d",&Month,&Day,&Year);
  445.  
  446.     StartTime = GetDayNum(Year,Month,Day);
  447.     printf("Starting Hour (UTC) :");
  448.     scanf("%lf",&hour);
  449.     StartTime += hour/24;
  450.  
  451.     printf("Duration (Days) :");
  452.     scanf("%lf",&duration);
  453.     EndTime = StartTime + duration;
  454.  
  455.     printf("Time Step (Minutes) :");
  456.     scanf("%lf",&StepTime);
  457.     StepTime /= MinutesPerDay;
  458. }
  459.  
  460. PrintMode(OutFile,Phase)
  461. FILE *OutFile;
  462. {
  463.     int CurMode;
  464.  
  465.     for (CurMode = 0; CurMode < NumModes; CurMode++)
  466.         if ((Phase >= Modes[CurMode].MinPhase
  467.                 && Phase < Modes[CurMode].MaxPhase)
  468.               || ((Modes[CurMode].MinPhase > Modes[CurMode].MaxPhase)
  469.                   && (Phase >= Modes[CurMode].MinPhase
  470.                         || Phase < Modes[CurMode].MaxPhase)))
  471.             {
  472.             fprintf(OutFile,"%s ",Modes[CurMode].ModeStr);
  473.             }
  474. }
  475.  
  476.  
  477. main()
  478. {
  479.     double ReferenceOrbit;      /* Floating point orbit # at epoch */
  480.     double CurrentTime;         /* In Days                         */
  481.     double CurrentOrbit;
  482.     double AverageMotion,       /* Corrected for drag              */
  483.         CurrentMotion;
  484.     double MeanAnomaly,TrueAnomaly;
  485.     double SemiMajorAxis;
  486.     double Radius;              /* From geocenter                  */
  487.     double SatX,SatY,SatZ;    /* In Right Ascension based system */
  488.     double SatVX,SatVY,SatVZ;   /* Kilometers/second           */
  489.     double SiteX,SiteY,SiteZ;
  490.     double SiteVX,SiteVY;
  491.     double SiteMatrix[3][3];
  492.     double Height;
  493.     double RAANPrecession,PerigeePrecession;
  494.     double SSPLat,SSPLong;
  495.     long OrbitNum,PrevOrbitNum;
  496.     long Day,PrevDay;
  497.     double Azimuth,Elevation,Range;
  498.     double RangeRate,Doppler;
  499.     int Phase;
  500.     char FileName[100];
  501.     FILE *OutFile;
  502.     int DidApogee;
  503.     double TmpTime,PrevTime;
  504.     int PrevVisible;
  505.  
  506.     printf("%s\n",VersionStr);
  507.  
  508.  
  509.     GetSatelliteParams();
  510.     GetSiteParams();
  511.     GetSimulationParams();
  512.  
  513.     InitOrbitRoutines((StartTime+EndTime)/2);
  514.  
  515.     printf("Output file (RETURN for TTY) :");
  516.     gets(FileName);     /* Skip previous RETURN */
  517.     gets(FileName);
  518.  
  519.  
  520.     if (strlen(FileName) > 0)
  521.         {
  522.         if ((OutFile = fopen(FileName,"w")) == 0)
  523.             {
  524.             printf("Can't write to %s\n",FileName);
  525.             exit(-1);
  526.             }
  527.         }
  528.       else OutFile = stdout;
  529.  
  530.     fprintf(OutFile,"%s Element Set %d\n",SatName,ElementSet);
  531.  
  532.     fprintf(OutFile,"%s\n",SiteName);
  533.  
  534.     fprintf(OutFile,"Doppler calculated for freq = %lf MHz\n",BeaconFreq);
  535.  
  536.     SemiMajorAxis = 331.25 * exp(2*log(MinutesPerDay/epochMeanMotion)/3);
  537.     GetPrecession(SemiMajorAxis,Eccentricity,Inclination,&RAANPrecession,
  538.                         &PerigeePrecession);
  539.  
  540.     ReferenceOrbit = EpochMeanAnomaly/PI2 + EpochOrbitNum;
  541.  
  542.     PrevDay = -10000; PrevOrbitNum = -10000;
  543.     PrevTime = StartTime-2*StepTime;
  544.  
  545.     BeaconFreq *= 1E6;          /* Convert to Hz */
  546.  
  547.     DidApogee = 0;
  548.  
  549.     for (CurrentTime = StartTime; CurrentTime <= EndTime;
  550.                 CurrentTime += StepTime)
  551.         {
  552.  
  553.         AverageMotion = epochMeanMotion
  554.        + (CurrentTime-EpochDay)*OrbitalDecay/2;
  555.         CurrentMotion = epochMeanMotion
  556.        + (CurrentTime-EpochDay)*OrbitalDecay;
  557.  
  558.         SemiMajorAxis = 331.25 * exp(2*log(MinutesPerDay/CurrentMotion)/3);
  559.  
  560.         CurrentOrbit = ReferenceOrbit +
  561.                         (CurrentTime-EpochDay)*AverageMotion;
  562.         OrbitNum = CurrentOrbit;
  563.  
  564.         MeanAnomaly = (CurrentOrbit-OrbitNum)*PI2;
  565.  
  566.         TmpTime = CurrentTime;
  567.         if (MeanAnomaly < PI)
  568.             DidApogee = 0;
  569.         if (PrintApogee && !DidApogee && MeanAnomaly > PI)
  570.             {                   /* Calculate Apogee */
  571.             TmpTime -= StepTime;   /* So we pick up later where we left off */
  572.             MeanAnomaly = PI;
  573.             CurrentTime=EpochDay+(OrbitNum-ReferenceOrbit+0.5)/AverageMotion;
  574.             }
  575.  
  576.         TrueAnomaly = Kepler(MeanAnomaly,Eccentricity);
  577.  
  578.     GetSatPosition(EpochDay,EpochRAAN,EpochArgPerigee,SemiMajorAxis,
  579.         Inclination,Eccentricity,RAANPrecession,PerigeePrecession,
  580.         CurrentTime,TrueAnomaly,&SatX,&SatY,&SatZ,&Radius,
  581.         &SatVX,&SatVY,&SatVZ);
  582.  
  583.     GetSitPosition(SiteLat,SiteLong,SiteAltitude,CurrentTime,
  584.         &SiteX,&SiteY,&SiteZ,&SiteVX,&SiteVY,SiteMatrix);
  585.  
  586.  
  587.     GetBearings(SatX,SatY,SatZ,SiteX,SiteY,SiteZ,SiteMatrix,
  588.         &Azimuth,&Elevation);
  589.  
  590.  
  591.         if (Elevation >= SiteMinElev && CurrentTime >= StartTime)
  592.             {
  593.  
  594.             Day = CurrentTime + HalfSecond;
  595.             if (((double) Day) > CurrentTime+HalfSecond)
  596.                 Day -= 1;    /* Correct for truncation of negative values */
  597.  
  598.         if (OrbitNum == PrevOrbitNum && Day == PrevDay && !PrevVisible)
  599.             fprintf(OutFile,"\n");    /* Dipped out of sight; print blank */
  600.  
  601.             if (OrbitNum != PrevOrbitNum || Day != PrevDay)
  602.                 {                       /* Print Header */
  603.         PrintDayOfWeek(OutFile,(long) Day);
  604.         fprintf(OutFile," ");
  605.                 PrintDate(OutFile,(long) Day);
  606.                 fprintf(OutFile,"  ----Orbit # %ld-----\n",OrbitNum);
  607.                 fprintf(OutFile," U.T.C.   Az  El ");
  608.         if (Flip)
  609.                     fprintf(OutFile," Az'  El' ");
  610.  
  611.         fprintf(OutFile,"Frequency  Range");
  612.                 fprintf(OutFile," Height  Lat  Long  Phase(%3.0lf)\n",
  613.                                 MaxPhase);
  614.                 }
  615.             PrevOrbitNum = OrbitNum; PrevDay = Day;
  616.             PrintTime(OutFile,CurrentTime + HalfSecond);
  617.  
  618.             fprintf(OutFile,"  %3.0lf %3.0lf",Azimuth*DegreesPerRadian,
  619.                 Elevation*DegreesPerRadian);
  620.         if (Flip)
  621.         {
  622.         Azimuth += PI;
  623.         if (Azimuth >= PI2)
  624.             Azimuth -= PI2;
  625.         Elevation = PI-Elevation;
  626.         fprintf(OutFile,"  %3.0lf  %3.0lf",Azimuth*DegreesPerRadian,
  627.             Elevation*DegreesPerRadian);
  628.         }
  629.  
  630.           GetRange(SiteX,SiteY,SiteZ,SiteVX,SiteVY,
  631.         SatX,SatY,SatZ,SatVX,SatVY,SatVZ,&Range,&RangeRate);
  632.             Doppler = -BeaconFreq*RangeRate/C;
  633.             fprintf(OutFile,"  %4.4lf %6.0lf",(BeaconFreq+Doppler)/1.E6,Range);
  634.  
  635.         GetSubSatPoint(SatX,SatY,SatZ,CurrentTime,
  636.             &SSPLat,&SSPLong,&Height);
  637.             fprintf(OutFile," %6.0lf  %3.0lf  %4.0lf",
  638.                 Height,SSPLat*DegreesPerRadian,
  639.                 SSPLong*DegreesPerRadian);
  640.  
  641.             Phase = (MeanAnomaly/PI2*MaxPhase + perigeePhase);
  642.             while (Phase < 0)
  643.                 Phase += MaxPhase;
  644.             while (Phase >= MaxPhase)
  645.                 Phase -= MaxPhase;
  646.  
  647.             fprintf(OutFile," %4d  ", Phase);
  648.             PrintMode(OutFile,Phase);
  649.  
  650.             if (PrintApogee && (MeanAnomaly == PI))
  651.                 fprintf(OutFile,"    Apogee");
  652.  
  653.         if (PrintEclipses)
  654.         if (Eclipsed(SatX,SatY,SatZ,Radius,CurrentTime))
  655.             fprintf(OutFile,"  Eclipse");
  656.  
  657.             fprintf(OutFile,"\n");
  658.         PrevVisible = 1;
  659.             }
  660.      else
  661.         PrevVisible = 0;    
  662.         if (PrintApogee && (MeanAnomaly == PI))
  663.             DidApogee = 1;
  664.  
  665.  
  666.         PrevTime = CurrentTime;
  667.         CurrentTime = TmpTime;
  668.         }
  669.     fclose(OutFile);
  670. }
  671.