// -----------------------------------------------------------------------------
// ----- concat_calib.C                                                    -----
// -----                                                                   -----
// ----- created by C. Simon on 2020-02-06                                 -----
// -----                                                                   -----
// -----------------------------------------------------------------------------

#include <iostream>
#include <vector>

#include "TROOT.h"
#include "TSystem.h"
#include "TFile.h"
#include "TH1.h"
#include "TMath.h"


void concat_calib()
{
  const Int_t kiNInputRuns = 3;//4;
  const Bool_t kbMergeFineTime = kTRUE;
  const Bool_t kbJointEdgeDetection = kTRUE;
  const Int_t kiNTDCBoards = 20;
  const Int_t kiNTDCChannels = 64;

  const char* kcAnalysisDirName = "/lustre/nyx/cbm/users/csimon/analysis/sps_feb15";
  const char* kcCalibFileName = "raw.calib.root";

//  const char* kcInputRunID[kiNInputRuns] = {"04Mar1923", "04Mar1937", "04Mar1945", "04Mar2002"};
  const char* kcInputRunID[kiNInputRuns] = {"05Mar0041", "05Mar0219", "05Mar0404"};

  const Int_t kiEventWindow = 75;
  const Int_t kiTriggerSet = 1;
  const Int_t kiTriggerMUL = 3;

  const Int_t kiReferenceRun = 0;//2;


  std::vector<std::vector<TH1I*>> tRefFineTimeHistosOUT(kiNInputRuns, std::vector<TH1I*>(kiNTDCBoards, NULL));
  std::vector<std::vector<std::vector<TH1I*>>> tFineTimeHistosOUT(kiNInputRuns, std::vector<std::vector<TH1I*>>(kiNTDCBoards, std::vector<TH1I*>(kiNTDCChannels, NULL)));
  std::vector<Int_t> iReferenceBoardOUT[kiNInputRuns];
  std::vector<UInt_t> uGlobalFineTimeLimitsOUT[kiNInputRuns];
  std::vector<Bool_t> bDisableBoardOUT[kiNInputRuns];
  std::vector<Double_t> dFinalOffsetOUT[kiNInputRuns];

  for(Int_t iInputRun = 0; iInputRun < kiNInputRuns; iInputRun++)
  {
    std::vector<TH1I*> tRefFineTimeHistosIN(kiNTDCBoards, NULL);
    std::vector<std::vector<TH1I*>> tFineTimeHistosIN(kiNTDCBoards, std::vector<TH1I*>(kiNTDCChannels, NULL));
    std::vector<Int_t>* iReferenceBoardIN = NULL;
    std::vector<UInt_t>* uGlobalFineTimeLimitsIN = NULL;
    std::vector<Bool_t>* bDisableBoardIN = NULL;
    std::vector<Double_t>* dFinalOffsetIN = NULL;

    TString tFileName = Form("%s/CbmTofAna_%s/UnpackDigis_EVT%dns_TRG%d_MUL%d/data/%s", kcAnalysisDirName, kcInputRunID[iInputRun],
                             kiEventWindow, kiTriggerSet, kiTriggerMUL, kcCalibFileName);

    TFile* tInputCalibFile = new TFile(tFileName, "READ");

    if(tInputCalibFile->IsZombie())
    {
      tInputCalibFile->Close();
      delete tInputCalibFile;

      cout<<Form("Could not open input calibration file of run %s.", kcInputRunID[iInputRun])<<endl;
      return;
    }

    if(!tInputCalibFile->cd("TofTrbCalib"))
    {
      cout<<Form("Could not find directory 'TofTrbCalib' in input calibration file of run %s.", kcInputRunID[iInputRun])<<endl;
      return;
    }

    TDirectory* tCalibParent = gDirectory;

    gDirectory->GetObject("iReferenceBoard", iReferenceBoardIN);
    if(!iReferenceBoardIN)
    {
      cout<<Form("Could not retrieve 'iReferenceBoard' from input calibration file of run %s.", kcInputRunID[iInputRun])<<endl;
      return;
    }
    else
    {
      iReferenceBoardOUT[iInputRun] = *iReferenceBoardIN;
    }

    gDirectory->GetObject("uGlobalFineTimeLimits", uGlobalFineTimeLimitsIN);
    if(!uGlobalFineTimeLimitsIN)
    {
      cout<<Form("Could not retrieve 'uGlobalFineTimeLimits' from input calibration file of run %s.", kcInputRunID[iInputRun])<<endl;
      return;
    }
    else
    {
      uGlobalFineTimeLimitsOUT[iInputRun] = *uGlobalFineTimeLimitsIN;
    }

    gDirectory->GetObject("bDisableBoard", bDisableBoardIN);
    if(!bDisableBoardIN)
    {
      cout<<Form("Could not retrieve 'bDisableBoard' from input calibration file of run %s.", kcInputRunID[iInputRun])<<endl;
      return;
    }
    else
    {
      bDisableBoardOUT[iInputRun] = *bDisableBoardIN;
    }

    gDirectory->GetObject("dFinalOffset", dFinalOffsetIN);
    if(!dFinalOffsetIN)
    {
      cout<<Form("Could not retrieve 'dFinalOffset' from input calibration file of run %s.", kcInputRunID[iInputRun])<<endl;
      return;
    }
    else
    {
      dFinalOffsetOUT[iInputRun] = *dFinalOffsetIN;
    }

    for(Int_t iTDCBoard = 0; iTDCBoard < kiNTDCBoards; iTDCBoard++)
    {
      if(!tCalibParent->cd(Form("tdc_b%03d", iTDCBoard)))
      {
        cout<<Form("Could not find directory 'tdc_b%03d' in input calibration file of run %s.", iTDCBoard, kcInputRunID[iInputRun])<<endl;
        return;
      }

      gDirectory->GetObject(Form("calib_trb_ft_b%03d_ref", iTDCBoard), tRefFineTimeHistosIN.at(iTDCBoard));

      if(!tRefFineTimeHistosIN.at(iTDCBoard))
      {
        cout<<Form("Could not retrieve the reference fine time histogram of TDC %d in run %s.", iTDCBoard, kcInputRunID[iInputRun])<<::endl;
        return;
      }
      else
      {
        tRefFineTimeHistosIN.at(iTDCBoard)->SetName(Form("tof_trb_ft_b%03d_ref", iTDCBoard));
        (tRefFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard) = dynamic_cast<TH1I*>(tRefFineTimeHistosIN.at(iTDCBoard)->Clone(Form("calib_trb_ft_r%03d_b%03d_ref",
                                                                                                                                 iInputRun, iTDCBoard)));
        (tRefFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)->SetDirectory(gROOT);
      }

      for(Int_t iTDCChannel = 0; iTDCChannel < kiNTDCChannels; iTDCChannel++)
      {
        gDirectory->GetObject(Form("calib_trb_ft_b%03d_ch%03d", iTDCBoard, iTDCChannel), (tFineTimeHistosIN.at(iTDCBoard)).at(iTDCChannel));

        if(!(tFineTimeHistosIN.at(iTDCBoard)).at(iTDCChannel))
        {
          cout<<Form("Could not retrieve the fine time histogram of channel %2d of TDC %2d in run %s.", iTDCChannel, iTDCBoard, kcInputRunID[iInputRun])<<endl;
          return;
        }
        else
        {
          (tFineTimeHistosIN.at(iTDCBoard)).at(iTDCChannel)->SetName(Form("tof_trb_ft_b%03d_ch%03d", iTDCBoard, iTDCChannel));
          ((tFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)).at(iTDCChannel) = dynamic_cast<TH1I*>((tFineTimeHistosIN.at(iTDCBoard)).at(iTDCChannel)->Clone(Form("calib_trb_ft_r%03d_b%03d_ch%03d",
                                                                                                                                                                 iInputRun, iTDCBoard, iTDCChannel)));
          ((tFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)).at(iTDCChannel)->SetDirectory(gROOT);
        }
      }
    }

    tInputCalibFile->Close();
    delete tInputCalibFile;

    gSystem->Exec(Form("mv -v %s %s.old", tFileName.Data(), tFileName.Data()));
  }


  if(kbMergeFineTime)
  {
    if(kbJointEdgeDetection)
    {
      for(Int_t iInputRun = 0; iInputRun < kiNInputRuns; iInputRun++)
      {
        for(Int_t iTDCBoard = 0; iTDCBoard < kiNTDCBoards; iTDCBoard++)
        {
          for(Int_t iTDCChannel = 0; iTDCChannel < kiNTDCChannels; iTDCChannel++)
          {
            if(0 == iTDCChannel%2)
            {
              ((tFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)).at(iTDCChannel)->Add(((tFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)).at(iTDCChannel + 1));
            }
            else
            {
              ((tFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)).at(iTDCChannel)->Reset("ICESM");
              ((tFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)).at(iTDCChannel)->Add(((tFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)).at(iTDCChannel - 1));
            }
          }
        }
      }
    }

    for(Int_t iTDCBoard = 0; iTDCBoard < kiNTDCBoards; iTDCBoard++)
    {
      // merge FT histograms for run 0
      for(Int_t iInputRun = 1; iInputRun < kiNInputRuns; iInputRun++)
      {
        (tRefFineTimeHistosOUT.at(0)).at(iTDCBoard)->Add((tRefFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard));
      }

      // copy merged FT histograms for run 0 to the other runs
      for(Int_t iInputRun = 1; iInputRun < kiNInputRuns; iInputRun++)
      {
        (tRefFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)->Reset("ICESM");
        (tRefFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)->Add((tRefFineTimeHistosOUT.at(0)).at(iTDCBoard));
      }

      for(Int_t iTDCChannel = 0; iTDCChannel < kiNTDCChannels; iTDCChannel++)
      {
        for(Int_t iInputRun = 1; iInputRun < kiNInputRuns; iInputRun++)
        {
          ((tFineTimeHistosOUT.at(0)).at(iTDCBoard)).at(iTDCChannel)->Add(((tFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)).at(iTDCChannel));
        }

        for(Int_t iInputRun = 1; iInputRun < kiNInputRuns; iInputRun++)
        {
          ((tFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)).at(iTDCChannel)->Reset("ICESM");
          ((tFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)).at(iTDCChannel)->Add(((tFineTimeHistosOUT.at(0)).at(iTDCBoard)).at(iTDCChannel));
        }
      }
    }
  }


  for(Int_t iInputRun = 0; iInputRun < kiNInputRuns; iInputRun++)
  {
    if(kiReferenceRun != iInputRun)
    {
      for(Int_t iTDCBoard = 0; iTDCBoard < kiNTDCBoards; iTDCBoard++)
      {
        if(0. != dFinalOffsetOUT[iInputRun].at(iTDCBoard))
        {
          dFinalOffsetOUT[iInputRun].at(iTDCBoard) -= dFinalOffsetOUT[kiReferenceRun].at(iTDCBoard);

          Double_t dTimeOffset = std::round(dFinalOffsetOUT[iInputRun].at(iTDCBoard));

          cout<<Form("TDC %03d clock cycle offset between runs %s and %s: %f -> %f", iTDCBoard, kcInputRunID[iInputRun], kcInputRunID[kiReferenceRun], dFinalOffsetOUT[iInputRun].at(iTDCBoard), dTimeOffset)<<endl;

          dFinalOffsetOUT[iInputRun].at(iTDCBoard) = dFinalOffsetOUT[kiReferenceRun].at(iTDCBoard) + dTimeOffset;
        }
      }

      cout<<endl;
    }
  }


  for(Int_t iInputRun = 0; iInputRun < kiNInputRuns; iInputRun++)
  {
    TString tFileName = Form("%s/CbmTofAna_%s/UnpackDigis_EVT%dns_TRG%d_MUL%d/data/%s", kcAnalysisDirName, kcInputRunID[iInputRun],
                             kiEventWindow, kiTriggerSet, kiTriggerMUL, kcCalibFileName);

    TFile* tOutputCalibFile = new TFile(tFileName, "RECREATE");

    tOutputCalibFile->mkdir("TofTrbCalib");
    tOutputCalibFile->cd("TofTrbCalib");

    for(Int_t iTDCBoard = 0; iTDCBoard < kiNTDCBoards; iTDCBoard++)
    {
      TDirectory* tBoardDirectory = gDirectory->mkdir(Form("tdc_b%03d", iTDCBoard));
      tBoardDirectory->cd();

      (tRefFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)->SetName(Form("calib_trb_ft_b%03d_ref", iTDCBoard));
      (tRefFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)->Write();
      (tRefFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)->SetName(Form("calib_trb_ft_r%03d_b%03d_ref", iInputRun, iTDCBoard));

      for(Int_t iTDCChannel = 0; iTDCChannel < kiNTDCChannels; iTDCChannel++)
      {
        ((tFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)).at(iTDCChannel)->SetName(Form("calib_trb_ft_b%03d_ch%03d", iTDCBoard, iTDCChannel));
        ((tFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)).at(iTDCChannel)->Write();
        ((tFineTimeHistosOUT.at(iInputRun)).at(iTDCBoard)).at(iTDCChannel)->SetName(Form("calib_trb_ft_r%03d_b%03d_ch%03d", iInputRun, iTDCBoard, iTDCChannel));
      }

      gDirectory->cd("../");
    }

    gDirectory->WriteObjectAny(&iReferenceBoardOUT[iInputRun], "vector<Int_t>", "iReferenceBoard");
    gDirectory->WriteObjectAny(&uGlobalFineTimeLimitsOUT[iInputRun], "vector<UInt_t>", "uGlobalFineTimeLimits");
    gDirectory->WriteObjectAny(&bDisableBoardOUT[iInputRun], "vector<Bool_t>", "bDisableBoard");
    gDirectory->WriteObjectAny(&dFinalOffsetOUT[iInputRun], "vector<Double_t>", "dFinalOffset");

    tOutputCalibFile->Close();
    delete tOutputCalibFile;
  }

}