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

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