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

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