#ifndef FILTER_H
#define FILTER_H

#include <vector>
#include <newmat.h>
#include "pppInclude.h"
#include "pppParlist.h"
#include "bnctime.h"
#include "t_prn.h"

namespace BNC_PPP {

class t_pppParlist;
class t_pppObsPool;
class t_pppSatObs;

class t_pppFilter {
 public:
  t_pppFilter(t_pppObsPool* obsPool);
  ~t_pppFilter();

  t_irc processEpoch(int num);

  void datumTransformation();

  const ColumnVector&    x() const {return _xFlt;}
  const SymmetricMatrix& Q() const {return _QFlt;}

  int    numSat() const {return _numSat;}
  double HDOP() const {return _dop.H;}
  double HDOV() const {return _dop.V;}
  double PDOP() const {return _dop.P;}
  double GDOP() const {return _dop.G;}
  double trp() const {
    const std::vector<t_pppParam*>& par = _parlist->params();
    for (unsigned ii = 0; ii < par.size(); ++ii) {
      if (par[ii]->type() == t_pppParam::trp) {
        return x()[ii];
      }
    }
    return 0.0;
  };
  double trpStdev() const {
    const std::vector<t_pppParam*>& par = _parlist->params();
    for (unsigned ii = 0; ii < par.size(); ++ii) {
      if (par[ii]->type() == t_pppParam::trp) {
        return sqrt(Q()[ii][ii]);
      }
    }
    return 0.0;
  };

 private:
  class t_slip {
   public:
    t_slip() {
      _slip            = false;
      _obsSlipCounter  = -1;
      _biasJumpCounter = -1;
    }
    bool _slip;
    int  _obsSlipCounter;
    int  _biasJumpCounter;
  };

  class t_dop {
   public:
    t_dop() {reset();}
    void reset() {H = V = P = T = G = 0.0;}
    double H;
    double V;
    double P;
    double T;
    double G;
  };

  class t_datumTrafo{
  public:
    t_datumTrafo () {initIndices();}
    void initIndices() {_firstRow = 1; _lastRow = 0;}
    void setFirstSystem(bool firstSys) {_firstSys = firstSys;}
    bool firstSystem() {return _firstSys;}
    void updateIndices(int maxObs) {
      if (_firstSys) {
        initIndices();
      }
      else {
        _firstRow += maxObs;
      }
      _lastRow += maxObs;
    };
    void initAA(int maxObs, int numPar) {
      _AA1.ReSize(maxObs, numPar); _AA1 = 0.0;
      _AA2.ReSize(maxObs, numPar); _AA2 = 0.0;
    }
    void prepareAA(Matrix& AA, int _numEpoProcessing, int nPar) {
      Matrix& Prep = _AA2;
      if (_numEpoProcessing == 1) {
        Prep = _AA1;
      }
      Prep.SubMatrix(_firstRow, _lastRow, 1, nPar) = AA;
    }
    Matrix varCov(const SymmetricMatrix& QFlt) {
      Matrix D21 = (_AA2.t() * _AA2).i() * _AA2.t() * _AA1;
      return D21 * QFlt * D21.t();
    }
    int     _firstRow;
    int     _lastRow;
    Matrix  _AA1;
    Matrix  _AA2;
    bool    _firstSys;
  };

  t_irc processSystem(const std::vector<t_lc::type>& LCs,
                      const std::vector<t_pppSatObs*>& obsVector,
                      const t_prn& refPrn,
                      bool pseudoObsIonoAvailable,
                      bool preProcessing);

  t_irc detectCycleSlips(const std::vector<t_lc::type>& LCs,
                         const std::vector<t_pppSatObs*>& obsVector,
                         const t_prn& refPrn,
                         bool preProcessing);

  t_irc resetAmb(t_prn prn, const std::vector<t_pppSatObs*>& obsVector,
                 SymmetricMatrix* QSav = 0, ColumnVector* xSav = 0);

  void cmpDOP(const std::vector<t_pppSatObs*>& obsVector);

  void predictCovCrdPart(const SymmetricMatrix& QFltOld);

  bncTime         _epoTime;
  t_pppParlist*   _parlist;
  t_pppObsPool*   _obsPool;
  t_datumTrafo*   _datumTrafo;
  SymmetricMatrix _QFlt;
  ColumnVector    _xFlt;
  ColumnVector    _x0;
  t_slip          _slips[t_prn::MAXPRN+1];
  int             _numSat;
  int             _numEpoProcessing;
  int             _numAllUsedLCs;
  t_dop           _dop;
  bncTime         _firstEpoTime;
  bncTime         _lastEpoTimeOK;
  t_prn           _refPrn;
};

}

#endif