SdSimulationCalibrationFillerOG/SdSimulationCalibrationFiller.cc

Go to the documentation of this file.

#include <det/Detector.h>

#include <sdet/SDetector.h>
#include <sdet/PMTConstants.h>

#include <evt/Event.h>

#include <fwk/CentralConfig.h>

#include <sevt/SEvent.h>
#include <sevt/Station.h>
#include <sevt/StationSimData.h>
#include <sevt/StationCalibData.h>
#include <sevt/PMTCalibData.h>
#include <sevt/PMTSimData.h>
#include <sevt/PMT.h>

#include <sdet/PMT.h>

#include <utl/Reader.h>

#include "SdSimulationCalibrationFiller.h"

using namespace std;
using namespace sevt;
using namespace utl;
using namespace fwk;
using namespace evt;

using namespace SdSimulationCalibrationFillerOG;


VModule::ResultFlag
SdSimulationCalibrationFiller::Init()
{
  auto topB = CentralConfig::GetInstance()->GetTopBranch("SdSimulationCalibrationFiller");

  const auto calibratorBranch = topB.GetChild("calibratorBranch").Get<string>();

  // get these from the config xml of SdCalibrator
  auto calibTopB = CentralConfig::GetInstance()->GetTopBranch(calibratorBranch);

  calibTopB.GetChild("peakOnlineToVEMFactor").GetData(fPeakOnlineToVEMFactor);
  calibTopB.GetChild("chargeOnlineToVEMFactor").GetData(fChargeOnlineToVEMFactor);

  ostringstream info;
  info << "from top branch <" << calibratorBranch << "> of the calibrator module got "
          "the peak online-to-vem factor " << fPeakOnlineToVEMFactor << " and "
          "the charge online-to-vem " << fChargeOnlineToVEMFactor;
  INFO(info);

  // continue with SdSimulationCalibrationFiller xml
  topB.GetChild("calibrationVersion").GetData(fCalibVersion);

  // SmallPMT calibration filling
  // TODO : temporarily here, to be moved to standard
  // XML sim. calib. consts configuration asap
  auto spmtCalibB = topB.GetChild("SmallPMTCalibration");
  if (spmtCalibB) {
    spmtCalibB.GetChild("version").GetData(fSPMTVersion);
    spmtCalibB.GetChild("beta").GetData(fBeta);
    spmtCalibB.GetChild("betaErr").GetData(fBetaError);
    spmtCalibB.GetChild("chi2").GetData(fChi2);
    info.str("");
    info << "Got calibration constants for SmallPMT -- version " << fSPMTVersion << "\n"
            "\tbeta = " << fBeta << " +- " << fBetaError << " (chi2 = " << fChi2 << ")";

    // Values with respect to each single LPMT
    const auto betas = spmtCalibB.GetChild("betasLPMT").Get<vector<double>>();
    const auto betaErrs = spmtCalibB.GetChild("betaErrsLPMT").Get<vector<double>>();
    const auto chi2s = spmtCalibB.GetChild("chi2sLPMT").Get<vector<double>>();

    // DV: these conditions should be encoded in the xsd!
    if (betas.size() != 3 || betaErrs.size() != 3 || chi2s.size() != 3) {
      ERROR("Wrong SmallPMT calibration values for single LPMTs");
      return eFailure;
    } else {
      info << "\tsingle LPMTs betas:";
      for (int i = 0; i < 3; ++i) {
        fBetaLPMT[i] = betas.at(i);
        fBetaErrorLPMT[i] = betaErrs[i];
        fChi2LPMT[i] = chi2s[i];
        info << "\n"
                "\tLPMT" << i+1
             << " beta = " << fBetaLPMT[i] << " +- " << fBetaErrorLPMT[i]
             << " (chi2 = " << fChi2LPMT[i] << ")";
      }
    }

    INFO(info);
  }

  return eSuccess;
}


VModule::ResultFlag
SdSimulationCalibrationFiller::Run(Event& event)
{
  if (!event.HasSEvent())
    event.MakeSEvent();

  auto& sEvent = event.GetSEvent();

  for (auto& station : sEvent.CandidateStationsRange()) {

    const bool isUUB = station.IsUUB();

    string signature;
    if (station.HasSimData())
      signature = station.GetSimData().GetSimulatorSignature();
    else
      ERROR("No SimData found, though it should exist at this point");

    if (!station.HasCalibData())
      station.MakeCalibData();

    auto& stationCalib = station.GetCalibData();

    stationCalib.SetVersion(fCalibVersion);
    stationCalib.SetStartSecond(1);
    stationCalib.SetEndSecond(3);
    stationCalib.SetNT1(100);
    stationCalib.SetNT2(20);
    stationCalib.SetNTot(1);

    for (auto& pmt : station.PMTsRange(sdet::PMTConstants::eAnyType)) {

      if (!pmt.HasCalibData())
        pmt.MakeCalibData();
      auto& pmtCalib = pmt.GetCalibData();

      if (!pmt.HasSimData())
        pmt.MakeSimData();
      auto& pmtSim = pmt.GetSimData();

      const auto& dStation = det::Detector::GetInstance().GetSDetector().GetStation(station);
      const auto pmtId = pmt.GetId();
      const auto& dPMT = dStation.GetPMT(pmtId);

      // apply calibration magic numbers from Laura (in reverse)
      // clarify correct numbers for the case of small PMT

      const double onlinePeak = dPMT.GetVEMPeak(signature) / fPeakOnlineToVEMFactor;
      pmtSim.SetTriggerVEMPeak(onlinePeak);  // trigger runs with biased online peak
      if (!isUUB)
        pmtCalib.SetOnlinePeak(onlinePeak);  // used by the calibrator
      else {
        // DV we assume here that compatibility VEM bias in UUB is the same as VEM bias in UB
        const double compatPeak = dPMT.GetCompatibilityVEMPeak(signature) / fPeakOnlineToVEMFactor;
        pmtCalib.SetOnlinePeak(compatPeak);  // UUB calibrator needs compat VEM
        pmtSim.SetCompatibilityTriggerVEMPeak(compatPeak);
      }
      const double onlineCharge = dPMT.GetVEMCharge(signature) / fChargeOnlineToVEMFactor;
      pmtCalib.SetOnlineCharge(onlineCharge);

      for (const auto gain : { sdet::PMTConstants::eHighGain, sdet::PMTConstants::eLowGain })
        pmtCalib.SetBaseline(dPMT.GetBaseline(gain), dPMT.GetBaselineRMS(gain), gain);

      pmtCalib.SetIsTubeOk(true);
      pmtCalib.SetIsLowGainOk(true);
      pmtCalib.SetRate(7e3);
      pmtCalib.SetGainRatio(dPMT.GetGainRatio(), dPMT.GetGainRatioRMS());

      pmtCalib.SetEvolution(0);
      pmtCalib.SetNumberTDA(0);
      pmtCalib.SetHighGainDelay(0, 0);
      pmtCalib.SetHighGainDelayChi2(0);

      if (pmt.GetType() == sdet::PMTConstants::eWaterCherenkovSmall) {
        auto& spmtData = station.GetSmallPMTData();
        if (!spmtData.HasCalibData())
          spmtData.MakeCalibData();
        auto& spmtCalib = spmtData.GetCalibData();
        spmtCalib.SetIsTubeOk(true);
        spmtCalib.SetVersion(fSPMTVersion);  // temporary
        spmtCalib.SetBeta(fBeta);  // temporary
        spmtCalib.SetBetaError(fBetaError);  // temporary
        spmtCalib.SetChi2(fChi2);  // temporary
        spmtCalib.SetCorrectionFactor(1);
        for (const auto& lpmt : station.PMTsRange(sdet::PMTConstants::eWaterCherenkovLarge)) {
          const auto lpmtId = lpmt.GetId();
          const auto i = lpmtId - dStation.GetFirstPMTId();
          spmtCalib.SetBeta(fBetaLPMT[i], lpmtId);  // temporary
          spmtCalib.SetBetaError(fBetaErrorLPMT[i], lpmtId);  // temporary
          spmtCalib.SetChi2(fChi2LPMT[i], lpmtId);  // temporary
          spmtCalib.SetCorrectionFactor(1, lpmtId);
        }
      }

    }
  }

  return eSuccess;
}