timetable.cpp

gpstk1.5的有关内容 对于刚刚接触gps有一定的帮助 很有用的啊

   if(!ttifs) {
      cerr << "Failed to open input time table file " << CI.TimeTableFile << endl;
      return -3;
   }

   //REF site site sat week use_first use_last data_start data_end
   do {
      getline(ttifs,line);
      StringUtils::stripTrailing(line,'\r');
      if(ttifs.eof() || !ttifs.good()) break;

      if(line.size() <= 0) continue;                              // skip blank lines
      if(line[0] == '#') continue;                                // skip comments
      if(StringUtils::numWords(line) < 9) continue; // TD msg?    // skip bad lines
      if(StringUtils::words(line,0,1) != string("REF")) continue; // only REF lines
      TTSegment ts;
      ts.site1 = StringUtils::words(line,1,1);
      ts.site2 = StringUtils::words(line,2,1);
      ts.sat.fromString(StringUtils::words(line,3,1));

      week = StringUtils::asInt(StringUtils::words(line,4,1));
      sow = StringUtils::asInt(StringUtils::words(line,5,1));
      tt.setGPSfullweek(week,sow);           // TD handle week rollover
      ts.first = int(0.5+(tt-FirstEpoch)/CI.DataInterval);

      sow = StringUtils::asInt(StringUtils::words(line,6,1));
      tt.setGPSfullweek(week,sow);
      ts.last = int(0.5+(tt-FirstEpoch)/CI.DataInterval);

      sow = StringUtils::asInt(StringUtils::words(line,7,1));
      tt.setGPSfullweek(week,sow);
      ts.start = int(0.5+(tt-FirstEpoch)/CI.DataInterval);

      sow = StringUtils::asInt(StringUtils::words(line,8,1));
      tt.setGPSfullweek(week,sow);
      ts.end = int(0.5+(tt-FirstEpoch)/CI.DataInterval);

      //ts.minelev = ts.maxelev = 0.0;
      ts.length = ts.end - ts.start + 1;
      ts.findElev();
      TimeTable.push_back(ts);

   } while(1);
   // close the output timetable file
   ttifs.close();

   oflog << "Read time table from file " << CI.TimeTableFile << endl;

   return 0;
}
catch(Exception& e) { GPSTK_RETHROW(e); }
catch(exception& e) { Exception E("std except: "+string(e.what())); GPSTK_THROW(E); }
catch(...) { Exception e("Unknown exception"); GPSTK_THROW(e); }
}

//------------------------------------------------------------------------------------
int ComputeBaselineTimeTable(const string& bl)
{
try {
   int i,j;
   map<SDid,SDData>::const_iterator it;
   list<TTSegment> SegList;

   SDmap.clear();
   i = TTComputeSingleDifferences(bl,40.0);
   if(i) return i;

   // now compute the timetable based on SDmap
   for(it=SDmap.begin(); it != SDmap.end(); it++) {
      TTSegment ts;
      ts.site1 = it->first.site1;
      ts.site2 = it->first.site2;
      ts.sat = it->first.sat;
      ts.start = it->second.count[0];
      ts.minelev = ts.maxelev = 0.0;
      for(i=0; i<it->second.count.size()-1; i++) {
         j = it->second.count.at(i+1) - it->second.count.at(i);
         if(j > 1) {
            TTSegment tts;
            tts = ts;
            tts.end = it->second.count.at(i);
            tts.length = tts.end - tts.start + 1;
            tts.findElev();
            SegList.push_back(tts);
            ts.start = it->second.count.at(i+1);
         }
      }
      ts.end = it->second.count[it->second.count.size()-1];
      ts.length = ts.end - ts.start + 1;
      ts.findElev();
      SegList.push_back(ts);
   }
   
   if(SegList.size() == 0) return -2;

   // figure out the time table from the list of segments
   list<TTSegment> TTable;
   i = TimeTableAlgorithm(SegList, TTable);
   if(i) return i;

   // add this timetable to the master timetable.
   list<TTSegment>::iterator ttit;
   for(ttit=TTable.begin(); ttit != TTable.end(); ttit++)
      TimeTable.push_back(*ttit);

   return 0;
}
catch(Exception& e) { GPSTK_RETHROW(e); }
catch(exception& e) { Exception E("std except: "+string(e.what())); GPSTK_THROW(E); }
catch(...) { Exception e("Unknown exception"); GPSTK_THROW(e); }
}

//------------------------------------------------------------------------------------
int TTComputeSingleDifferences(const string& bl, const double ElevLimit)
{
try {
   int i,j,k;
   double elevi,elevj;
   GSatID sat;
   map<GSatID,RawData>::const_iterator it,jt;
   map<SDid,SDData>::const_iterator kt;
   ofstream rawofs;
   format f62(6,2),f133(13,3);

   // loop over buffered raw data of sats common to both
   string est=StringUtils::word(bl,0,'-');
   string fix=StringUtils::word(bl,1,'-');

   for(it=Stations[est].RawDataBuffers.begin();
       it != Stations[est].RawDataBuffers.end(); it++) {

      // consider each satellite
      sat = it->first;
      if(CI.Verbose) oflog << "Single difference " << est << " " << fix << " " << sat;

      // is sat also found at fixed site?
      jt = Stations[fix].RawDataBuffers.find(sat);
      if(jt == Stations[fix].RawDataBuffers.end()) {
         if(CI.Verbose) oflog << " not found on both sites" << endl;
         continue;
      }

      if(CI.Verbose) oflog << " (raw buffers size: " << it->second.count.size()
         << " " << jt->second.count.size() << ")";

      // compute continuous segments of SD data
      // sdd.count is the intersection of the two count vectors
      SDid sdid(fix,est,sat);
      SDData sdd;
      sdd.elevmin = 100.0;
      sdd.elevmax = -1.0;
      i = j = 0;
      do {
         if(it->second.count[i] == jt->second.count[j]) {
            elevi = RotatedAntennaElevation(it->second.elev[i],it->second.az[i]);
            elevj = RotatedAntennaElevation(jt->second.elev[j],jt->second.az[j]);
            if(elevi >= ElevLimit && elevj >= ElevLimit) {
               sdd.count.push_back(it->second.count[i]);
               if(elevi < sdd.elevmin) sdd.elevmin = elevi;
               if(elevi > sdd.elevmax) sdd.elevmax = elevi;
            }
            i++;
            j++;
         }
         else if(it->second.count[i] < jt->second.count[j]) {
            i++;
         }
         else {
            j++;
         }
      } while( i < it->second.count.size() &&
               j < jt->second.count.size() );

      // TD ?
      if(sdd.count.size() < 10) {
         if(CI.Verbose) oflog << " size is too small ("
            << sdd.count.size() << ")" << endl;
         continue;
      }

      // save it in the map
      SDmap[sdid] = sdd;

      if(CI.Verbose) oflog << endl;

   }  // end loop over raw buffered data common to both sites

   // write out a summary of single differences
   oflog << "Single differences summary :" << endl;
   for(k=1,kt=SDmap.begin(); kt != SDmap.end(); kt++,k++) {
      oflog << " " << setw(2) << k << " " << kt->first
         << " " << setw(5) << kt->second.count.size()
         << " " << setw(5) << kt->second.count.at(0)
         << " - " << setw(5) << kt->second.count.at(kt->second.count.size()-1);
         // elevation
      oflog << " elev: "
         << fixed << setw(4) << setprecision(1) << kt->second.elevmin
         << " - " << setw(4) << setprecision(1) << kt->second.elevmax;
         // gaps - (count : number of pts)
      for(i=0; i<kt->second.count.size()-1; i++) {
         j = kt->second.count.at(i+1) - kt->second.count.at(i);
         if(j > 1)
            oflog << " (" << kt->second.count.at(i)+1 << ":" << j-1 << ")";
      }
      oflog << endl;
   }

   if(SDmap.size() == 0) {
      oflog << "Returning error code -1 from TTComputeSingleDifferences()" << endl;
      return -1;
   }

   return 0;
}
catch(Exception& e) { GPSTK_RETHROW(e); }
catch(exception& e) { Exception E("std except: "+string(e.what())); GPSTK_THROW(E); }
catch(...) { Exception e("Unknown exception"); GPSTK_THROW(e); }
}

//------------------------------------------------------------------------------------
// preprocess the segment list - get first and last counts
int TimeTableAlgorithm(list<TTSegment>& TTS, list<TTSegment>& TTab)
{
try {
   bool keep;
   int i,beg,end;
   int begcount,endcount;
   list<TTSegment>::iterator ttit,ttjt;
   typedef pair<int,int> seg;
   list<seg> Segs;
   list<seg>::iterator lit,ljt;

   // 0 sort list in order of decreasing metric()
   // 0.5 delete very small segments
   // 0.6 run through the list, finding the smallest begin, and largest end, counts;
   // call these the begin and end counts.
   // 1. run through the sorted list, starting at largest metric(), and pick out
   // first segments which have the begin and/or end counts; add these to TTab
   // and erase from TTS.
   // 2. run through TTS again, starting at the largest metric(); if a segment
   // includes some counts that have not been covered before, than add this
   // to TTab and erase from TTS. Quit when either the entire range from
   // begcount to endcount is covered, a given minimum metric() is reached,
   // or when the end of TTS is reached.
   // 3. if gaps remain in the coverage, these are real gaps in the data and
   // the estimation will have to reset.
   // 4. sort TTab in increasing order. Run through TTab looking for
   // segments which can be removed without generating gaps; remove these.

   // 0. sort in reverse order (largest metric() first)
   // [ list has its own sort sort(TTS.rbegin(),TTS.rend()); ]
   TTS.sort(decreasingMetricSort);

   // 0.5 delete very small segments and output the sorted list
   // 0.6 find begcount and endcount (after deletion)
   begcount = endcount = -1;
   oflog << "Here is the sorted list of segments:" << endl;
   for(i=1,ttit=TTS.begin(); ttit != TTS.end(); i++) {
      oflog << " " << setw(4) << i << *ttit;
      if(ttit->length < 10) {
         oflog << " -- delete this segment: too small";
         ttit = TTS.erase(ttit); // ttit now points to next seg
      }
      else {
         if(begcount < 0 || ttit->start < begcount)
            begcount = ttit->start;
         if(endcount < 0 || ttit->end > endcount)
            endcount = ttit->end;
         ttit++;
      }
      oflog << endl;
   }
   oflog << "End the sorted list; limits : " << begcount << " - " << endcount << endl;

   for(ttit=TTS.begin(); ttit != TTS.end(); ttit++) {
      if(ttit->length < 10) {
      }
   }
   
   // 1.find the begin point
   for(ttit=TTS.begin(); ttit != TTS.end(); ttit++) {
      if(ttit->start == begcount) {