/
line2.cpp
165 lines (148 loc) · 5.66 KB
/
line2.cpp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
/* See LICENSE.
The Minor Planet Center accepts astrometry from spacecraft using
a modification of their usual 80-column "punched-card" format. A
second line is used to tell you where the spacecraft was, relative
to the geocenter.
https://minorplanetcenter.net/iau/info/SatelliteObs.html
This code can reset those positions for Earth-orbiting spacecraft using
TLEs ('two-line elements'; https://www.projectpluto.com/tle_info.htm.)
This helps when accurate positions are not provided or completely trusted
or appear to be bad.
This code will read the 80-column astrometry and, when it finds a
spacecraft position report (the "second line"), look through the TLE
file for a matching TLE. If it finds one, it'll compute the spacecraft
location for that time and modify the position accordingly.
This was written specifically to handle a problem with a satellite in
LEO. TLEs don't exist for heliocentric objects, but I may revise this
code eventually to add positions for them (this could be done by, for
example, requesting them from JPL Horizons). TLEs exist for TESS,
computed by me, but some work would be required to make those actually
function properly (the numbers are larger and the format is somewhat
different as a result). So this really works, at present, only for
(C51) WISE, (C52) Swift, and (C53) NEOSSat. */
#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <assert.h>
#include <errno.h>
#include <math.h>
#include "norad.h"
#include "watdefs.h"
#include "afuncs.h"
#include "mpc_func.h"
#include "stringex.h"
#define PI \
3.1415926535897932384626433832795028841971693993751058209749445923
int find_tle( tle_t *tle, const char *filename, const int norad_no)
{
FILE *ifile = fopen( filename, "rb");
char line0[100], line1[100], line2[100];
int rval = -1;
if( !ifile)
{
fprintf( stderr, "Couldn't open TLE file '%s'\n", filename);
exit( -1);
}
*line0 = *line1 = '\0';
while( rval && fgets( line2, sizeof( line2), ifile))
{
if( *line1 == '1' && *line2 == '2' && atoi( line1 + 2) == norad_no)
if( parse_elements( line1, line2, tle) >= 0)
rval = 0;
strlcpy_error( line0, line1);
strlcpy_error( line1, line2);
}
fclose( ifile);
if( rval)
{
fprintf( stderr, "Couldn't find TLEs for %5d in TLE file '%s'\n",
norad_no, filename);
exit( -1);
}
return( rval);
}
/*
C49 = 29510 = 2006-047A = STEREO-A
C50 = 29511 = 2006-047B = STEREO-B
C51 = 36119 = 2009-071A = WISE *
C52 = 28485 = 2004-047A = Swift *
C53 = 39089 = 2013-009D = NEOSSat *
C54 = 28928 = 2006-001A = New Horizons
C55 = 34380 = 2009-011A = Kepler
C56 = 41043 = 2015-070A = LISA-Pathfinder
C57 = 43435 = 2018-038A = TESS *
Note that only the asterisked objects are actually in earth orbit
and therefore have TLEs. But I may try to figure out some way to add
in positions for the heliocentric spacecraft later. */
static int mpc_code_to_norad_number( const char *mpc_code)
{
const char *codes = "C49 C50 C51 C52 C53 C54 C55 C56 C57 ";
const int norad_numbers[] = { 29510, 29511, 36119, 28485, 39089,
28928, 34380, 41043, 43435 };
int i;
for( i = 0; codes[i * 4]; i++)
if( !memcmp( codes + i * 4, mpc_code, 3))
return( norad_numbers[i]);
assert( 1); /* not supposed to ever happen, unless a new */
return( 0); /* satellite has been added */
}
int main( const int argc, const char **argv)
{
const char *tle_filename = (argc > 2 ? argv[2] : "all_tle.txt");
const char *astrometry_filename = argv[1];
FILE *ifile = (argc > 1 ? fopen( astrometry_filename, "rb") : NULL);
int is_deep = 0, curr_norad = 0;
char buff[200];
double params[N_SAT_PARAMS];
const double J2000 = 2451545.;
tle_t tle;
if( argc > 1 && !ifile)
fprintf( stderr, "'%s' not found : %s\n", astrometry_filename, strerror( errno));
if( !ifile)
{
fprintf( stderr, "'line2' takes as a command-line argument the name of the input\n"
"astrometry file. It sends astrometry with corrected/added\n"
"spacecraft locations to stdout.\n");
exit( -1);
}
while( fgets( buff, sizeof( buff), ifile))
{
double jd;
if( strlen( buff) > 80 && buff[14] == 's'
&& (jd = extract_date_from_mpc_report( buff, NULL)) > 2400000.)
{
double state[6], state_j2000[6], precess_matrix[9];
double t_since;
const int norad_number = mpc_code_to_norad_number( buff + 77);
int i;
if( curr_norad != norad_number)
{
curr_norad = norad_number;
find_tle( &tle, tle_filename, norad_number);
is_deep = select_ephemeris( &tle);
if( is_deep)
SDP4_init( params, &tle);
else
SGP4_init( params, &tle);
}
t_since = (jd - tle.epoch) * minutes_per_day;
if( is_deep)
SDP4( t_since, &tle, params, state, state + 3);
else
SGP4( t_since, &tle, params, state, state + 3);
setup_precession( precess_matrix, 2000. + (jd - J2000) / 365.25, 2000.);
precess_vector( precess_matrix, state, state_j2000);
precess_vector( precess_matrix, state + 3, state_j2000 + 3);
for( i = 0; i < 3; i++)
{
char *tptr = buff + 34 + i * 12;
snprintf_err( tptr, 12, "%11.4f", fabs( state_j2000[i]));
*tptr = (state_j2000[i] > 0. ? '+' : '-');
tptr[11] = ' ';
}
}
printf( "%s", buff);
}
fclose( ifile);
return( 0);
}