-
Notifications
You must be signed in to change notification settings - Fork 1
Expand file tree
/
Copy pathGNSSdataFromGRD.cpp
More file actions
2373 lines (2317 loc) · 128 KB
/
GNSSdataFromGRD.cpp
File metadata and controls
2373 lines (2317 loc) · 128 KB
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
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
/** @file GNSSdataFromGRD.cpp
* Contains the implementation of the GNSSdataFromGRD class for ORD and NRD raw data files.
*
*/
#include "GNSSdataFromGRD.h"
#define LOG_MSG_COUNT (" @" + to_string(msgCount))
/**Constructs a GNSSdataFromGRD object using parameter passed.
*
*@param pl a pointer to the Logger to be used to record logging messages
*/
GNSSdataFromGRD::GNSSdataFromGRD(Logger* pl) {
plog = pl;
dynamicLog = false;
setInitValues();
}
/**Constructs a GNSSdataFromGRD object logging data into the stderr.
*
*/
GNSSdataFromGRD::GNSSdataFromGRD() {
plog = new Logger();
dynamicLog = true;
setInitValues();
}
/**Destroys a GNSSdataFromGRD object
*/
GNSSdataFromGRD::~GNSSdataFromGRD(void) {
if (dynamicLog) delete plog;
}
/**openInputGRD opens the GRD input file with the name and in the path given.
* @param inputFilePath the full path to the file to be open
* @param inputFileName the name of the file to be open
* @return true if input file succesfully open, false otherwhise
*/
bool GNSSdataFromGRD::openInputGRD(string inputFilePath, string inputFileName) {
msgCount = 0;
//open input raw data file
bool retVal = true;
string inFileName = inputFilePath + inputFileName;
if ((grdFile = fopen(inFileName.c_str(), "r")) == NULL) {
plog->warning(LOG_MSG_ERROPEN + inFileName);
return false;
}
return true;
}
/**rewindInputGRD rewinds the GRD input file already open
*/
void GNSSdataFromGRD::rewindInputGRD() {
msgCount = 0;
rewind(grdFile);
}
/**closeInputGRD closes the currently open input GRD file
*/
void GNSSdataFromGRD::closeInputGRD() {
fclose(grdFile);
}
/**collectHeaderData extracts data from the current ORD or NRD file for the RINEX file header.
* Also other parameters useful for processing observation or navigation data are collected here.
* <p>To collect these data the whole file is parsed from begin to end, and lines with message types not containing
* data useful for header are skipped.
* <p>Most data are collected only from the first file. Data from MT_SATOBS useful to identify systems and signals
* being tracked, and from MT_SATNAV_GLONASS_L1_CA related to slot number and carrier frequency are collected from all files.
* <p>Data collected are saved in the RinexData object passed.
*
* @param rinex the RinexData object where header data will be saved
* @param inFileNum the number of the current input raw data file from a total of inFileTotal. First value = 0
* @param inFileLast the last number of input file to be processed
* @return true if header data have been extracted, false otherwise (file cannot be processed)
*/
bool GNSSdataFromGRD::collectHeaderData(RinexData &rinex, int inFileNum = 0, int inFileLast = 0) {
char msgBuffer[100];
char smallBuffer[10];
double dvoid, dvoid2;
long long llvoid;
int ivoid;
string svoid;
string pgm, runby, date;
string observer, agency;
string rcvNum, rcvType, rcvVer;
vector <string> aVectorStr;
string logMsg;
//message parameters in ORD or NRD files
int msgType;
char constId;
int satNum, trackState, phaseState;
double carrierFrequencyMHz;
//for data extracted from the navigation message
bool hasGLOsats;
bool tofoUnset = true; //time of first observation not set
int weekNumber; //GPS week number without roll over
string msgEpoch = "Fist epoch";
rewind(grdFile);
while ((feof(grdFile) == 0) && (fscanf(grdFile, "%d;", &msgType) == 1)) {
//there are messages in the raw data file
msgCount++;
logMsg = getMsgDescription(msgType);
switch(msgType) {
case MT_GRDVER:
if ((fgets(msgBuffer, sizeof msgBuffer, grdFile) != msgBuffer)
|| !trimBuffer(msgBuffer, "\r \t\f\v\n")
|| !processHdData(rinex, msgType, string(msgBuffer))) {
plog->severe(logMsg + "CANNOT process this file (.type;version): " + string(msgBuffer));
return false;
}
//fgets already reads the EOL (skipToEOM not neccesary)!
continue;
case MT_DATE:
case MT_RINEXVER:
case MT_PGM:
case MT_RUN_BY:
case MT_INTERVALMS:
case MT_SIGU:
case MT_MARKER_NAME:
case MT_MARKER_TYPE:
case MT_OBSERVER:
case MT_AGENCY:
case MT_RECNUM:
case MT_DVTYPE:
case MT_DVVER:
case MT_LLA:
case MT_FIT:
//they include data directly used for RINEX header lines
if (fgets(msgBuffer, sizeof msgBuffer, grdFile) == msgBuffer) {
trimBuffer(msgBuffer, "\r \t\f\v\n");
if (inFileNum == 0) processHdData(rinex, msgType, string(msgBuffer));
} else plog->warning(logMsg + LOG_MSG_PARERR);
//fgets already reads the EOL (skipToEOM not neccesary)!
continue;
case MT_SATOBS:
//it includes data used to identify systems and signals being tracked
memset(smallBuffer, 0, sizeof smallBuffer); //signal identification to be stored here
if (fscanf(grdFile, "%c%d;%c%c;%d;%*lld;%*lf;%d;%*lf;%*lf;%lf", &constId, &satNum, smallBuffer, smallBuffer+1, &trackState, &phaseState, &carrierFrequencyMHz) == 7) {
if (constId == 'R') satNum = gloOSN(satNum, *smallBuffer, carrierFrequencyMHz, true);
//ignore unknown measurements or not having at least a valid pseudorrange or carrier phase
if (isKnownMeasur(constId, satNum, *smallBuffer, *(smallBuffer+1))) {
if (!isPsAmbiguous(constId, smallBuffer, trackState, dvoid, dvoid2, llvoid) || !isCarrierPhInvalid(constId, smallBuffer, phaseState)) {
if (addSignal(constId, string(smallBuffer)))
plog->config(logMsg + " added signal " + string(1, constId) + MSG_SPACE + string(smallBuffer));
}
}
} else plog->warning(logMsg + LOG_MSG_PARERR);
break;
case MT_EPOCH:
//it includes data used in time related header lines
collectAndSetEpochTime(rinex, dvoid, ivoid, logMsg + msgEpoch);
if (tofoUnset && inFileNum == 0) {
//set Time of Firts and Last Observation
rinex.setHdLnData(rinex.TOFO);
tofoUnset = false;
rinex.setHdLnData(rinex.TOLO);
//compute number of roll overs occurred at current epoch time for different constellations
rinex.getEpochTime(weekNumber, dvoid, dvoid2, ivoid);
nGPSrollOver = weekNumber / 1024;
nGALrollOver = (weekNumber - 1024) / 4096; //GALT starts at first GPST roll over (22/8/1999 00:00:00 GPST)
nBDSrollOver = (weekNumber - 1356) / 8192; //1356: week number of Jan 1, 2006 (BDST start) = 1/1/2006 00:00:14 GPST
msgEpoch = "Epoch ";
} else {
rinex.setHdLnData(rinex.TOLO);
}
break;
case MT_SATNAV_GPS_L1_CA:
//it includes data used in corrections (ionospheric, time, etc.) header lines
if (collectGPSL1CACorrections(rinex, msgType)) {
rinex.setHdLnData(RinexData::SYS, 'G', aVectorStr);
}
break;
case MT_SATNAV_GPS_L5_C:
case MT_SATNAV_GPS_C2:
case MT_SATNAV_GPS_L2_C:
//TODO extract data for header lines from these navigation message signals
rinex.setHdLnData(RinexData::SYS, 'G', aVectorStr);
break;
case MT_SATNAV_GLONASS_L1_CA:
//it includes data used in corrections (ionospheric, time, etc.) header lines, and identify satellite OSN
if (collectGLOL1CACorrections(rinex, msgType)) {
rinex.setHdLnData(RinexData::SYS, 'R', aVectorStr);
}
break;
case MT_SATNAV_GALILEO_INAV:
//it includes data used in corrections (ionospheric, time, etc.) header lines
if (collectGALINCorrections(rinex, msgType)) {
rinex.setHdLnData(RinexData::SYS, 'E', aVectorStr);
}
break;
case MT_SATNAV_GALILEO_FNAV:
//TODO extract data for header lines from this navigation message signal
rinex.setHdLnData(RinexData::SYS, 'E', aVectorStr);
break;
case MT_SATNAV_BEIDOU_D1:
//it includes data used in corrections (ionospheric, time, etc.) header lines
if (collectBDSD1Corrections(rinex, msgType)) {
rinex.setHdLnData(RinexData::SYS, 'C', aVectorStr);
}
break;
case MT_SATNAV_BEIDOU_D2:
//TODO extract data for header lines from this navigation message signal
rinex.setHdLnData(RinexData::SYS, 'C', aVectorStr);
break;
default:
plog->warning(logMsg + to_string(msgType));
break;
}
skipToEOM();
}
if (inFileNum == inFileLast) {
setHdSys(rinex);
processFilterData(rinex);
hasGLOsats = false;
for (int i=0; rinex.getHdLnData(RinexData::SYS, constId, aVectorStr, i); i++) {
//set empty PHSH records because Android does not provide specific data on this subject
aVectorStr.clear();
rinex.setHdLnData(RinexData::PHSH, constId, string(), 0.0, aVectorStr);
if (constId == 'R') hasGLOsats = true;
}
if (hasGLOsats) {
//set empty GLPHS records because Android does not provides data on this subject
rinex.setHdLnData(RinexData::GLPHS, "C1C", 0.0);
rinex.setHdLnData(RinexData::GLPHS, "C1P", 0.0);
rinex.setHdLnData(RinexData::GLPHS, "C2C", 0.0);
rinex.setHdLnData(RinexData::GLPHS, "C2P", 0.0);
//set GLSLT data for the existing slots
//?for (int i=0; i<GLO_MAXSATELLITES; i++) {
for (int i=0; i<GLO_MAXOSN; i++) {
if ((glonassOSN_FCN[i].osn != 0) && glonassOSN_FCN[i].fcnSet) {
rinex.setHdLnData(RinexData::GLSLT, glonassOSN_FCN[i].osn, glonassOSN_FCN[i].fcn);
}
}
//log the OSN - FCN table
plog->config("Table from GLONASS almanacs [Sat, nA(OSN), HnA(FCN)]:");
ivoid = 1; //offset to obtain satellite number from table index
for (int i=0; i<GLO_MAXSATELLITES; i++) {
if (i == GLO_MAXOSN) ivoid = GLO_FCN2OSN;
logMsg = "R" + to_string(i+ivoid) + MSG_COMMA;
if (glonassOSN_FCN[i].osn != 0) logMsg += to_string(glonassOSN_FCN[i].osn);
logMsg += MSG_COMMA;
if (glonassOSN_FCN[i].fcnSet) logMsg += to_string(glonassOSN_FCN[i].fcn);
plog->config(logMsg);
}
}
}
return true;
}
/**collectEpochObsData extracts observation and time data from observation raw data (ORD) file messages for one epoch.
*<p>Epoch RINEX data are contained in a sequence of messages in the form: MT_EPOCH {MT_SATOBS}.
*<p>The epoch starts with a MT_EPOCH message. It contains clock data for the current epoch, including its time,
* and the number of MT_SATOBS that would follow.
*<p>Each MT_SATOBS message contains observables for a given satellite (in each epoch the receiver sends a MT_SATOBS
* message for each satellite being tracked).
* An unexpected end of epoch occurs when the number of MT_SATOBS messages does not correspond with what has been
* stated in MT_EPOCH. This event would be logged at warning level.
*<p>This method acquires data reading messages recorded in the raw data file according to the above sequence and
* computes the obsevables from the acquired raw data. These data are saved in the RINEX data structure for further
* generation/printing of RINEX observation records in the observation file.
* The method ends when it is read the last message of the current epoch.
*<p>Other messages in the input file different from MT_EPOCH ot MT_SATOBS are ignored.
*
* @param rinex the RinexData object where data got from receiver will be placed
* @return true when data from all observation messages of an epoch have been acquired, false otherwise (End Of File reached)
*/
bool GNSSdataFromGRD::collectEpochObsData(RinexData &rinex) {
//variables to get data from ORD observation records
int msgType;
char constellId;
int satNum;
char signalId[4] = {0};
int synchState;
long long tTx = 0; //the satellite transmitted clock in nanosec. (from receivedSatTimeNanos given by Android I/F)
long long tTxUncert;
double timeOffsetNanos;
int carrierPhaseState;
double carrierPhase, cn0db, carrierFrequencyMHz;
double psRangeRate, psRangeRateUncert;
//variables to obtain observables
//This app version assumes for tRxGNSS the GPS time reference (tRxGNSS is here tRxGPS)
double tRx = 0.0; //the receiver clock in nanosececonds from the beginning of the current GPS week
double tRxGNSS = 0; //the receiver clock for the GNSS system of the satellite being processed
double tow = 0; //time of week in seconds from the beginning of the current GPS week
double pseudorange = 0.0;
double dopplerShift = 0.0;
int lli = 0; //loss o lock indicator as per RINEX V3.04
int sn_rnx = 0; //signal strength as per RINEX V3.04
int numMeasur = 0; //number of satellite measurements in current epoch
bool psAmbiguous = false;
bool phInvalid = false;
while (fscanf(grdFile, "%d;", &msgType) == 1) {
msgCount++;
switch(msgType) {
case MT_EPOCH:
if (numMeasur > 0) plog->warning(getMsgDescription(msgType) + "Few MT_SATOBS in epoch");
tRx = collectAndSetEpochTime(rinex, tow, numMeasur, getMsgDescription(msgType) + "Epoch");
break;
case MT_SATOBS:
if (numMeasur <= 0) {
plog->warning(getMsgDescription(msgType) + "MT_SATOBS before MT_EPOCH");
break;
}
numMeasur--; //an MT_SATOBS has been read, get its parameters
if (fscanf(grdFile, "%c%d;%c%c;%d;%lld;%lf;%d;%lf;%lf;%lf;%lf;%lf;%lld",
&constellId, &satNum, signalId+1, signalId+2,
&synchState, &tTx, &timeOffsetNanos,
&carrierPhaseState, &carrierPhase,
&cn0db, &carrierFrequencyMHz,
&psRangeRate, &psRangeRateUncert, &tTxUncert) != 14) {
plog->warning(getMsgDescription(msgType) + "MT_SATOBS params");
break;
}
//solve the issue of Glonass satellite number
if (constellId == 'R') satNum = gloOSN(satNum);
if (isKnownMeasur(constellId, satNum, signalId[1], signalId[2])) {
psAmbiguous = isPsAmbiguous(constellId, signalId+1, synchState, tRx, tRxGNSS, tTx);
phInvalid = isCarrierPhInvalid(constellId, signalId, carrierPhaseState);
if (!psAmbiguous || !phInvalid) {
//from data available compute signal to noise RINEX index
sn_rnx = (int) (cn0db / 6);
if (sn_rnx < 1) sn_rnx = 1;
else if (sn_rnx > 9) sn_rnx = 1;
//set and comupute pseudorrange values. Computation depends on tracking state and constellation time frame
signalId[0] = 'C';
pseudorange = (tRxGNSS - (double) tTx - timeOffsetNanos) * SPEED_OF_LIGTH_MxNS;
if (psAmbiguous || pseudorange < 0) pseudorange = 0.0;
rinex.saveObsData(constellId, satNum, string(signalId), pseudorange, 0, sn_rnx, tow);
//set carrier phase values and LLI
signalId[0] = 'L';
lli = 0; //valid or unkown by default
if (phInvalid) {
carrierPhase = 0.0; //invalid carrier phase
} else {
if ((carrierPhaseState & ADR_ST_CYCLE_SLIP) != 0) lli |= 0x01; //cycle slip detected
if ((carrierPhaseState & ADR_ST_RESET) != 0) lli |= 0x01; //reset detected
if ((carrierPhaseState & ADR_ST_HALF_CYCLE_RESOLVED) != 0) lli |= 0x01;
}
//phase, given in meters, shall be converted to full cycles
//TODO to analyse taking into account apply bias and half cycle
carrierPhase *= carrierFrequencyMHz * WLFACTOR;
rinex.saveObsData(constellId, satNum, string(signalId), carrierPhase, lli, sn_rnx, tow);
//set doppler values
signalId[0] = 'D';
dopplerShift = - psRangeRate * carrierFrequencyMHz * DOPPLER_FACTOR;
rinex.saveObsData(constellId, satNum, string(signalId), dopplerShift, 0, sn_rnx, tow);
//set signal to noise values
signalId[0] = 'S';
rinex.saveObsData(constellId, satNum, string(signalId), cn0db, 0, sn_rnx, tow);
plog->finer(getMsgDescription(msgType) + string(1, constellId) + to_string(satNum) + MSG_SPACE +
string(signalId+1) + MSG_SPACE +
to_string(pseudorange) + MSG_SPACE + to_string(carrierPhase) + MSG_SPACE +
to_string(dopplerShift) + MSG_SPACE + to_string(cn0db));
} else {
plog->fine(getMsgDescription(msgType) + string(1, constellId) + to_string(satNum) + MSG_SPACE +
string(signalId+1) + LOG_MSG_INVM);
}
} else plog->warning(getMsgDescription(msgType) + string(1, constellId) + to_string(satNum) + MSG_SPACE +
string(signalId+1) + MSG_SPACE + LOG_MSG_UNK);
if (numMeasur <= 0) {
skipToEOM();
return true;
}
break;
case MT_SATNAV_GPS_L1_CA:
case MT_SATNAV_GLONASS_L1_CA:
case MT_SATNAV_GALILEO_FNAV:
case MT_SATNAV_BEIDOU_D1:
case MT_SATNAV_GPS_L5_C:
plog->warning(getMsgDescription(msgType) + LOG_MSG_NINO);
default:
break;
}
skipToEOM();
}
return false;
}
/**collectNavData iterates over the input raw data file extracting navigation messages to process their data.
* Data extracted are saved in a RinexData object for futher printing in a navigation RINEX file.
* It is assumed that grdFile file type and version are the correct ones.
*
* @param rinex the RinexData object where data got from receiver will be placed
* @return true when navigation data from at least one epoch messages have been acquired, false otherwise
*/
bool GNSSdataFromGRD::collectNavData(RinexData &rinex) {
int msgType;
bool acquiredNavData = false;
msgCount = 0;
while (fscanf(grdFile, "%d;", &msgType) == 1) {
msgCount++;
switch(msgType) {
case MT_SATNAV_GPS_L1_CA:
acquiredNavData |= collectGPSL1CAEphemeris(rinex, msgType);
break;
case MT_SATNAV_GLONASS_L1_CA:
acquiredNavData |= collectGLOL1CAEphemeris(rinex, msgType);
break;
case MT_SATNAV_GALILEO_INAV:
acquiredNavData |= collectGALINEphemeris(rinex, msgType);
break;
case MT_SATNAV_BEIDOU_D1:
acquiredNavData |= collectBDSD1Ephemeris(rinex, msgType);
break;
case MT_SATNAV_GPS_L5_C:
case MT_SATNAV_GPS_C2:
case MT_SATNAV_GPS_L2_C:
case MT_SATNAV_GALILEO_FNAV:
case MT_SATNAV_BEIDOU_D2:
plog->warning(getMsgDescription(msgType) + MSG_NOT_IMPL + LOG_MSG_NAVIG);
break;
case MT_EPOCH:
break;
case MT_SATOBS:
plog->warning(getMsgDescription(msgType) + LOG_MSG_NONI);
default:
break;
}
skipToEOM();
}
return acquiredNavData;
}
/**processHdData sets the Rinex header record data from the contents of the given message.
* Each message type contains data to be processed and stored in the corresponding header record.
*
* @param rinex the RinexData object which header data will be updated
* @param msgType the raw data file message type
* @param msgContent the content of the message with data to be extracted
* @return true when data has been processed without error, false otherwise
*/
bool GNSSdataFromGRD::processHdData(RinexData &rinex, int msgType, string msgContent) {
double da, db, dc;
double x, y, z;
string msgError;
vector<string> selElements;
size_t position;
double dvoid = 0.0;
string svoid = string();
plog->config(getMsgDescription(msgType) + msgContent);
try {
switch(msgType) {
case MT_GRDVER:
position = msgContent.find(';');
if ((position != string::npos)
&& isGoodGRDver(msgContent.substr(0, position), stoi(msgContent.substr(position+1, string::npos)))) return true;
break;
case MT_SITE:
plog->finer(getMsgDescription(msgType) + " currently ignored");
return true;
case MT_PGM:
return rinex.setHdLnData(rinex.RUNBY, msgContent, svoid, svoid);
case MT_DVTYPE:
return rinex.setHdLnData(rinex.RECEIVER, svoid, msgContent, svoid);
case MT_DVVER:
return rinex.setHdLnData(rinex.RECEIVER, svoid, svoid, msgContent);
case MT_LLA:
if (sscanf(msgContent.c_str(), "%lf;%lf;%lf", &da, &db, &dc) == 3) {
llaTOxyz(da * dgrToRads, db * dgrToRads, dc, x, y, z);
return rinex.setHdLnData(rinex.APPXYZ, x, y, z);
}
plog->warning(getMsgDescription(msgType) + LOG_MSG_PARERR);
return false;
case MT_DATE: //Date of the file
return rinex.setHdLnData(rinex.RUNBY, svoid, svoid, msgContent);
case MT_INTERVALMS:
return rinex.setHdLnData(rinex.INT, stod(msgContent) / 1000., dvoid, dvoid);
case MT_SIGU:
return rinex.setHdLnData(rinex.SIGU, msgContent, svoid, svoid);
case MT_RINEXVER:
return rinex.setHdLnData(rinex.VERSION, stod(msgContent), dvoid, dvoid);
case MT_RUN_BY:
return rinex.setHdLnData(rinex.RUNBY, svoid, msgContent, svoid);
case MT_MARKER_NAME:
return rinex.setHdLnData(rinex.MRKNAME, msgContent, svoid, svoid);
case MT_MARKER_TYPE:
return rinex.setHdLnData(rinex.MRKTYPE, msgContent, svoid, svoid);
case MT_OBSERVER:
return rinex.setHdLnData(rinex.AGENCY, msgContent, svoid, svoid);
case MT_AGENCY:
return rinex.setHdLnData(rinex.AGENCY, svoid, msgContent, svoid);
case MT_RECNUM:
return rinex.setHdLnData(rinex.RECEIVER, msgContent, svoid, svoid);
case MT_COMMENT:
return rinex.setHdLnData(rinex.COMM, rinex.RUNBY, msgContent);
case MT_MARKER_NUM:
return rinex.setHdLnData(rinex.MRKNUMBER, msgContent, svoid, svoid);
case MT_CLKOFFS:
clkoffset = stoi(msgContent);
rinex.setHdLnData(rinex.CLKOFFS, clkoffset);
applyBias = clkoffset == 1;
plog->config(getMsgDescription(msgType) + " applyBias:" + (applyBias?string("TRUE"):string("FALSE")));
return true;
case MT_FIT:
if (msgContent.find("TRUE") != string::npos) fitInterval = true;
else fitInterval = false;
return true;
case MT_LOGLEVEL:
plog->setLevel(msgContent);
return true;
case MT_CONSTELLATIONS:
selElements = getElements(msgContent, "[], ");
for (vector<string>::iterator it = selElements.begin(); it != selElements.end(); ++it) {
if ((*it).compare("GPS") == 0) *it = "G";
else if ((*it).compare("GLONASS") == 0) *it = "R";
else if ((*it).compare("GALILEO") == 0) *it = "E";
else if ((*it).compare("BEIDOU") == 0) *it = "C";
else if ((*it).compare("SBAS") == 0) *it = "S";
else if ((*it).compare("QZSS") == 0) *it = "J";
else {
plog->warning(getMsgDescription(msgType) + LOG_MSG_UNKSELSYS + (*it));
selElements.erase(it);
}
}
selSatellites.insert(selSatellites.end(), selElements.begin(), selElements.end());
return true;
case MT_SATELLITES:
selElements = getElements(msgContent, "[],;.:- ");
selSatellites.insert(selSatellites.end(), selElements.begin(), selElements.end());
return true;
case MT_OBSERVABLES:
selObservables = getElements(msgContent, "[], ");
return true;
default:
plog->warning(getMsgDescription(msgType) + to_string(msgType));
break;
}
} catch (string error) {
plog->severe(error + LOG_MSG_COUNT);
}
return false;
}
/**processFilterData sets filtering data (constellations, satellites, observations, if any) to be used when printing
* RINEX files.
*
* @param rinex the RinexData object where filtering data will be set
*/
void GNSSdataFromGRD::processFilterData(RinexData &rinex) {
rinex.setFilter(selSatellites, selObservables);
}
//PRIVATE METHODS
//===============
/**setInitValues set initial values in class variables and tables containig f.e. conversion parameters
* used to translate scaled normalized GPS message data to values in actual units.
* <b>Called by the construtors
*/
void GNSSdataFromGRD::setInitValues() {
ordVersion = 0;
nrdVersion = 0;
msgCount = 0;
clkoffset = 0;
applyBias = false;
fitInterval = false;
//default values for roll overs
nGPSrollOver = 2;
nGALrollOver = 0;
nBDSrollOver = 0;
//scale factors to apply to ephemeris
//set default values
double COMMON_SCALEFACTOR[BO_LINSTOTAL][BO_MAXCOLS]; //the scale factors to apply to GPS, GAL and BDS broadcast orbit data to obtain ephemeris
for (int i = 0; i < BO_LINSTOTAL; i++) {
for (int j = 0; j < BO_MAXCOLS; j++) {
COMMON_SCALEFACTOR[i][j] = 1.0;
}
}
//factors common to GPS, Glonass, Galileo or Beidou?
//SV clock data have specific factors for each constellation
//broadcast orbit 1
COMMON_SCALEFACTOR[1][1] = pow(2.0, -5.0); //Crs
COMMON_SCALEFACTOR[1][2] = pow(2.0, -43.0) * ThisPI; //Delta N
COMMON_SCALEFACTOR[1][3] = pow(2.0, -31.0) * ThisPI; //M0
//broadcast orbit 2
COMMON_SCALEFACTOR[2][0] = pow(2.0, -29.0); //Cuc
COMMON_SCALEFACTOR[2][1] = pow(2.0, -33.0); //e
COMMON_SCALEFACTOR[2][2] = pow(2.0, -29.0); //Cus
COMMON_SCALEFACTOR[2][3] = pow(2.0, -19.0); //sqrt(A) ????
//broadcast orbit 3
COMMON_SCALEFACTOR[3][1] = pow(2.0, -29.0); //Cic
COMMON_SCALEFACTOR[3][2] = pow(2.0, -31.0) * ThisPI; //Omega0
COMMON_SCALEFACTOR[3][3] = pow(2.0, -29.0); //Cis
//broadcast orbit 4
COMMON_SCALEFACTOR[4][0] = pow(2.0, -31.0) * ThisPI; //i0
COMMON_SCALEFACTOR[4][1] = pow(2.0, -5.0); //Crc
COMMON_SCALEFACTOR[4][2] = pow(2.0, -31.0) * ThisPI; //w
COMMON_SCALEFACTOR[4][3] = pow(2.0, -43.0) * ThisPI; //w dot
//broadcast orbit 5
COMMON_SCALEFACTOR[5][0] = pow(2.0, -43.0) * ThisPI; //Idot
//Iono scale factors
COMMON_SCALEFACTOR[BO_LIN_IONOA][0] = pow(2.0, -30); //Iono alfa0
COMMON_SCALEFACTOR[BO_LIN_IONOA][1] = pow(2.0, -27); //Iono alfa1
COMMON_SCALEFACTOR[BO_LIN_IONOA][2] = pow(2.0, -24); //Iono alfa2
COMMON_SCALEFACTOR[BO_LIN_IONOA][3] = pow(2.0, -24); //Iono alfa3
COMMON_SCALEFACTOR[BO_LIN_IONOB][0] = pow(2.0, 11); //Iono beta0
COMMON_SCALEFACTOR[BO_LIN_IONOB][1] = pow(2.0, 14); //Iono beta1
COMMON_SCALEFACTOR[BO_LIN_IONOB][2] = pow(2.0, 16); //Iono beta2
COMMON_SCALEFACTOR[BO_LIN_IONOB][3] = pow(2.0, 16); //Iono beta3
//Time corrections common scale factors
COMMON_SCALEFACTOR[BO_LIN_TIMEU][0] = pow(2.0, -30); //TCORR A0UTC
COMMON_SCALEFACTOR[BO_LIN_TIMEU][1] = pow(2.0, -50); //TCORR A1UTC
COMMON_SCALEFACTOR[BO_LIN_TIMEU][3] = 1.0; //WNt scale factor
COMMON_SCALEFACTOR[BO_LIN_TIMEG][2] = 1.0;
COMMON_SCALEFACTOR[BO_LIN_TIMEG][3] = 1.0;
//set scale factors for GPS
memcpy(GPS_SCALEFACTOR, COMMON_SCALEFACTOR, sizeof(GPS_SCALEFACTOR));
//SV clock data
GPS_SCALEFACTOR[0][0] = pow(2.0, 4.0); //T0c
GPS_SCALEFACTOR[0][1] = pow(2.0, -31.0); //Af0: SV clock bias
GPS_SCALEFACTOR[0][2] = pow(2.0, -43.0); //Af1: SV clock drift
GPS_SCALEFACTOR[0][3] = pow(2.0, -55.0); //Af2: SV clock drift rate GAL
//broadcast orbit 3
GPS_SCALEFACTOR[3][0] = pow(2.0, 4.0); //TOE
//broadcast orbit 6
GPS_SCALEFACTOR[6][2] = pow(2.0, -31.0); //TGD
//broadcast orbit 7
GPS_SCALEFACTOR[7][0] = 0.01; //Transmission time of message in sec x 100
GPS_SCALEFACTOR[7][2] = 0.0; //Spare
GPS_SCALEFACTOR[7][3] = 0.0; //Spare
//GPST - UTC factors
GPS_SCALEFACTOR[BO_LIN_TIMEU][2] = pow(2.0, 12); //t0t scale factor
//GPS_URA table to obtain in meters the value associated to the satellite accuracy index
//fill GPS_URA table as per GPS ICD 20.3.3.3.1.3
GPS_URA[0] = 2.0; //0
GPS_URA[1] = 2.8; //1
GPS_URA[2] = 4.0; //2
GPS_URA[3] = 5.7; //3
GPS_URA[4] = 8.0; //4
GPS_URA[5] = 11.3; //5
GPS_URA[6] = pow(2.0, 4.0); //<=6
GPS_URA[7] = pow(2.0, 5.0);
GPS_URA[8] = pow(2.0, 6.0);
GPS_URA[9] = pow(2.0, 7.0);
GPS_URA[10] = pow(2.0, 8.0);
GPS_URA[11] = pow(2.0, 9.0);
GPS_URA[12] = pow(2.0, 10.0);
GPS_URA[13] = pow(2.0, 11.0);
GPS_URA[14] = pow(2.0, 12.0);
GPS_URA[15] = 6144.00;
//set factors for Galileo
memcpy(GAL_SCALEFACTOR, COMMON_SCALEFACTOR, sizeof(GAL_SCALEFACTOR));
//SV clock data
GAL_SCALEFACTOR[0][0] = 60.; //T0c
GAL_SCALEFACTOR[0][1] = pow(2.0, -34.0); //Af0: SV clock bias
GAL_SCALEFACTOR[0][2] = pow(2.0, -46.0); //Af1: SV clock drift
GAL_SCALEFACTOR[0][3] = pow(2.0, -59.0); //Af2: SV clock drift rate
//broadcast orbit 3
GAL_SCALEFACTOR[3][0] = 60.; //TOE
//broadcast orbit 5
GAL_SCALEFACTOR[5][3] = 0.0; //spare
//broadcast orbit 6
GAL_SCALEFACTOR[6][2] = pow(2.0, -32.0); //BGD E5a / E1
GAL_SCALEFACTOR[6][3] = pow(2.0, -32.0); //BGD E5b / E1
//broadcast orbit 7
//GAL_SCALEFACTOR[7][0] = 1.0; //Transmission time of message
//GAL_SCALEFACTOR[7][1] = 0.0; //Spare
//GAL_SCALEFACTOR[7][2] = 0.0; //Spare
//GAL_SCALEFACTOR[7][3] = 0.0; //Spare
//other data stored in bom
GAL_SCALEFACTOR[BO_LIN_IONOA][0] = pow(2.0, -2); //Iono Ai0
GAL_SCALEFACTOR[BO_LIN_IONOA][1] = pow(2.0, -8); //Iono Ai1
GAL_SCALEFACTOR[BO_LIN_IONOA][2] = pow(2.0, -15); //Iono Ai2
GAL_SCALEFACTOR[BO_LIN_IONOA][3] = 0.0; //Iono blank
GAL_SCALEFACTOR[BO_LIN_TIMEU][2] = 3600; //GST-UTC t0t
GAL_SCALEFACTOR[BO_LIN_TIMEG][0] = pow(2.0, -35); //GST - GPS time corr. A0GPS
GAL_SCALEFACTOR[BO_LIN_TIMEG][1] = pow(2.0, -51); //GST - GPS time corr. A1GPS
GAL_SCALEFACTOR[BO_LIN_TIMEG][2] = 3600; //GST-GPS t0G
//Glonass scale factors
//SV clock data
GLO_SCALEFACTOR[0][0] = 1.; //T0c
GLO_SCALEFACTOR[0][1] = pow(2.0, -30.0); //-TauN
GLO_SCALEFACTOR[0][2] = pow(2.0, -40.0); //+GammaN
GLO_SCALEFACTOR[0][3] = 1.; //Message frame time
//broadcast orbit 1
GLO_SCALEFACTOR[1][0] = pow(2.0, -11.0); //X
GLO_SCALEFACTOR[1][1] = pow(2.0, -20.0); //Vel X
GLO_SCALEFACTOR[1][2] = pow(2.0, -30.0); //Accel X
GLO_SCALEFACTOR[1][3] = 1.; //SV Health
//broadcast orbit 2
GLO_SCALEFACTOR[2][0] = pow(2.0, -11.0); //Y
GLO_SCALEFACTOR[2][1] = pow(2.0, -20.0); //Vel Y
GLO_SCALEFACTOR[2][2] = pow(2.0, -30.0); //Accel Y
GLO_SCALEFACTOR[2][3] = 1.; //Frequency number
//broadcast orbit 3
GLO_SCALEFACTOR[3][0] = pow(2.0, -11.0); //Z
GLO_SCALEFACTOR[3][1] = pow(2.0, -20.0); //Vel Z
GLO_SCALEFACTOR[3][2] = pow(2.0, -30.0); //Accel Z
GLO_SCALEFACTOR[3][3] = 1.; //Age of oper.
//other data stored in bom
GLO_SCALEFACTOR[BO_LIN_TIMEU][0] = pow(2.0, -31); //TauC
GLO_SCALEFACTOR[BO_LIN_TIMEU][1] = 0.0;
GLO_SCALEFACTOR[BO_LIN_TIMEU][2] = 0.0;
GLO_SCALEFACTOR[BO_LIN_TIMEU][3] = 0.0;
GLO_SCALEFACTOR[BO_LIN_TIMEG][0] = pow(2.0, -30); //TauGPS
GLO_SCALEFACTOR[BO_LIN_TIMEG][1] = 0.0;
GLO_SCALEFACTOR[BO_LIN_TIMEG][2] = 0.0;
GLO_SCALEFACTOR[BO_LIN_TIMEG][3] = 0.0;
//set scale factors for BDS
memcpy(BDS_SCALEFACTOR, COMMON_SCALEFACTOR, sizeof(BDS_SCALEFACTOR));
//SV clock data
BDS_SCALEFACTOR[0][0] = pow(2.0, 3.0); //T0c
BDS_SCALEFACTOR[0][1] = pow(2.0, -33.0); //Af0: SV clock bias
BDS_SCALEFACTOR[0][2] = pow(2.0, -50.0); //Af1: SV clock drift
BDS_SCALEFACTOR[0][3] = pow(2.0, -66.0); //Af2: SV clock drift rate GAL
//broadcast orbit 1
BDS_SCALEFACTOR[1][1] = pow(2.0, -6.0); //Crs
//broadcast orbit 2
BDS_SCALEFACTOR[2][0] = pow(2.0, -31.0); //Cuc
BDS_SCALEFACTOR[2][2] = pow(2.0, -31.0); //Cus
//broadcast orbit 3
BDS_SCALEFACTOR[3][0] = pow(2.0, 3.0); //T0e
BDS_SCALEFACTOR[3][1] = pow(2.0, -31.0); //Cic
BDS_SCALEFACTOR[3][3] = pow(2.0, -31.0); //Cis
//broadcast orbit 4
BDS_SCALEFACTOR[4][1] = pow(2.0, -6.0); //Crc
//broadcast orbit 5
BDS_SCALEFACTOR[5][1] = 0.0; //Spare
BDS_SCALEFACTOR[5][3] = 0.0; //Spare
//broadcast orbit 6
//broadcast orbit 7
BDS_SCALEFACTOR[7][2] = 0.0; //Spare
BDS_SCALEFACTOR[7][3] = 0.0; //Spare
//other data stored in bom
BDS_SCALEFACTOR[BO_LIN_TIMEG][0] = BDS_SCALEFACTOR[BO_LIN_TIMEG][1] = 1.0E-10; //A0GPS, A1GPS
//BDS_URA table to obtain in meters the value associated to the satellite accuracy index
BDS_URA[0] = 2.4;
BDS_URA[1] = 3.4;
BDS_URA[2] = 4.85;
BDS_URA[3] = 6.85;
BDS_URA[4] = 9.65;
BDS_URA[5] = 13.65;
BDS_URA[6] = 24.0;
BDS_URA[7] = 48.0;
BDS_URA[8] = 96.0;
BDS_URA[9] = 192,0;
BDS_URA[10] = 384,0;
BDS_URA[11] = 768,0;
BDS_URA[12] = 1536.0;
BDS_URA[13] = 3072.0;
BDS_URA[14] = 6144.0;
BDS_URA[15] = 6144.0;
//set tables to 0
memset(gpsSatFrame, 0, sizeof(gpsSatFrame));
memset(gloSatFrame, 0, sizeof(gloSatFrame));
memset(glonassOSN_FCN, 0, sizeof(glonassOSN_FCN));
for(int i=0; i<GLO_MAXSATELLITES; i++) glonassOSN_FCN[i].fcnSet = false;
memset(nAhnA, 0, sizeof(nAhnA));
memset(galInavSatFrame, 0, sizeof(galInavSatFrame));
}
//methods for GPS L1 CA message processing
//A macro to get bit position in a subframe from the bit position (BITNUMBER) in the message subframe stated in the GPS ICD
//Each GPS 30 bits message word is stored in an uint32_t (aligned to the rigth)
//Each GPS subframe, comprising 10 words, is tored in an stream of 10 uint32_t
//BITNUMBER is the position in the subframe, GPSL1CA_BIT gets the position in the stream
#define GPSL1CA_BIT(BITNUMBER) ((BITNUMBER-1)/30*32 + (BITNUMBER-1)%30 + 2)
/**collectGPSL1CAEphemeris gets GPS L1 C/A navigation data from a raw message and store them to generate RINEX navigation files.
* Satellite navigation messages are temporary stored until all needed subframes for the current satellite and epoch have been received.
* When completed, satellite ephemeris are extracted from subframes, scaled to working units, and stored into "broadcast orbit lines"
* of the RinexData object.
* For GPS L1 C/A each subframe of the navigation message contains 10 30-bit words.
* Each word (30 bits) should be fit into the last 30 bits in a 4-byte word (skip B31 and B32), with MSB first,
* for a total of 40 bytes, covering a time period of 6, 6, and 0.6 seconds, respectively.
* Only have interest subframes: 1,2,3 & page 18 of subframe 4 (pgID = 56 in GPS ICD Table 20-V)
*
* @param rinex the class instance where data are stored
* @param msgType the true type of the message
* @return true if data properly extracted (data correctly read, correct parity status, relevant subframe), false otherwise
*/
bool GNSSdataFromGRD::collectGPSL1CAEphemeris(RinexData &rinex, int msgType) {
char constId; //the constellation identifier this satellite belongs
int satNum; //the satellite number this navigation message belongssatt
int sfrmNum; //navigation message subframe number
int pageNum; //navigation message page number
int bom[BO_LINSTOTAL][BO_MAXCOLS]; //a RINEX broadcats orbit like arrangement for satellite ephemeris mantissa
double bo[BO_LINSTOTAL][BO_MAXCOLS]; //the RINEX broadcats orbit arrangement for satellite ephemeris
double tTag; //the time tag for ephemeris data
GPSFrameData *pframe;
string logMsg = getMsgDescription(msgType);
//read MT_SATNAV_GPS_L1_CA message data
if (!readGPSL1CANavMsg(constId, satNum, sfrmNum, pageNum, logMsg)) return false;
pframe = &gpsSatFrame[satNum - 1];
//check if all ephemeris have been received, that is:
//-subframes 1, 2, and 3 have data
//-and their data belong belong to the same Issue Of Data (IOD)
bool allRec = pframe->hasData;
for (int i = 0; allRec && i < 3; ++i) allRec = allRec && pframe->gpsSatSubframes[i].hasData;
if (allRec) {
logMsg += LOG_MSG_FRM;
//IODC (8LSB in subframe 1) must be equal to IODE in subframe 2 and IODE in subframe 3
uint32_t iodcLSB = getBits(pframe->gpsSatSubframes[0].words, GPSL1CA_BIT(211), 8);
uint32_t iode2 = getBits(pframe->gpsSatSubframes[1].words, GPSL1CA_BIT(61), 8);
uint32_t iode3 = getBits(pframe->gpsSatSubframes[2].words, GPSL1CA_BIT(271), 8);
if ((iodcLSB == iode2) && (iodcLSB == iode3)) {
//all ephemerides have been already received; extract and store them into the RINEX object
logMsg += LOG_MSG_IOD;
plog->fine(logMsg);
extractGPSL1CAEphemeris(satNum - 1, bom);
tTag = scaleGPSEphemeris(bom, bo);
rinex.saveNavData('G', satNum, bo, tTag);
//clear satellite frame storage
pframe->hasData = false;
for (int i = 0; i < GPS_MAXSUBFRS; ++i) pframe->gpsSatSubframes[i].hasData = false;
} else plog->fine(logMsg + " and IODs different.");
} else plog->finer(logMsg);
return true;
}
/**collectGPSL1CACorrections gets from GPS L1 C/A navigation raw data message the ionospheric, clock and leap corrections data
* needed to generate some RINEX navigation files header records.
* Satellite navigation messages are temporary stored until subframes 1 and 4 of a given satellite and epoch have been received.
* When completed, corrections are extracted from subframes, scaled to working units, and stored into header records of the RinexData object.
*
* @param rinex the class instance where data are stored
* @param msgType the true type of the message
* @return true if data properly extracted (data correctly read, correct parity status, relevant subframe), false otherwise
*/
bool GNSSdataFromGRD::collectGPSL1CACorrections(RinexData &rinex, int msgType) {
char constId; //the constellation identifier this satellite belongs
int satNum; //the satellite number this navigation message belongssatt
int sfrmNum; //navigation message subframe number
int pageNum; //navigation message page number
int bom[BO_LINSTOTAL][BO_MAXCOLS]; //a RINEX broadcats orbit like arrangement for satellite ephemeris mantissa
double bo[BO_LINSTOTAL][BO_MAXCOLS]; //the RINEX broadcats orbit arrangement for satellite ephemeris
double tTag; //the time tag for ephemeris data
string logMsg = getMsgDescription(msgType);
//read MT_SATNAV_GPS_L1_CA message data
if (!readGPSL1CANavMsg(constId, satNum, sfrmNum, pageNum, logMsg)) return false;
GPSFrameData* pframe = &gpsSatFrame[satNum - 1];
//check if corrections have been received, that is, subframes 1 and 4 have data
if (pframe->hasData && pframe->gpsSatSubframes[0].hasData && pframe->gpsSatSubframes[3].hasData) {
logMsg += LOG_MSG_CORR;
extractGPSL1CAEphemeris(satNum - 1, bom);
tTag = scaleGPSEphemeris(bom, bo);
rinex.setHdLnData(rinex.IONC, rinex.IONC_GPSA, bo[BO_LIN_IONOA], bom[BO_LIN_TIMEG][2], satNum);
rinex.setHdLnData(rinex.IONC, rinex.IONC_GPSB, bo[BO_LIN_IONOB], bom[BO_LIN_TIMEG][2], satNum);
rinex.setHdLnData(rinex.TIMC, rinex.TIMC_GPUT, bo[BO_LIN_TIMEU], 0, satNum);
rinex.setHdLnData(rinex.LEAP, (int) bo[BO_LIN_LEAPS][0], (int) bo[BO_LIN_LEAPS][1], (int) bo[BO_LIN_LEAPS][2], (int) bo[BO_LIN_LEAPS][3], 'G');
logMsg += "IONA&B TIMEG TIMEU LEAPS";
plog->fine(logMsg);
//clear satellite frame storage
pframe->hasData = false;
for (int i = 0; i < GPS_MAXSUBFRS; ++i) pframe->gpsSatSubframes[i].hasData = false;
} else plog->finer(logMsg);
return true;
}
/**readGPSL1CANavMsg read a GPS L1 C/A navigation message data from the navigation raw data file.
* Data are stored in gpsSatFrame of the corresponding satNum.
* For GPS L1 C/A each subframe of the navigation message contains 10 30-bit words.
* Each word (30 bits) should be fit into the last 30 bits in a 4-byte word (skip B31 and B32), with MSB first,
* for a total of 40 bytes, covering a time period of 6, 6, and 0.6 seconds, respectively.
* Only have interest subframes: 1,2,3 & page 18 of subframe 4 (pgID = 56 in GPS ICD Table 20-V)
*
* @param constId constellation identification (G, R, E, ...)
* @param satNum satellite number in the constellation
* @param sfrmNum subframe number
* @param pageNum page number
* @param logMsg is the loggin message
* @return true if data succesful read, false otherwise
*/
bool GNSSdataFromGRD::readGPSL1CANavMsg(char &constId, int &satNum, int &sfrmNum, int &pageNum, string &logMsg) {
int status; //the status of the navigation message:0=UNKNOWN, 1=PARITY PASSED, 2=REBUILT
int msgSize; //the message size in bytes
unsigned int navMsg[GPS_L1_CA_MSGSIZE]; //to store GPS message bytes from receiver
GPSSubframeData *psubframe;
GPSFrameData *pframe;
try {
//read MT_SATNAV_GPS_L1_CA message data
if (fscanf(grdFile, "%d;%c%d;%d;%d;%d", &status, &constId, &satNum, &sfrmNum, &pageNum, &msgSize) != 6
|| status < 1
|| constId != 'G'
|| satNum < GPS_MINPRN
|| satNum > GPS_MAXPRN
|| msgSize != GPS_L1_CA_MSGSIZE) {
plog->warning(logMsg + LOG_MSG_INMP + LOG_MSG_OSIZ);
return false;
}
logMsg += " sat:" + to_string(satNum) + " subfr:" + to_string(sfrmNum) + " pg:" + to_string(pageNum);
for (int i = 0; i < GPS_L1_CA_MSGSIZE; ++i) {
if (fscanf(grdFile,";%X", navMsg+i) != 1) {
plog->warning(logMsg + LOG_MSG_ERRO + LOG_MSG_INMP);
return false;
}
}
if ((sfrmNum < 1 || sfrmNum > GPS_MAXSUBFRS) || (sfrmNum == 4 && pageNum != 18)) {
plog->finer(logMsg + LOG_MSG_NAVIG);
return false;
}
pframe = &gpsSatFrame[satNum - 1];
psubframe = &pframe->gpsSatSubframes[sfrmNum - 1];
//pack message bytes into the ten words with navigation data
for (int i = 0, n=0; i < GPS_SUBFRWORDS; ++i, n+=4) {
psubframe->words[i] = navMsg[n]<<24 | navMsg[n+1]<<16 | navMsg[n+2]<<8 | navMsg[n+3];
}
psubframe->hasData = true;
pframe->hasData = true;
logMsg += LOG_MSG_SFR;
/*TODO to analyze including the code for checking satellite health
//check for SVhealth
uint32_t svHealth = (navMsq[4]>>10) & 0x3F;
if ((svHealth & 0x20) != 0) {
sprintf(errorBuff, "Wrong SV %2u health <%02hx>", sv, svHealth);
plog->info(string(errorBuff));
return false;
}
*/
} catch (int error) {
plog->severe(logMsg + LOG_MSG_ERRO + to_string((long long) error));
return false;
}
return true;
}
/**extractGPSL1CAEphemeris extract satellite ephemeris from the stored navigation frames transmitted for a given GPS satellite.
* <p>The navigation message data of interest here have been stored in GPS frame (gpsSatFrame) data words (without parity).
* Ephemeris are extracted from data words and stored into a RINEX Broadcast Orbit lines like arrangement.
* <p>This method stores the MANTISSA of each satellite ephemeris into a broadcast orbit array, without applying any scale factor.
* <p>In addition to the 8 lines x 4 ephemeris per satellite needed for RINEX, they are extracted data for iono and time corrections,
* and the UTC leap seconds available in the navigation message and needed for the header of the file. They are stored in the
* rows BO_LIN_IONOA, BO_LIN_IONOB, BO_LIN_TIMEx, and BO_LIN_LEAP respectively.
*
* @param satIdx the satellite index in the gpsSatFrame
* @param bom an array of broadcats orbit data containing the mantissa of each satellite ephemeris
*/
void GNSSdataFromGRD::extractGPSL1CAEphemeris(int satIdx, int (&bom)[BO_LINSTOTAL][BO_MAXCOLS]) {
uint32_t *psubfr1data = gpsSatFrame[satIdx].gpsSatSubframes[0].words;
uint32_t *psubfr2data = gpsSatFrame[satIdx].gpsSatSubframes[1].words;
uint32_t *psubfr3data = gpsSatFrame[satIdx].gpsSatSubframes[2].words;
uint32_t *psubfr4data = gpsSatFrame[satIdx].gpsSatSubframes[3].words;
for (int i = 0; i < BO_LINSTOTAL; ++i) {
for (int j = 0; j < BO_MAXCOLS; ++j) {
bom[i][j] = 0;
}
}
//broadcast line 0
bom[0][0] = getBits(psubfr1data, GPSL1CA_BIT(219), 16); //T0C
bom[0][1] = getTwosComplement(getBits(psubfr1data, GPSL1CA_BIT(271), 22), 22); //Af0
bom[0][2] = getTwosComplement(getBits(psubfr1data, GPSL1CA_BIT(249), 16), 16); //Af1
bom[0][3] = getTwosComplement(getBits(psubfr1data, GPSL1CA_BIT(241), 8), 8); //Af2
//broadcast line 1
bom[1][0] = getBits(psubfr2data, GPSL1CA_BIT(61), 8); //IODE
bom[1][1] = getTwosComplement(getBits(psubfr2data, GPSL1CA_BIT(69), 16), 16); //Crs
bom[1][2] = getTwosComplement(getBits(psubfr2data, GPSL1CA_BIT(91), 16), 16); //Delta n
bom[1][3] = (getBits(psubfr2data, GPSL1CA_BIT(107), 8) << 24) | getBits(psubfr2data, GPSL1CA_BIT(121), 24);//M0
//broadcast line 2
bom[2][0] = getTwosComplement(getBits(psubfr2data, GPSL1CA_BIT(151), 16), 16); //Cuc
bom[2][1] = (getBits(psubfr2data, GPSL1CA_BIT(167), 8) << 24) | getBits(psubfr2data, GPSL1CA_BIT(181), 24); //e
bom[2][2] = getTwosComplement(getBits(psubfr2data, GPSL1CA_BIT(211), 16), 16); //Cus
bom[2][3] = (getBits(psubfr2data, GPSL1CA_BIT(227), 8) << 24) | getBits(psubfr2data, GPSL1CA_BIT(241), 24); //sqrt(A)
//broadcast line 3
bom[3][0] = getBits(psubfr2data, GPSL1CA_BIT(271), 16); //Toe
bom[3][1] = getTwosComplement(getBits(psubfr3data, GPSL1CA_BIT(61), 16), 16); //Cic
bom[3][2] = (getBits(psubfr3data, GPSL1CA_BIT(77), 8) << 24) | getBits(psubfr3data, GPSL1CA_BIT(91), 24); //OMEGA0
bom[3][3] = getTwosComplement(getBits(psubfr3data, GPSL1CA_BIT(121), 16), 16); //CIS
//broadcast line 4
bom[4][0] = (getBits(psubfr3data, GPSL1CA_BIT(137), 8) << 24) | getBits(psubfr3data, GPSL1CA_BIT(151), 24); //i0
bom[4][1] = getTwosComplement(getBits(psubfr3data, GPSL1CA_BIT(181), 16), 16); //Crc
bom[4][2] = (getBits(psubfr3data, GPSL1CA_BIT(197), 8) << 24) | getBits(psubfr3data, GPSL1CA_BIT(211), 24); //w (omega)
bom[4][3] = getTwosComplement(getBits(psubfr3data, GPSL1CA_BIT(241), 24), 24); //w dot
//broadcast line 5
bom[5][0] = getTwosComplement(getBits(psubfr3data, GPSL1CA_BIT(279), 14), 14); //IDOT
bom[5][1] = getBits(psubfr1data, GPSL1CA_BIT(71), 2); //Codes on L2
bom[5][2] = getBits(psubfr1data, GPSL1CA_BIT(61), 10) + nGPSrollOver * 1024; //GPS week#
bom[5][3] = getBits(psubfr1data, GPSL1CA_BIT(91), 1); //L2P data flag
//broadcast line 6
bom[6][0] = getBits(psubfr1data, GPSL1CA_BIT(73), 4); //URA index
bom[6][1] = getBits(psubfr1data, GPSL1CA_BIT(77), 6); //SV health
bom[6][2] = getTwosComplement(getBits(psubfr1data, GPSL1CA_BIT(197), 8), 8); //TGD
bom[6][3] = (getBits(psubfr1data, GPSL1CA_BIT(83), 2) << 8) | getBits(psubfr1data, GPSL1CA_BIT(211), 8); //IODC
//broadcast line 7
bom[7][0] = getBits(psubfr1data, GPSL1CA_BIT(31), 17) *6 *100; //Transmission time of message:
// the 17 MSB of the Zcount in HOW converted to sec and scaled by 100
bom[7][1] = getBits(psubfr2data, GPSL1CA_BIT(287), 1); //Fit interval flag
bom[BO_LIN_IONOA][0] = getTwosComplement(getBits(psubfr4data, GPSL1CA_BIT(69), 8), 8); //alfa0
bom[BO_LIN_IONOA][1] = getTwosComplement(getBits(psubfr4data, GPSL1CA_BIT(77), 8), 8); //alfa1
bom[BO_LIN_IONOA][2] = getTwosComplement(getBits(psubfr4data, GPSL1CA_BIT(91), 8), 8); //alfa2
bom[BO_LIN_IONOA][3] = getTwosComplement(getBits(psubfr4data, GPSL1CA_BIT(99), 8), 8); //alfa3
bom[BO_LIN_IONOB][0] = getTwosComplement(getBits(psubfr4data, GPSL1CA_BIT(107), 8), 8); //beta0
bom[BO_LIN_IONOB][1] = getTwosComplement(getBits(psubfr4data, GPSL1CA_BIT(121), 8), 8); //beta1
bom[BO_LIN_IONOB][2] = getTwosComplement(getBits(psubfr4data, GPSL1CA_BIT(129), 8), 8); //beta2
bom[BO_LIN_IONOB][3] = getTwosComplement(getBits(psubfr4data, GPSL1CA_BIT(137), 8), 8); //beta3
bom[BO_LIN_TIMEU][0] = (getBits(psubfr4data, GPSL1CA_BIT(181), 24) << 8) | getBits(psubfr4data, GPSL1CA_BIT(211), 8); //Time correction A0
bom[BO_LIN_TIMEU][1] = getTwosComplement(getBits(psubfr4data, GPSL1CA_BIT(151), 24), 24); //Time correction A1
bom[BO_LIN_TIMEU][2] = getBits(psubfr4data, GPSL1CA_BIT(219), 8); //t0t, ref time for for UTC data
bom[BO_LIN_TIMEU][3] = getBits(psubfr4data, GPSL1CA_BIT(227), 8); //WNt, UTC ref week number
bom[BO_LIN_TIMEU][3] |= bom[5][2] & (~MASK8b); //put WNt as a continuous week number
bom[BO_LIN_TIMEG][2] = getBits(psubfr4data, GPSL1CA_BIT(31), 17) * 6; //Transmission time of message: the 17 MSB of the Zcount in HOw
bom[BO_LIN_TIMEG][3] = bom[5][2]; //the week number when message was transmitted
bom[BO_LIN_LEAPS][0] = getTwosComplement(getBits(psubfr4data, GPSL1CA_BIT(241), 8), 8); //delta tLS (leap seconds)
bom[BO_LIN_LEAPS][1] = getTwosComplement(getBits(psubfr4data, GPSL1CA_BIT(271), 8), 8); //delta tLSF
bom[BO_LIN_LEAPS][2] = getBits(psubfr4data, GPSL1CA_BIT(249), 8); //WN_LSF
bom[BO_LIN_LEAPS][2] |= bom[5][2] & (~MASK8b); //put WN_LSF as a continuous week number
bom[BO_LIN_LEAPS][3] = getBits(psubfr4data, GPSL1CA_BIT(257), 8); //DN_LSF
}