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

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