1 //============================================================================== 2 // 3 // This file is part of GPSTk, the GPS Toolkit. 4 // 5 // The GPSTk is free software; you can redistribute it and/or modify 6 // it under the terms of the GNU Lesser General Public License as published 7 // by the Free Software Foundation; either version 3.0 of the License, or 8 // any later version. 9 // 10 // The GPSTk is distributed in the hope that it will be useful, 11 // but WITHOUT ANY WARRANTY; without even the implied warranty of 12 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 13 // GNU Lesser General Public License for more details. 14 // 15 // You should have received a copy of the GNU Lesser General Public 16 // License along with GPSTk; if not, write to the Free Software Foundation, 17 // Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110, USA 18 // 19 // This software was developed by Applied Research Laboratories at the 20 // University of Texas at Austin. 21 // Copyright 2004-2020, The Board of Regents of The University of Texas System 22 // 23 //============================================================================== 24 25 //============================================================================== 26 // 27 // This software was developed by Applied Research Laboratories at the 28 // University of Texas at Austin, under contract to an agency or agencies 29 // within the U.S. Department of Defense. The U.S. Government retains all 30 // rights to use, duplicate, distribute, disclose, or release this software. 31 // 32 // Pursuant to DoD Directive 523024 33 // 34 // DISTRIBUTION STATEMENT A: This software has been approved for public 35 // release, distribution is unlimited. 36 // 37 //============================================================================== 38 39 /** 40 * @file Rinex3ObsHeader.cpp 41 * Encapsulate header of Rinex observation file, including I/O 42 */ 43 44 #include <sstream> 45 #include <algorithm> 46 #include <set> 47 #include <string.h> 48 49 #include "StringUtils.hpp" 50 #include "SystemTime.hpp" 51 #include "TimeString.hpp" 52 #include "Rinex3ObsStream.hpp" 53 #include "Rinex3ObsHeader.hpp" 54 55 using namespace std; 56 using namespace gpstk::StringUtils; 57 58 namespace gpstk 59 { 60 const string Rinex3ObsHeader::hsVersion = "RINEX VERSION / TYPE"; 61 const string Rinex3ObsHeader::hsRunBy = "PGM / RUN BY / DATE"; 62 const string Rinex3ObsHeader::hsComment = "COMMENT"; 63 const string Rinex3ObsHeader::hsMarkerName = "MARKER NAME"; 64 const string Rinex3ObsHeader::hsMarkerNumber = "MARKER NUMBER"; 65 const string Rinex3ObsHeader::hsMarkerType = "MARKER TYPE"; 66 const string Rinex3ObsHeader::hsObserver = "OBSERVER / AGENCY"; 67 const string Rinex3ObsHeader::hsReceiver = "REC # / TYPE / VERS"; 68 const string Rinex3ObsHeader::hsAntennaType = "ANT # / TYPE"; 69 const string Rinex3ObsHeader::hsAntennaPosition = "APPROX POSITION XYZ"; 70 const string Rinex3ObsHeader::hsAntennaDeltaHEN = "ANTENNA: DELTA H/E/N"; 71 const string Rinex3ObsHeader::hsAntennaDeltaXYZ = "ANTENNA: DELTA X/Y/Z"; 72 const string Rinex3ObsHeader::hsAntennaPhaseCtr = "ANTENNA: PHASECENTER"; 73 const string Rinex3ObsHeader::hsAntennaBsightXYZ = "ANTENNA: B.SIGHT XYZ"; 74 const string Rinex3ObsHeader::hsAntennaZeroDirAzi = "ANTENNA: ZERODIR AZI"; 75 const string Rinex3ObsHeader::hsAntennaZeroDirXYZ = "ANTENNA: ZERODIR XYZ"; 76 const string Rinex3ObsHeader::hsCenterOfMass = "CENTER OF MASS: XYZ"; 77 const string Rinex3ObsHeader::hsNumObs = "# / TYPES OF OBSERV"; 78 const string Rinex3ObsHeader::hsSystemNumObs = "SYS / # / OBS TYPES"; 79 const string Rinex3ObsHeader::hsWaveFact = "WAVELENGTH FACT L1/2"; 80 const string Rinex3ObsHeader::hsSigStrengthUnit = "SIGNAL STRENGTH UNIT"; 81 const string Rinex3ObsHeader::hsInterval = "INTERVAL"; 82 const string Rinex3ObsHeader::hsFirstTime = "TIME OF FIRST OBS"; 83 const string Rinex3ObsHeader::hsLastTime = "TIME OF LAST OBS"; 84 const string Rinex3ObsHeader::hsReceiverOffset = "RCV CLOCK OFFS APPL"; 85 const string Rinex3ObsHeader::hsSystemDCBSapplied = "SYS / DCBS APPLIED"; 86 const string Rinex3ObsHeader::hsSystemPCVSapplied = "SYS / PCVS APPLIED"; 87 const string Rinex3ObsHeader::hsSystemScaleFac = "SYS / SCALE FACTOR"; 88 const string Rinex3ObsHeader::hsSystemPhaseShift = "SYS / PHASE SHIFT"; 89 const string Rinex3ObsHeader::hsGlonassSlotFreqNo = "GLONASS SLOT / FRQ #"; 90 const string Rinex3ObsHeader::hsGlonassCodPhsBias = "GLONASS COD/PHS/BIS"; 91 const string Rinex3ObsHeader::hsLeapSeconds = "LEAP SECONDS"; 92 const string Rinex3ObsHeader::hsNumSats = "# OF SATELLITES"; 93 const string Rinex3ObsHeader::hsPrnObs = "PRN / # OF OBS"; 94 const string Rinex3ObsHeader::hsEoH = "END OF HEADER"; 95 96 int Rinex3ObsHeader::debug = 0; 97 98 const Rinex3ObsHeader::Fields Rinex3ObsHeader::allValid2({ 99 Rinex3ObsHeader::validVersion, 100 Rinex3ObsHeader::validRunBy, 101 Rinex3ObsHeader::validMarkerName, 102 Rinex3ObsHeader::validObserver, 103 Rinex3ObsHeader::validReceiver, 104 Rinex3ObsHeader::validAntennaType, 105 Rinex3ObsHeader::validAntennaPosition, 106 Rinex3ObsHeader::validAntennaDeltaHEN, 107 Rinex3ObsHeader::validNumObs, 108 Rinex3ObsHeader::validFirstTime 109 }); 110 const Rinex3ObsHeader::Fields Rinex3ObsHeader::allValid30({ 111 Rinex3ObsHeader::validVersion, 112 Rinex3ObsHeader::validRunBy, 113 Rinex3ObsHeader::validMarkerName, 114 /** @note MARKER TYPE is actually required but in the 115 * interest of supporting the invalid IGS data files, 116 * we make it optional */ 117 //Rinex3ObsHeader::validMarkerType, 118 Rinex3ObsHeader::validObserver, 119 Rinex3ObsHeader::validReceiver, 120 Rinex3ObsHeader::validAntennaType, 121 Rinex3ObsHeader::validAntennaPosition, 122 Rinex3ObsHeader::validAntennaDeltaHEN, 123 Rinex3ObsHeader::validSystemNumObs, 124 Rinex3ObsHeader::validFirstTime 125 }); 126 const Rinex3ObsHeader::Fields Rinex3ObsHeader::allValid301({ 127 Rinex3ObsHeader::validVersion, 128 Rinex3ObsHeader::validRunBy, 129 Rinex3ObsHeader::validMarkerName, 130 //Rinex3ObsHeader::validMarkerType, 131 Rinex3ObsHeader::validObserver, 132 Rinex3ObsHeader::validReceiver, 133 Rinex3ObsHeader::validAntennaType, 134 Rinex3ObsHeader::validAntennaPosition, 135 Rinex3ObsHeader::validAntennaDeltaHEN, 136 Rinex3ObsHeader::validSystemNumObs, 137 Rinex3ObsHeader::validFirstTime, 138 Rinex3ObsHeader::validSystemPhaseShift 139 }); 140 const Rinex3ObsHeader::Fields Rinex3ObsHeader::allValid302({ 141 Rinex3ObsHeader::validVersion, 142 Rinex3ObsHeader::validRunBy, 143 Rinex3ObsHeader::validMarkerName, 144 //Rinex3ObsHeader::validMarkerType, 145 Rinex3ObsHeader::validObserver, 146 Rinex3ObsHeader::validReceiver, 147 Rinex3ObsHeader::validAntennaType, 148 Rinex3ObsHeader::validAntennaPosition, 149 Rinex3ObsHeader::validAntennaDeltaHEN, 150 Rinex3ObsHeader::validSystemNumObs, 151 Rinex3ObsHeader::validFirstTime, 152 Rinex3ObsHeader::validSystemPhaseShift 153 }); 154 const Rinex3ObsHeader::Fields Rinex3ObsHeader::allValid303({ 155 Rinex3ObsHeader::validVersion, 156 Rinex3ObsHeader::validRunBy, 157 Rinex3ObsHeader::validMarkerName, 158 //Rinex3ObsHeader::validMarkerType, 159 Rinex3ObsHeader::validObserver, 160 Rinex3ObsHeader::validReceiver, 161 Rinex3ObsHeader::validAntennaType, 162 Rinex3ObsHeader::validAntennaPosition, 163 Rinex3ObsHeader::validAntennaDeltaHEN, 164 Rinex3ObsHeader::validSystemNumObs, 165 Rinex3ObsHeader::validFirstTime, 166 Rinex3ObsHeader::validSystemPhaseShift 167 }); 168 Rinex3ObsHeader()169 Rinex3ObsHeader::Rinex3ObsHeader() 170 : PisY(false) 171 { 172 clear(); 173 } 174 175 clear()176 void Rinex3ObsHeader::clear() 177 { 178 R2ObsTypes.clear(); 179 mapSysR2toR3ObsID.clear(); 180 version = currentVersion; 181 fileType = "O"; // observation data 182 fileSys = "G"; // GPS only by default 183 preserveVerType = false; // let the write methods chose the above 184 fileSysSat = SatID(-1,SatelliteSystem::GPS); 185 fileProgram.clear(); 186 fileAgency.clear(); 187 date.clear(); 188 preserveDate = false; 189 commentList.clear(); 190 markerName.clear(); 191 markerNumber.clear(); 192 markerType.clear(); 193 observer.clear(); 194 agency.clear(); 195 recNo.clear(); 196 recType.clear(); 197 recVers.clear(); 198 antNo.clear(); 199 antType.clear(); 200 antennaPosition = Triple(); 201 antennaDeltaHEN = Triple(); 202 antennaDeltaXYZ = Triple(); 203 antennaSatSys.clear(); 204 antennaObsCode.clear(); 205 antennaPhaseCtr = Triple(); 206 antennaBsightXYZ = Triple(); 207 antennaZeroDirAzi = 0.; 208 antennaZeroDirXYZ = Triple(); 209 centerOfMass = Triple(); 210 mapObsTypes.clear(); 211 wavelengthFactor[0] = wavelengthFactor[1] = 1; 212 extraWaveFactList.clear(); 213 sigStrengthUnit.clear(); 214 interval = 0.; 215 firstObs = CivilTime(); 216 lastObs = CivilTime(); 217 receiverOffset = 0; 218 infoDCBS.clear(); 219 infoPCVS.clear(); 220 sysSfacMap.clear(); 221 sysPhaseShift.clear(); 222 glonassFreqNo.clear(); 223 glonassCodPhsBias.clear(); 224 leapSeconds = 0; 225 numSVs = 0; 226 numObsForSat.clear(); 227 obsTypeList.clear(); 228 valid.clear(); 229 validEoH = false; 230 // Only do this in the constructor so the desired handling of 231 // "P" code in RINEX 2 stays the same. 232 //PisY = false; 233 sysPhaseShiftObsID = RinexObsID(); 234 satSysTemp.clear(); 235 satSysPrev.clear(); 236 numObs = 0; 237 numObsPrev = 0; 238 lastPRN = RinexSatID(); 239 factor = 0; 240 factorPrev = 0; 241 } 242 243 reallyPutRecord(FFStream & ffs) const244 void Rinex3ObsHeader::reallyPutRecord(FFStream& ffs) const 245 { 246 using gpstk::StringUtils::asString; 247 Rinex3ObsStream& strm = dynamic_cast<Rinex3ObsStream&>(ffs); 248 249 strm.header = *this; 250 251 Fields allValid = Fields::getRequired(version); 252 if (allValid.empty()) 253 { 254 FFStreamError err("Unknown RINEX version: " + 255 gpstk::StringUtils::asString(version,2)); 256 err.addText("Make sure to set the version correctly."); 257 GPSTK_THROW(err); 258 } 259 260 if((valid & allValid) != allValid) 261 { 262 FFStreamError err("Incomplete or invalid header."); 263 err.addText("Make sure you set all header valid bits for all of the" 264 " available data."); 265 allValid.describeMissing(valid, err); 266 GPSTK_THROW(err); 267 } 268 269 try 270 { 271 writeHeaderRecords(strm); 272 } 273 catch(FFStreamError& e) 274 { 275 GPSTK_RETHROW(e); 276 } 277 catch(StringException& e) 278 { 279 GPSTK_RETHROW(e); 280 } 281 282 } // end reallyPutRecord 283 284 285 // This function computes the number of valid header records 286 // which writeHeaderRecords will write. 287 // NB not used in Rinex3Obs.... numberHeaderRecordsToBeWritten(void) const288 int Rinex3ObsHeader::numberHeaderRecordsToBeWritten(void) const throw() 289 { 290 int n = 0; 291 292 if(valid & validVersion ) n++; 293 if(valid & validRunBy ) n++; 294 if(valid & validComment ) n += commentList.size(); 295 if(valid & validMarkerName ) n++; 296 if(valid & validMarkerNumber ) n++; 297 if(version >= 3 && (valid & validMarkerType)) n++; 298 if(valid & validObserver ) n++; 299 if(valid & validReceiver ) n++; 300 if(valid & validAntennaType ) n++; 301 if(valid & validAntennaPosition ) n++; 302 if(valid & validAntennaDeltaHEN ) n++; 303 if(version >= 3 && (valid & validAntennaDeltaXYZ)) n++; 304 if(version >= 3 && (valid & validAntennaPhaseCtr)) n++; 305 if(version >= 3 && (valid & validAntennaBsightXYZ)) n++; 306 if(version >= 3 && (valid & validAntennaZeroDirAzi)) n++; 307 if(version >= 3 && (valid & validAntennaZeroDirXYZ)) n++; 308 if(version >= 3 && (valid & validCenterOfMass)) n++; 309 if(version < 3 && (valid & validNumObs) && R2ObsTypes.size() != 0) 310 n += 1 + (R2ObsTypes.size()-1)/9; 311 if(version >= 3 && (valid & validSystemNumObs) && numObs != 0) 312 n += 1 + (numObs-1)/9; 313 if(version < 3 && (valid & validWaveFact)) 314 { 315 n++; 316 if(extraWaveFactList.size() != 0) n += extraWaveFactList.size(); 317 } 318 if(version >= 3 && (valid & validSigStrengthUnit)) n++; 319 if(valid & validInterval ) n++; 320 if(valid & validFirstTime ) n++; 321 if(valid & validLastTime ) n++; 322 if(valid & validReceiverOffset ) n++; 323 if(version >= 3 && (valid & validSystemDCBSapplied)) n++; 324 if(version >= 3 && (valid & validSystemPCVSapplied)) n++; 325 if(version >= 3 && (valid & validSystemScaleFac)) n++; 326 if(version >= 3.01 && (valid & validSystemPhaseShift)) n++; // one per system at least 327 if(version >= 3.01 && (valid & validGlonassSlotFreqNo)) n++; // TODO: continuation lines... 328 if(version >= 3.02 && (valid & validGlonassCodPhsBias)) n++; 329 if(valid & validLeapSeconds ) n++; 330 if(valid & validNumSats ) n++; 331 if(valid & validPrnObs ) 332 n += numObsForSat.size() * (1+numObsForSat.begin()->second.size()/9); 333 if(validEoH ) n++; 334 335 return n; 336 } // end numberHeaderRecordsToBeWritten 337 338 339 // This function writes all valid header records. writeHeaderRecords(FFStream & ffs) const340 void Rinex3ObsHeader::writeHeaderRecords(FFStream& ffs) const 341 { 342 using gpstk::StringUtils::asString; 343 Rinex3ObsStream& strm = dynamic_cast<Rinex3ObsStream&>(ffs); 344 string line; 345 346 if(valid & validVersion) 347 { 348 line = rightJustify(asString(version,2), 9); 349 line += string(11, ' '); 350 351 if((fileType[0] != 'O') && (fileType[0] != 'o')) 352 { 353 FFStreamError err("File type is not Observation: " + fileType); 354 GPSTK_THROW(err); 355 } 356 357 if (preserveVerType) 358 { 359 line += leftJustify(fileType, 20); 360 line += leftJustify(fileSys, 20); 361 } 362 else 363 { 364 if(fileSysSat.system == SatelliteSystem::Unknown) 365 { 366 FFStreamError err("Invalid satellite system"); 367 GPSTK_THROW(err); 368 } 369 370 line += leftJustify(string("OBSERVATION DATA"), 20); 371 string str; 372 if(fileSysSat.system == SatelliteSystem::Mixed) 373 str = "MIXED"; 374 else 375 { 376 RinexSatID sat(fileSysSat); 377 str = sat.systemChar(); 378 str = str + " (" + sat.systemString() + ")"; 379 } 380 line += leftJustify(str, 20); 381 } 382 line += hsVersion; 383 strm << line << endl; 384 strm.lineNumber++; 385 } 386 if(valid & validRunBy) 387 { 388 line = leftJustify(fileProgram, 20); 389 line += leftJustify(fileAgency , 20); 390 if (preserveDate) 391 { 392 line += leftJustify(date, 20); 393 } 394 else 395 { 396 SystemTime sysTime; 397 string curDate; 398 curDate = printTime(sysTime,"%04Y%02m%02d %02H%02M%02S %P"); 399 line += leftJustify(curDate, 20); 400 } 401 line += hsRunBy; 402 strm << line << endl; 403 strm.lineNumber++; 404 } 405 if(valid & validComment) 406 { 407 vector<string>::const_iterator itr = commentList.begin(); 408 while (itr != commentList.end()) 409 { 410 line = leftJustify((*itr), 60); 411 line += hsComment; 412 strm << line << endl; 413 strm.lineNumber++; 414 itr++; 415 } 416 } 417 if(!R2DisambiguityMap.empty()) 418 { 419 DisAmbMap::const_iterator itr = R2DisambiguityMap.begin(); 420 while (itr != R2DisambiguityMap.end()) 421 { 422 line = leftJustify(itr->first.substr(0,1) + " Obs Type " + itr->first.substr(1) + " originated from " + itr->second, 60); 423 line += hsComment; 424 strm << line << endl; 425 strm.lineNumber++; 426 itr++; 427 } 428 } 429 if(valid & validMarkerName) 430 { 431 line = leftJustify(markerName, 60); 432 line += hsMarkerName; 433 strm << line << endl; 434 strm.lineNumber++; 435 } 436 if(valid & validMarkerNumber) 437 { 438 line = leftJustify(markerNumber, 20); 439 line += string(40, ' '); 440 line += hsMarkerNumber; 441 strm << line << endl; 442 strm.lineNumber++; 443 } 444 if(version >= 3 && (valid & validMarkerType)) 445 { 446 line = leftJustify(markerType, 20); 447 line += string(40, ' '); 448 line += hsMarkerType; 449 strm << line << endl; 450 strm.lineNumber++; 451 } 452 if(valid & validObserver) 453 { 454 line = leftJustify(observer, 20); 455 line += leftJustify(agency , 40); 456 line += hsObserver; 457 strm << line << endl; 458 strm.lineNumber++; 459 } 460 if(valid & validReceiver) 461 { 462 line = leftJustify(recNo , 20); 463 line += leftJustify(recType, 20); 464 line += leftJustify(recVers, 20); 465 line += hsReceiver; 466 strm << line << endl; 467 strm.lineNumber++; 468 } 469 if(valid & validAntennaType) 470 { 471 line = leftJustify(antNo , 20); 472 line += leftJustify(antType, 20); 473 line += string(20, ' '); 474 line += hsAntennaType; 475 strm << line << endl; 476 strm.lineNumber++; 477 } 478 if(valid & validAntennaPosition) 479 { 480 line = rightJustify(asString(antennaPosition[0], 4), 14); 481 line += rightJustify(asString(antennaPosition[1], 4), 14); 482 line += rightJustify(asString(antennaPosition[2], 4), 14); 483 line += string(18, ' '); 484 line += hsAntennaPosition; 485 strm << line << endl; 486 strm.lineNumber++; 487 } 488 if(valid & validAntennaDeltaHEN) 489 { 490 line = rightJustify(asString(antennaDeltaHEN[0], 4), 14); 491 line += rightJustify(asString(antennaDeltaHEN[1], 4), 14); 492 line += rightJustify(asString(antennaDeltaHEN[2], 4), 14); 493 line += string(18, ' '); 494 line += hsAntennaDeltaHEN; 495 strm << line << endl; 496 strm.lineNumber++; 497 } 498 if(version >= 3 && (valid & validAntennaDeltaXYZ)) 499 { 500 line = rightJustify(asString(antennaDeltaXYZ[0], 4), 14); 501 line += rightJustify(asString(antennaDeltaXYZ[1], 4), 14); 502 line += rightJustify(asString(antennaDeltaXYZ[2], 4), 14); 503 line += string(18, ' '); 504 line += hsAntennaDeltaXYZ; 505 strm << line << endl; 506 strm.lineNumber++; 507 } 508 if(version >= 3 && (valid & validAntennaPhaseCtr)) 509 { 510 line = leftJustify(antennaSatSys , 1); 511 line += string(1, ' '); 512 line += rightJustify(antennaObsCode, 3); 513 line += rightJustify(asString(antennaPhaseCtr[0], 4), 9); 514 line += rightJustify(asString(antennaPhaseCtr[1], 4), 14); 515 line += rightJustify(asString(antennaPhaseCtr[2], 4), 14); 516 line += string(18, ' '); 517 line += hsAntennaPhaseCtr; 518 strm << line << endl; 519 strm.lineNumber++; 520 } 521 if(version >= 3 && (valid & validAntennaBsightXYZ)) 522 { 523 line = rightJustify(asString(antennaBsightXYZ[0], 4), 14); 524 line += rightJustify(asString(antennaBsightXYZ[1], 4), 14); 525 line += rightJustify(asString(antennaBsightXYZ[2], 4), 14); 526 line += string(18, ' '); 527 line += hsAntennaBsightXYZ; 528 strm << line << endl; 529 strm.lineNumber++; 530 } 531 if(version >= 3 && (valid & validAntennaZeroDirAzi)) 532 { 533 line = rightJustify(asString(antennaZeroDirAzi, 4), 14); 534 line += string(46, ' '); 535 line += hsAntennaZeroDirAzi; 536 strm << line << endl; 537 strm.lineNumber++; 538 } 539 if(version >= 3 && (valid & validAntennaZeroDirXYZ)) 540 { 541 line = rightJustify(asString(antennaZeroDirXYZ[0], 4), 14); 542 line += rightJustify(asString(antennaZeroDirXYZ[1], 4), 14); 543 line += rightJustify(asString(antennaZeroDirXYZ[2], 4), 14); 544 line += string(18, ' '); 545 line += hsAntennaZeroDirXYZ; 546 strm << line << endl; 547 strm.lineNumber++; 548 } 549 if(version >= 3 && (valid & validCenterOfMass)) 550 { 551 line = rightJustify(asString(centerOfMass[0], 4), 14); 552 line += rightJustify(asString(centerOfMass[1], 4), 14); 553 line += rightJustify(asString(centerOfMass[2], 4), 14); 554 line += string(18, ' '); 555 line += hsCenterOfMass; 556 strm << line << endl; 557 strm.lineNumber++; 558 } 559 if(version < 3 && (valid & validNumObs)) // R2 only 560 { 561 if(R2ObsTypes.empty()) 562 { 563 InvalidRequest er("Header contains no R2ObsTypes. " 564 "You must run prepareVer2Write before outputting an R2 file"); 565 GPSTK_THROW(er); 566 } 567 // write out RinexObsTypes 568 const int maxObsPerLine = 9; 569 int obsWritten = 0; 570 line = ""; // make sure the line contents are reset. 571 572 for(size_t i=0; i<R2ObsTypes.size(); i++) 573 { 574 string val; 575 // the first line needs to have the # of obs 576 if(obsWritten == 0) 577 line = rightJustify(asString(R2ObsTypes.size()), 6); 578 // if you hit 9, write out the line and start a new one 579 else if((obsWritten % maxObsPerLine) == 0) 580 { 581 line += hsNumObs; 582 strm << line << endl; 583 strm.lineNumber++; 584 line = string(6, ' '); 585 } 586 val = R2ObsTypes[i]; 587 line += rightJustify(val, 6); 588 obsWritten++; 589 } 590 591 line += string(60 - line.size(), ' '); 592 line += hsNumObs; 593 strm << line << endl; 594 strm.lineNumber++; 595 } 596 if(version >= 3 && (valid & validSystemNumObs)) 597 { 598 static const int maxObsPerLine = 13; 599 600 map<string,vector<RinexObsID> >::const_iterator mapIter; 601 602 // because of the way pseudo-observables are mapped from 603 // ObsID/RinexObsID to the RINEX file, we have to manually 604 // count the observables first. 605 map<string,unsigned> obsCount; 606 RinexObsMap dummy; 607 remapObsTypes(dummy, obsCount); 608 // now do the actual writing 609 for(mapIter = mapObsTypes.begin(); mapIter != mapObsTypes.end(); 610 mapIter++) 611 { 612 bool addedChannel = false; 613 std::set<CarrierBand> addedIono; 614 int obsWritten = 0; 615 line = ""; // make sure the line contents are reset 616 617 vector<RinexObsID> ObsTypeList = mapIter->second; 618 619 for(size_t i = 0; i < ObsTypeList.size(); i++) 620 { 621 if (ObsTypeList[i].type == ObservationType::Iono) 622 { 623 if (addedIono.count(ObsTypeList[i].band) > 0) 624 continue; // only write this pseudo-obs once 625 addedIono.insert(ObsTypeList[i].band); 626 } 627 else if (ObsTypeList[i].type == ObservationType::Channel) 628 { 629 if (addedChannel) 630 continue; // only write this pseudo-obs once 631 addedChannel = true; 632 } 633 // the first line needs to have the GNSS type and # of obs 634 if(obsWritten == 0) 635 { 636 line = leftJustify(mapIter->first, 1); 637 line += string(2, ' '); 638 line += rightJustify(asString(obsCount[mapIter->first]), 3); 639 } 640 // if you hit 13, write out the line and start a new one 641 else if((obsWritten % maxObsPerLine) == 0) 642 { 643 line += string(2, ' '); 644 line += hsSystemNumObs; 645 strm << line << endl; 646 strm.lineNumber++; 647 line = string(6, ' '); 648 } 649 line += string(1, ' '); 650 line += rightJustify(ObsTypeList[i].asString(version), 3); 651 obsWritten++; 652 } 653 line += string(60 - line.size(), ' '); 654 line += hsSystemNumObs; 655 strm << line << endl; 656 strm.lineNumber++; 657 } 658 659 } 660 if(version < 3 && (valid & validWaveFact)) 661 { 662 line = rightJustify(asString<short>(wavelengthFactor[0]),6); 663 line += rightJustify(asString<short>(wavelengthFactor[1]),6); 664 line += string(48, ' '); 665 line += hsWaveFact; 666 strm << line << endl; 667 strm.lineNumber++; 668 669 // handle continuation lines 670 if(!extraWaveFactList.empty()) 671 { 672 vector<ExtraWaveFact>::const_iterator itr = extraWaveFactList.begin(); 673 674 while (itr != extraWaveFactList.end()) 675 { 676 const int maxSatsPerLine = 7; 677 short satsWritten = 0, satsLeft = (*itr).satList.size(), satsThisLine; 678 vector<SatID>::const_iterator vecItr = (*itr).satList.begin(); 679 680 while ((vecItr != (*itr).satList.end())) 681 { 682 if(satsWritten == 0) 683 { 684 line = rightJustify(asString<short>((*itr).wavelengthFactor[0]),6); 685 line += rightJustify(asString<short>((*itr).wavelengthFactor[1]),6); 686 satsThisLine = (satsLeft > maxSatsPerLine ? maxSatsPerLine : satsLeft); 687 line += rightJustify(asString<short>(satsThisLine),6); 688 } 689 try 690 { 691 line += string(3, ' ') + RinexSatID(*vecItr).toString(); 692 } 693 catch (Exception& e) 694 { 695 FFStreamError ffse(e); 696 GPSTK_THROW(ffse); 697 } 698 satsWritten++; 699 satsLeft--; 700 if(satsWritten==maxSatsPerLine || satsLeft==0) 701 { 702 // output a complete line 703 line += string(60 - line.size(), ' '); 704 line += hsWaveFact; 705 strm << line << endl; 706 strm.lineNumber++; 707 satsWritten = 0; 708 } 709 vecItr++; 710 } 711 itr++; 712 } 713 } 714 } 715 if(version >= 3 && valid & validSigStrengthUnit) 716 { 717 line = leftJustify(sigStrengthUnit, 20); 718 line += string(40, ' '); 719 line += hsSigStrengthUnit; 720 strm << line << endl; 721 strm.lineNumber++; 722 } 723 if(valid & validInterval) 724 { 725 line = rightJustify(asString(interval, 3), 10); 726 line += string(50, ' '); 727 line += hsInterval; 728 strm << line << endl; 729 strm.lineNumber++; 730 } 731 if(valid & validFirstTime) 732 { 733 line = writeTime(firstObs); 734 line += string(60 - line.size(), ' '); 735 line += hsFirstTime; 736 strm << line << endl; 737 strm.lineNumber++; 738 } 739 if(valid & validLastTime) 740 { 741 line = writeTime(lastObs); 742 line += string(60 - line.size(), ' '); 743 line += hsLastTime; 744 strm << line << endl; 745 strm.lineNumber++; 746 } 747 if(valid & validReceiverOffset) 748 { 749 line = rightJustify(asString(receiverOffset), 6); 750 line += string(54, ' '); 751 line += hsReceiverOffset; 752 strm << line << endl; 753 strm.lineNumber++; 754 } 755 if(version >= 3 && (valid & validSystemDCBSapplied)) 756 { 757 for(size_t i = 0; i < infoDCBS.size(); i++) 758 { 759 line = leftJustify(infoDCBS[i].satSys, 1); 760 line += string(1, ' '); 761 line += leftJustify(infoDCBS[i].name , 17); 762 line += string(1, ' '); 763 line += leftJustify(infoDCBS[i].source, 40); 764 line += hsSystemDCBSapplied; 765 strm << line << endl; 766 strm.lineNumber++; 767 } 768 } 769 if(version >= 3 && (valid & validSystemPCVSapplied)) 770 { 771 for(size_t i = 0; i < infoPCVS.size(); i++) 772 { 773 line = leftJustify(infoPCVS[i].satSys, 1); 774 line += string(1, ' '); 775 line += leftJustify(infoPCVS[i].name , 17); 776 line += string(1, ' '); 777 line += leftJustify(infoPCVS[i].source, 40); 778 line += hsSystemPCVSapplied; 779 strm << line << endl; 780 strm.lineNumber++; 781 } 782 } 783 if(version >= 3 && (valid & validSystemScaleFac)) 784 { 785 static const int maxObsPerLine = 12; 786 787 static const int size = 4; 788 static const int factors[size] = {1,10,100,1000}; 789 vector<string> obsTypes; 790 791 // loop over GNSSes 792 map<string, ScaleFacMap>::const_iterator mapIter; 793 for(mapIter = sysSfacMap.begin(); mapIter != sysSfacMap.end(); mapIter++) 794 { 795 map<RinexObsID, int>::const_iterator iter; 796 797 for(int i = 0; i < size; i++) // loop over possible factors (above) 798 { 799 int count = 0; 800 obsTypes.clear(); // clear the list of Obs Types we're going to make 801 802 for(iter = mapIter->second.begin(); // loop over scale factor map 803 iter != mapIter->second.end(); iter++) 804 { 805 if(iter->second == factors[i] ) 806 { 807 count++; 808 obsTypes.push_back(iter->first.asString(version)); 809 } 810 } 811 812 if(count == 0 ) continue; 813 814 line = leftJustify(mapIter->first , 1); 815 line += string(1, ' '); 816 line += rightJustify(asString(factors[i]), 4); 817 line += string(2, ' '); 818 line += rightJustify(asString(count ), 2); 819 820 for(int j = 0; j < count; j++) 821 { 822 if(j > maxObsPerLine-1 && (j % maxObsPerLine) == 0 ) 823 { 824 // need continuation; end current line 825 line += string(2, ' '); 826 line += hsSystemScaleFac; 827 strm << line << endl; 828 strm.lineNumber++; 829 line = string(10, ' '); 830 } 831 line += string(1, ' '); 832 line += rightJustify(obsTypes[j], 3); 833 } 834 int space = 60 - 10 - 4*(count % maxObsPerLine); 835 line += string(space, ' '); 836 line += hsSystemScaleFac; 837 strm << line << endl; 838 strm.lineNumber++; 839 } 840 } 841 } 842 if(version >= 3.01 && (valid & validSystemPhaseShift)) 843 { 844 //map<string, map<RinexObsID, map<RinexSatID,double> > > sysPhaseShift; 845 map<string, map<RinexObsID, map<RinexSatID,double> > >::const_iterator it; 846 for(it=sysPhaseShift.begin(); it!=sysPhaseShift.end(); ++it) 847 { 848 string sys(it->first); 849 map<RinexObsID, map<RinexSatID,double> >::const_iterator jt(it->second.begin()); 850 if(jt == it->second.end()) 851 { 852 line = sys; 853 line += string(60-line.length(), ' '); 854 line += hsSystemPhaseShift; 855 strm << line << endl; 856 strm.lineNumber++; 857 } 858 else 859 { 860 for( ; jt!=it->second.end(); ++jt) 861 { 862 RinexSatID sat(jt->second.begin()->first); 863 double corr(jt->second.begin()->second); 864 if (jt->first.type != ObservationType::Phase) 865 { 866 // Phase shift only makes sense for phase measurements. 867 continue; 868 } 869 line = sys + " "; 870 line += leftJustify(jt->first.asString(version),3) + " "; 871 line += rightJustify(asString(corr,5),8); 872 if(sat.id == -1) 873 { 874 line += string(60-line.length(), ' '); 875 line += hsSystemPhaseShift; 876 strm << line << endl; 877 strm.lineNumber++; 878 } 879 else 880 { 881 // list of sats 882 setfill('0'); 883 line += string(" ") + rightJustify(asString(jt->second.size()),2); 884 setfill(' '); 885 886 int n(0); 887 map<RinexSatID,double>::const_iterator kt,lt; 888 for(kt=jt->second.begin(); kt!=jt->second.end(); ++kt) 889 { 890 line += string(" ") + kt->first.toString(); 891 if(++n == 10 || ++(lt=kt) == jt->second.end()) 892 { 893 // end this line 894 line += string(60-line.length(), ' '); 895 line += hsSystemPhaseShift; 896 strm << line << endl; 897 strm.lineNumber++; 898 n = 0; 899 // are there more for a continuation line? 900 if(lt != jt->second.end()) 901 line = string(18,' '); 902 } // if(++n == 10 || ++(lt=kt) == jt->second.end()) 903 } // for(kt=jt->second.begin(); kt!=jt->second.end(); ++kt) 904 } // else 905 } // for( ; jt!=it->second.end(); ++jt) 906 } // else 907 } // for(it=sysPhaseShift.begin(); it!=sysPhaseShift.end(); ++it) 908 } // if(version >= 3.01 && (valid & validSystemPhaseShift)) 909 if(version >= 3.01) 910 { 911 if ((valid & validGlonassSlotFreqNo)) 912 { 913 size_t n(0), nsat(glonassFreqNo.size()); 914 line = rightJustify(asString(nsat), 3) + " "; 915 GLOFreqNumMap::const_iterator it, kt; 916 for (it = glonassFreqNo.begin(); it != glonassFreqNo.end(); ++it) { 917 line += it->first.toString(); 918 line += rightJustify(asString(it->second), 3) + " "; 919 if (++n == 8 || ++(kt = it) == glonassFreqNo.end()) { 920 // write it 921 line += string(60 - line.length(), ' '); 922 line += hsGlonassSlotFreqNo; 923 strm << line << endl; 924 strm.lineNumber++; 925 n = 0; 926 // are there more for a continuation line? 927 if (kt != glonassFreqNo.end()) 928 line = string(4, ' '); 929 } 930 } 931 } 932 else if(mapObsTypes.find("R") != mapObsTypes.end()) 933 { 934 FFStreamError err("Glonass Slot Freq No required for files containing Glonass Observations "); 935 GPSTK_THROW(err); 936 } 937 } 938 if(version >= 3.02) 939 { 940 if ((valid & validGlonassCodPhsBias)) 941 { 942 line.clear(); 943 GLOCodPhsBias::const_iterator it; 944 const string labs[4] = {"C1C", "C1P", "C2C", "C2P"}; 945 for (int i = 0; i < 4; i++) { 946 RinexObsID obsid(RinexObsID("R" + labs[i], version)); 947 it = glonassCodPhsBias.find(obsid); 948 double bias = 0.0; 949 if (it != glonassCodPhsBias.end()) 950 bias = it->second; 951 line += " " + labs[i] + " " + rightJustify(asString(bias, 3), 8); 952 } 953 line += string(60 - line.length(), ' '); 954 line += hsGlonassCodPhsBias; 955 strm << line << endl; 956 strm.lineNumber++; 957 } 958 else if(mapObsTypes.find("R") != mapObsTypes.end()) 959 { 960 FFStreamError err("Glonass Code Phase Bias required for files containing Glonass Observations "); 961 GPSTK_THROW(err); 962 } 963 } 964 if(valid & validLeapSeconds) 965 { 966 line = rightJustify(asString(leapSeconds), 6); 967 line += string(54, ' '); 968 line += hsLeapSeconds; 969 strm << line << endl; 970 strm.lineNumber++; 971 } 972 if(valid & validNumSats) 973 { 974 line = rightJustify(asString(numSVs), 6); 975 line += string(54, ' '); 976 line += hsNumSats; 977 strm << line << endl; 978 strm.lineNumber++; 979 } 980 if(valid & validPrnObs) 981 { 982 static const int maxObsPerLine = 9; 983 map<RinexSatID, vector<int> >::const_iterator itr(numObsForSat.begin()); 984 // loop over satellites 985 while(itr != numObsForSat.end()) 986 { 987 int numObsWritten = 0; // # of counts written for this sat 988 RinexSatID sat(itr->first); // the sat 989 const vector<int>& numObs(itr->second); // the vector of ints stored 990 vector<int> vec; // the vector of ints to write 991 992 if(version >= 3) 993 vec = numObs; 994 else 995 { 996 // fill in zeros for version 2 997 int j; 998 size_t i; 999 string sys(string(1,sat.systemChar())); 1000 map<string, map<string, RinexObsID> >::const_iterator jt(mapSysR2toR3ObsID.find(sys)); 1001 const map<string, RinexObsID> mapVec(jt->second); 1002 map<string, RinexObsID>::const_iterator kt; 1003 for(i=0,j=0; i<R2ObsTypes.size(); i++) 1004 { 1005 kt = mapVec.find(R2ObsTypes[i]); 1006 string obsid(kt->second.asString(version)); 1007 if(obsid == string(" ")) vec.push_back(0.0); 1008 else vec.push_back(numObs[j++]); 1009 } 1010 } 1011 1012 vector<int>::const_iterator vecItr(vec.begin()); 1013 while (vecItr != vec.end()) 1014 { 1015 if(numObsWritten == 0) 1016 { 1017 // start of line 1018 try 1019 { 1020 line = string(3, ' ') + sat.toString(); // ' G01' 1021 } 1022 catch (Exception& e) 1023 { 1024 FFStreamError ffse(e); 1025 GPSTK_RETHROW(ffse); 1026 } 1027 } 1028 else if((numObsWritten % maxObsPerLine) == 0) 1029 { 1030 // end of line 1031 line += hsPrnObs; 1032 strm << line << endl; 1033 strm.lineNumber++; 1034 line = string(6, ' '); 1035 } 1036 1037 line += rightJustify(asString(*vecItr), 6); // add num obs to line 1038 ++vecItr; 1039 ++numObsWritten; 1040 } 1041 1042 // finish last line 1043 line += string(60 - line.size(), ' '); 1044 line += hsPrnObs; 1045 strm << line << endl; 1046 strm.lineNumber++; 1047 itr++; 1048 } 1049 } 1050 if(validEoH) 1051 { 1052 line = string(60, ' '); 1053 line += hsEoH; 1054 strm << line << endl; 1055 strm.lineNumber++; 1056 } 1057 } // end writeHeaderRecords 1058 1059 1060 // This function parses a single header record. parseHeaderRecord(string & line)1061 void Rinex3ObsHeader::parseHeaderRecord(string& line) 1062 { 1063 int i; 1064 string label(line, 60, 20); 1065 1066 if(label == hsVersion) 1067 { 1068 version = asDouble(line.substr( 0,20)); 1069 fileType = strip( line.substr(20,20)); 1070 fileSys = strip( line.substr(40,20)); 1071 1072 if(fileSys[0] != 'M' && fileSys[0] != 'm') 1073 { 1074 RinexSatID sat; 1075 sat.fromString(fileSys); 1076 fileSysSat = SatID(sat); 1077 } 1078 else 1079 fileSysSat = SatID(-1,SatelliteSystem::Mixed); 1080 1081 if(fileType[0] != 'O' && fileType[0] != 'o') 1082 { 1083 FFStreamError e("This isn't a RINEX 3 Obs file."); 1084 GPSTK_THROW(e); 1085 } 1086 1087 valid |= validVersion; 1088 } 1089 else if(label == hsRunBy) 1090 { 1091 fileProgram = strip(line.substr( 0,20)); 1092 fileAgency = strip(line.substr(20,20)); 1093 date = strip(line.substr(40,20)); 1094 valid |= validRunBy; 1095 } 1096 else if(label == hsComment) 1097 { 1098 commentList.push_back(strip(line.substr(0,60))); 1099 valid |= validComment; 1100 } 1101 else if(label == hsMarkerName) 1102 { 1103 markerName = strip(line.substr(0,60)); 1104 valid |= validMarkerName; 1105 } 1106 else if(label == hsMarkerNumber) 1107 { 1108 markerNumber = strip(line.substr(0,20)); 1109 valid |= validMarkerNumber; 1110 } 1111 else if(label == hsMarkerType) 1112 { 1113 markerType = strip(line.substr(0,20)); 1114 valid |= validMarkerType; 1115 } 1116 else if(label == hsObserver) 1117 { 1118 observer = strip(line.substr( 0,20)); 1119 agency = strip(line.substr(20,40)); 1120 valid |= validObserver; 1121 } 1122 else if(label == hsReceiver) 1123 { 1124 recNo = strip(line.substr( 0,20)); 1125 recType = strip(line.substr(20,20)); 1126 recVers = strip(line.substr(40,20)); 1127 valid |= validReceiver; 1128 } 1129 else if(label ==hsAntennaType) 1130 { 1131 antNo = strip(line.substr( 0,20)); 1132 antType = strip(line.substr(20,20)); 1133 valid |= validAntennaType; 1134 } 1135 else if(label == hsAntennaPosition) 1136 { 1137 antennaPosition[0] = asDouble(line.substr( 0,14)); 1138 antennaPosition[1] = asDouble(line.substr(14,14)); 1139 antennaPosition[2] = asDouble(line.substr(28,14)); 1140 valid |= validAntennaPosition; 1141 } 1142 else if(label == hsAntennaDeltaHEN) 1143 { 1144 antennaDeltaHEN[0] = asDouble(line.substr( 0,14)); 1145 antennaDeltaHEN[1] = asDouble(line.substr(14,14)); 1146 antennaDeltaHEN[2] = asDouble(line.substr(28,14)); 1147 valid |= validAntennaDeltaHEN; 1148 } 1149 else if(label == hsAntennaDeltaXYZ) 1150 { 1151 antennaDeltaXYZ[0] = asDouble(line.substr( 0,14)); 1152 antennaDeltaXYZ[1] = asDouble(line.substr(14,14)); 1153 antennaDeltaXYZ[2] = asDouble(line.substr(28,14)); 1154 valid |= validAntennaDeltaXYZ; 1155 } 1156 else if(label == hsAntennaPhaseCtr) 1157 { 1158 antennaSatSys = strip(line.substr(0,2)); 1159 antennaObsCode = strip(line.substr(2,3)); 1160 antennaPhaseCtr[0] = asDouble(line.substr( 5, 9)); 1161 antennaPhaseCtr[1] = asDouble(line.substr(14,14)); 1162 antennaPhaseCtr[2] = asDouble(line.substr(28,14)); 1163 valid |= validAntennaPhaseCtr; 1164 } 1165 else if(label == hsAntennaBsightXYZ) 1166 { 1167 antennaBsightXYZ[0] = asDouble(line.substr( 0,14)); 1168 antennaBsightXYZ[1] = asDouble(line.substr(14,14)); 1169 antennaBsightXYZ[2] = asDouble(line.substr(28,14)); 1170 valid |= validAntennaBsightXYZ; 1171 } 1172 else if(label == hsAntennaZeroDirAzi) 1173 { 1174 antennaZeroDirAzi = asDouble(line.substr(0,14)); 1175 valid |= validAntennaBsightXYZ; 1176 } 1177 else if(label == hsAntennaZeroDirXYZ) 1178 { 1179 antennaZeroDirXYZ[0] = asDouble(line.substr( 0,14)); 1180 antennaZeroDirXYZ[1] = asDouble(line.substr(14,14)); 1181 antennaZeroDirXYZ[2] = asDouble(line.substr(28,14)); 1182 valid |= validAntennaBsightXYZ; 1183 } 1184 else if(label == hsCenterOfMass) 1185 { 1186 centerOfMass[0] = asDouble(line.substr( 0,14)); 1187 centerOfMass[1] = asDouble(line.substr(14,14)); 1188 centerOfMass[2] = asDouble(line.substr(28,14)); 1189 valid |= validCenterOfMass; 1190 } 1191 else if(label == hsNumObs) // R2 only 1192 { 1193 if(version >= 3) 1194 { 1195 FFStreamError e("RINEX 2 record in RINEX 3 file: " + label); 1196 GPSTK_THROW(e); 1197 } 1198 1199 // process the first line 1200 if(!(valid & validNumObs)) 1201 { 1202 numObs = asInt(line.substr(0,6)); 1203 valid |= validNumObs; 1204 } 1205 1206 const int maxObsPerLine = 9; 1207 for(i = 0; (R2ObsTypes.size() < numObs) && (i < maxObsPerLine); i++) 1208 R2ObsTypes.push_back(line.substr(i * 6 + 10, 2)); 1209 } 1210 else if(label == hsSystemNumObs) 1211 { 1212 if(version < 3) 1213 { 1214 FFStreamError e("RINEX 3 record in RINEX 2 file: " + label); 1215 GPSTK_THROW(e); 1216 } 1217 1218 string satSys = strip(line.substr(0,1)); 1219 if (satSys != "") 1220 { 1221 numObs = asInt(line.substr(3,3)); 1222 valid |= validSystemNumObs; 1223 satSysPrev = satSys; 1224 } 1225 else 1226 satSys = satSysPrev; 1227 1228 try 1229 { 1230 const int maxObsPerLine = 13; 1231 for (int i=0; 1232 i < maxObsPerLine && mapObsTypes[satSys].size() < numObs; i++) 1233 { 1234 string obstype(line.substr(4 * i + 7, 3)); 1235 mapObsTypes[satSys].push_back( 1236 RinexObsID(satSys+obstype,version)); 1237 } 1238 } 1239 catch(InvalidParameter& ip) 1240 { 1241 FFStreamError fse("InvalidParameter: "+ip.what()); 1242 GPSTK_THROW(fse); 1243 } 1244 } 1245 else if(label == hsWaveFact) // R2 only 1246 { 1247 // first time reading this 1248 if(!(valid & validWaveFact)) 1249 { 1250 wavelengthFactor[0] = asInt(line.substr(0,6)); 1251 wavelengthFactor[1] = asInt(line.substr(6,6)); 1252 valid |= validWaveFact; 1253 } 1254 else 1255 { 1256 // additional wave fact lines 1257 const int maxSatsPerLine = 7; 1258 int Nsats; 1259 ExtraWaveFact ewf; 1260 ewf.wavelengthFactor[0] = asInt(line.substr(0,6)); 1261 ewf.wavelengthFactor[1] = asInt(line.substr(6,6)); 1262 Nsats = asInt(line.substr(12,6)); 1263 1264 if(Nsats > maxSatsPerLine) // > not >= 1265 { 1266 FFStreamError e("Invalid number of Sats for " + hsWaveFact); 1267 GPSTK_THROW(e); 1268 } 1269 1270 for(i = 0; i < Nsats; i++) 1271 { 1272 try 1273 { 1274 RinexSatID prn(line.substr(21+i*6,3)); 1275 ewf.satList.push_back(prn); 1276 } 1277 catch (Exception& e) 1278 { 1279 FFStreamError ffse(e); 1280 GPSTK_RETHROW(ffse); 1281 } 1282 } 1283 1284 extraWaveFactList.push_back(ewf); 1285 } 1286 } 1287 else if(label == hsSigStrengthUnit) 1288 { 1289 sigStrengthUnit = strip(line.substr(0,20)); 1290 valid |= validSigStrengthUnit; 1291 } 1292 else if(label == hsInterval) 1293 { 1294 interval = asDouble(line.substr(0,10)); 1295 valid |= validInterval; 1296 } 1297 else if(label == hsFirstTime) 1298 { 1299 firstObs = parseTime(line); 1300 valid |= validFirstTime; 1301 } 1302 else if(label == hsLastTime) 1303 { 1304 lastObs = parseTime(line); 1305 valid |= validLastTime; 1306 } 1307 else if(label == hsReceiverOffset) 1308 { 1309 receiverOffset = asInt(line.substr(0,6)); 1310 valid |= validReceiverOffset; 1311 } 1312 1313 else if(label == hsSystemDCBSapplied) 1314 { 1315 Rinex3CorrInfo tempInfo; 1316 tempInfo.satSys = strip(line.substr( 0, 1)); 1317 tempInfo.name = strip(line.substr( 2,17)); 1318 tempInfo.source = strip(line.substr(20,40)); 1319 infoDCBS.push_back(tempInfo); 1320 valid |= validSystemDCBSapplied; 1321 } 1322 else if(label == hsSystemPCVSapplied) 1323 { 1324 Rinex3CorrInfo tempInfo; 1325 tempInfo.satSys = strip(line.substr( 0, 1)); 1326 tempInfo.name = strip(line.substr( 2,17)); 1327 tempInfo.source = strip(line.substr(20,40)); 1328 infoPCVS.push_back(tempInfo); 1329 valid |= validSystemPCVSapplied; 1330 } 1331 else if(label == hsSystemScaleFac) 1332 { 1333 static const int maxObsPerLine = 12; 1334 1335 satSysTemp = strip(line.substr(0,1)); 1336 factor = asInt(line.substr(2,4)); 1337 numObs = asInt(line.substr(8,2)); 1338 1339 int startPosition = 0; 1340 1341 if(satSysTemp == "" ) 1342 { 1343 // it's a continuation line; use prev. info., end pt. to start 1344 satSysTemp = satSysPrev; 1345 factor = factorPrev; 1346 numObs = numObsPrev; 1347 1348 startPosition = sysSfacMap[satSysTemp].size(); 1349 } 1350 1351 // 0/blank numObs means factor applies to all obs types 1352 // in appropriate obsTypeList 1353 if(numObs == 0) 1354 numObs = mapObsTypes[satSysTemp].size(); 1355 1356 ScaleFacMap tempSfacMap = sysSfacMap[satSysTemp]; 1357 for(i = startPosition; 1358 (i < numObs) && ((i % maxObsPerLine) < maxObsPerLine); i++) 1359 { 1360 int position = 4*(i % maxObsPerLine) + 10 + 1; 1361 RinexObsID tempType(satSysTemp+strip(line.substr(position,3)), 1362 version); 1363 tempSfacMap.insert(make_pair(tempType,factor)); 1364 } 1365 sysSfacMap[satSysTemp] = tempSfacMap; 1366 1367 ScaleFacMap::const_iterator iter; 1368 ScaleFacMap tempmap; 1369 tempmap = sysSfacMap[satSysTemp]; 1370 1371 // save values in case next line is a continuation line 1372 satSysPrev = satSysTemp; 1373 factorPrev = factor; 1374 numObsPrev = numObs; 1375 1376 valid |= validSystemScaleFac; 1377 } 1378 else if(label == hsSystemPhaseShift || label == hsSystemPhaseShift+"S") ///< "SYS / PHASE SHIFT" R3.01 1379 { 1380 // The above +"S" is that some files pluralize this... 1381 RinexSatID sat; 1382 // system 1383 satSysTemp = strip(line.substr(0,1)); 1384 1385 if(satSysTemp.empty()) 1386 { 1387 // continuation line 1388 satSysTemp = satSysPrev; 1389 1390 if(sysPhaseShift[satSysTemp].find(sysPhaseShiftObsID) 1391 == sysPhaseShift[satSysTemp].end()) 1392 { 1393 //FFStreamError e("SYS / PHASE SHIFT: unexpected continuation line"); 1394 //GPSTK_THROW(e); 1395 // lenient - some writers have only a single blank record 1396 return; 1397 } 1398 1399 map<RinexSatID,double>& satcorrmap(sysPhaseShift[satSysTemp][sysPhaseShiftObsID]); 1400 double cor(sysPhaseShift[satSysTemp][sysPhaseShiftObsID].begin()->second); 1401 for(i=0; i<10; i++) 1402 { 1403 string str = strip(line.substr(19+4*i,3)); 1404 if(str.empty()) break; 1405 sat = RinexSatID(str); 1406 satcorrmap.insert(make_pair(sat,cor)); 1407 } 1408 } 1409 else 1410 { 1411 // not a cont. line 1412 sat.fromString(satSysTemp); 1413 if(sysPhaseShift.find(satSysTemp) == sysPhaseShift.end()) 1414 { 1415 map<RinexObsID, map<RinexSatID, double> > obssatcormap; 1416 sysPhaseShift.insert(make_pair(satSysTemp,obssatcormap)); 1417 } 1418 1419 // obs id 1420 string str = strip(line.substr(2,3)); 1421 1422 // obsid and correction may be blank <=> unknown: ignore this 1423 if(!str.empty()) 1424 { 1425 RinexObsID obsid(satSysTemp+str, version); 1426 double cor(asDouble(strip(line.substr(6,8)))); 1427 int nsat(asInt(strip(line.substr(16,2)))); 1428 if(nsat > 0) 1429 { 1430 // list of sats 1431 map<RinexSatID,double> satcorrmap; 1432 for(i=0; i<(nsat < 10 ? nsat : 10); i++) 1433 { 1434 sat = RinexSatID(strip(line.substr(19+4*i,3))); 1435 satcorrmap.insert(make_pair(sat,cor)); 1436 } 1437 sysPhaseShift[satSysTemp].insert(make_pair(obsid,satcorrmap)); 1438 if(nsat > 10) // expect continuation 1439 sysPhaseShiftObsID = obsid; 1440 } 1441 else 1442 { 1443 // no sat, just system 1444 map<RinexSatID,double> satcorrmap; 1445 satcorrmap.insert(make_pair(sat,cor)); 1446 sysPhaseShift[satSysTemp].insert(make_pair(obsid,satcorrmap)); 1447 } 1448 } 1449 1450 // save for continuation lines 1451 satSysPrev = satSysTemp; 1452 1453 valid |= validSystemPhaseShift; 1454 } 1455 } 1456 else if(label == hsGlonassSlotFreqNo) 1457 { 1458 //map<RinexSatID,int> glonassFreqNo; 1459 int tmp; 1460 RinexSatID sat; 1461 string str(strip(line.substr(0,3))); 1462 1463 for(i=0; i<8; i++) 1464 { 1465 str = strip(line.substr(4+i*7,3)); 1466 if(str.empty()) break; 1467 sat = RinexSatID(str); 1468 str = strip(line.substr(8+i*7,2)); 1469 tmp = asInt(str); 1470 glonassFreqNo.insert(make_pair(sat,tmp)); 1471 } 1472 1473 valid |= validGlonassSlotFreqNo; 1474 } 1475 else if(label == hsGlonassCodPhsBias || label == hsGlonassCodPhsBias+"#") 1476 { 1477 // several IGS MGEX 3.02 files do the extra "#" 1478 //std::map<RinexObsID,double> glonassCodPhsBias; ///< "GLONASS COD/PHS/BIS" R3.02 1479 for(i=0; i<4; i++) 1480 { 1481 string str(strip(line.substr(i*13+1,3))); 1482 if(str.empty()) continue; 1483 RinexObsID obsid("R"+str, version); 1484 double bias(asDouble(strip(line.substr(i*13+5,8)))); 1485 glonassCodPhsBias[obsid] = bias; 1486 } 1487 valid |= validGlonassCodPhsBias; 1488 } 1489 else if(label == hsLeapSeconds) 1490 { 1491 leapSeconds = asInt(line.substr(0,6)); 1492 valid |= validLeapSeconds; 1493 } 1494 else if(label == hsNumSats) 1495 { 1496 numSVs = asInt(line.substr(0,6)) ; 1497 valid |= validNumSats; 1498 } 1499 else if(label == hsPrnObs) 1500 { 1501 // this assumes 'PRN / # OF OBS' comes after '# / TYPES OF OBSERV' or 'SYS / # / OBS TYPES' 1502 // NOT a good assumption for auxiliary header... ignore in that case 1503 if(version >= 3.0 && mapObsTypes.size() == 0) 1504 { 1505 commentList.push_back( 1506 string("Warning - can't read PRN/OBS in auxHeader: no" 1507 " mapObsTypes")); 1508 return; 1509 } 1510 1511 static const int maxObsPerLine = 9; 1512 1513 int j,otmax; 1514 RinexSatID PRN; 1515 string prn, GNSS; 1516 vector<int> numObsList; 1517 1518 prn = strip(line.substr(3,3)); 1519 1520 if(prn == "" ) // this is a continuation line; use last PRN 1521 { 1522 PRN = lastPRN; 1523 GNSS = PRN.systemChar(); 1524 if(version < 3) 1525 otmax = R2ObsTypes.size(); 1526 else 1527 { 1528 if(mapObsTypes.find(GNSS) == mapObsTypes.end()) 1529 { 1530 Exception e("PRN/#OBS for system "+PRN.toString()+" not found in SYS/#/OBS"); 1531 GPSTK_THROW(e); 1532 } 1533 otmax = mapObsTypes[GNSS].size(); 1534 } 1535 1536 numObsList = numObsForSat[PRN]; // grab the existing list 1537 1538 for(j=0,i=numObsList.size(); j<maxObsPerLine && i<otmax; i++,j++) 1539 numObsList.push_back(asInt(line.substr(6*j+6,6))); 1540 1541 numObsForSat[PRN] = numObsList; 1542 } 1543 else // this is a new PRN line 1544 { 1545 PRN = RinexSatID(prn); 1546 GNSS = PRN.systemChar(); 1547 if(version < 3) 1548 otmax = R2ObsTypes.size(); 1549 else 1550 { 1551 if(mapObsTypes.find(GNSS) == mapObsTypes.end()) 1552 { 1553 Exception e("PRN/#OBS for system "+PRN.toString()+" not found in SYS/#/OBS"); 1554 GPSTK_THROW(e); 1555 } 1556 otmax = mapObsTypes[GNSS].size(); 1557 } 1558 1559 for(i=0; i<maxObsPerLine && i<otmax; i++) 1560 numObsList.push_back(asInt(line.substr(6*i+6,6))); 1561 1562 numObsForSat[PRN] = numObsList; 1563 1564 lastPRN = PRN; 1565 } 1566 1567 valid |= validPrnObs; 1568 } 1569 else if(label == hsEoH) 1570 { 1571 validEoH = true; 1572 } 1573 else 1574 { 1575 FFStreamError e("Unidentified label: >" + label + "<"); 1576 GPSTK_THROW(e); 1577 } 1578 } // end of parseHeaderRecord 1579 1580 1581 // This function parses the entire header from the given stream reallyGetRecord(FFStream & ffs)1582 void Rinex3ObsHeader::reallyGetRecord(FFStream& ffs) 1583 { 1584 using gpstk::StringUtils::asString; 1585 Rinex3ObsStream& strm = dynamic_cast<Rinex3ObsStream&>(ffs); 1586 1587 // If already read, just return. 1588 if(strm.headerRead == true) return; 1589 1590 // Since we're reading a new header, we need to reinitialize 1591 // all our list structures. All the other objects should be 1592 // ok. This also applies if we threw an exception the first 1593 // time we read the header and are now re-reading it. Some 1594 // of these could be full and we need to empty them. 1595 clear(); 1596 1597 string line; 1598 1599 while (!validEoH) 1600 { 1601 strm.formattedGetLine(line); 1602 StringUtils::stripTrailing(line); 1603 1604 if(line.length() == 0) 1605 { 1606 FFStreamError e("No data read"); 1607 GPSTK_THROW(e); 1608 } 1609 else if(line.length() < 60 || line.length() > 80) 1610 { 1611 FFStreamError e("Invalid line length"); 1612 GPSTK_THROW(e); 1613 } 1614 1615 try 1616 { 1617 parseHeaderRecord(line); 1618 } 1619 catch(FFStreamError& e) 1620 { 1621 GPSTK_RETHROW(e); 1622 } 1623 catch(Exception& e) 1624 { 1625 FFStreamError fse("Exception: "+e.what()); 1626 GPSTK_THROW(fse); 1627 } 1628 1629 } // end while(not end of header) 1630 1631 // if RINEX 2, define mapObsTypes from R2ObsTypes and 1632 // system(s) this may have to be corrected later using 1633 // wavelengthFactor also define mapSysR2toR3ObsID in case 1634 // version 2 is written out later 1635 if(version < 3) 1636 { 1637 // try to determine systems included in the file 1638 vector<string> syss; // 1-char strings "G" "R" "E" ... 1639 if(numObsForSat.size() > 0) 1640 { 1641 // get syss from PRN/#OBS 1642 map<RinexSatID, vector<int> >::const_iterator it; 1643 for(it=numObsForSat.begin(); it != numObsForSat.end(); ++it) 1644 { 1645 string sys(string(1,(it->first).systemChar())); 1646 if(find(syss.begin(),syss.end(),sys) == syss.end()) 1647 syss.push_back(sys); 1648 } 1649 } 1650 else if(fileSysSat.system != SatelliteSystem::Mixed) 1651 { 1652 // only one system in this file 1653 syss.push_back(string(1,RinexSatID(fileSysSat).systemChar())); 1654 } 1655 else 1656 { 1657 // have to replicate obs type list for all RINEX2 systems 1658 syss.push_back("G"); 1659 syss.push_back("R"); 1660 syss.push_back("S"); // ?? 1661 syss.push_back("E"); 1662 } 1663 1664 // given systems and list of R2ObsTypes, compute 1665 // mapObsTypes and mapSysR2toR3ObsID 1666 mapSysR2toR3ObsID.clear(); 1667 for(size_t i=0; i<syss.size(); i++) 1668 { 1669 const string s(syss[i]); 1670 vector<RinexObsID> obsids; 1671 1672 try 1673 { 1674 if (s=="G") 1675 obsids = mapR2ObsToR3Obs_G(); 1676 else if (s=="R") 1677 obsids = mapR2ObsToR3Obs_R(); 1678 else if (s=="E") 1679 obsids = mapR2ObsToR3Obs_E(); 1680 else if (s=="S") 1681 obsids = mapR2ObsToR3Obs_S(); 1682 } 1683 catch(FFStreamError fse) 1684 { 1685 GPSTK_RETHROW(fse); 1686 } 1687 1688 // TD if GPS and have wavelengthFactors, add more 1689 // ObsIDs with tc=N 1690 1691 mapObsTypes[syss[i]] = obsids; 1692 } 1693 1694 // modify numObsForSat if necessary 1695 map<RinexSatID, vector<int> >::const_iterator it(numObsForSat.begin()); 1696 for( ; it != numObsForSat.end(); ++it) 1697 { 1698 RinexSatID sat(it->first); 1699 string sys; 1700 sys = sat.systemChar(); 1701 vector<int> vec; 1702 for(size_t i=0; i<R2ObsTypes.size(); i++) 1703 { 1704 if(mapSysR2toR3ObsID[sys][R2ObsTypes[i]].asString(version) == string(" ")) 1705 ; 1706 else 1707 vec.push_back(it->second[i]); 1708 } 1709 numObsForSat[sat] = vec; 1710 } 1711 } 1712 1713 // Since technically the Phase Shift record is required in ver 3.01, 1714 // create SystemPhaseShift record(s) if not present. 1715 //map<string, map<RinexObsID, map<RinexSatID,double> > > sysPhaseShift; 1716 if(version >= 3.01 && (valid & validSystemNumObs) 1717 && !(valid & validSystemPhaseShift)) 1718 { 1719 // loop over obs types to get systems 1720 map<string,vector<RinexObsID> >::const_iterator iter; 1721 for(iter=mapObsTypes.begin(); iter != mapObsTypes.end(); iter++) 1722 { 1723 string sys(iter->first); 1724 if(sysPhaseShift.find(sys) == sysPhaseShift.end()) 1725 { 1726 map<RinexObsID, map<RinexSatID, double> > dummy; 1727 sysPhaseShift.insert(make_pair(sys,dummy)); 1728 } 1729 } 1730 valid |= validSystemPhaseShift; 1731 } 1732 1733 // is the header valid? 1734 Fields allValid = Fields::getRequired(version); 1735 if (allValid.empty()) 1736 { 1737 FFStreamError e("Unknown or unsupported RINEX version " + 1738 asString(version,2)); 1739 GPSTK_THROW(e); 1740 } 1741 1742 if((valid & allValid) != allValid) 1743 { 1744 FFStreamError e("Incomplete or invalid header"); 1745 allValid.describeMissing(valid, e); 1746 GPSTK_THROW(e); 1747 } 1748 1749 // If we get here, we should have reached the end of header line. 1750 strm.header = *this; 1751 strm.headerRead = true; 1752 1753 // determine the time system of epochs in this file; cf. R3.02 Table A2 1754 // 1.determine time system from time tag in TIME OF FIRST OBS record 1755 // 2.if not given, determine from type in RINEX VERSION / TYPE record 1756 // 3.(if the type is MIXED, the time system in firstObs is required by RINEX) 1757 strm.timesystem = firstObs.getTimeSystem(); 1758 if(strm.timesystem == TimeSystem::Any || 1759 strm.timesystem == TimeSystem::Unknown) 1760 { 1761 if(fileSysSat.system == SatelliteSystem::GPS) 1762 { 1763 strm.timesystem = TimeSystem::GPS; 1764 firstObs.setTimeSystem(TimeSystem::GPS); 1765 } 1766 else if(fileSysSat.system == SatelliteSystem::Glonass) 1767 { 1768 strm.timesystem = TimeSystem::UTC; 1769 firstObs.setTimeSystem(TimeSystem::UTC); 1770 } 1771 else if(fileSysSat.system == SatelliteSystem::Galileo) 1772 { 1773 strm.timesystem = TimeSystem::GAL; 1774 firstObs.setTimeSystem(TimeSystem::GAL); 1775 } 1776 else if(fileSysSat.system == SatelliteSystem::QZSS) 1777 { 1778 strm.timesystem = TimeSystem::QZS; 1779 firstObs.setTimeSystem(TimeSystem::QZS); 1780 } 1781 else if(fileSysSat.system == SatelliteSystem::BeiDou) 1782 { 1783 strm.timesystem = TimeSystem::BDT; 1784 firstObs.setTimeSystem(TimeSystem::BDT); 1785 } 1786 else if(fileSysSat.system == SatelliteSystem::Mixed) 1787 { 1788 // lenient 1789 strm.timesystem = TimeSystem::GPS; 1790 firstObs.setTimeSystem(TimeSystem::GPS); 1791 // RINEX 3 requires this 1792 //FFStreamError e("TimeSystem in MIXED files must be given by first obs"); 1793 //GPSTK_THROW(e); 1794 } 1795 else 1796 { 1797 FFStreamError e("Unknown file system type"); 1798 GPSTK_THROW(e); 1799 } 1800 } 1801 1802 } // end reallyGetRecord 1803 1804 // This method maps v2.11 GPS observation types to the v3 equivalent. 1805 // Since only GPS and only v2.11 are of interest, only L1/L2/L5 1806 // are considered. mapR2ObsToR3Obs_G()1807 vector<RinexObsID> Rinex3ObsHeader::mapR2ObsToR3Obs_G() 1808 { 1809 vector<RinexObsID> obsids; 1810 1811 // Assume D1, S1, and L1 come from C/A unless P is being treated as Y and P1 is present 1812 // Furthermore, if P1 is present and P is NOT being treated as Y, assume that P1 1813 // is some Z-mode or equivalent "smart" codeless process. 1814 // 1815 // Condition Result 1816 // PisY P1? 1817 // N Y L1,D1,S1 considered C, P1 becomes C1W 1818 // N N L1,D1,S1 considered C 1819 // Y Y L1,D1,S1 considered Y, P1 becomes C1Y 1820 // Y N L1,D1,S1 considered C 1821 // 1822 bool hasL1P = find(R2ObsTypes.begin(),R2ObsTypes.end(), string("P1")) != R2ObsTypes.end(); 1823 string code1 = "C"; 1824 string code1P = "W"; 1825 if (PisY && hasL1P) 1826 { 1827 code1 = "Y"; 1828 code1P = "Y"; 1829 } 1830 1831 // Assume D2, S2, and L2 come from Y if P is being treated as Y and P2 is present 1832 // codeless unless L2C is tracked. 1833 // If BOTH C2 and P2 are present, and P is NOT being treated as Y, assume C2 1834 // is code tracking the open signal and that P2 is codelessly tracking an 1835 // authorized signal. 1836 // 1837 // Condition Result 1838 // PisY C2? P2? 1839 // N Y N L2,D2,S2 considered X, 1840 // N Y Y L2,D2,S2 considered W, P2 becomes C2W** 1841 // N N Y L2,D2,S2 considered W, P2 becomes C2W 1842 // N N N L2,D2,S2 considered X* 1843 // Y Y N L2,D2,S2 considered X 1844 // Y Y Y L2,D2,S2 considered Y, P2 becomes C2Y 1845 // Y N Y L2,D2,S2 considered Y, P2 becomes C2Y 1846 // Y N N L2,D2,S2 considered X* 1847 // * - Probably not a reasonable set of conditions. It implies no L2 pseudoranges 1848 // were collected on any tracking code. 1849 // **- Interesting case. Currently presence of C2 in the header means 1850 // that the data MAY be present. However, since only some of the GPS 1851 // SVs have L2C, the C2 data field will frequently be empty. 1852 // Therefore, we'll go with "W" if P2 is present. The other option 1853 // would be to add smarts to the SV-by-SV record reading process to 1854 // coerce this to X if there are actually data in the C2 field at 1855 // the time the observations are read. That would really do violence 1856 // to the existing logic. Better to hope for a transition to Rinex 3 1857 // before this becomes a real issue. 1858 // 1859 // N.B.: This logic (both for P1 and P2) assumes P is NEVER P. If we want to allow for 1860 // live sky (or simulator capture) P code, we'll have to add more logic 1861 // to differentate between PisY, PisW, and PisP. That will have to be 1862 // "beyond RINEX v2.11" extra-special handling. 1863 // 1864 bool hasL2P = find(R2ObsTypes.begin(),R2ObsTypes.end(), string("P2")) != R2ObsTypes.end(); 1865 bool hasL2C = find(R2ObsTypes.begin(),R2ObsTypes.end(), string("C2")) != R2ObsTypes.end(); 1866 1867 string code2 = "X"; // Correct condition as long as P2 is not in the list 1868 string code2P = "X"; // Condition is irrelvant unless P2 is in the list 1869 if (hasL2P) 1870 { 1871 if (PisY) 1872 { 1873 code2 = "Y"; 1874 code2P = "Y"; 1875 } 1876 else 1877 { 1878 code2 = "W"; 1879 code2P = "W"; 1880 } 1881 } 1882 string syss("G"); 1883 for(size_t j=0; j<R2ObsTypes.size(); ++j) 1884 { 1885 string ot(R2ObsTypes[j]); 1886 string obsid(syss); 1887 1888 if (ot == "C1") obsid += "C1C"; 1889 else if (ot == "P1") obsid += "C1" + code1P; 1890 else if (ot == "L1") obsid += "L1" + code1; 1891 else if (ot == "D1") obsid += "D1" + code1; 1892 else if (ot == "S1") obsid += "S1" + code1; 1893 1894 else if (ot == "C2") obsid += "C2X"; 1895 else if (ot == "P2") obsid += "C2" + code2P; 1896 else if (ot == "L2") obsid += "L2" + code2; 1897 else if (ot == "D2") obsid += "D2" + code2; 1898 else if (ot == "S2") obsid += "S2" + code2; 1899 1900 else if (ot == "C5") obsid += "C5X"; 1901 else if (ot == "L5") obsid += "L5X"; 1902 else if (ot == "D5") obsid += "D5X"; 1903 else if (ot == "S5") obsid += "S5X"; 1904 1905 // If the obs type isn't valid for GPS, skip it. 1906 else continue; 1907 1908 try 1909 { 1910 RinexObsID OT(obsid, version); 1911 obsids.push_back(OT); 1912 mapSysR2toR3ObsID[syss][ot] = OT; //map<string, map<string, RinexObsID> > 1913 } 1914 catch(InvalidParameter& ip) 1915 { 1916 FFStreamError fse("InvalidParameter: "+ip.what()); 1917 GPSTK_THROW(fse); 1918 } 1919 } 1920 return obsids; 1921 } 1922 1923 // This method maps v2.11 GLONASS observation types to the v3 equivalent. 1924 // Since only GLONASS and only v2.11 are of interest, only L1/L2 1925 // are considered. mapR2ObsToR3Obs_R()1926 vector<RinexObsID> Rinex3ObsHeader::mapR2ObsToR3Obs_R( ) 1927 { 1928 vector<RinexObsID> obsids; 1929 1930 // Assume D1, S1, and L1 come from C/A 1931 // This assumes that any files claiming to track GLONASS P1 is 1932 // actually doing so with a codeless technique. There is no RINEX V3 1933 // "C1W" for GLONASS, so we'll leave P1 as C1P as the closest approximation. 1934 bool hasL1P = find(R2ObsTypes.begin(),R2ObsTypes.end(), string("P1")) != R2ObsTypes.end(); 1935 string code1 = "C"; 1936 1937 // Assume D2, S2, and L2 come from C/A. Same logic as above. 1938 bool hasL2P = find(R2ObsTypes.begin(),R2ObsTypes.end(), string("P2")) != R2ObsTypes.end(); 1939 string code2 = "C"; 1940 1941 string syss("R"); 1942 for(size_t j=0; j<R2ObsTypes.size(); ++j) 1943 { 1944 string ot(R2ObsTypes[j]); 1945 string obsid(syss); 1946 1947 if (ot == "C1") obsid += "C1C"; 1948 else if (ot == "P1") obsid += "C1P"; 1949 else if (ot == "L1") obsid += "L1" + code1; 1950 else if (ot == "D1") obsid += "D1" + code1; 1951 else if (ot == "S1") obsid += "S1" + code1; 1952 1953 else if (ot == "C2") obsid += "C2C"; 1954 else if (ot == "P2") obsid += "C2P"; 1955 else if (ot == "L2") obsid += "L2" + code2; 1956 else if (ot == "D2") obsid += "D2" + code2; 1957 else if (ot == "S2") obsid += "S2" + code2; 1958 1959 // If the obs type isn't valid for GLONASS, skip it. 1960 else continue; 1961 1962 try 1963 { 1964 RinexObsID OT(obsid, version); 1965 obsids.push_back(OT); 1966 mapSysR2toR3ObsID[syss][ot] = OT; //map<string, map<string, RinexObsID> > 1967 } 1968 catch(InvalidParameter& ip) 1969 { 1970 FFStreamError fse("InvalidParameter: "+ip.what()); 1971 GPSTK_THROW(fse); 1972 } 1973 } 1974 return obsids; 1975 } 1976 1977 // This method maps v2.11 Galileo observation types to the v3 equivalent. 1978 // Since only Galileo and only v2.11 are of interest no L2 types 1979 // are considered. Furthermore, Rinex v2.11 states that there is no 1980 // P for Galileo. (Where that leaves the PRS is a good question.) 1981 // 1982 // In RINEX v3, there are 3-5 tracking codes defined for each carrier. 1983 // Given the current lack of experience, the code makes some 1984 // guesses on what the v2.11 translations should mean. mapR2ObsToR3Obs_E()1985 vector<RinexObsID> Rinex3ObsHeader::mapR2ObsToR3Obs_E() 1986 { 1987 vector<RinexObsID> obsids; 1988 1989 string code1 = "B"; // Corresponds to the open service 1990 string code5 = "I"; // Corresponds to the open service 1991 string code7 = "X"; // Corresponds to I + Q tracking 1992 string code8 = "X"; // Corresponds to I + Q tracking 1993 string code6 = "X"; // Corresponds to B + C tracking 1994 1995 string syss("E"); 1996 for(size_t j=0; j<R2ObsTypes.size(); ++j) 1997 { 1998 string ot(R2ObsTypes[j]); 1999 string obsid(syss); 2000 if (ot == "C1") obsid += "C1" + code1; 2001 else if (ot == "L1") obsid += "L1" + code1; 2002 else if (ot == "D1") obsid += "D1" + code1; 2003 else if (ot == "S1") obsid += "S1" + code1; 2004 2005 else if (ot == "C5") obsid += "C5" + code5; 2006 else if (ot == "L5") obsid += "L5" + code5; 2007 else if (ot == "D5") obsid += "D5" + code5; 2008 else if (ot == "S5") obsid += "S5" + code5; 2009 2010 else if (ot == "C6") obsid += "C6" + code6; 2011 else if (ot == "L6") obsid += "L6" + code6; 2012 else if (ot == "D6") obsid += "D6" + code6; 2013 else if (ot == "S6") obsid += "S6" + code6; 2014 2015 else if (ot == "C7") obsid += "C7" + code7; 2016 else if (ot == "L7") obsid += "L7" + code7; 2017 else if (ot == "D7") obsid += "D7" + code7; 2018 else if (ot == "S7") obsid += "S7" + code7; 2019 2020 else if (ot == "C8") obsid += "C8" + code8; 2021 else if (ot == "L8") obsid += "L8" + code8; 2022 else if (ot == "D8") obsid += "D8" + code8; 2023 else if (ot == "S8") obsid += "S8" + code8; 2024 2025 // If the obs type isn't valid for Galileo, skip it. 2026 else continue; 2027 2028 try 2029 { 2030 RinexObsID OT(obsid, version); 2031 obsids.push_back(OT); 2032 mapSysR2toR3ObsID[syss][ot] = OT; //map<string, map<string, RinexObsID> > 2033 } 2034 catch(InvalidParameter& ip) 2035 { 2036 FFStreamError fse("InvalidParameter: "+ip.what()); 2037 GPSTK_THROW(fse); 2038 } 2039 } 2040 return obsids; 2041 } 2042 2043 2044 // This method maps v2.11 SBAS observation types to the v3 equivalent. 2045 // Since only SBAS and only v2.11 are of interest only L1/L5 2046 // are considered. mapR2ObsToR3Obs_S()2047 vector<RinexObsID> Rinex3ObsHeader::mapR2ObsToR3Obs_S() 2048 { 2049 vector<RinexObsID> obsids; 2050 2051 string code1 = "C"; // Only option 2052 string code5 = "X"; // Corresponds to I + Q tracking 2053 2054 string syss("S"); 2055 for(size_t j=0; j<R2ObsTypes.size(); ++j) 2056 { 2057 string ot(R2ObsTypes[j]); 2058 string obsid(syss); 2059 if (ot == "C1") obsid += "C1" + code1; 2060 else if (ot == "L1") obsid += "L1" + code1; 2061 else if (ot == "D1") obsid += "D1" + code1; 2062 else if (ot == "S1") obsid += "S1" + code1; 2063 2064 else if (ot == "C5") obsid += "C5" + code5; 2065 else if (ot == "L5") obsid += "L5" + code5; 2066 else if (ot == "D5") obsid += "D5" + code5; 2067 else if (ot == "S5") obsid += "S5" + code5; 2068 2069 // If the obs type isn't valid for SBAS, skip it. 2070 else continue; 2071 2072 try 2073 { 2074 RinexObsID OT(obsid, version); 2075 obsids.push_back(OT); 2076 mapSysR2toR3ObsID[syss][ot] = OT; //map<string, map<string, RinexObsID> > 2077 } 2078 catch(InvalidParameter& ip) 2079 { 2080 FFStreamError fse("InvalidParameter: "+ip.what()); 2081 GPSTK_THROW(fse); 2082 } 2083 } 2084 return obsids; 2085 } 2086 2087 parseTime(const string & line) const2088 CivilTime Rinex3ObsHeader::parseTime(const string& line) const 2089 { 2090 int year, month, day, hour, min; 2091 double sec; 2092 string tsys; 2093 TimeSystem ts; 2094 2095 year = asInt( line.substr(0, 6)); 2096 month = asInt( line.substr(6, 6)); 2097 day = asInt( line.substr(12, 6)); 2098 hour = asInt( line.substr(18, 6)); 2099 min = asInt( line.substr(24, 6)); 2100 sec = asDouble(line.substr(30, 13)); 2101 tsys = line.substr(48, 3) ; 2102 2103 ts = gpstk::StringUtils::asTimeSystem(tsys); 2104 2105 return CivilTime(year, month, day, hour, min, sec, ts); 2106 } // end parseTime 2107 2108 writeTime(const CivilTime & civtime) const2109 string Rinex3ObsHeader::writeTime(const CivilTime& civtime) const 2110 { 2111 using gpstk::StringUtils::asString; 2112 string line, tsStr(gpstk::StringUtils::asString(civtime.getTimeSystem())); 2113 2114 line = rightJustify(asString<short>(civtime.year ) , 6); 2115 line += rightJustify(asString<short>(civtime.month ) , 6); 2116 line += rightJustify(asString<short>(civtime.day ) , 6); 2117 line += rightJustify(asString<short>(civtime.hour ) , 6); 2118 line += rightJustify(asString<short>(civtime.minute ) , 6); 2119 line += rightJustify(asString( civtime.second,7) , 13); 2120 line += rightJustify(tsStr, 8); 2121 2122 return line; 2123 } // end writeTime 2124 2125 // Compute map of obs types for use in writing version 2 header 2126 // and data, call before writing prepareVer2Write(void)2127 void Rinex3ObsHeader::prepareVer2Write(void) 2128 { 2129 if(version > 3) 2130 { 2131 version = 2.11; 2132 } 2133 valid |= Rinex3ObsHeader::validWaveFact; 2134 if (valid.isSet(Rinex3ObsHeader::validSystemNumObs)) 2135 { 2136 // change from RINEX 3 header to RINEX 2 equivalent 2137 valid.clear(Rinex3ObsHeader::validSystemNumObs); 2138 valid.set(Rinex3ObsHeader::validNumObs); 2139 } 2140 /// @todo unset R3-specific header members? 2141 2142 // make a list of R2 obstype strings, and a map R3ObsIDs <= 2143 // R2 obstypes for each system 2144 R2ObsTypes.clear(); 2145 map<string,vector<RinexObsID> >::const_iterator mit; 2146 for (mit = mapObsTypes.begin(); mit != mapObsTypes.end(); mit++) 2147 { 2148 string sysString = mit->first; 2149 // Compass/BeiDou and QZSS and IRNSS are unsupported by RINEX 2.11. 2150 // Transit was deprecated in 2.11. 2151 // Skip them. 2152 if ((sysString == "I") || (sysString == "J") || (sysString == "C") || 2153 (sysString == "T")) 2154 { 2155 continue; 2156 } 2157 // 2.11 supports only GPS, GLONASS, Geo/SBAS and Galileo 2158 if(sysString!="G" && sysString!="R" && sysString!="E" && 2159 sysString!="S") 2160 { 2161 FFStreamError er( 2162 "Invalid system char string in header.mapObsTypes: "+sysString); 2163 GPSTK_THROW(er); 2164 } 2165 // mit->first is system char as a 1-char string 2166 map<string, RinexObsID> mapR2toR3ObsID; 2167 2168 vector<string> uniqueCheck; 2169 // loop over all ObsIDs for this system 2170 for(size_t i=0; i<mit->second.size(); i++) 2171 { 2172 string R2ot, lab(mit->second[i].asString(version)); 2173 // the list of all tracking code characters for this sys, freq 2174 string allCodes( 2175 RinexObsID::validRinexTrackingCodes[mit->first[0]][lab[1]]); 2176 2177 if (lab == string("C1C")) 2178 R2ot = string("C1"); 2179 else if (lab == string("C2X") && mit->first == "G") 2180 R2ot = string("C2"); 2181 else if (lab == string("C2C") && mit->first == "R") 2182 R2ot = string("C2"); 2183 // R2 has C5 but not P5 2184 else if (lab.substr(0,2) == "C5") 2185 R2ot = string("C5"); 2186 else if (lab[0] == 'C') 2187 R2ot = string("P")+string(1,lab[1]); 2188 else 2189 R2ot = lab.substr(0,2); 2190 // add to list, if not already there 2191 vector<string>::iterator it; 2192 it = find(R2ObsTypes.begin(),R2ObsTypes.end(),R2ot); 2193 if (it == R2ObsTypes.end()) 2194 { 2195 // its not there - add it 2196 R2ObsTypes.push_back(R2ot); 2197 mapR2toR3ObsID[R2ot] = mit->second[i]; 2198 } 2199 else 2200 { 2201 // its already there - in list of R2 ots 2202 if (mapR2toR3ObsID.find(R2ot) == mapR2toR3ObsID.end()) 2203 { 2204 // must also add to sys map 2205 mapR2toR3ObsID[R2ot] = mit->second[i]; 2206 } 2207 else 2208 { 2209 // its already in sys map ... 2210 // .. but is the new tc 'better'? 2211 string::size_type posold,posnew; 2212 posold = allCodes.find((mapR2toR3ObsID[R2ot].asString(version))[2]); 2213 posnew = allCodes.find(lab[2]); 2214 if(posnew < posold) 2215 { 2216 // replace the R3ObsID in the map 2217 mapR2toR3ObsID[R2ot] = mit->second[i]; 2218 } 2219 if (R2DisambiguityMap.find(sysString + R2ot) == 2220 R2DisambiguityMap.end()) 2221 { 2222 R2DisambiguityMap.insert( 2223 std::pair<string,string>( 2224 sysString + R2ot, 2225 mapR2toR3ObsID[R2ot].asString(version))); 2226 } 2227 else 2228 { 2229 R2DisambiguityMap[sysString + R2ot] = lab; 2230 } 2231 } 2232 } 2233 } 2234 // save for this system 2235 mapSysR2toR3ObsID[mit->first] = mapR2toR3ObsID; 2236 } 2237 } // end prepareVer2Write() 2238 dump(ostream & s,double dumpVersion) const2239 void Rinex3ObsHeader::dump(ostream& s, double dumpVersion) const 2240 { 2241 using gpstk::StringUtils::asString; 2242 size_t i; 2243 2244 string str; 2245 if(fileSysSat.system == SatelliteSystem::Mixed) 2246 str = "MIXED"; 2247 else 2248 { 2249 RinexSatID sat(fileSysSat); 2250 str = sat.systemChar(); 2251 str = str + " (" + sat.systemString() + ")"; 2252 } 2253 2254 s << "---------------------------------- REQUIRED " 2255 << "----------------------------------" << endl; 2256 s << "Rinex Version " << fixed << setw(5) << setprecision(2) << dumpVersion 2257 << ", File type " << fileType << ", System " << str << "." << endl; 2258 s << "Prgm: " << fileProgram << ", Run: " << date 2259 << ", By: " << fileAgency << endl; 2260 s << "Marker name: " << markerName << ", "; 2261 s << "Marker type: " << markerType << "." << endl; 2262 s << "Observer : " << observer << ", Agency: " << agency << endl; 2263 s << "Rec#: " << recNo << ", Type: " << recType 2264 << ", Vers: " << recVers << endl; 2265 s << "Antenna # : " << antNo << ", Type : " << antType << endl; 2266 s << "Position (XYZ,m) : " << setprecision(4) << antennaPosition 2267 << "." << endl; 2268 s << "Antenna Delta (HEN,m) : " << setprecision(4) << antennaDeltaHEN 2269 << "." << endl; 2270 map<string,vector<RinexObsID> >::const_iterator iter; 2271 for(iter = mapObsTypes.begin(); iter != mapObsTypes.end(); iter++) 2272 { 2273 RinexSatID rsid; 2274 rsid.fromString(iter->first); 2275 s << rsid.systemString() << " Observation types (" 2276 << iter->second.size() << "):" << endl; 2277 for(i = 0; i < iter->second.size(); i++) 2278 { 2279 s << " Type #" << setw(2) << setfill('0') << i+1 << setfill(' ') 2280 << " (" << iter->second[i].asString(dumpVersion) << ") " 2281 << asString(static_cast<ObsID>(iter->second[i])) << endl; 2282 } 2283 } 2284 2285 s << "R2ObsTypes: "; 2286 for (StringVec::const_iterator i=R2ObsTypes.begin(); i != R2ObsTypes.end(); i++) 2287 s << *i << " "; 2288 s << endl; 2289 for (VersionObsMap::const_iterator i = mapSysR2toR3ObsID.begin(); i != mapSysR2toR3ObsID.end(); i++) 2290 { 2291 s << "mapSysR2toR3ObsID[" << i->first << "] "; 2292 for (ObsIDMap::const_iterator j = i->second.begin(); j != i->second.end(); j++) 2293 s << j->first << ":" << j->second.asString(dumpVersion) << " "; 2294 s << endl; 2295 } 2296 2297 s << "Time of first obs " 2298 << printTime(firstObs,"%04Y/%02m/%02d %02H:%02M:%06.3f %P") << endl; 2299 2300 s << "(This header is "; 2301 if (Fields::isValid(dumpVersion, valid)) 2302 s << "VALID)" << endl; 2303 else 2304 { 2305 s << "NOT VALID"; 2306 s << " RINEX " << setprecision(2) << dumpVersion << ")" << endl; 2307 s << "valid = " << valid << endl; 2308 Fields required = Fields::getRequired(dumpVersion); 2309 s << "allValid = " << required << endl; 2310 2311 s << "Invalid or missing header records:" << endl; 2312 // print all header fields in required not in valid 2313 for (const auto& reqi : required.fieldsSet) 2314 { 2315 if (valid.isSet(reqi) == 0) 2316 { 2317 s << " " << setw(2) << reqi << " " 2318 << Rinex3ObsHeader::asString(reqi) << endl; 2319 } 2320 } 2321 s << "END Invalid header records." << endl; 2322 } 2323 2324 s << "---------------------------------- OPTIONAL " 2325 << "----------------------------------" << endl; 2326 if(valid & validMarkerNumber ) 2327 s << "Marker number : " << markerNumber << endl; 2328 if(valid & validMarkerType ) 2329 s << "Marker type : " << markerType << endl; 2330 if(valid & validAntennaDeltaXYZ ) 2331 s << "Antenna Delta (XYZ,m) : " 2332 << setprecision(4) << antennaDeltaXYZ << endl; 2333 if(valid & validAntennaPhaseCtr ) 2334 s << "Antenna PhaseCtr (XYZ,m) : " 2335 << setprecision(4) << antennaPhaseCtr << endl; 2336 if(valid & validAntennaBsightXYZ ) 2337 s << "Antenna B.sight (XYZ,m) : " 2338 << setprecision(4) << antennaBsightXYZ << endl; 2339 if(valid & validAntennaZeroDirAzi) 2340 s << "Antenna ZeroDir (deg) : " 2341 << setprecision(4) << antennaZeroDirAzi << endl; 2342 if(valid & validAntennaZeroDirXYZ) 2343 s << "Antenna ZeroDir (XYZ,m) : " 2344 << setprecision(4) << antennaZeroDirXYZ << endl; 2345 if(valid & validCenterOfMass ) 2346 s << "Center of Mass (XYZ,m) : " 2347 << setprecision(4) << antennaPhaseCtr << endl; 2348 if(valid & validSigStrengthUnit ) 2349 s << "Signal Strength Unit = " << sigStrengthUnit << endl; 2350 if(valid & validInterval ) 2351 s << "Interval = " 2352 << fixed << setw(7) << setprecision(3) << interval << endl; 2353 if(valid & validLastTime ) 2354 s << "Time of Last Obs " 2355 << printTime(lastObs,"%04Y/%02m/%02d %02H:%02M:%06.3f %P") << endl; 2356 if(valid & validReceiverOffset ) 2357 s << "Clock offset record is present and offsets " 2358 << (receiverOffset ? "ARE" : "are NOT") << " applied." << endl; 2359 if(dumpVersion < 3 && (valid & validWaveFact)) // TD extraWaveFactList 2360 s << "Wavelength factor L1: " << wavelengthFactor[0] 2361 << " L2: " << wavelengthFactor[1] << endl; 2362 if(valid & validSystemDCBSapplied) 2363 { 2364 for(i = 0; i < infoDCBS.size(); i++) 2365 { 2366 RinexSatID rsid; 2367 rsid.fromString(infoDCBS[i].satSys); 2368 s << "System DCBS Correction Applied to " << rsid.systemString() 2369 << " data using program " << infoDCBS[i].name << endl; 2370 s << " from source " << infoDCBS[i].source << "." << endl; 2371 } 2372 } 2373 if(valid & validSystemPCVSapplied) 2374 { 2375 for(i = 0; i < infoPCVS.size(); i++) 2376 { 2377 RinexSatID rsid; 2378 rsid.fromString(infoPCVS[i].satSys); 2379 s << "System PCVS Correction Applied to " << rsid.systemString() 2380 << " data using program " << infoPCVS[i].name << endl; 2381 s << " from source " << infoPCVS[i].source << "." << endl; 2382 } 2383 } 2384 if(valid & validSystemScaleFac ) 2385 { 2386 map<string, ScaleFacMap>::const_iterator mapIter; 2387 // loop over GNSSes 2388 for(mapIter = sysSfacMap.begin(); mapIter != sysSfacMap.end(); mapIter++) 2389 { 2390 RinexSatID rsid; 2391 rsid.fromString(mapIter->first); 2392 s << rsid.systemString() << " scale factors applied:" << endl; 2393 map<RinexObsID,int>::const_iterator iter; 2394 // loop over scale factor map 2395 for(iter = mapIter->second.begin(); iter != mapIter->second.end(); iter++) 2396 s << " " << iter->first.asString(dumpVersion) << " " << iter->second << endl; 2397 } 2398 } 2399 if(valid & validSystemPhaseShift ) 2400 { 2401 map<string, map<RinexObsID, map<RinexSatID,double> > >::const_iterator it; 2402 for(it=sysPhaseShift.begin(); it!=sysPhaseShift.end(); ++it) 2403 { 2404 string sys(it->first); 2405 map<RinexObsID, map<RinexSatID, double> >::const_iterator jt; 2406 jt = it->second.begin(); 2407 if(jt == it->second.end()) 2408 s << "Phase shift correction for system " << sys << " is empty." << endl; 2409 for( ; jt!=it->second.end(); ++jt) 2410 { 2411 map<RinexSatID,double>::const_iterator kt; 2412 for(kt=jt->second.begin(); kt!=jt->second.end(); ++kt) 2413 s << "Phase shift correction for system " << sys << ": " 2414 << fixed << setprecision(5) 2415 << setw(8) << kt->second << " cycles applied to obs type " 2416 << jt->first.asString(dumpVersion) << " " 2417 << RinexSatID(sys).systemString() << endl; 2418 } 2419 } 2420 } 2421 if(valid & validGlonassSlotFreqNo) 2422 { 2423 int n(0); 2424 map<RinexSatID,int>::const_iterator it; 2425 s << "GLONASS frequency channels:\n"; 2426 for(it=glonassFreqNo.begin(); it!=glonassFreqNo.end(); ++it) 2427 { 2428 s << " " << it->first.toString() << " " << setw(2) << it->second; 2429 if(++n > 1 && (n%8)==0) s << endl; 2430 } 2431 if((n%8) != 0) s << endl; 2432 } 2433 if(valid & validGlonassCodPhsBias) 2434 { 2435 map<RinexObsID,double>::const_iterator it; 2436 s << "GLONASS Code-phase biases:\n" << fixed << setprecision(3); 2437 for(it=glonassCodPhsBias.begin(); it!=glonassCodPhsBias.end(); ++it) 2438 s << " " << it->first.asString(dumpVersion) << " " << setw(8) << it->second; 2439 s << endl; 2440 } 2441 if(valid & validLeapSeconds) 2442 s << "Leap seconds: " << leapSeconds << endl; 2443 if(valid & validNumSats) 2444 s << "Number of Satellites with data : " << numSVs << endl; 2445 if(valid & validPrnObs) 2446 { 2447 RinexSatID sat, sys(-1,SatelliteSystem::Unknown); 2448 s << " PRN and number of observations for each obs type:" << endl; 2449 map<RinexSatID, vector<int> >::const_iterator it = numObsForSat.begin(); 2450 while (it != numObsForSat.end()) 2451 { 2452 sat = it->first; 2453 if(sat.system != sys.system) 2454 { 2455 // print a header: SYS OT OT OT ... 2456 s << " " << sat.systemString3() << " "; 2457 iter = mapObsTypes.find(string(1,sat.systemChar())); 2458 const vector<RinexObsID>& vec(iter->second); 2459 for(i=0; i<vec.size(); i++) 2460 s << setw(7) << vec[i].asString(dumpVersion); 2461 s << endl; 2462 sys = sat; 2463 } 2464 vector<int> obsvec = it->second; 2465 s << " " << sat.toString() << " "; 2466 for(i = 0; i < obsvec.size(); i++) // print the numbers of obss 2467 s << " " << setw(6) << obsvec[i]; 2468 s << endl; 2469 it++; 2470 } 2471 } 2472 if(commentList.size()) 2473 { 2474 if(!(valid & validComment)) s << " Comment list is NOT valid" << endl; 2475 s << "Comments (" << commentList.size() << ") :" << endl; 2476 for(i=0; i<commentList.size(); i++) s << commentList[i] << endl; 2477 } 2478 2479 s << "-------------------------------- END OF HEADER " 2480 << "--------------------------------" << endl; 2481 } // end dump 2482 2483 2484 /* This method returns the numerical index of a given observation 2485 * 2486 * @param type String representing the observation type. 2487 */ getObsIndex(const string & type) const2488 size_t Rinex3ObsHeader::getObsIndex( const string& type ) const 2489 { 2490 string newType(type); 2491 2492 // 'old-style' type: Let's change it to 'new style'. 2493 if( newType.size() == 2 ) 2494 { 2495 if( newType == "C1" ) newType = "C1C"; 2496 else if( newType == "P1" ) newType = "C1P"; 2497 else if( newType == "L1" ) newType = "L1P"; 2498 else if( newType == "D1" ) newType = "D1P"; 2499 else if( newType == "S1" ) newType = "S1P"; 2500 else if( newType == "C2" ) newType = "C2C"; 2501 else if( newType == "P2" ) newType = "C2P"; 2502 else if( newType == "L2" ) newType = "L2P"; 2503 else if( newType == "D2" ) newType = "D2P"; 2504 else if( newType == "S2" ) newType = "S2P"; 2505 else 2506 { 2507 InvalidRequest exc("Invalid type."); 2508 GPSTK_THROW(exc); 2509 } 2510 } 2511 2512 // Add GNSS code. By default the system is GPS 2513 if( newType.size() == 3 ) 2514 { 2515 newType = "G" + newType; 2516 } 2517 2518 // Check if resulting 'newType' is valid 2519 if( !isValidRinexObsID(newType) ) 2520 { 2521 InvalidRequest ir(newType + " is not a valid RinexObsID!."); 2522 GPSTK_THROW(ir); 2523 } 2524 2525 // Extract the GNSS from the newType 2526 string sys( newType, 0, 1 ); 2527 return getObsIndex(sys, RinexObsID(newType, version)); 2528 } 2529 2530 getObsIndex(const string & sys,const RinexObsID & obsID) const2531 size_t Rinex3ObsHeader::getObsIndex(const string& sys, 2532 const RinexObsID& obsID ) const 2533 { 2534 /// typedef std::vector<RinexObsID> RinexObsVec; 2535 /// typedef std::map<std::string, RinexObsVec> RinexObsMap; 2536 /// RinexObsMap mapObsTypes; ///< SYS / # / OBS TYPES 2537 2538 // find the GNSS in the map 2539 RinexObsMap remapped; 2540 std::map<std::string,unsigned> obsCount; 2541 remapObsTypes(remapped, obsCount); 2542 RinexObsMap::const_iterator it = remapped.find(sys); 2543 2544 if (it == remapped.end()) 2545 { 2546 InvalidRequest ir("GNSS system " + sys + " not stored."); 2547 GPSTK_THROW(ir); 2548 } 2549 2550 const RinexObsVec& rov = it->second; 2551 for (size_t i=0; i<rov.size(); i++) 2552 { 2553 if (rov[i].equalIndex(obsID)) 2554 return i; 2555 } 2556 2557 InvalidRequest ir(obsID.asString(version) + " is not stored in system " + sys + "."); 2558 GPSTK_THROW(ir); 2559 return 0; 2560 } 2561 2562 compare(const Rinex3ObsHeader & right,std::vector<std::string> & diffs,const std::vector<std::string> & inclExclList,bool incl)2563 bool Rinex3ObsHeader::compare(const Rinex3ObsHeader& right, 2564 std::vector<std::string>& diffs, 2565 const std::vector<std::string>& inclExclList, 2566 bool incl) 2567 { 2568 // map header token to comparison result 2569 std::map<std::string,bool> lineMap; 2570 std::map<std::string,bool>::const_iterator lmi; 2571 // Put the comments in a sorted set, we don't really care 2572 // about the ordering. 2573 std::set<std::string> 2574 lcomments(commentList.begin(), commentList.end()), 2575 rcomments(right.commentList.begin(), right.commentList.end()); 2576 std::set<RinexObsID> 2577 lobs(obsTypeList.begin(), obsTypeList.end()), 2578 robs(right.obsTypeList.begin(), right.obsTypeList.end()); 2579 // Compare everything first... 2580 // deliberately ignoring valid flags 2581 2582 // only comparing first character of file type because that's 2583 // all that matters according to RINEX 2584 lineMap[hsVersion] = 2585 ((version == right.version) && 2586 (fileType[0] == right.fileType[0]) && 2587 (fileSysSat.system == right.fileSysSat.system)); 2588 lineMap[hsRunBy] = 2589 ((fileProgram == right.fileProgram) && 2590 (fileAgency == right.fileAgency) && 2591 (date == right.date)); 2592 lineMap[hsComment] = (lcomments == rcomments); 2593 lineMap[hsMarkerName] = (markerName == right.markerName); 2594 lineMap[hsMarkerNumber] = (markerNumber == right.markerNumber); 2595 lineMap[hsMarkerType] = (markerType == right.markerType); 2596 lineMap[hsObserver] = 2597 ((observer == right.observer) && 2598 (agency == right.agency)); 2599 lineMap[hsReceiver] = 2600 ((recNo == right.recNo) && 2601 (recType == right.recType) && 2602 (recVers == right.recVers)); 2603 lineMap[hsAntennaType] = 2604 ((antNo == right.antNo) && 2605 (antType == right.antType)); 2606 lineMap[hsAntennaPosition] = 2607 (antennaPosition == right.antennaPosition); 2608 lineMap[hsAntennaDeltaHEN] = 2609 (antennaDeltaHEN == right.antennaDeltaHEN); 2610 lineMap[hsAntennaDeltaXYZ] = 2611 (antennaDeltaXYZ == right.antennaDeltaXYZ); 2612 lineMap[hsAntennaPhaseCtr] = 2613 (antennaPhaseCtr == right.antennaPhaseCtr); 2614 lineMap[hsAntennaBsightXYZ] = 2615 (antennaBsightXYZ == right.antennaBsightXYZ); 2616 lineMap[hsAntennaZeroDirAzi] = 2617 (antennaZeroDirAzi == right.antennaZeroDirAzi); 2618 lineMap[hsAntennaZeroDirXYZ] = 2619 (antennaZeroDirXYZ == right.antennaZeroDirXYZ); 2620 lineMap[hsCenterOfMass] = (centerOfMass == right.centerOfMass); 2621 lineMap[hsNumObs] = (lobs == robs); 2622 lineMap[hsSystemNumObs] = true; 2623 lineMap[hsWaveFact] = 2624 (memcmp(wavelengthFactor, right.wavelengthFactor, 2625 sizeof(wavelengthFactor)) == 0); 2626 lineMap[hsSigStrengthUnit] = 2627 (sigStrengthUnit == right.sigStrengthUnit); 2628 lineMap[hsInterval] = (interval == right.interval); 2629 lineMap[hsFirstTime] = (firstObs == right.firstObs); 2630 lineMap[hsLastTime] = (lastObs == right.lastObs); 2631 lineMap[hsReceiverOffset] = (receiverOffset == right.receiverOffset); 2632 lineMap[hsSystemDCBSapplied] = true; 2633 lineMap[hsSystemPCVSapplied] = true; 2634 lineMap[hsSystemScaleFac] = true; 2635 lineMap[hsSystemPhaseShift] = true; 2636 lineMap[hsGlonassSlotFreqNo] = true; 2637 lineMap[hsGlonassCodPhsBias] = true; 2638 lineMap[hsLeapSeconds] = (leapSeconds == right.leapSeconds); 2639 lineMap[hsNumSats] = (numSVs == right.numSVs); 2640 lineMap[hsPrnObs] = true; 2641 // ...then filter by inclExclList later 2642 if (incl) 2643 { 2644 std::map<std::string,bool> oldLineMap(lineMap); 2645 std::map<std::string,bool>::const_iterator olmi; 2646 lineMap.clear(); 2647 for (unsigned i = 0; i < inclExclList.size(); i++) 2648 { 2649 if ((olmi = oldLineMap.find(inclExclList[i])) != oldLineMap.end()) 2650 { 2651 lineMap[olmi->first] = olmi->second; 2652 } 2653 } 2654 } 2655 else 2656 { 2657 // exclude, remove items in inclExclList 2658 for (unsigned i = 0; i < inclExclList.size(); i++) 2659 { 2660 lineMap.erase(inclExclList[i]); 2661 } 2662 } 2663 // check the equality of the final remaining set of header lines 2664 bool rv = true; 2665 for (lmi = lineMap.begin(); lmi != lineMap.end(); lmi++) 2666 { 2667 if (!lmi->second) 2668 { 2669 diffs.push_back(lmi->first); 2670 rv = false; 2671 } 2672 } 2673 return rv; 2674 } // bool Rinex3ObsHeader::compare 2675 2676 2677 Rinex3ObsHeader::Fields Rinex3ObsHeader::Fields :: getRequired(double version)2678 getRequired(double version) 2679 { 2680 if (version < 3.00) return allValid2; 2681 else if(version < 3.01) return allValid30; 2682 else if(version < 3.02) return allValid301; 2683 else if(version < 3.03) return allValid302; 2684 else if(version < 3.04) return allValid303; 2685 else if(version < 3.05) return allValid303; 2686 return Fields(); 2687 } 2688 2689 2690 Rinex3ObsHeader::Fields Rinex3ObsHeader::Fields :: operator &(const Rinex3ObsHeader::Fields & rhs) const2691 operator&(const Rinex3ObsHeader::Fields& rhs) const 2692 { 2693 FieldSet results; 2694 set_intersection(fieldsSet.begin(), fieldsSet.end(), 2695 rhs.fieldsSet.begin(), rhs.fieldsSet.end(), 2696 inserter(results, results.begin())); 2697 return Rinex3ObsHeader::Fields(results); 2698 } 2699 2700 Rinex3ObsHeader::Fields Rinex3ObsHeader::Fields :: operator |(const Rinex3ObsHeader::Fields & rhs) const2701 operator|(const Rinex3ObsHeader::Fields& rhs) const 2702 { 2703 FieldSet results; 2704 set_union(fieldsSet.begin(), fieldsSet.end(), 2705 rhs.fieldsSet.begin(), rhs.fieldsSet.end(), 2706 inserter(results, results.begin())); 2707 return Rinex3ObsHeader::Fields(results); 2708 } 2709 2710 2711 Rinex3ObsHeader::Field Rinex3ObsHeader::Fields :: operator &(Rinex3ObsHeader::Field rhs) const2712 operator&(Rinex3ObsHeader::Field rhs) const 2713 { 2714 if (fieldsSet.count(rhs)) 2715 return rhs; 2716 return validInvalid; 2717 } 2718 2719 2720 bool Rinex3ObsHeader::Fields :: isValid(const Rinex3ObsHeader::Fields & present) const2721 isValid(const Rinex3ObsHeader::Fields& present) const 2722 { 2723 FieldSet results; 2724 set_difference(fieldsSet.begin(), fieldsSet.end(), 2725 present.fieldsSet.begin(), present.fieldsSet.end(), 2726 inserter(results, results.begin())); 2727 return results.empty(); 2728 } 2729 2730 2731 void Rinex3ObsHeader::Fields :: describeMissing(const Rinex3ObsHeader::Fields & valid,Exception & exc)2732 describeMissing(const Rinex3ObsHeader::Fields& valid, 2733 Exception& exc) 2734 { 2735 for (const auto& f : fieldsSet) 2736 { 2737 if (!valid.isSet(f)) 2738 { 2739 exc.addText("Missing required header field: " + asString(f)); 2740 } 2741 } 2742 } 2743 2744 2745 std::string Rinex3ObsHeader :: asString(Rinex3ObsHeader::Field b)2746 asString(Rinex3ObsHeader::Field b) 2747 { 2748 switch (b) 2749 { 2750 case validVersion: return hsVersion; 2751 case validRunBy: return hsRunBy; 2752 case validComment: return hsComment; 2753 case validMarkerName: return hsMarkerName; 2754 case validMarkerNumber: return hsMarkerNumber; 2755 case validMarkerType: return hsMarkerType; 2756 case validObserver: return hsObserver; 2757 case validReceiver: return hsReceiver; 2758 case validAntennaType: return hsAntennaType; 2759 case validAntennaPosition: return hsAntennaPosition; 2760 case validAntennaDeltaHEN: return hsAntennaDeltaHEN; 2761 case validAntennaDeltaXYZ: return hsAntennaDeltaXYZ; 2762 case validAntennaPhaseCtr: return hsAntennaPhaseCtr; 2763 case validAntennaBsightXYZ: return hsAntennaBsightXYZ; 2764 case validAntennaZeroDirAzi: return hsAntennaZeroDirAzi; 2765 case validAntennaZeroDirXYZ: return hsAntennaZeroDirXYZ; 2766 case validCenterOfMass: return hsCenterOfMass; 2767 case validNumObs: return hsNumObs; 2768 case validSystemNumObs: return hsSystemNumObs; 2769 case validWaveFact: return hsWaveFact; 2770 case validSigStrengthUnit: return hsSigStrengthUnit; 2771 case validInterval: return hsInterval; 2772 case validFirstTime: return hsFirstTime; 2773 case validLastTime: return hsLastTime; 2774 case validReceiverOffset: return hsReceiverOffset; 2775 case validSystemDCBSapplied: return hsSystemDCBSapplied; 2776 case validSystemPCVSapplied: return hsSystemPCVSapplied; 2777 case validSystemScaleFac: return hsSystemScaleFac; 2778 case validSystemPhaseShift: return hsSystemPhaseShift; 2779 case validGlonassSlotFreqNo: return hsGlonassSlotFreqNo; 2780 case validGlonassCodPhsBias: return hsGlonassCodPhsBias; 2781 case validLeapSeconds: return hsLeapSeconds; 2782 case validNumSats: return hsNumSats; 2783 case validPrnObs: return hsPrnObs; 2784 default: return "???"; 2785 } 2786 } 2787 2788 2789 Rinex3ObsHeader::Field Rinex3ObsHeader :: asField(const std::string & s)2790 asField(const std::string& s) 2791 { 2792 if (s == hsVersion) return validVersion; 2793 if (s == hsRunBy) return validRunBy; 2794 if (s == hsComment) return validComment; 2795 if (s == hsMarkerName) return validMarkerName; 2796 if (s == hsMarkerNumber) return validMarkerNumber; 2797 if (s == hsMarkerType) return validMarkerType; 2798 if (s == hsObserver) return validObserver; 2799 if (s == hsReceiver) return validReceiver; 2800 if (s == hsAntennaType) return validAntennaType; 2801 if (s == hsAntennaPosition) return validAntennaPosition; 2802 if (s == hsAntennaDeltaHEN) return validAntennaDeltaHEN; 2803 if (s == hsAntennaDeltaXYZ) return validAntennaDeltaXYZ; 2804 if (s == hsAntennaPhaseCtr) return validAntennaPhaseCtr; 2805 if (s == hsAntennaBsightXYZ) return validAntennaBsightXYZ; 2806 if (s == hsAntennaZeroDirAzi) return validAntennaZeroDirAzi; 2807 if (s == hsAntennaZeroDirXYZ) return validAntennaZeroDirXYZ; 2808 if (s == hsCenterOfMass) return validCenterOfMass; 2809 if (s == hsNumObs) return validNumObs; 2810 if (s == hsSystemNumObs) return validSystemNumObs; 2811 if (s == hsWaveFact) return validWaveFact; 2812 if (s == hsSigStrengthUnit) return validSigStrengthUnit; 2813 if (s == hsInterval) return validInterval; 2814 if (s == hsFirstTime) return validFirstTime; 2815 if (s == hsLastTime) return validLastTime; 2816 if (s == hsReceiverOffset) return validReceiverOffset; 2817 if (s == hsSystemDCBSapplied) return validSystemDCBSapplied; 2818 if (s == hsSystemPCVSapplied) return validSystemPCVSapplied; 2819 if (s == hsSystemScaleFac) return validSystemScaleFac; 2820 if (s == hsSystemPhaseShift) return validSystemPhaseShift; 2821 if (s == hsGlonassSlotFreqNo) return validGlonassSlotFreqNo; 2822 if (s == hsGlonassCodPhsBias) return validGlonassCodPhsBias; 2823 if (s == hsLeapSeconds) return validLeapSeconds; 2824 if (s == hsNumSats) return validNumSats; 2825 if (s == hsPrnObs) return validPrnObs; 2826 return validInvalid; 2827 } 2828 operator <<(std::ostream & s,const Rinex3ObsHeader::Fields & v)2829 std::ostream& operator<<(std::ostream& s, const Rinex3ObsHeader::Fields& v) 2830 { 2831 Rinex3ObsHeader::FieldSet::const_iterator i; 2832 for (i = v.fieldsSet.begin(); i != v.fieldsSet.end(); i++) 2833 { 2834 if (i != v.fieldsSet.begin()) 2835 s << ","; 2836 s << *i; 2837 } 2838 return s; 2839 } 2840 2841 2842 void Rinex3ObsHeader :: remapObsTypes(RinexObsMap & remapped,map<string,unsigned> & obsCount) const2843 remapObsTypes(RinexObsMap& remapped, map<string,unsigned>& obsCount) 2844 const 2845 { 2846 for (const auto& mapIter : mapObsTypes) 2847 { 2848 // I and X pseudo-observables are special cases, and can 2849 // only be listed once (or once per band in the case of 2850 // the ionospheric delay) in any sane manner. 2851 bool addedChannel = false; 2852 std::set<CarrierBand> addedIono; 2853 for(size_t i = 0; i < mapIter.second.size(); i++) 2854 { 2855 if (mapIter.second[i].type == ObservationType::Iono) 2856 { 2857 if (addedIono.count(mapIter.second[i].band) > 0) 2858 continue; // only write this pseudo-obs once 2859 addedIono.insert(mapIter.second[i].band); 2860 } 2861 else if (mapIter.second[i].type == ObservationType::Channel) 2862 { 2863 if (addedChannel) 2864 continue; // only write this pseudo-obs once 2865 addedChannel = true; 2866 } 2867 remapped[mapIter.first].push_back(mapIter.second[i]); 2868 obsCount[mapIter.first]++; 2869 } 2870 } 2871 } 2872 2873 } // namespace gpstk 2874