// PCABits.cc
// Contact person: Freija Descamps <fbdescamps@lbl.gov>
// See PCABits.hh for more details
//———————————————————————//

// GEANT stuff
#include <G4UImessenger.hh>

// RAT stuff
#include <RAT/PCABits.hh>
#include <RAT/DB.hh>
#include <RAT/Log.hh>
#include <RAT/BitManip.hh>
#include <RAT/DS/Run.hh>

// C++ stuff
#include <iostream>
#include <sstream>
#include <cmath>
#include <vector>

using namespace std;
namespace RAT {

  PCABits::PCABits()
  {
    // Empty for now
  }

  PCABits::~PCABits()
  {
    // Empty for now
  }

  void PCABits::BeginOfRun( DS::Run& )
  {
    // Get PCA user-defined input concerning how to
    // squeeze the data in an integer
    DBLinkPtr fPCAbank;
    fPCAbank = DB::Get()->GetLink("PCA_GENERATION");
    fQbits=fPCAbank->GetI("qbits");
    fTbits=fPCAbank->GetI("tbits");
    fIbits=fPCAbank->GetI("ibits");
    fRbits=fPCAbank->GetI("rbits");
    fQfact=fPCAbank->GetI("qfact");
    fTfact=fPCAbank->GetI("tfact");
    fIfact=fPCAbank->GetI("ifact");
    fRfact=fPCAbank->GetI("rfact");
  }

  int PCABits::PCAPacker(Direction direction, DataType datatype,
                         int &packeddata,
                         double &point1, double &point2)
  {
    if( direction==1 )
      {  // pack!
        int packedint=0;
        if( datatype==1 )
          { // (q,t) interpolation point
            double thisq=static_cast<int>(point1*fQfact);
            double thist=static_cast<int>(point2*fTfact);
            // Set the status
            if( point1<-9000 && point2<-9000 )
              {
                packedint= BitManip::SetBit(packedint, 0);
                return -1;
              }
            else
              {
                // Set the signs
                if( thisq>0 )
                  {
                    packedint=BitManip::SetBit(packedint, 1);
                  }
                else{
                  thisq=thisq*-1;
                }
                if( thist>0 )
                  {
                    packedint=BitManip::SetBit(packedint, fQbits+2);
                  }
                else
                  {
                    thist=thist*-1;
                  }
                // Fill the data
                packedint=BitManip::SetBits(packedint,2,thisq);
                packedint=BitManip::SetBits(packedint,fQbits+3,thist);}
          }
        if( datatype==2 )
          { // fitted intercept and gradient
            double thisi=static_cast<int>(point1*fIfact);
            double thisr=static_cast<int>(point2*fRfact);
            if( point1<-9000 && point2<-9000 )
              {
                packedint= BitManip::SetBit(packedint,0);
                return -1;
              }
            else
              {
                //Check the sign
                if( thisi>0 )
                  {
                    packedint=BitManip::SetBit(packedint,1);
                  }
                else
                  {
                    thisi=thisi*-1;
                  }
                if( thisr>0 )
                  {
                    packedint=BitManip::SetBit(packedint,fIbits+2);
                  }
                else
                  {
                    thisr=thisr*-1;
                  }
                // Fill the data
                packedint=BitManip::SetBits(packedint,2,thisi);
                packedint=BitManip::SetBits(packedint,fIbits+3,thisr);
              }
          }
        packeddata=packedint;
      }
    if( direction==2 )
      { // unpack!
        if( BitManip::GetBits(packeddata,0,1)==1 )
          {
            // this means that there were no datapoints
            return -1;
          }
        else
          {
            if( datatype==1 )
              { // (q,t) interpolation point
                int signofq=0; // 0=negative
                int signoft=0; // 0=negative
                int thisqint=0;
                int thistint=0;
                double thisq=0.;
                double thist=0.;
                thisqint=BitManip::GetBits(packeddata,2,fQbits);
                thistint=BitManip::GetBits(packeddata,fQbits+3,fTbits);
                signofq=BitManip::GetBits(packeddata,1,1);
                signoft=BitManip::GetBits(packeddata,fQbits+2,1);

                if( signofq==0 )
                  {
                    thisq=thisqint*-1.0;
                  }
                else
                  {
                    thisq=thisqint*1.0;
                  }
                if( signoft==0 )
                  {
                    thist=thistint*-1.0;
                  }
                else
                  {
                    thist=thistint*1.0;
                  }
                // Fill the data
                point1=thisq/(1.0*fQfact);
                point2=thist/(1.0*fTfact);
              }
            if( datatype==2 )
              { // fitted intercept and gradient
                int signofr=0; // 0=negative
                int signofi=0; // 0=negative
                int thisrint=0;
                int thisiint=0;
                double thisr=0.;
                double thisi=0.;
                thisiint=BitManip::GetBits(packeddata,2,fIbits);
                thisrint=BitManip::GetBits(packeddata,fIbits+3,fRbits);
                signofi=BitManip::GetBits(packeddata,1,1);
                signofr=BitManip::GetBits(packeddata,fIbits+2,1);
                if( signofi==0 )
                  {
                    thisi=thisiint*-1.0;
                  }
                else
                  {
                    thisi=thisiint*1.0;
                  }
                if( signofr==0 )
                  {
                    thisr=thisrint*-1.0;
                  }
                else
                  {
                    thisr=thisrint*1.0;
                  }
                // Fill the data
                point1=thisi/(1.0*fIfact);
                point2=thisr/(1.0*fRfact);
              }
          }
      }
    return 0;
  }


  void PCABits::InitialiseStatusWords()
  {
    fOverallMask.resize(32);
    fPMTTWMask.resize(32);
    fPMTGFMask.resize(32);
    // OR of bits 1-31
    fOverallMask[0]  = "status PCA";
    // OR of bits 8-19
    fOverallMask[1]  = "status TW";
    // OR of bits 20-31
    fOverallMask[2]  = "status GF";
    fOverallMask[3]  = "too many channels offline";
    fOverallMask[4]  = "too many online channels that have 0 occupancy";
    fOverallMask[5]  = "too many channels that have low occupancy";
    fOverallMask[6]  = "overall high occupancy";
    fOverallMask[7]  = "did not reach end of run";

    fOverallMask[8]  = "N channels with very high RMS > max";
    fOverallMask[9]  = "N channels with too high RMS > max";
    fOverallMask[10] = "N channels with no fit to high Q tail > max";
    fOverallMask[11] = "N channels with Tstep warning > max";
    fOverallMask[12] = "N channels with Flate warning > max";
    fOverallMask[13] = "N channels with Fout warning > max";
    fOverallMask[14] = "N channels with missing points warning > max";
    fOverallMask[15] = "N channels with gradient too big > max";
    fOverallMask[16] = "N channels with gradient too small > max";
    fOverallMask[17] = "Bootstrap: correlation warning";
    fOverallMask[18] = "Bootstrap failed";
    fOverallMask[19] = "TW spare";

    fOverallMask[20] = "too many (#PMTs with d(QHS TH) too large)";
    fOverallMask[21] = "too many (#PMTs with d(QHS PK) too large)";
    fOverallMask[22] = "too many (#PMTs with d(QHS HP) too large)";
    fOverallMask[23] = "too many (#PMTs with d(QHL TH) too large)";
    fOverallMask[24] = "too many (#PMTs with d(QHL PK) too large)";
    fOverallMask[25] = "too many (#PMTs with d(QHL HP) too large)";
    fOverallMask[26] = "too many (#PMTs with QHS TH too high)";
    fOverallMask[27] = "too many (#PMTs with QHL TH too high)";
    fOverallMask[28] = "too many (#PMTs with peak finder)";
    fOverallMask[29] = "too many (#PMTs with QHS TH too low)";
    fOverallMask[30] = "too many (#PMTs with QHL TH too low)";
    fOverallMask[31] = "GF spare";

    // OR of bits 1-31
    fPMTTWMask[0]  = "status";
    // same as GF
    fPMTTWMask[1]  = "channel offline";
    // no hits on this PMT, same as GF
    fPMTTWMask[2]  = "zero occupancy";
    //Less than min hits on this PMT, could be different than GF
    fPMTTWMask[3]  = "low occupancy";
    fPMTTWMask[4]  = "spare";
    fPMTTWMask[5]  = "spare";
    fPMTTWMask[6]  = "spare";
    fPMTTWMask[7]  = "RMS is very high";
    fPMTTWMask[8]  = "RMS is too high";
    fPMTTWMask[9]  = "high Q tail was not fitted";
    // Possible bad ADC charge conversion. This will affect a complete card.
    fPMTTWMask[10] = "Tstep warning";
    // The fraction of late light is too large
    fPMTTWMask[11] = "Flate warning";
    // The fraction of hits <QHSmin and < Tpeak-10 is too large
    fPMTTWMask[12] = "Fout warning";
    // The gradient of the fit to the high Q tail is too large
    fPMTTWMask[13] = "Gradient too big";
    // The gradient of the fit to the high Q tail is too small
    fPMTTWMask[14] = "Gradient too small";
    fPMTTWMask[15] = "spare";
    fPMTTWMask[16] = "spare";
    fPMTTWMask[17] = "spare";
    fPMTTWMask[18] = "spare";
    fPMTTWMask[19] = "spare";
    fPMTTWMask[20] = "spare";
    fPMTTWMask[21] = "spare";
    fPMTTWMask[22] = "spare";
    fPMTTWMask[23] = "spare";
    fPMTTWMask[24] = "spare";
    fPMTTWMask[25] = "spare";
    fPMTTWMask[26] = "spare";
    fPMTTWMask[27] = "spare";
    fPMTTWMask[28] = "spare";
    fPMTTWMask[29] = "spare";
    fPMTTWMask[30] = "spare";
    fPMTTWMask[31] = "N Bins failed interpolation point calculation > Max";

    // OR of bits 1-31
    fPMTGFMask[0]  = "status";
    // FAIL tube is marked as off in DQXX
    fPMTGFMask[1]  = "channel offline";
    // FAIL tube is marked as on in DQXX but did not see any hits
    fPMTGFMask[2]  = "zero occupancy QHS";
    // FAIL tube saw less than min hits
    fPMTGFMask[3]  = "low occupancy";
    fPMTGFMask[4]  = "< 100 hits in 100-bin window for QHS";
    fPMTGFMask[5]  = "< 100 hits in 100-bin window for QHL";
    fPMTGFMask[6]  = "spare";
    fPMTGFMask[7]  = "spare";
    fPMTGFMask[8]  = "spare";
    fPMTGFMask[9]  = "spare";

    fPMTGFMask[10] = "QHS TH diff too large"; // WARN
    fPMTGFMask[11] = "QHS PK diff too large"; // WARN
    fPMTGFMask[12] = "QHS HP diff too large"; // WARN
    fPMTGFMask[13] = "QHL TH diff too large"; // WARN
    fPMTGFMask[14] = "QHL PK diff too large"; // WARN
    fPMTGFMask[15] = "QHL HP diff too large"; // WARN
    fPMTGFMask[16] = "spare";
    fPMTGFMask[17] = "spare";
    fPMTGFMask[18] = "spare";
    fPMTGFMask[19] = "spare";

    fPMTGFMask[20] = "QHS TH too high"; // WARN, bad pedestal?
    fPMTGFMask[21] = "QHL TH too high"; // WARN, bad pedestal?
    // WARN, peakfinder was called and used to determine position of second peak
    fPMTGFMask[22] = "Peakfinder called and used";
    // WARN, Possible noise peak fitted for threshold
    fPMTGFMask[23] = "QHS TH too low";
    fPMTGFMask[24] = "QHL TH too low";
    fPMTGFMask[25] = "spare";
    fPMTGFMask[26] = "spare";
    fPMTGFMask[27] = "spare";
    fPMTGFMask[28] = "spare";
    fPMTGFMask[29] = "spare";
    fPMTGFMask[30] = "spare";
    fPMTGFMask[31] = "spare";

  }
} // namespace RAT