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

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