source: ntrip/trunk/BNC/src/rinex/reqcanalyze.cpp@ 6532

Last change on this file since 6532 was 6532, checked in by mervart, 9 years ago
File size: 33.5 KB
Line 
1// Part of BNC, a utility for retrieving decoding and
2// converting GNSS data streams from NTRIP broadcasters.
3//
4// Copyright (C) 2007
5// German Federal Agency for Cartography and Geodesy (BKG)
6// http://www.bkg.bund.de
7// Czech Technical University Prague, Department of Geodesy
8// http://www.fsv.cvut.cz
9//
10// Email: euref-ip@bkg.bund.de
11//
12// This program is free software; you can redistribute it and/or
13// modify it under the terms of the GNU General Public License
14// as published by the Free Software Foundation, version 2.
15//
16// This program is distributed in the hope that it will be useful,
17// but WITHOUT ANY WARRANTY; without even the implied warranty of
18// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19// GNU General Public License for more details.
20//
21// You should have received a copy of the GNU General Public License
22// along with this program; if not, write to the Free Software
23// Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
24
25/* -------------------------------------------------------------------------
26 * BKG NTRIP Client
27 * -------------------------------------------------------------------------
28 *
29 * Class: t_reqcAnalyze
30 *
31 * Purpose: Analyze RINEX Files
32 *
33 * Author: L. Mervart
34 *
35 * Created: 11-Apr-2012
36 *
37 * Changes:
38 *
39 * -----------------------------------------------------------------------*/
40
41#include <iostream>
42#include <iomanip>
43#include <qwt_plot_renderer.h>
44
45#include "reqcanalyze.h"
46#include "bnccore.h"
47#include "bncsettings.h"
48#include "reqcedit.h"
49#include "bncutils.h"
50#include "graphwin.h"
51#include "polarplot.h"
52#include "availplot.h"
53#include "eleplot.h"
54#include "dopplot.h"
55#include "bncephuser.h"
56
57using namespace std;
58
59// Constructor
60////////////////////////////////////////////////////////////////////////////
61t_reqcAnalyze::t_reqcAnalyze(QObject* parent) : QThread(parent) {
62
63 bncSettings settings;
64
65 _logFileName = settings.value("reqcOutLogFile").toString(); expandEnvVar(_logFileName);
66 _logFile = 0;
67 _log = 0;
68 _currEpo = 0;
69 _obsFileNames = settings.value("reqcObsFile").toString().split(",", QString::SkipEmptyParts);
70 _navFileNames = settings.value("reqcNavFile").toString().split(",", QString::SkipEmptyParts);
71
72 connect(this, SIGNAL(dspSkyPlot(const QString&, const QString&, QVector<t_polarPoint*>*,
73 const QString&, QVector<t_polarPoint*>*,
74 const QByteArray&, double)),
75 this, SLOT(slotDspSkyPlot(const QString&, const QString&, QVector<t_polarPoint*>*,
76 const QString&, QVector<t_polarPoint*>*,
77 const QByteArray&, double)));
78
79 connect(this, SIGNAL(dspAvailPlot(const QString&, const QByteArray&)),
80 this, SLOT(slotDspAvailPlot(const QString&, const QByteArray&)));
81}
82
83// Destructor
84////////////////////////////////////////////////////////////////////////////
85t_reqcAnalyze::~t_reqcAnalyze() {
86 for (int ii = 0; ii < _rnxObsFiles.size(); ii++) {
87 delete _rnxObsFiles[ii];
88 }
89 for (int ii = 0; ii < _ephs.size(); ii++) {
90 delete _ephs[ii];
91 }
92 delete _log; _log = 0;
93 delete _logFile; _logFile = 0;
94 if (BNC_CORE->mode() != t_bncCore::interactive) {
95 qApp->exit(0);
96 }
97}
98
99//
100////////////////////////////////////////////////////////////////////////////
101void t_reqcAnalyze::run() {
102
103 static const double QC_FORMAT_VERSION = 1.0;
104
105 // Open Log File
106 // -------------
107 _logFile = new QFile(_logFileName);
108 if (_logFile->open(QIODevice::WriteOnly | QIODevice::Text)) {
109 _log = new QTextStream();
110 _log->setDevice(_logFile);
111 }
112 *_log << "QC Format Version : " << QString("%1").arg(QC_FORMAT_VERSION,3,'f',1) << endl << endl;
113
114 // Check Ephemerides
115 // -----------------
116 checkEphemerides();
117
118 // Initialize RINEX Observation Files
119 // ----------------------------------
120 t_reqcEdit::initRnxObsFiles(_obsFileNames, _rnxObsFiles, _log);
121
122 // Read Ephemerides
123 // ----------------
124 t_reqcEdit::readEphemerides(_navFileNames, _ephs);
125
126 // Loop over all RINEX Files
127 // -------------------------
128 for (int ii = 0; ii < _rnxObsFiles.size(); ii++) {
129 analyzeFile(_rnxObsFiles[ii]);
130 }
131
132 // Exit
133 // ----
134 emit finished();
135 deleteLater();
136}
137
138//
139////////////////////////////////////////////////////////////////////////////
140void t_reqcAnalyze::analyzeFile(t_rnxObsFile* obsFile) {
141
142 _qcFile.clear();
143
144 // A priori Coordinates
145 // --------------------
146 ColumnVector xyzSta = obsFile->xyz();
147
148 // Loop over all Epochs
149 // --------------------
150 try {
151 QMap<QString, bncTime> lastObsTime;
152 bool firstEpo = true;
153 while ( (_currEpo = obsFile->nextEpoch()) != 0) {
154 if (firstEpo) {
155 firstEpo = false;
156 _qcFile._startTime = _currEpo->tt;
157 _qcFile._antennaName = obsFile->antennaName();
158 _qcFile._markerName = obsFile->markerName();
159 _qcFile._receiverType = obsFile->receiverType();
160 _qcFile._interval = obsFile->interval();
161 }
162 _qcFile._endTime = _currEpo->tt;
163
164 t_qcEpo qcEpo;
165 qcEpo._epoTime = _currEpo->tt;
166 qcEpo._PDOP = cmpDOP(xyzSta);
167
168 // Loop over all satellites
169 // ------------------------
170 for (unsigned iObs = 0; iObs < _currEpo->rnxSat.size(); iObs++) {
171 const t_rnxObsFile::t_rnxSat& rnxSat = _currEpo->rnxSat[iObs];
172 t_satObs satObs;
173 t_rnxObsFile::setObsFromRnx(obsFile, _currEpo, rnxSat, satObs);
174 t_qcSat& qcSat = qcEpo._qcSat[satObs._prn];
175 setQcObs(qcEpo._epoTime, xyzSta, satObs, lastObsTime, qcSat);
176 updateQcSat(qcSat, _qcFile._qcSatSum[satObs._prn]);
177 }
178 _qcFile._qcEpo.push_back(qcEpo);
179 }
180
181 analyzeMultipath();
182
183 preparePlotData(obsFile);
184
185 printReport(obsFile);
186 }
187 catch (QString str) {
188 if (_log) {
189 *_log << "Exception " << str << endl;
190 }
191 else {
192 qDebug() << str;
193 }
194 }
195}
196
197// Compute Dilution of Precision
198////////////////////////////////////////////////////////////////////////////
199double t_reqcAnalyze::cmpDOP(const ColumnVector& xyzSta) const {
200
201 if ( xyzSta.size() != 3 || (xyzSta[0] == 0.0 && xyzSta[1] == 0.0 && xyzSta[2] == 0.0) ) {
202 return 0.0;
203 }
204
205 unsigned nSat = _currEpo->rnxSat.size();
206
207 if (nSat < 4) {
208 return 0.0;
209 }
210
211 Matrix AA(nSat, 4);
212
213 unsigned nSatUsed = 0;
214 for (unsigned iSat = 0; iSat < nSat; iSat++) {
215
216 const t_rnxObsFile::t_rnxSat& rnxSat = _currEpo->rnxSat[iSat];
217 const t_prn& prn = rnxSat.prn;
218
219 t_eph* eph = 0;
220 for (int ie = 0; ie < _ephs.size(); ie++) {
221 if (_ephs[ie]->prn() == prn) {
222 eph = _ephs[ie];
223 break;
224 }
225 }
226 if (eph) {
227 ColumnVector xSat(4);
228 ColumnVector vv(3);
229 if (eph->getCrd(_currEpo->tt, xSat, vv, false) == success) {
230 ++nSatUsed;
231 ColumnVector dx = xSat.Rows(1,3) - xyzSta;
232 double rho = dx.norm_Frobenius();
233 AA(nSatUsed,1) = dx(1) / rho;
234 AA(nSatUsed,2) = dx(2) / rho;
235 AA(nSatUsed,3) = dx(3) / rho;
236 AA(nSatUsed,4) = 1.0;
237 }
238 }
239 }
240
241 if (nSatUsed < 4) {
242 return 0.0;
243 }
244
245 AA = AA.Rows(1, nSatUsed);
246
247 SymmetricMatrix QQ;
248 QQ << AA.t() * AA;
249 QQ = QQ.i();
250
251 return sqrt(QQ.trace());
252}
253
254//
255////////////////////////////////////////////////////////////////////////////
256void t_reqcAnalyze::updateQcSat(const t_qcSat& qcSat, t_qcSatSum& qcSatSum) {
257
258 for (int ii = 0; ii < qcSat._qcFrq.size(); ii++) {
259 const t_qcFrq& qcFrq = qcSat._qcFrq[ii];
260 t_qcFrqSum& qcFrqSum = qcSatSum._qcFrqSum[qcFrq._rnxType2ch];
261 qcFrqSum._numObs += 1;
262 if (qcFrq._slip) {
263 qcFrqSum._numSlipsFlagged += 1;
264 }
265 if (qcFrq._gap) {
266 qcFrqSum._numGaps += 1;
267 }
268 if (qcFrq._SNR > 0.0) {
269 qcFrqSum._numSNR += 1;
270 qcFrqSum._sumSNR += qcFrq._SNR;
271 }
272 }
273}
274
275//
276////////////////////////////////////////////////////////////////////////////
277void t_reqcAnalyze::setQcObs(const bncTime& epoTime, const ColumnVector& xyzSta,
278 const t_satObs& satObs, QMap<QString, bncTime>& lastObsTime,
279 t_qcSat& qcSat) {
280
281 t_eph* eph = 0;
282 for (int ie = 0; ie < _ephs.size(); ie++) {
283 if (_ephs[ie]->prn() == satObs._prn) {
284 eph = _ephs[ie];
285 break;
286 }
287 }
288 if (eph) {
289 ColumnVector xc(4);
290 ColumnVector vv(3);
291 if ( xyzSta.size() == 3 && (xyzSta[0] != 0.0 || xyzSta[1] != 0.0 || xyzSta[2] != 0.0) &&
292 eph->getCrd(epoTime, xc, vv, false) == success) {
293 double rho, eleSat, azSat;
294 topos(xyzSta(1), xyzSta(2), xyzSta(3), xc(1), xc(2), xc(3), rho, eleSat, azSat);
295 qcSat._eleSet = true;
296 qcSat._azDeg = azSat * 180.0/M_PI;
297 qcSat._eleDeg = eleSat * 180.0/M_PI;
298 }
299 if (satObs._prn.system() == 'R') {
300 qcSat._slotSet = true;
301 qcSat._slotNum = eph->slotNum();
302 }
303 }
304
305 // Availability and Slip Flags
306 // ---------------------------
307 for (unsigned ii = 0; ii < satObs._obs.size(); ii++) {
308 const t_frqObs* frqObs = satObs._obs[ii];
309
310 qcSat._qcFrq.push_back(t_qcFrq());
311 t_qcFrq& qcFrq = qcSat._qcFrq.back();
312
313 qcFrq._rnxType2ch = QString(frqObs->_rnxType2ch.c_str());
314 qcFrq._SNR = frqObs->_snr;
315 qcFrq._slip = frqObs->_slip;
316 qcFrq._phaseValid = frqObs->_phaseValid;
317 qcFrq._codeValid = frqObs->_codeValid;
318 // Check Gaps
319 // ----------
320 QString key = QString(satObs._prn.toString().c_str()) + qcFrq._rnxType2ch;
321 if (lastObsTime[key].valid()) {
322 double dt = epoTime - lastObsTime[key];
323 if (dt > 1.5 * _qcFile._interval) {
324 qcFrq._gap = true;
325 }
326 }
327 lastObsTime[key] = epoTime;
328
329 // Compute the Multipath Linear Combination
330 // ----------------------------------------
331 if (frqObs->_codeValid) {
332 t_frequency::type fA;
333 t_frequency::type fB;
334 if (satObs._prn.system() == 'G') {
335 if (frqObs->_rnxType2ch[0] == '1') {
336 fA = t_frequency::G1;
337 fB = t_frequency::G2;
338 }
339 else if (frqObs->_rnxType2ch[0] == '2') {
340 fA = t_frequency::G2;
341 fB = t_frequency::G1;
342 }
343 }
344 else if (satObs._prn.system() == 'R') {
345 if (frqObs->_rnxType2ch[0] == '1') {
346 fA = t_frequency::R1;
347 fB = t_frequency::R2;
348 }
349 else if (frqObs->_rnxType2ch[0] == '2') {
350 fA = t_frequency::R2;
351 fB = t_frequency::R1;
352 }
353 }
354 else if (satObs._prn.system() == 'E') {
355 if (frqObs->_rnxType2ch[0] == '1') {
356 fA = t_frequency::E1;
357 fB = t_frequency::E5;
358 }
359 else if (frqObs->_rnxType2ch[0] == '5') {
360 fA = t_frequency::E5;
361 fB = t_frequency::E1;
362 }
363 }
364 else if (satObs._prn.system() == 'J') {
365 if (frqObs->_rnxType2ch[0] == '1') {
366 fA = t_frequency::J1;
367 fB = t_frequency::J2;
368 }
369 else if (frqObs->_rnxType2ch[0] == '2') {
370 fA = t_frequency::J2;
371 fB = t_frequency::J1;
372 }
373 }
374 else if (satObs._prn.system() == 'S') {
375 if (frqObs->_rnxType2ch[0] == '1') {
376 fA = t_frequency::S1;
377 fB = t_frequency::S5;
378 }
379 else if (frqObs->_rnxType2ch[0] == '5') {
380 fA = t_frequency::S5;
381 fB = t_frequency::S1;
382 }
383 }
384 else if (satObs._prn.system() == 'C') {
385 if (frqObs->_rnxType2ch[0] == '1') {
386 fA = t_frequency::C1;
387 fB = t_frequency::C7;
388 }
389 else if (frqObs->_rnxType2ch[0] == '7') {
390 fA = t_frequency::C7;
391 fB = t_frequency::C1;
392 }
393 }
394 if (fA != t_frequency::dummy && fB != t_frequency::dummy) {
395 double f_a = t_CST::freq(fA, qcSat._slotNum);
396 double f_b = t_CST::freq(fB, qcSat._slotNum);
397 double C_a = frqObs->_code;
398
399 bool foundA = false;
400 double L_a = 0.0;
401 bool foundB = false;
402 double L_b = 0.0;
403 for (unsigned jj = 0; jj < satObs._obs.size(); jj++) {
404 const t_frqObs* frqObsHlp = satObs._obs[jj];
405 if (frqObsHlp->_rnxType2ch[0] == t_frequency::toString(fA)[1] &&
406 frqObsHlp->_phaseValid) {
407 foundA = true;
408 L_a = frqObsHlp->_phase * t_CST::c / f_a;
409 }
410 else if (frqObsHlp->_rnxType2ch[0] == t_frequency::toString(fB)[1] &&
411 frqObsHlp->_phaseValid) {
412 foundB = true;
413 L_b = frqObsHlp->_phase * t_CST::c / f_b;
414 }
415 }
416 if (foundA && foundB) {
417 qcFrq._setMP = true;
418 qcFrq._rawMP = C_a - L_a - 2.0*f_b*f_b/(f_a*f_a-f_b*f_b) * (L_a - L_b);
419 }
420 }
421 }
422 } // satObs loop
423}
424
425//
426////////////////////////////////////////////////////////////////////////////
427void t_reqcAnalyze::analyzeMultipath() {
428
429 const double SLIPTRESH = 10.0; // cycle-slip threshold (meters)
430 const double chunkStep = 600.0; // 10 minutes
431
432 // Loop over all satellites available
433 // ----------------------------------
434 QMutableMapIterator<t_prn, t_qcSatSum> itSat(_qcFile._qcSatSum);
435 while (itSat.hasNext()) {
436 itSat.next();
437 const t_prn& prn = itSat.key();
438 t_qcSatSum& qcSatSum = itSat.value();
439
440 // Loop over all frequencies available
441 // -----------------------------------
442 QMutableMapIterator<QString, t_qcFrqSum> itFrq(qcSatSum._qcFrqSum);
443 while (itFrq.hasNext()) {
444 itFrq.next();
445 const QString& frqType = itFrq.key();
446 t_qcFrqSum& qcFrqSum = itFrq.value();
447
448 // Loop over all Chunks of Data
449 // ----------------------------
450 for (bncTime chunkStart = _qcFile._startTime;
451 chunkStart < _qcFile._endTime; chunkStart += chunkStep) {
452
453 bncTime chunkEnd = chunkStart + chunkStep;
454
455 QVector<t_qcFrq*> frqVec;
456 QVector<double> MP;
457
458 // Loop over all Epochs within one Chunk of Data
459 // ---------------------------------------------
460 for (int iEpo = 0; iEpo < _qcFile._qcEpo.size(); iEpo++) {
461 t_qcEpo& qcEpo = _qcFile._qcEpo[iEpo];
462 if (chunkStart <= qcEpo._epoTime && qcEpo._epoTime < chunkEnd) {
463 if (qcEpo._qcSat.contains(prn)) {
464 t_qcSat& qcSat = qcEpo._qcSat[prn];
465 for (int iFrq = 0; iFrq < qcSat._qcFrq.size(); iFrq++) {
466 t_qcFrq& qcFrq = qcSat._qcFrq[iFrq];
467 if (qcFrq._rnxType2ch == frqType) {
468 frqVec << &qcFrq;
469 if (qcFrq._setMP) {
470 MP << qcFrq._rawMP;
471 }
472 }
473 }
474 }
475 }
476 }
477
478 // Compute the multipath mean and standard deviation
479 // -------------------------------------------------
480 if (MP.size() > 1) {
481 double meanMP = 0.0;
482 for (int ii = 0; ii < MP.size(); ii++) {
483 meanMP += MP[ii];
484 }
485 meanMP /= MP.size();
486
487 bool slipMP = false;
488
489 double stdMP = 0.0;
490 for (int ii = 0; ii < MP.size(); ii++) {
491 double diff = MP[ii] - meanMP;
492 if (fabs(diff) > SLIPTRESH) {
493 slipMP = true;
494 break;
495 }
496 stdMP += diff * diff;
497 }
498
499 if (slipMP) {
500 stdMP = 0.0;
501 stdMP = 0.0;
502 qcFrqSum._numSlipsFound += 1;
503 }
504 else {
505 stdMP = sqrt(stdMP / (MP.size()-1));
506 qcFrqSum._numMP += 1;
507 qcFrqSum._sumMP += stdMP;
508 }
509
510 for (int ii = 0; ii < frqVec.size(); ii++) {
511 t_qcFrq* qcFrq = frqVec[ii];
512 if (slipMP) {
513 qcFrq->_slip = true;
514 }
515 else {
516 qcFrq->_stdMP = stdMP;
517 }
518 }
519 }
520 } // chunk loop
521 } // frq loop
522 } // sat loop
523}
524
525//
526////////////////////////////////////////////////////////////////////////////
527void t_reqcAnalyze::preparePlotData(const t_rnxObsFile* obsFile) {
528
529 QString mp1Title = "Multipath\n";
530 QString mp2Title = "Multipath\n";
531 QString sn1Title = "Signal-to-Noise Ratio\n";
532 QString sn2Title = "Signal-to-Noise Ratio\n";
533
534 bncSettings settings;
535
536 QString reqcSkyPlotSignals = settings.value("reqcSkyPlotSignals").toString();
537 if (reqcSkyPlotSignals.isEmpty()) {
538 reqcSkyPlotSignals = "G:1&2 R:1&2 J:1&2 E:1&7 C:1&6";
539 }
540
541 QStringList signalsOpt = reqcSkyPlotSignals.split(" ", QString::SkipEmptyParts);
542 QMap<char, QString> signalTypes1;
543 QMap<char, QString> signalTypes2;
544 for (int ii = 0; ii < signalsOpt.size(); ii++) {
545 QStringList hlp = signalsOpt.at(ii).split(QRegExp("[:&]"), QString::SkipEmptyParts);
546 if (hlp.size() > 1 && hlp[0].length() == 1) {
547 char system = hlp[0].toAscii().constData()[0];
548 signalTypes1[system] = hlp[1];
549 mp1Title += hlp[0] + ":" + hlp[1] + " ";
550 sn1Title += hlp[0] + ":" + hlp[1] + " ";
551 if (hlp.size() > 2) {
552 signalTypes2[system] = hlp[2];
553 mp2Title += hlp[0] + ":" + hlp[2] + " ";
554 sn2Title += hlp[0] + ":" + hlp[2] + " ";
555 }
556 }
557 }
558
559 QVector<t_polarPoint*>* dataMP1 = new QVector<t_polarPoint*>;
560 QVector<t_polarPoint*>* dataMP2 = new QVector<t_polarPoint*>;
561 QVector<t_polarPoint*>* dataSNR1 = new QVector<t_polarPoint*>;
562 QVector<t_polarPoint*>* dataSNR2 = new QVector<t_polarPoint*>;
563
564 // Loop over all observations
565 // --------------------------
566 for (int iEpo = 0; iEpo < _qcFile._qcEpo.size(); iEpo++) {
567 t_qcEpo& qcEpo = _qcFile._qcEpo[iEpo];
568 QMapIterator<t_prn, t_qcSat> it(qcEpo._qcSat);
569 while (it.hasNext()) {
570 it.next();
571 const t_prn& prn = it.key();
572 const t_qcSat& qcSat = it.value();
573 if (qcSat._eleSet) {
574
575 QString frqType[2];
576
577 for (int iFrq = 0; iFrq < qcSat._qcFrq.size(); iFrq++) {
578 const t_qcFrq& qcFrq = qcSat._qcFrq[iFrq];
579
580 for (int ii = 0; ii < 2; ii++) {
581 const QMap<char, QString>& signalTypes = (ii == 0 ? signalTypes1 : signalTypes2);
582 if (frqType[ii].isEmpty()) {
583 QMapIterator<char, QString> it(signalTypes);
584 while (it.hasNext()) {
585 it.next();
586 if (it.key() == prn.system()) {
587 if (it.value() == qcFrq._rnxType2ch || it.value() == qcFrq._rnxType2ch.left(1)) {
588 frqType[ii] = qcFrq._rnxType2ch;
589 break;
590 }
591 }
592 }
593 }
594 }
595 if (qcFrq._rnxType2ch == frqType[0]) {
596 (*dataSNR1) << (new t_polarPoint(qcSat._azDeg, 90.0 - qcSat._eleDeg, qcFrq._SNR));
597 (*dataMP1) << (new t_polarPoint(qcSat._azDeg, 90.0 - qcSat._eleDeg, qcFrq._stdMP));
598 }
599 else if (qcFrq._rnxType2ch == frqType[1]) {
600 (*dataSNR2) << (new t_polarPoint(qcSat._azDeg, 90.0 - qcSat._eleDeg, qcFrq._SNR));
601 (*dataMP2) << (new t_polarPoint(qcSat._azDeg, 90.0 - qcSat._eleDeg, qcFrq._stdMP));
602 }
603 }
604 }
605 }
606 }
607
608 // Show the plots
609 // --------------
610 if (BNC_CORE->GUIenabled()) {
611 QFileInfo fileInfo(obsFile->fileName());
612 QByteArray title = fileInfo.fileName().toAscii();
613 emit dspSkyPlot(obsFile->fileName(), mp1Title, dataMP1, mp2Title, dataMP2, "Meters", 2.0);
614 emit dspSkyPlot(obsFile->fileName(), sn1Title, dataSNR1, sn2Title, dataSNR2, "dbHz", 54.0);
615 emit dspAvailPlot(obsFile->fileName(), title);
616 }
617 else {
618 for (int ii = 0; ii < dataMP1->size(); ii++) {
619 delete dataMP1->at(ii);
620 }
621 delete dataMP1;
622 for (int ii = 0; ii < dataMP2->size(); ii++) {
623 delete dataMP2->at(ii);
624 }
625 delete dataMP2;
626 for (int ii = 0; ii < dataSNR1->size(); ii++) {
627 delete dataSNR1->at(ii);
628 }
629 delete dataSNR1;
630 for (int ii = 0; ii < dataSNR2->size(); ii++) {
631 delete dataSNR2->at(ii);
632 }
633 delete dataSNR2;
634 }
635}
636
637//
638////////////////////////////////////////////////////////////////////////////
639void t_reqcAnalyze::slotDspSkyPlot(const QString& fileName, const QString& title1,
640 QVector<t_polarPoint*>* data1, const QString& title2,
641 QVector<t_polarPoint*>* data2, const QByteArray& scaleTitle,
642 double maxValue) {
643
644 if (BNC_CORE->GUIenabled()) {
645
646 if (maxValue == 0.0) {
647 if (data1) {
648 for (int ii = 0; ii < data1->size(); ii++) {
649 double val = data1->at(ii)->_value;
650 if (maxValue < val) {
651 maxValue = val;
652 }
653 }
654 }
655 if (data2) {
656 for (int ii = 0; ii < data2->size(); ii++) {
657 double val = data2->at(ii)->_value;
658 if (maxValue < val) {
659 maxValue = val;
660 }
661 }
662 }
663 }
664
665 QwtInterval scaleInterval(0.0, maxValue);
666
667 QVector<QWidget*> plots;
668 if (data1) {
669 QwtText title(title1);
670 QFont font = title.font(); font.setPointSize(font.pointSize()-1); title.setFont(font);
671 t_polarPlot* plot1 = new t_polarPlot(title, scaleInterval, BNC_CORE->mainWindow());
672 plot1->addCurve(data1);
673 plots << plot1;
674 }
675 if (data2) {
676 QwtText title(title2);
677 QFont font = title.font(); font.setPointSize(font.pointSize()-1); title.setFont(font);
678 t_polarPlot* plot2 = new t_polarPlot(title, scaleInterval, BNC_CORE->mainWindow());
679 plot2->addCurve(data2);
680 plots << plot2;
681 }
682
683 t_graphWin* graphWin = new t_graphWin(0, fileName, plots,
684 &scaleTitle, &scaleInterval);
685
686 graphWin->show();
687
688 bncSettings settings;
689 QString dirName = settings.value("reqcPlotDir").toString();
690 if (!dirName.isEmpty()) {
691 QByteArray ext = (scaleTitle == "Meters") ? "_M.png" : "_S.png";
692 graphWin->savePNG(dirName, ext);
693 }
694 }
695}
696
697//
698////////////////////////////////////////////////////////////////////////////
699void t_reqcAnalyze::slotDspAvailPlot(const QString& fileName, const QByteArray& title) {
700
701 t_plotData plotData;
702 QMap<t_prn, t_plotData> plotDataMap;
703
704 for (int ii = 0; ii < _qcFile._qcEpo.size(); ii++) {
705 const t_qcEpo& qcEpo = _qcFile._qcEpo[ii];
706 double mjdX24 = qcEpo._epoTime.mjddec() * 24.0;
707
708 plotData._mjdX24 << mjdX24;
709 plotData._PDOP << qcEpo._PDOP;
710 plotData._numSat << qcEpo._qcSat.size();
711
712 QMapIterator<t_prn, t_qcSat> it(qcEpo._qcSat);
713 while (it.hasNext()) {
714 it.next();
715 const t_prn& prn = it.key();
716 const t_qcSat& qcSat = it.value();
717 t_plotData& data = plotDataMap[prn];
718
719 if (qcSat._eleSet) {
720 data._mjdX24 << mjdX24;
721 data._eleDeg << qcSat._eleDeg;
722 }
723
724 char frqChar1 = '1';
725 char frqChar2 = '2';
726 if (prn.system() == 'E' || prn.system() == 'S') {
727 frqChar2 = '5';
728 }
729 else if (prn.system() == 'C') {
730 frqChar2 = '7';
731 }
732 QString frqType1;
733 QString frqType2;
734 for (int iFrq = 0; iFrq < qcSat._qcFrq.size(); iFrq++) {
735 const t_qcFrq& qcFrq = qcSat._qcFrq[iFrq];
736 if (qcFrq._rnxType2ch[0] == frqChar1 && frqType1.isEmpty()) {
737 frqType1 = qcFrq._rnxType2ch;
738 }
739 if (qcFrq._rnxType2ch[0] == frqChar2 && frqType2.isEmpty()) {
740 frqType2 = qcFrq._rnxType2ch;
741 }
742 if (qcFrq._rnxType2ch == frqType1) {
743 if (qcFrq._slip) {
744 data._L1slip << mjdX24;
745 }
746 else if (qcFrq._gap) {
747 data._L1gap << mjdX24;
748 }
749 else {
750 data._L1ok << mjdX24;
751 }
752 }
753 else if (qcFrq._rnxType2ch == frqType2) {
754 if (qcFrq._slip) {
755 data._L2slip << mjdX24;
756 }
757 else if (qcFrq._gap) {
758 data._L2gap << mjdX24;
759 }
760 else {
761 data._L2ok << mjdX24;
762 }
763 }
764 }
765 }
766 }
767
768 if (BNC_CORE->GUIenabled()) {
769 t_availPlot* plotA = new t_availPlot(0, plotDataMap);
770 plotA->setTitle(title);
771
772 t_elePlot* plotZ = new t_elePlot(0, plotDataMap);
773
774 t_dopPlot* plotD = new t_dopPlot(0, plotData);
775
776 QVector<QWidget*> plots;
777 plots << plotA << plotZ << plotD;
778 t_graphWin* graphWin = new t_graphWin(0, fileName, plots, 0, 0);
779
780 int ww = QFontMetrics(graphWin->font()).width('w');
781 graphWin->setMinimumSize(120*ww, 40*ww);
782
783 graphWin->show();
784
785 bncSettings settings;
786 QString dirName = settings.value("reqcPlotDir").toString();
787 if (!dirName.isEmpty()) {
788 QByteArray ext = "_A.png";
789 graphWin->savePNG(dirName, ext);
790 }
791 }
792}
793
794// Finish the report
795////////////////////////////////////////////////////////////////////////////
796void t_reqcAnalyze::printReport(const t_rnxObsFile* obsFile) {
797
798 if (!_log) {
799 return;
800 }
801
802 QFileInfo obsFi(obsFile->fileName());
803 QString obsFileName = obsFi.fileName();
804
805 // Summary
806 // -------
807 *_log << "Observation File : " << obsFileName << endl
808 << "RINEX Version : " << QString("%1").arg(obsFile->version(),4,'f',2) << endl
809 << "Marker Name : " << _qcFile._markerName << endl
810 << "Marker Number : " << obsFile->markerNumber() << endl
811 << "Receiver : " << _qcFile._receiverType << endl
812 << "Antenna : " << _qcFile._antennaName << endl
813 << "Position XYZ : " << QString("%1 %2 %3").arg(obsFile->xyz()(1), 14, 'f', 4)
814 .arg(obsFile->xyz()(2), 14, 'f', 4)
815 .arg(obsFile->xyz()(3), 14, 'f', 4) << endl
816 << "Antenna dH/dE/dN : " << QString("%1 %2 %3").arg(obsFile->antNEU()(3), 8, 'f', 4)
817 .arg(obsFile->antNEU()(2), 8, 'f', 4)
818 .arg(obsFile->antNEU()(1), 8, 'f', 4) << endl
819 << "Start Time : " << _qcFile._startTime.datestr().c_str() << ' '
820 << _qcFile._startTime.timestr(1,'.').c_str() << endl
821 << "End Time : " << _qcFile._endTime.datestr().c_str() << ' '
822 << _qcFile._endTime.timestr(1,'.').c_str() << endl
823 << "Interval : " << _qcFile._interval << endl;
824
825 // Number of systems
826 // -----------------
827 QMap<QChar, QVector<const t_qcSatSum*> > systemMap;
828 QMapIterator<t_prn, t_qcSatSum> itSat(_qcFile._qcSatSum);
829 while (itSat.hasNext()) {
830 itSat.next();
831 const t_prn& prn = itSat.key();
832 const t_qcSatSum& qcSatSum = itSat.value();
833 systemMap[prn.system()].push_back(&qcSatSum);
834 }
835
836 *_log << "Navigation Systems: " << systemMap.size() << " ";
837 QMapIterator<QChar, QVector<const t_qcSatSum*> > itSys(systemMap);
838 while (itSys.hasNext()) {
839 itSys.next();
840 *_log << ' ' << itSys.key();
841 }
842 *_log << endl;
843
844 itSys.toFront();
845 while (itSys.hasNext()) {
846 itSys.next();
847 const QChar& sys = itSys.key();
848 const QVector<const t_qcSatSum*>& qcSatVec = itSys.value();
849 QString prefixSys = QString(" ") + sys + QString(": ");
850 QMap<QString, QVector<const t_qcFrqSum*> > frqMap;
851 for (int ii = 0; ii < qcSatVec.size(); ii++) {
852 const t_qcSatSum* qcSatSum = qcSatVec[ii];
853 QMapIterator<QString, t_qcFrqSum> itFrq(qcSatSum->_qcFrqSum);
854 while (itFrq.hasNext()) {
855 itFrq.next();
856 QString frqType = itFrq.key(); if (frqType.length() < 2) frqType += '?';
857 const t_qcFrqSum& qcFrqSum = itFrq.value();
858 frqMap[frqType].push_back(&qcFrqSum);
859 }
860 }
861 *_log << endl
862 << prefixSys << "Satellites: " << qcSatVec.size() << endl
863 << prefixSys << "Signals : " << frqMap.size() << " ";
864 QMapIterator<QString, QVector<const t_qcFrqSum*> > itFrq(frqMap);
865 while (itFrq.hasNext()) {
866 itFrq.next();
867 QString frqType = itFrq.key(); if (frqType.length() < 2) frqType += '?';
868 *_log << ' ' << frqType;
869 }
870 *_log << endl;
871 QString prefixSys2 = " " + prefixSys;
872 itFrq.toFront();
873 while (itFrq.hasNext()) {
874 itFrq.next();
875 QString frqType = itFrq.key(); if (frqType.length() < 2) frqType += '?';
876 const QVector<const t_qcFrqSum*> qcFrqVec = itFrq.value();
877 QString prefixFrq = QString(" ") + frqType + QString(": ");
878
879 int numObs = 0;
880 int numSlipsFlagged = 0;
881 int numSlipsFound = 0;
882 int numGaps = 0;
883 int numSNR = 0;
884 int numMP = 0;
885 double sumSNR = 0.0;
886 double sumMP = 0.0;
887 for (int ii = 0; ii < qcFrqVec.size(); ii++) {
888 const t_qcFrqSum* qcFrqSum = qcFrqVec[ii];
889 numObs += qcFrqSum->_numObs ;
890 numSlipsFlagged += qcFrqSum->_numSlipsFlagged;
891 numSlipsFound += qcFrqSum->_numSlipsFound ;
892 numGaps += qcFrqSum->_numGaps ;
893 numSNR += qcFrqSum->_numSNR;
894 numMP += qcFrqSum->_numMP;
895 sumSNR += qcFrqSum->_sumSNR;
896 sumMP += qcFrqSum->_sumMP;
897 }
898 if (numSNR > 0) {
899 sumSNR /= numSNR;
900 }
901 if (numMP > 0) {
902 sumMP /= numMP;
903 }
904 *_log << endl
905 << prefixSys2 << prefixFrq << "Observations : " << QString("%1\n").arg(numObs, 6)
906 << prefixSys2 << prefixFrq << "Slips (file+found): " << QString("%1 +").arg(numSlipsFlagged, 6)
907 << QString("%1\n").arg(numSlipsFound, 6)
908 << prefixSys2 << prefixFrq << "Gaps : " << QString("%1\n").arg(numGaps, 6)
909 << prefixSys2 << prefixFrq << "Mean SNR : " << QString("%1\n").arg(sumSNR, 6, 'f', 1)
910 << prefixSys2 << prefixFrq << "Mean Multipath : " << QString("%1\n").arg(sumMP, 6, 'f', 2);
911 }
912 }
913
914 // Epoch-Specific Output
915 // ---------------------
916 bncSettings settings;
917 if (Qt::CheckState(settings.value("reqcLogSummaryOnly").toInt()) == Qt::Checked) {
918 return;
919 }
920 *_log << endl;
921 for (int iEpo = 0; iEpo < _qcFile._qcEpo.size(); iEpo++) {
922 const t_qcEpo& qcEpo = _qcFile._qcEpo[iEpo];
923
924 unsigned year, month, day, hour, min;
925 double sec;
926 qcEpo._epoTime.civil_date(year, month, day);
927 qcEpo._epoTime.civil_time(hour, min, sec);
928
929 QString dateStr;
930 QTextStream(&dateStr) << QString("> %1 %2 %3 %4 %5%6")
931 .arg(year, 4)
932 .arg(month, 2, 10, QChar('0'))
933 .arg(day, 2, 10, QChar('0'))
934 .arg(hour, 2, 10, QChar('0'))
935 .arg(min, 2, 10, QChar('0'))
936 .arg(sec, 11, 'f', 7);
937
938 *_log << dateStr << QString(" %1").arg(qcEpo._qcSat.size(), 2)
939 << QString(" %1").arg(qcEpo._PDOP, 4, 'f', 1)
940 << endl;
941
942 QMapIterator<t_prn, t_qcSat> itSat(qcEpo._qcSat);
943 while (itSat.hasNext()) {
944 itSat.next();
945 const t_prn& prn = itSat.key();
946 const t_qcSat& qcSat = itSat.value();
947
948 *_log << prn.toString().c_str()
949 << QString(" %1 %2").arg(qcSat._eleDeg, 6, 'f', 2).arg(qcSat._azDeg, 7, 'f', 2);
950
951 int numObsTypes = 0;
952 for (int iFrq = 0; iFrq < qcSat._qcFrq.size(); iFrq++) {
953 const t_qcFrq& qcFrq = qcSat._qcFrq[iFrq];
954 if (qcFrq._phaseValid) {
955 numObsTypes += 1;
956 }
957 if (qcFrq._codeValid) {
958 numObsTypes += 1;
959 }
960 }
961 *_log << QString(" %1").arg(numObsTypes, 2);
962
963 for (int iFrq = 0; iFrq < qcSat._qcFrq.size(); iFrq++) {
964 const t_qcFrq& qcFrq = qcSat._qcFrq[iFrq];
965 if (qcFrq._phaseValid) {
966 *_log << " L" << qcFrq._rnxType2ch << ' ';
967 if (qcFrq._slip) {
968 *_log << 's';
969 }
970 else {
971 *_log << '.';
972 }
973 if (qcFrq._gap) {
974 *_log << 'g';
975 }
976 else {
977 *_log << '.';
978 }
979 *_log << QString(" %1").arg(qcFrq._SNR, 4, 'f', 1);
980 }
981 if (qcFrq._codeValid) {
982 *_log << " C" << qcFrq._rnxType2ch << ' ';
983 if (qcFrq._gap) {
984 *_log << " g";
985 }
986 else {
987 *_log << " .";
988 }
989 *_log << QString(" %1").arg(qcFrq._stdMP, 3, 'f', 2);
990 }
991 }
992 *_log << endl;
993 }
994 }
995 _log->flush();
996}
997
998//
999////////////////////////////////////////////////////////////////////////////
1000void t_reqcAnalyze::checkEphemerides() {
1001
1002 QString navFileName;
1003 QStringListIterator namIt(_navFileNames);
1004 bool firstName = true;
1005 while (namIt.hasNext()) {
1006 QFileInfo navFi(namIt.next());
1007 if (firstName) {
1008 firstName = false;
1009 navFileName += navFi.fileName();
1010 }
1011 else {
1012 navFileName += ", " + navFi.fileName();
1013 }
1014 }
1015
1016 *_log << "Navigation File(s): " << navFileName << endl;
1017
1018 QStringListIterator it(_navFileNames);
1019 while (it.hasNext()) {
1020 const QString& fileName = it.next();
1021 unsigned numOK = 0;
1022 unsigned numBad = 0;
1023 bncEphUser ephUser(false);
1024 t_rnxNavFile rnxNavFile(fileName, t_rnxNavFile::input);
1025 for (unsigned ii = 0; ii < rnxNavFile.ephs().size(); ii++) {
1026 t_eph* eph = rnxNavFile.ephs()[ii];
1027 ephUser.putNewEph(eph, true);
1028 if (eph->checkState() == t_eph::bad) {
1029 ++numBad;
1030 }
1031 else {
1032 ++numOK;
1033 }
1034 }
1035 *_log << "Ephemeris : " << numOK << " OK " << numBad << " BAD" << endl;
1036 if (numBad > 0) {
1037 for (unsigned ii = 0; ii < rnxNavFile.ephs().size(); ii++) {
1038 t_eph* eph = rnxNavFile.ephs()[ii];
1039 if (eph->checkState() == t_eph::bad) {
1040 QFileInfo navFi(fileName);
1041 *_log << " Bad Ephemeris : " << navFi.fileName() << ' '
1042 << eph->prn().toString().c_str() << ' '
1043 << eph->TOC().datestr(' ').c_str() << ' '
1044 << eph->TOC().timestr(0, ' ').c_str() << endl;
1045 }
1046 }
1047 }
1048 }
1049
1050 *_log << endl;
1051}
Note: See TracBrowser for help on using the repository browser.