// Part of BNC, a utility for retrieving decoding and
// converting GNSS data streams from NTRIP broadcasters.
//
// Copyright (C) 2007
// German Federal Agency for Cartography and Geodesy (BKG)
// http://bkg.bund.de
// Czech Technical University Prague, Department of Geodesy
// http://www.fsv.cvut.cz
//
// Email: euref-ip@bkg.bund.de
//
// This program is free software; you can redistribute it and/or
// modify it under the terms of the GNU General Public License
// as published by the Free Software Foundation, version 2.
//
// This program 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 General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.

/* -------------------------------------------------------------------------
 * BKG NTRIP Client
 * -------------------------------------------------------------------------
 *
 * Class:      bncBiasSnx
 *
 * Purpose:    Biases from Bias SINEX File
 *
 * Author:     A. Stuerze
 *
 * Created:    15-Mar-2022
 *
 * Changes:
 *
 * -----------------------------------------------------------------------*/

#include "bncbiassnx.h"
#include "ssrRefSig.h"

#include <iostream>
#include <newmatio.h>

using namespace std;

// Constructor
////////////////////////////////////////////////////////////////////////////
bncBiasSnx::bncBiasSnx(bool allGnss) {
  bncSettings settings;
  if (allGnss) {
    for (char sys : QString("GRECJSI").toLatin1()) {
      _useGnss[sys] = true;
    }
    return;
  }
  _useGnss['G'] = (Qt::CheckState(settings.value("cmbGps").toInt()) == Qt::Checked) ? true : false;
  _useGnss['R'] = (Qt::CheckState(settings.value("cmbGlo").toInt()) == Qt::Checked) ? true : false;
  _useGnss['E'] = (Qt::CheckState(settings.value("cmbGal").toInt()) == Qt::Checked) ? true : false;
  _useGnss['C'] = (Qt::CheckState(settings.value("cmbBds").toInt()) == Qt::Checked) ? true : false;
  _useGnss['J'] = (Qt::CheckState(settings.value("cmbQzss").toInt()) == Qt::Checked) ? true : false;
  _useGnss['S'] = (Qt::CheckState(settings.value("cmbSbas").toInt()) == Qt::Checked) ? true : false;
  _useGnss['I'] = (Qt::CheckState(settings.value("cmbNavic").toInt()) == Qt::Checked) ? true : false;
}

// Destructor
////////////////////////////////////////////////////////////////////////////
bncBiasSnx::~bncBiasSnx() {
  clear();
}

// Clear Satellite Map
////////////////////////////////////////////////////////////////////////////
void bncBiasSnx::clear() {
  QMapIterator<QString,snxSatCodeBias*>it(_snxSatCodeBiasMap);
  while (it.hasNext()) {
    it.next();
    delete it.value();
  }
  _snxSatCodeBiasMap.clear();
}

// Print
////////////////////////////////////////////////////////////////////////////
void bncBiasSnx::print() const {
  QMapIterator<QString,snxSatCodeBias*> it(_snxSatCodeBiasMap);
  cout << "============= bncBiasSnx::print() ==========" << endl;
  while (it.hasNext()) {
    it.next();
    cout << it.key().toStdString() << " ";
    snxSatCodeBias* satCbMap = it.value();
    cout << satCbMap->_startTime.datestr() << "_" << satCbMap->_startTime.timestr() << " ";
    QMapIterator <QString, double> itB(satCbMap->_biasMap);
    cout << satCbMap->_biasMap.size() << " ";
    while (itB.hasNext()) {
      itB.next();
      QString osb  = itB.key();
      double  value = itB.value();
      cout << osb.toStdString() << " " << value << " ";
    }
    cout << endl;
  }
}


// Read Bias SINEX File
////////////////////////////////////////////////////////////////////////////
t_irc bncBiasSnx::readFile(const QString& fileName) {
  clear();
  QFile inFile(fileName);
  if (!inFile.open(QIODevice::ReadOnly | QIODevice::Text)) {
    return failure;
  }

  QTextStream in(&inFile);

  while ( !in.atEnd() ) {
    QString biasType, biasAtrib, svn, prn, station, obs1, obs2, startTimeStr, endTimeStr, unit;
    double value;
    bncTime startTime;

    QString line = in.readLine();//cout << line.toStdString() << endl;
    if      (line.indexOf("+BIAS/SOLUTION") == 14) {
      clear();
    }
    else if (line.indexOf("BIAS_MODE") == 1) {
       if      (line.contains("ABSOLUTE")) {
         _type = ABS;
       }
       else if (line.contains("RELATIVE")) {
         _type = REL;
       }
    }
    else if (line.indexOf("OSB") == 1 || line.indexOf("DSB") == 1) {
      QTextStream inLine(&line, QIODevice::ReadOnly);
      inLine >> biasType;
      if      (biasType == "OSB") {
        inLine >> svn >> prn >> obs1 >> startTimeStr >> endTimeStr >> unit >> value;
        biasAtrib = obs1;
      }
      else if (biasType == "DSB") {
        inLine >> svn;
        if (svn.size() == 1) { // station data
          continue;
        } else {
          inLine >> prn >> obs1 >> obs2 >> startTimeStr >> endTimeStr >> unit >> value;
          biasAtrib = QString("%1-%2").arg(obs1).arg(obs2);
        }
      }
      // bncTime from snx time
      // ---------------------
      if (bncTimeFromSinex(startTimeStr, startTime) != success) {
        continue;
      }
      // change of sign
      // --------------
      // post-processing convention is to subtract biases
      // real-time convention is to add biases.
      value *= -1.0;

      // conversion of Code Biases
      // -------------------------
      if (unit.contains("ns")) {
        // from nano seconds into meters
        value /= 1.e9;
        value *= t_CST::c;

      }

      // delete existing entries if there are newer once available
      // ---------------------------------------------------------
      if (_snxSatCodeBiasMap.contains(prn)) {
        if (_snxSatCodeBiasMap[prn]->_startTime < startTime) {
          delete _snxSatCodeBiasMap[prn];
          _snxSatCodeBiasMap[prn] = new snxSatCodeBias(startTime);
        }
      }
      else {
        _snxSatCodeBiasMap[prn] = new snxSatCodeBias(startTime);
      }
      // save biases
      // -----------
      QMap<QString, double>& biasMap = _snxSatCodeBiasMap[prn]->_biasMap;
      if (!biasMap.contains(biasAtrib)) {
        biasMap.insert(biasAtrib, value);
      }
    }
  }
  inFile.close();

  if (_type == ABS) {
    getDsbFromOsb();
  }
  cleanDsb();
  //print();
  return success;
}

// bncTime from SINEX Time String
////////////////////////////////////////////////////////////////////////////
t_irc bncBiasSnx::bncTimeFromSinex(QString snxStr, bncTime& time) {
  bool ok = false;
  QStringList hlp = snxStr.split(":");// yyyy:doy:dssec

  int year = hlp[0].toInt(&ok, 10);
  if (!ok) {
    return failure;
  }
  time.set(year, 1, 1, 0, 0, 0.0);

  int doy  = hlp[1].toInt(&ok, 10);
  if (!ok) {
    return failure;
  }
  time += (doy * 24.0 * 3600.0);

  int daysec = hlp[2].toInt(&ok, 10);
  if (!ok) {
    return failure;
  }

  time += daysec;

  return success;
}


// get DSBs from OSBs
////////////////////////////////////////////////////////////////////////////
t_irc bncBiasSnx::getDsbFromOsb(void) {
  QMapIterator<QString,snxSatCodeBias*>itSat(_snxSatCodeBiasMap);
  while (itSat.hasNext()) {
    itSat.next();
    char sys = itSat.key().at(0).toLatin1();
    QMap<QString, double>& cb = itSat.value()->_biasMap;
    switch (sys) {
      case 'G':
        if  (cb.contains("C1W") && cb.contains("C2W")) {cb.insert("C1W-C2W", cb["C1W"] - cb["C2W"]);}
        if  (cb.contains("C1C") && cb.contains("C1W")) {cb.insert("C1C-C1W", cb["C1C"] - cb["C1W"]);}
        if  (cb.contains("C2C") && cb.contains("C2W")) {cb.insert("C2C-C2W", cb["C2C"] - cb["C2W"]);}
        if  (cb.contains("C2W") && cb.contains("C2S")) {cb.insert("C2W-C2S", cb["C2W"] - cb["C2S"]);}
        if  (cb.contains("C2W") && cb.contains("C2L")) {cb.insert("C2W-C2L", cb["C2W"] - cb["C2L"]);}
        if  (cb.contains("C2W") && cb.contains("C2X")) {cb.insert("C2W-C2X", cb["C2W"] - cb["C2X"]);}
        if  (cb.contains("C1C") && cb.contains("C2W")) {cb.insert("C1C-C2W", cb["C1C"] - cb["C2W"]);}
        if  (cb.contains("C1C") && cb.contains("C5Q")) {cb.insert("C1C-C5Q", cb["C1C"] - cb["C5Q"]);}
        if  (cb.contains("C1C") && cb.contains("C5X")) {cb.insert("C1C-C5X", cb["C1C"] - cb["C5X"]);}
        break;
      case 'R':
        if  (cb.contains("C1P") && cb.contains("C2P")) {cb.insert("C1P-C2P", cb["C1P"] - cb["C2P"]);}
        if  (cb.contains("C1C") && cb.contains("C1P")) {cb.insert("C1C-C1P", cb["C1C"] - cb["C1P"]);}
        if  (cb.contains("C2C") && cb.contains("C2P")) {cb.insert("C2C-C2P", cb["C2C"] - cb["C2P"]);}
        if  (cb.contains("C1C") && cb.contains("C2C")) {cb.insert("C1C-C2C", cb["C1C"] - cb["C2C"]);}
        if  (cb.contains("C1C") && cb.contains("C2P")) {cb.insert("C1C-C2P", cb["C1C"] - cb["C2P"]);}
        break;
      case 'E':
        if  (cb.contains("C1C") && cb.contains("C5Q")) {cb.insert("C1C-C5Q", cb["C1C"] - cb["C5Q"]);}
        if  (cb.contains("C1C") && cb.contains("C6C")) {cb.insert("C1C-C6C", cb["C1C"] - cb["C6C"]);}
        if  (cb.contains("C1C") && cb.contains("C7Q")) {cb.insert("C1C-C7Q", cb["C1C"] - cb["C7Q"]);}
        if  (cb.contains("C1C") && cb.contains("C8Q")) {cb.insert("C1C-C8Q", cb["C1C"] - cb["C8Q"]);}
        break;
      case 'C':
        if  (cb.contains("C2I") && cb.contains("C6I")) {cb.insert("C2I-C6I", cb["C2I"] - cb["C6I"]);}
        if  (cb.contains("C2I") && cb.contains("C7I")) {cb.insert("C2I-C7I", cb["C2I"] - cb["C7I"]);}
        if  (cb.contains("C1X") && cb.contains("C6I")) {cb.insert("C1X-C6I", cb["C1X"] - cb["C6I"]);}
        if  (cb.contains("C1P") && cb.contains("C6I")) {cb.insert("C1P-C6I", cb["C1P"] - cb["C6I"]);}
        if  (cb.contains("C1D") && cb.contains("C6I")) {cb.insert("C1D-C6I", cb["C1D"] - cb["C6I"]);}
        if  (cb.contains("C1X") && cb.contains("C7Z")) {cb.insert("C1X-C7Z", cb["C1X"] - cb["C7Z"]);}
        if  (cb.contains("C1X") && cb.contains("C8X")) {cb.insert("C1X-C8X", cb["C1X"] - cb["C8X"]);}
        if  (cb.contains("C1X") && cb.contains("C5X")) {cb.insert("C1X-C5X", cb["C1X"] - cb["C5X"]);}
        if  (cb.contains("C1P") && cb.contains("C5P")) {cb.insert("C1P-C5P", cb["C1P"] - cb["C5P"]);}
        if  (cb.contains("C1D") && cb.contains("C5D")) {cb.insert("C1D-C5D", cb["C1D"] - cb["C5D"]);}
        break;
      case 'J':
        if  (cb.contains("C1C") && cb.contains("C1X")) {cb.insert("C1C-C1X", cb["C1C"] - cb["C1X"]);}
        if  (cb.contains("C1C") && cb.contains("C2L")) {cb.insert("C1C-C2L", cb["C1C"] - cb["C2L"]);}
        if  (cb.contains("C1C") && cb.contains("C5X")) {cb.insert("C1C-C5X", cb["C1C"] - cb["C5X"]);}
        if  (cb.contains("C1C") && cb.contains("C5Q")) {cb.insert("C1C-C5Q", cb["C1C"] - cb["C5Q"]);}
        if  (cb.contains("C1X") && cb.contains("C2X")) {cb.insert("C1X-C2X", cb["C1X"] - cb["C2X"]);}
        if  (cb.contains("C1X") && cb.contains("C5X")) {cb.insert("C1X-C5X", cb["C1X"] - cb["C5X"]);}
        break;
      default:
        break;
    }
  }
  return success;
}

//
////////////////////////////////////////////////////////////////////////////
t_irc bncBiasSnx::cleanDsb() {

  for (auto itS = _snxSatCodeBiasMap.begin(); itS != _snxSatCodeBiasMap.end();) {
    char sys = itS.key().at(0).toLatin1();
    if (!_useGnss[sys]) {
      delete itS.value();
      itS = _snxSatCodeBiasMap.erase(itS);
    }
    else {
      QMap<QString, double>& cb = itS.value()->_biasMap;
      for (auto itB = cb.begin(); itB != cb.end();){
        if  (itB.key().size() == 3) {
          itB = cb.erase(itB);
        }
        else {
          ++itB;
        }
      }
      ++itS;
    }
  }
  return success;
}

// DCBs (OSB_A - OSB_B) used to derive the code biases from the default
// reference signals (see t_ssrRefSig), in the order of evaluation
////////////////////////////////////////////////////////////////////////////
static const QVector<QPair<QString, QString> >& defaultDcbs(char sys) {
  static QMap<char, QVector<QPair<QString, QString> > > dcbs;
  static QMutex mutex;
  QMutexLocker locker(&mutex);
  if (dcbs.isEmpty()) {
    auto add = [](QVector<QPair<QString, QString> >& vv, const char* list) {
      for (const QString& dcb : QString(list).split(' ', Qt::SkipEmptyParts)) {
        vv.append(qMakePair(dcb.section('-', 0, 0), dcb.section('-', 1, 1)));
      }
    };
    add(dcbs['G'], "1C-1W 1C-5Q 1C-5X 1C-2W 2C-2W 2W-2S 2W-2L 2W-2X");
    add(dcbs['R'], "1C-1P 2C-2P");
    add(dcbs['E'], "1C-6C 1C-7Q 1C-8Q");
    add(dcbs['C'], "2I-7I 1X-6I 1X-7Z 1X-8X 1X-5X 1P-6I 1P-5P 1D-6I 1D-5D");
    add(dcbs['J'], "1C-1X 1X-2X 1C-5X");
  }
  static const QVector<QPair<QString, QString> > none;
  return dcbs.contains(sys) ? dcbs[sys] : none;
}

// Satellite code biases (OSB) from the DCBs, referred to the reference
// signals of the clocks: aIF1 * OSB_ref1 + aIF2 * OSB_ref2 = 0
////////////////////////////////////////////////////////////////////////////
void bncBiasSnx::determineSsrSatCodeBiases(QString prn, double aIF1, double aIF2,
                                           t_satCodeBias& satCodeBias) {
  if (!_snxSatCodeBiasMap.contains(prn)) {
    return;
  }
  char sys = prn.at(0).toLatin1();
  t_ssrRefSig::t_pair pair;
  if (!t_ssrRefSig::get(sys, pair)) {
    return;
  }
  const QString ref1 = QString::fromStdString(pair.sig1());
  const QString ref2 = QString::fromStdString(pair.sig2());
  const QMap<QString, double>& cb = _snxSatCodeBiasMap[prn]->_biasMap;
  auto dcb = [&cb](const QString& sig1, const QString& sig2, double& value) {
    QString key = "C" + sig1 + "-C" + sig2;
    if (cb.contains(key)) {
      value = cb[key];
      return true;
    }
    return false;
  };

  QMap<QString, double> osb;
  QStringList           order;

  // Default reference signals: the DCB of the reference signals and the
  // DCBs of the IGS conventions, evaluated in one pass as before
  if (t_ssrRefSig::isDefault(sys)) {
    const QVector<QPair<QString, QString> >& dcbs = defaultDcbs(sys);
    double dcbRef;
    if (dcbs.isEmpty() || !dcb(ref1, ref2, dcbRef)) {
      return;
    }
    osb[ref1] =  aIF2 * dcbRef;
    osb[ref2] = -aIF1 * dcbRef;
    if (sys == 'E') {   // order of the output as before
      order << ref2 << ref1;
    }
    else {
      order << ref1 << ref2;
    }
    for (const QPair<QString, QString>& dd : dcbs) {
      double value;
      if (!dcb(dd.first, dd.second, value)) {
        continue;
      }
      if (osb.contains(dd.first) && !osb.contains(dd.second)) {
        osb[dd.second] = osb[dd.first] - value;
        order << dd.second;
      }
      else if (osb.contains(dd.second) && !osb.contains(dd.first)) {
        osb[dd.first] = value + osb[dd.second];
        order << dd.first;
      }
    }
  }

  // Other reference signals: all DCBs of the satellite, the DCB of the
  // reference signals may also be derived via other signals
  else {
    QMap<QString, double> rel;   // OSB - OSB_ref1
    QStringList           relOrder;
    rel[ref1] = 0.0;
    relOrder << ref1;
    bool changed = true;
    while (changed) {
      changed = false;
      QMapIterator<QString, double> it(cb);
      while (it.hasNext()) {
        it.next();
        QStringList hlp = it.key().split('-');
        if (hlp.size() != 2 || hlp[0].size() != 3 || hlp[1].size() != 3) {
          continue;
        }
        QString sigA = hlp[0].mid(1);
        QString sigB = hlp[1].mid(1);
        if (rel.contains(sigA) && !rel.contains(sigB)) {
          rel[sigB] = rel[sigA] - it.value();
          relOrder << sigB;
          changed = true;
        }
        else if (rel.contains(sigB) && !rel.contains(sigA)) {
          rel[sigA] = rel[sigB] + it.value();
          relOrder << sigA;
          changed = true;
        }
      }
    }
    if (!rel.contains(ref2)) {
      return;
    }
    double osbRef1 = -aIF2 * rel[ref2];
    order << ref1 << ref2;
    for (const QString& sig : relOrder) {
      osb[sig] = osbRef1 + rel[sig];
      if (sig != ref1 && sig != ref2) {
        order << sig;
      }
    }
  }

  for (const QString& sig : order) {
    setSsrSatCodeBias(sig.toStdString(), osb[sig], satCodeBias);
  }
}


void bncBiasSnx::setSsrSatCodeBias(std::string rnxType2ch, double value,
                                    t_satCodeBias& satCodeBias) {
  t_frqCodeBias frqCodeBias;

  frqCodeBias._rnxType2ch = rnxType2ch;
  frqCodeBias._value = value;

  satCodeBias._bias.push_back(frqCodeBias);
  return;
}