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

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