SBNDFlashFinder_module.cc

Go to the documentation of this file.

////////////////////////////////////////////////////////////////////////
//// Class:       SBNDFlashFinder
//// Module Type: producer
//// File:        SBNDFlashFinder_module.cc
////
//// Adapted form ICARUSFlashFinder by Kazuhiro Terao
//////////////////////////////////////////////////////////////////////////

#include "art/Framework/Core/EDProducer.h"
#include "art/Framework/Core/ModuleMacros.h"
#include "art/Framework/Principal/Event.h"
#include "art/Framework/Principal/Handle.h"
#include "art/Framework/Principal/Run.h"
#include "art/Framework/Principal/SubRun.h"
#include "canvas/Utilities/InputTag.h"
#include "fhiclcpp/ParameterSet.h"
#include "messagefacility/MessageLogger/MessageLogger.h"
#include "art/Utilities/make_tool.h"

#include "lardataobj/RecoBase/OpHit.h"
#include "lardataobj/RecoBase/OpFlash.h"
#include "lardata/Utilities/AssociationUtil.h"

#include <memory>
#include <string>
#include "sbndcode/OpDetReco/OpFlash/FlashFinder/FlashFinderManager.h"
#include "sbndcode/OpDetReco/OpFlash/FlashFinder/FlashFinderFMWKInterface.h"
#include "sbndcode/OpDetReco/OpFlash/FlashFinder/PECalib.h"
#include "sbndcode/OpDetReco/OpFlash/FlashTools/FlashGeoBase.hh"
#include "sbndcode/OpDetReco/OpFlash/FlashTools/FlashT0Base.hh"
#include "sbndcode/OpDetReco/OpFlash/FlashTools/DriftEstimatorBase.hh"

namespace opdet{

  class SBNDFlashFinder;

  class SBNDFlashFinder : public art::EDProducer {
  public:
    explicit SBNDFlashFinder(fhicl::ParameterSet const & p);
    // The destructor generated by the compiler is fine for classes
    // without bare pointers or other resource use.

    // Plugins should not be copied or assigned.
    SBNDFlashFinder(SBNDFlashFinder const &) = delete;
    SBNDFlashFinder(SBNDFlashFinder &&) = delete;
    SBNDFlashFinder & operator = (SBNDFlashFinder const &) = delete;
    SBNDFlashFinder & operator = (SBNDFlashFinder &&) = delete;

    // Required functions.
    void produce(art::Event & e) override;

  private:

    ::lightana::FlashFinderManager _mgr;
    ::lightana::PECalib _pecalib;
    std::vector<std::string> _hit_producers;
    std::string _ophit_input_time;
    bool _use_t0tool;
    bool _correct_light_propagation;

    // Tool for calculating the OpFlash Y and Z centers
    std::unique_ptr<lightana::FlashGeoBase> _flashgeo;

    // Tool for calculating the OpFlash t0
    std::unique_ptr<lightana::FlashT0Base> _flasht0calculator;

    // Tool for light propagation correction
    std::unique_ptr<lightana::DriftEstimatorBase> _driftestimator;

    ::lightana::LiteOpHitArray_t GetAssociatedLiteHits(::lightana::LiteOpFlash_t lite_flash,
                                                       ::lightana::LiteOpHitArray_t lite_hits_v);

  };

  SBNDFlashFinder::SBNDFlashFinder(lightana::Config_t const & p)
  : EDProducer{p}
  // Initialize member data here.
  {
    _hit_producers = p.get<std::vector<std::string>>("OpHitProducers");

    auto const flash_algo  = p.get<std::string>("FlashFinderAlgo");
    auto const flash_pset = p.get<lightana::Config_t>("AlgoConfig");
    auto algo_ptr = ::lightana::FlashAlgoFactory::get().create(flash_algo,flash_algo);
    algo_ptr->Configure(flash_pset);
    _mgr.SetFlashAlgo(algo_ptr);
    _pecalib.Configure(p.get<lightana::Config_t>("PECalib"));
    _ophit_input_time = p.get<std::string>("OpHitInputTime", "PeakTime");
    _use_t0tool = p.get<bool>("UseT0Tool", false);
    _correct_light_propagation = p.get<bool>("CorrectLightPropagation", false);

    auto const flashgeo_pset = p.get<lightana::Config_t>("FlashGeoConfig");
    _flashgeo = art::make_tool<lightana::FlashGeoBase>(flashgeo_pset);

    if(_use_t0tool){
      auto const flasht0_pset = p.get<lightana::Config_t>("FlashT0Config");
      _flasht0calculator = art::make_tool<lightana::FlashT0Base>(flasht0_pset);
    }

    if(_correct_light_propagation){
      auto const driftestimator_pset = p.get<lightana::Config_t>("DriftEstimatorConfig");
      _driftestimator = art::make_tool<lightana::DriftEstimatorBase>(driftestimator_pset);
    }

    produces< std::vector<recob::OpFlash>   >();
    produces< art::Assns <recob::OpHit, recob::OpFlash> >();
  }

  void SBNDFlashFinder::produce(art::Event & e)
  {

    // produce OpFlash data-product to be filled within module
    std::unique_ptr< std::vector<recob::OpFlash> > opflashes(new std::vector<recob::OpFlash>);
    std::unique_ptr< art::Assns <recob::OpHit, recob::OpFlash> > flash2hit_assn_v
      (new art::Assns<recob::OpHit, recob::OpFlash>);

    std::vector<art::Ptr<recob::OpHit>> ophit_v;

    ::lightana::LiteOpHitArray_t ophits;
    double trigger_time=1.1e20;

    for (auto producer : _hit_producers) {

      // load OpHits previously created
      art::Handle<std::vector<recob::OpHit> > ophit_h;
      e.getByLabel(producer, ophit_h);

      // make sure hits look good
      if(!ophit_h.isValid()) {
        std::cerr << "\033[93m[ERROR]\033[00m ... could not locate OpHit!" << std::endl;
        throw std::exception();
      }

      std::vector<art::Ptr<recob::OpHit>> temp_v;
      art::fill_ptr_vector(temp_v, ophit_h);
      ophit_v.insert(ophit_v.end(), temp_v.begin(), temp_v.end());
    }

    for(auto const oph : ophit_v) {
      ::lightana::LiteOpHit_t loph;
      if(trigger_time > 1.e20) trigger_time = oph->PeakTimeAbs() - oph->PeakTime();

      if(_ophit_input_time=="RiseTime") loph.peak_time = oph->StartTime()+oph->RiseTime();
      else if(_ophit_input_time=="StartTime") loph.peak_time = oph->StartTime();
      else loph.peak_time = oph->PeakTime();

      size_t opdet = ::lightana::OpDetFromOpChannel(oph->OpChannel());
      loph.pe = _pecalib.Calibrate(opdet,oph->Area());
      loph.channel = oph->OpChannel();
      ophits.emplace_back(std::move(loph));
    }

    auto const flash_v = _mgr.RecoFlash(ophits);

    for(const auto& lflash :  flash_v) {

      // Get Flash Barycenter
      double Ycenter, Zcenter, Ywidth, Zwidth;
      _flashgeo->GetFlashLocation(lflash.channel_pe, Ycenter, Zcenter, Ywidth, Zwidth);

      // Get flasht0
      double flasht0 = lflash.time;

      // Refine t0 calculation
      if(_use_t0tool)
        flasht0 = _flasht0calculator->GetFlashT0(lflash.time, GetAssociatedLiteHits(lflash, ophits));

      // Estimate drift location of the interaction and
      // make t0 unbias (lght propagation time correction)
      if(_correct_light_propagation){
        
        double drift_distance = _driftestimator->GetDriftPosition( lflash.channel_pe );
        double propagation_time = _driftestimator->GetPropagationTime( drift_distance );
        flasht0 = flasht0-propagation_time * 1e-3;

        drift_distance = (lflash.tpc==0 ? -drift_distance : drift_distance);

        recob::OpFlash flash(flasht0, lflash.time_err, trigger_time + flasht0,
                           (trigger_time + flasht0) / 1600., lflash.channel_pe,
                           0, 0, 1, // this are just default values
                           drift_distance, -1, Ycenter, Ywidth, Zcenter, Zwidth);
        opflashes->emplace_back(std::move(flash));
      }
      else{
        recob::OpFlash flash(flasht0, lflash.time_err, trigger_time + flasht0,
                           (trigger_time + flasht0) / 1600., lflash.channel_pe,
                           0, 0, 1, // this are just default values
                           100., -1., Ycenter, Ywidth, Zcenter, Zwidth);
        opflashes->emplace_back(std::move(flash));
      }

      // Create OpHit association
      for(auto const& hitidx : lflash.asshit_idx) {
        const art::Ptr<recob::OpHit> hit_ptr(ophit_v.at(hitidx));
        util::CreateAssn(*this, e, *opflashes, hit_ptr, *flash2hit_assn_v);
      }
    }

    e.put(std::move(opflashes));
    e.put(std::move(flash2hit_assn_v));
  }

  ::lightana::LiteOpHitArray_t SBNDFlashFinder::GetAssociatedLiteHits(::lightana::LiteOpFlash_t lite_flash,
                                                                      ::lightana::LiteOpHitArray_t lite_hits_v)
  {

    ::lightana::LiteOpHitArray_t flash_hits_v;

    for(auto const& hitidx : lite_flash.asshit_idx) {
      flash_hits_v.emplace_back(std::move(lite_hits_v.at(hitidx)));
    }

    return flash_hits_v;
  }

  DEFINE_ART_MODULE(SBNDFlashFinder)

}//end namespace