TutAnalyzer_module.cc

Go to the documentation of this file.

////////////////////////////////////////////////////////////////////////
/// \brief   TODO
/// \author  bckhouse@hep.caltech.edu
////////////////////////////////////////////////////////////////////////

// Framework includes
#include "art/Framework/Core/EDAnalyzer.h"
#include "art/Framework/Core/ModuleMacros.h"
#include "art/Framework/Principal/Event.h"
#include "art/Framework/Principal/Handle.h"
#include "art/Framework/Services/Optional/TFileService.h"
#include "art/Framework/Services/Registry/ServiceHandle.h"

// ROOT includes
#include "TH1F.h"

// NOvASoft includes
#include "RecoBase/Prong.h"
#include "RecoBase/RecoHit.h"

namespace tut
{
  class TutAnalyzer : public art::EDAnalyzer
  {
  public:
    explicit TutAnalyzer(const fhicl::ParameterSet& pset);
    ~TutAnalyzer();

    void analyze(const art::Event& evt);

    void reconfigure(const fhicl::ParameterSet& pset);

    void beginJob();

  protected:
    std::string fProngLabel;

    TH1* fMassPeak;
  };
}


////////////////////////////////////////////////////////////////////////
namespace tut
{
  //.......................................................................
  TutAnalyzer::TutAnalyzer(const fhicl::ParameterSet& pset)
    : EDAnalyzer(pset)
  {
    reconfigure(pset);
  }

  //......................................................................
  TutAnalyzer::~TutAnalyzer()
  {
  }

  //......................................................................
  void TutAnalyzer::reconfigure(const fhicl::ParameterSet& pset)
  {
    fProngLabel = pset.get<std::string>("ProngLabel");
  }

  //......................................................................
  void TutAnalyzer::beginJob()
  {
    art::ServiceHandle<art::TFileService> tfs;
    fMassPeak = tfs->make<TH1F>("massPeak", ";Reconstructed mass (MeV)", 125, 0, 250);
  }

  //......................................................................
  void TutAnalyzer::analyze(const art::Event& evt)
  {
    art::Handle<std::vector<rb::Prong> > prongs;
    evt.getByLabel(fProngLabel, prongs);

    if(prongs->size() != 2) return;

    double E[2] = {0,};
    for(int photIdx = 0; photIdx < 2; ++photIdx){
      const rb::Prong& phot = (*prongs)[photIdx];
      for(unsigned int cellIdx = 0; cellIdx < phot.NCell(); ++cellIdx){
        // The prong knows how to use its 3D information to
        // attenuation-calibrate this hit.
        const rb::RecoHit rhit = phot.RecoHit(phot.Cell(cellIdx));
        if(rhit.IsCalibrated()) E[photIdx] += rhit.GeV();
      }
    }

    // Prong directions are guaranteed to be unit vectors
    const double dot = (*prongs)[0].Dir().Dot((*prongs)[1].Dir());

    const double m = sqrt(2*E[0]*E[1]*(1-dot));

    fMassPeak->Fill(1000*m);
  }

} // end namespace tut

////////////////////////////////////////////////////////////////////////
namespace tut
{
  DEFINE_ART_MODULE(TutAnalyzer)
}