gdc.hpp

Go to the documentation of this file.

//==============================================================================
//
//  This file is part of GNSSTk, the ARL:UT GNSS Toolkit.
//
//  The GNSSTk is free software; you can redistribute it and/or modify
//  it under the terms of the GNU Lesser General Public License as published
//  by the Free Software Foundation; either version 3.0 of the License, or
//  any later version.
//
//  The GNSSTk is distributed in the hope that it will be useful,
//  but WITHOUT ANY WARRANTY; without even the implied warranty of
//  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
//  GNU Lesser General Public License for more details.
//
//  You should have received a copy of the GNU Lesser General Public
//  License along with GNSSTk; if not, write to the Free Software Foundation,
//  Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110, USA
//
//  This software was developed by Applied Research Laboratories at the
//  University of Texas at Austin.
//  Copyright 2004-2022, The Board of Regents of The University of Texas System
//
//==============================================================================
 
//==============================================================================
//
//  This software was developed by Applied Research Laboratories at the
//  University of Texas at Austin, under contract to an agency or agencies
//  within the U.S. Department of Defense. The U.S. Government retains all
//  rights to use, duplicate, distribute, disclose, or release this software.
//
//  Pursuant to DoD Directive 523024
//
//  DISTRIBUTION STATEMENT A: This software has been approved for public
//                            release, distribution is unlimited.
//
//==============================================================================
 
#ifndef GNSSTK_DISCONTINUITY_CORRECTOR_INCLUDE
#define GNSSTK_DISCONTINUITY_CORRECTOR_INCLUDE
 
#include <fstream>
#include <iostream>
#include <map>
#include <sstream>
#include <string>
#include <vector>
 
#include "gnsstk_export.h"
#include "Exception.hpp"
#include "Epoch.hpp"
#include "Exception.hpp"
#include "FDiffFilter.hpp"
#include "FirstDiffFilter.hpp"
#include "RinexSatID.hpp"
#include "SatPass.hpp"
#include "WindowFilter.hpp"
 
namespace gnsstk
{
   //---------------------------------------------------------------------------
   class Arc
   {
      friend class gdc;
 
   private:
      static std::map<unsigned, std::string> create_mark_string_map()
      {
         std::map<unsigned, std::string> m;
         m[BEG]    = std::string("Begin");
         m[WLSLIP] = std::string("WLslip");
         m[GFSLIP] = std::string("GFslip");
         m[WLFIX]  = std::string("WLfix");
         m[GFFIX]  = std::string("GFfix");
         m[WLMARK] = std::string("WLmark");
         m[GFMARK] = std::string("GFmark");
         m[REJ]    = std::string("Reject");
         return m;
      }
 
   public:
      // NB cannot use markStr[...] since its const -- use iterators
      static const std::map<unsigned, std::string> markStr;
 
      static const unsigned
         BEG; 
      static const unsigned WLSLIP; 
      static const unsigned GFSLIP; 
      static const unsigned WLFIX;  
      static const unsigned GFFIX;  
      static const unsigned GFMARK; 
      static const unsigned WLMARK; 
      static const unsigned
         REJ; 
 
      // -------------------------------------------------------------------
      class Arcinfo
      {
       public:
         Arcinfo() : Nslip(0), step(0.0), n(0), ave(0.0), sig(0.0) {}
 
         // member data
         double
            step; 
         double sigma; 
         // TD half-cycles?
         int Nslip; 
         int n; 
         double ave; 
         double sig; 
      };             // end class Arcinfo
 
      // -------------------------------------------------------------------
      // member data
 
      // these are used more for "break" information
      unsigned
         mark;   
      int index; 
      unsigned npts;  
      unsigned ngood; 
 
      Arcinfo WLinfo, GFinfo; // others? array? IF?
 
      std::string message; 
 
      // -------------------------------------------------------------------
      // member functions
 
      Arc() { init(); }
 
      Arc(int ind, unsigned n, unsigned ng, unsigned m)
      {
         init();
         mark  = m;
         index = ind;
         npts  = n;
         ngood = ng;
      }
 
      void init()
      {
         mark  = 0;
         index = -1;
         npts = ngood = 0;
         WLinfo       = Arcinfo();
         GFinfo       = Arcinfo();
         message      = std::string();
      }
 
      Arc(const FilterHit<double>& seg)
      {
         init();
         index = seg.index;
         if (seg.type == FilterHit<double>::BOD)
         {
            mark = BEG; // NB BEG != FilterHit::BOD
         }
         else
         {
            mark = 0;
         }
         npts  = seg.npts;
         ngood = seg.ngood;
      }
 
      Arc(const Arc& right)
      {
         init();
         *this = right;
      }
 
      Arc& operator=(const Arc& right)
      {
         index   = right.index;
         npts    = right.npts;
         ngood   = right.ngood;
         mark    = right.mark;
         message = right.message;
         WLinfo  = right.WLinfo;
         GFinfo  = right.GFinfo;
         return *this;
      }
 
      static std::string markString(const unsigned mark)
      {
         std::map<unsigned, std::string>::const_iterator it;
         it = markStr.find(mark);
         if (it != markStr.end())
         {
            return it->second;
         }
         else
         {
            return std::string("Unknown");
         }
      }
 
      std::string markString() const
      {
         std::string msg;
         std::map<unsigned, std::string>::const_iterator it;
         //  NB don't use markStr[it->first]; b/c its const
         for (it = markStr.begin(); it != markStr.end(); ++it)
            if (mark & it->first)
            {
               msg += std::string(msg.empty() ? "" : "/") + it->second;
            }
         return msg;
      }
 
      std::string brkString() const
      {
         std::ostringstream os;
         os << markString() << " good/tot=" << ngood << "/" << npts << "="
            << std::fixed << std::setprecision(2) << 100. * ngood / double(npts)
            << "%";
         return os.str();
      }
 
      std::string asString(const int prec = 3) const
      {
         std::ostringstream os;
         os << std::fixed << std::setprecision(prec);
         os << brkString();
 
         if (markString() != "BEG" && markString() != "REJ")
         {
            os << " NWL=" << WLinfo.Nslip << " dWL=" << WLinfo.step << "wl";
         }
         else
         {
            os << " WL";
         }
         if (WLinfo.n > 0)
         {
            os << " " << WLinfo.n << " " << WLinfo.ave << " +- " << WLinfo.sig;
         }
         else
         {
            os << " NoAn";
         }
 
         if (markString() != "BEG" && markString() != "REJ")
         {
            os << " NGF=" << GFinfo.Nslip << " dGF=" << GFinfo.step << "wl";
         }
         else
         {
            os << " GF";
         }
         if (GFinfo.n > 0)
         {
            os << " " << GFinfo.n << " " << GFinfo.ave << " +- " << GFinfo.sig;
         }
         else
         {
            os << " NoAn";
         }
 
         return os.str();
      }
 
   }; // end class Arc
 
   //---------------------------------------------------------------------------------
   //---------------------------------------------------------------------------------
   class gdc
   {
   public:
      gdc() { init(); }
 
      ~gdc() {}
 
      bool setParameter(std::string cmd);
 
      bool setParameter(std::string label, double value);
 
      double getParameter(std::string label) { return CFG[label]; }
 
      void DisplayParameterUsage(std::ostream& os,
                                 std::string tag = std::string(),
                                 bool advanced   = false);
 
      double cfg_func(std::string a)
      {
         if (CFGdesc[a] == std::string())
         {
            Exception e("cfg(UNKNOWN LABEL) : " + a);
            GNSSTK_THROW(e);
         }
         return CFG[a];
      }
 
      int getUniqueNumber() { return unique; }
 
      void ForceUniqueNumber(int n) { unique = n; }
 
      //---------------------------------------------------------------------------
      int DiscontinuityCorrector(SatPass& SP, std::string& retMsg,
                                 std::vector<std::string>& cmds,
                                 int GLOn = -99);
 
      //---------------------------------------------------------------------------
      int DiscontinuityCorrector(
         const RinexSatID& sat, const double& nominalDT, const Epoch& beginTime,
         std::vector<double> dataL1, std::vector<double> dataL2,
         std::vector<double> dataP1, std::vector<double> dataP2,
         std::vector<double> dt, std::vector<int> flags, std::string& retMsg,
         std::vector<std::string>& cmds, int GLOn = -99,
         std::string outfmt = std::string("%4F %10.3g"));
 
   private:
      static std::vector<unsigned> create_vector_SLIP()
      {
         std::vector<unsigned> v;
         v.push_back(Arc::WLSLIP);
         v.push_back(Arc::GFSLIP);
         return v;
      }
      static std::vector<unsigned> create_vector_FIX()
      {
         std::vector<unsigned> v;
         v.push_back(Arc::WLFIX);
         v.push_back(Arc::GFFIX);
         return v;
      }
      static std::vector<std::string> create_vector_LAB()
      {
         std::vector<std::string> v;
         v.push_back("WL");
         v.push_back("GF");
         return v;
      }
 
   protected:
      GNSSTK_EXPORT
      static const std::string GDCVersion;
 
      // NB flags[] is either good (0) or not (!0);
      // flags have NOTHING to do with either SatPass or Arc::marks
      // Arcs handles BEG, SLIP/FIX, and GAP(new BEG); flags handle bad data
      // w/in Arc NB BAD != SatPass::BAD but ONLY SatPass::BAD on input => BAD
      static const unsigned OK;  
      static const unsigned BAD; 
      static const unsigned
         WLOUTLIER; 
      static const unsigned
         GFOUTLIER; 
      static const unsigned WLSHORT; 
      static const unsigned GFSHORT; 
      static const unsigned
         ISOLATED; 
      // others?
 
      static const unsigned WL; 
      static const unsigned GF; 
      static const std::vector<unsigned> SLIP, FIX; 
      static const std::vector<std::string> LAB; 
      static std::vector<double> wl;             
 
      //---------------------------------------------------------------------------
      // configuration
      std::map<std::string, double> CFG;
 
      std::map<std::string, std::string> CFGdesc;
 
      int CFGindex;
      std::map<int, std::string> CFGlist;
 
      //---------------------------------------------------------------------------
      int unique;
      std::string
         tag; 
      std::string SPSstr; 
 
      RinexSatID sat; 
      bool isGLO;     
      double dt;      
      Epoch beginT;   
 
      std::string outfmt; 
      int oswidth;        
      int osprec;         
 
      int GLOchan; 
      // just to keep numerical range down
      double WLbias; 
      double GFbias; 
      // approximate; for output of editing commands only
      long long
         N1bias; 
      long long
         N2bias; 
 
      double wl1;   
      double wl2;   
      double alpha; 
      double beta;  
      double wlWL;  
      double wlGF;  
      double wlNL;  
 
      std::vector<double> dataWL, dataGF;
 
      std::vector<double> xdata;
 
      // NB flags must be int, not unsigned, for StatsFilter processing
      std::vector<int> flags;
 
      std::map<int, Arc> Arcs;
 
      // member functions
      //---------------------------------------------------------------------------
 
      void init();
 
      int ProcessOneCombo(const unsigned which);
 
      int GrossProcessing(const unsigned which);
 
      int FineProcessing(const unsigned which);
 
      int filterFirstDiff(const unsigned which, const std::string label,
                          double limit, std::vector<FilterHit<double>>& hits);
 
      int filterWindow(const unsigned which, const std::string label,
                       double limit, std::vector<FilterHit<double>>& hits);
 
      int mergeFilterResultsIntoArcs(std::vector<FilterHit<double>>& hits,
                                     const unsigned which);
 
      // the next three are called by mergeFilterResultsIntoArcs()
 
      void flagBadData(const FilterHit<double>& hit, const unsigned flagvalue);
 
      bool addArc(const int index, const unsigned mark);
 
      void fixUpArcs();
 
      void recomputeArcs();
 
      inline unsigned computeNgood(Arc& arc)
      {
         arc.ngood = 0;
         for (unsigned int i = arc.index; i < arc.index + arc.npts; ++i)
         {
            if (flags[i] == OK)
            {
               arc.ngood++;
            }
         }
         return arc.ngood;
      }
 
      void getArcStats(const unsigned which)
      {
         try
         {
            for (std::map<int, Arc>::iterator ait = Arcs.begin();
                 ait != Arcs.end(); ++ait)
               getArcStats(ait, which);
         }
         catch (Exception& e)
         {
            GNSSTK_RETHROW(e);
         }
      }
 
      void getArcStats(std::map<int, Arc>::iterator& ait, const unsigned which);
 
      void findLargeGaps();
 
      int fixSlips(const unsigned which);
 
      std::map<int, int> findGaps(const Arc& arc);
 
      void dumpData(std::ostream& os, const std::string msg = std::string());
 
      Epoch xtime(const int& i)
      {
         if (i >= dataWL.size())
         {
            GNSSTK_THROW(Exception("Index out of range in xtime"));
         }
         Epoch tt(beginT);
         tt += xdata[i];
         return tt;
      }
 
      void findArc(const unsigned int ind, std::map<int, Arc>::iterator& ait)
      {
         if (int(ind) < ait->second.index)
         {
            GNSSTK_THROW(Exception("index before given Arc"));
         }
 
         while (ait != Arcs.end() &&
                ind >= ait->second.index + ait->second.npts)
            ++ait;
 
         if (ait == Arcs.end()) // index not found
         {
            GNSSTK_THROW(Exception("index after end of Arcs"));
         }
      }
 
      void DumpArcs(const std::string& tag, const std::string& label,
                    int prec = -1);
 
      void DumpHits(const std::vector<FilterHit<double>>& filterResults,
                    const std::string& tag, const std::string& label,
                    int prec = -1);
 
      std::string returnMessage(int prec = -1, int wid = -1);
 
      void applyFixesToSatPass(SatPass& SP);
 
      void generateCmds(std::vector<std::string>& cmds);
 
      int FinalCheck();
 
   }; // end class gdc
 
} // end namespace gnsstk
 
//------------------------------------------------------------------------------------
#endif // GNSSTK_DISCONTINUITY_CORRECTOR_INCLUDE