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

Last change on this file since 7068 was 6867, checked in by stuerze, 10 years ago

minor changes to prevent crash, if no logfile is specified for qc

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