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

Last change on this file since 7548 was 7548, checked in by stuerze, 9 years ago

minor format changes
File size: 33.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#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 // Observation types per system
833 // -----------------------------
834 for (int iSys = 0; iSys < obsFile->numSys(); iSys++) {
835 char sys = obsFile->system(iSys);
836 if (sys != ' ') {
837 *_log << "Observation Types " << sys << ":";
838 for (int iType = 0; iType < obsFile->nTypes(sys); iType++) {
839 QString type = obsFile->obsType(sys, iType);
840 *_log << " " << type;
841 }
842 *_log << endl;
843 }
844 }
845
846 itSys.toFront();
847 while (itSys.hasNext()) {
848 itSys.next();
849 const QChar& sys = itSys.key();
850 const QVector<const t_qcSatSum*>& qcSatVec = itSys.value();
851 QString prefixSys = QString(" ") + sys + QString(": ");
852 QMap<QString, QVector<const t_qcFrqSum*> > frqMap;
853 for (int ii = 0; ii < qcSatVec.size(); ii++) {
854 const t_qcSatSum* qcSatSum = qcSatVec[ii];
855 QMapIterator<QString, t_qcFrqSum> itFrq(qcSatSum->_qcFrqSum);
856 while (itFrq.hasNext()) {
857 itFrq.next();
858 QString frqType = itFrq.key(); if (frqType.length() < 2) frqType += '?';
859 const t_qcFrqSum& qcFrqSum = itFrq.value();
860 frqMap[frqType].push_back(&qcFrqSum);
861 }
862 }
863 *_log << endl
864 << prefixSys << "Satellites: " << qcSatVec.size() << endl
865 << prefixSys << "Signals : " << frqMap.size() << " ";
866 QMapIterator<QString, QVector<const t_qcFrqSum*> > itFrq(frqMap);
867 while (itFrq.hasNext()) {
868 itFrq.next();
869 QString frqType = itFrq.key(); if (frqType.length() < 2) frqType += '?';
870 *_log << ' ' << frqType;
871 }
872 *_log << endl;
873 QString prefixSys2 = " " + prefixSys;
874 itFrq.toFront();
875 while (itFrq.hasNext()) {
876 itFrq.next();
877 QString frqType = itFrq.key(); if (frqType.length() < 2) frqType += '?';
878 const QVector<const t_qcFrqSum*> qcFrqVec = itFrq.value();
879 QString prefixFrq = QString(" ") + frqType + QString(": ");
880
881 int numObs = 0;
882 int numSlipsFlagged = 0;
883 int numSlipsFound = 0;
884 int numGaps = 0;
885 int numSNR = 0;
886 int numMP = 0;
887 double sumSNR = 0.0;
888 double sumMP = 0.0;
889 for (int ii = 0; ii < qcFrqVec.size(); ii++) {
890 const t_qcFrqSum* qcFrqSum = qcFrqVec[ii];
891 numObs += qcFrqSum->_numObs ;
892 numSlipsFlagged += qcFrqSum->_numSlipsFlagged;
893 numSlipsFound += qcFrqSum->_numSlipsFound ;
894 numGaps += qcFrqSum->_numGaps ;
895 numSNR += qcFrqSum->_numSNR;
896 numMP += qcFrqSum->_numMP;
897 sumSNR += qcFrqSum->_sumSNR;
898 sumMP += qcFrqSum->_sumMP;
899 }
900 if (numSNR > 0) {
901 sumSNR /= numSNR;
902 }
903 if (numMP > 0) {
904 sumMP /= numMP;
905 }
906 *_log << endl
907 << prefixSys2 << prefixFrq << "Observations : " << QString("%1\n").arg(numObs, 6)
908 << prefixSys2 << prefixFrq << "Slips (file+found): " << QString("%1 +").arg(numSlipsFlagged, 6)
909 << QString("%1\n").arg(numSlipsFound, 6)
910 << prefixSys2 << prefixFrq << "Gaps : " << QString("%1\n").arg(numGaps, 6)
911 << prefixSys2 << prefixFrq << "Mean SNR : " << QString("%1\n").arg(sumSNR, 6, 'f', 1)
912 << prefixSys2 << prefixFrq << "Mean Multipath : " << QString("%1\n").arg(sumMP, 6, 'f', 2);
913 }
914 }
915
916 // Epoch-Specific Output
917 // ---------------------
918 bncSettings settings;
919 if (Qt::CheckState(settings.value("reqcLogSummaryOnly").toInt()) == Qt::Checked) {
920 return;
921 }
922 *_log << endl;
923 for (int iEpo = 0; iEpo < _qcFile._qcEpo.size(); iEpo++) {
924 const t_qcEpo& qcEpo = _qcFile._qcEpo[iEpo];
925
926 unsigned year, month, day, hour, min;
927 double sec;
928 qcEpo._epoTime.civil_date(year, month, day);
929 qcEpo._epoTime.civil_time(hour, min, sec);
930
931 QString dateStr;
932 QTextStream(&dateStr) << QString("> %1 %2 %3 %4 %5%6")
933 .arg(year, 4)
934 .arg(month, 2, 10, QChar('0'))
935 .arg(day, 2, 10, QChar('0'))
936 .arg(hour, 2, 10, QChar('0'))
937 .arg(min, 2, 10, QChar('0'))
938 .arg(sec, 11, 'f', 7);
939
940 *_log << dateStr << QString(" %1").arg(qcEpo._qcSat.size(), 2)
941 << QString(" %1").arg(qcEpo._PDOP, 4, 'f', 1)
942 << endl;
943
944 QMapIterator<t_prn, t_qcSat> itSat(qcEpo._qcSat);
945 while (itSat.hasNext()) {
946 itSat.next();
947 const t_prn& prn = itSat.key();
948 const t_qcSat& qcSat = itSat.value();
949
950 *_log << prn.toString().c_str()
951 << QString(" %1 %2").arg(qcSat._eleDeg, 6, 'f', 2).arg(qcSat._azDeg, 7, 'f', 2);
952
953 int numObsTypes = 0;
954 for (int iFrq = 0; iFrq < qcSat._qcFrq.size(); iFrq++) {
955 const t_qcFrq& qcFrq = qcSat._qcFrq[iFrq];
956 if (qcFrq._phaseValid) {
957 numObsTypes += 1;
958 }
959 if (qcFrq._codeValid) {
960 numObsTypes += 1;
961 }
962 }
963 *_log << QString(" %1").arg(numObsTypes, 2);
964
965 for (int iFrq = 0; iFrq < qcSat._qcFrq.size(); iFrq++) {
966 const t_qcFrq& qcFrq = qcSat._qcFrq[iFrq];
967 if (qcFrq._phaseValid) {
968 *_log << " L" << qcFrq._rnxType2ch << ' ';
969 if (qcFrq._slip) {
970 *_log << 's';
971 }
972 else {
973 *_log << '.';
974 }
975 if (qcFrq._gap) {
976 *_log << 'g';
977 }
978 else {
979 *_log << '.';
980 }
981 *_log << QString(" %1").arg(qcFrq._SNR, 4, 'f', 1);
982 }
983 if (qcFrq._codeValid) {
984 *_log << " C" << qcFrq._rnxType2ch << ' ';
985 if (qcFrq._gap) {
986 *_log << " g";
987 }
988 else {
989 *_log << " .";
990 }
991 *_log << QString(" %1").arg(qcFrq._stdMP, 3, 'f', 2);
992 }
993 }
994 *_log << endl;
995 }
996 }
997 _log->flush();
998}
999
1000//
1001////////////////////////////////////////////////////////////////////////////
1002void t_reqcAnalyze::checkEphemerides() {
1003
1004 QString navFileName;
1005 QStringListIterator namIt(_navFileNames);
1006 bool firstName = true;
1007 while (namIt.hasNext()) {
1008 QFileInfo navFi(namIt.next());
1009 if (firstName) {
1010 firstName = false;
1011 navFileName += navFi.fileName();
1012 }
1013 else {
1014 navFileName += ", " + navFi.fileName();
1015 }
1016 }
1017 if (_log) {
1018 *_log << "Navigation File(s) : " << navFileName << endl;
1019 }
1020 QStringListIterator it(_navFileNames);
1021 while (it.hasNext()) {
1022 const QString& fileName = it.next();
1023 unsigned numOK = 0;
1024 unsigned numBad = 0;
1025 bncEphUser ephUser(false);
1026 t_rnxNavFile rnxNavFile(fileName, t_rnxNavFile::input);
1027 for (unsigned ii = 0; ii < rnxNavFile.ephs().size(); ii++) {
1028 t_eph* eph = rnxNavFile.ephs()[ii];
1029 ephUser.putNewEph(eph, true);
1030 if (eph->checkState() == t_eph::bad) {
1031 ++numBad;
1032 }
1033 else {
1034 ++numOK;
1035 }
1036 }
1037 if (_log) {
1038 *_log << "Ephemeris : " << numOK << " OK " << numBad << " BAD" << endl;
1039 }
1040 if (numBad > 0) {
1041 for (unsigned ii = 0; ii < rnxNavFile.ephs().size(); ii++) {
1042 t_eph* eph = rnxNavFile.ephs()[ii];
1043 if (eph->checkState() == t_eph::bad) {
1044 QFileInfo navFi(fileName);
1045 if (_log) {
1046 *_log << " Bad Ephemeris : " << navFi.fileName() << ' '
1047 << eph->toString(3.0).left(24) << endl;
1048 }
1049 }
1050 }
1051 }
1052 }
1053 if (_log) {
1054 *_log << endl;
1055 }
1056}
Note: See TracBrowser for help on using the repository browser.