evgen::NDKGen Class Reference

A module to check the results from the Monte Carlo generator. More...

Inheritance diagram for evgen::NDKGen:

Public Types
using	ModuleType = EDProducer
using	WorkerType = WorkerT< EDProducer >
template<typename UserConfig , typename KeysToIgnore = void>
using	Table = ProducerBase::Table< UserConfig, KeysToIgnore >

Public Member Functions
	NDKGen (fhicl::ParameterSet const &pset)
virtual	~NDKGen ()
void	produce (art::Event &evt)
void	beginJob ()
void	beginRun (art::Run &run)
void	reconfigure (fhicl::ParameterSet const &p)
void	endJob ()
template<typename PROD , BranchType B = InEvent>
ProductID	getProductID (std::string const &instanceName={}) const
template<typename PROD , BranchType B>
ProductID	getProductID (ModuleDescription const &moduleDescription, std::string const &instanceName) const
bool	modifiesEvent () const
template<typename T , BranchType = InEvent>
ProductToken< T >	consumes (InputTag const &)
template<typename T , art::BranchType BT>
art::ProductToken< T >	consumes (InputTag const &it)
template<typename T , BranchType = InEvent>
void	consumesMany ()
template<typename Element , BranchType = InEvent>
ViewToken< Element >	consumesView (InputTag const &)
template<typename T , art::BranchType BT>
art::ViewToken< T >	consumesView (InputTag const &it)
template<typename T , BranchType = InEvent>
ProductToken< T >	mayConsume (InputTag const &)
template<typename T , art::BranchType BT>
art::ProductToken< T >	mayConsume (InputTag const &it)
template<typename T , BranchType = InEvent>
void	mayConsumeMany ()
template<typename Element , BranchType = InEvent>
ViewToken< Element >	mayConsumeView (InputTag const &)
template<typename T , art::BranchType BT>
art::ViewToken< T >	mayConsumeView (InputTag const &it)
base_engine_t &	createEngine (seed_t seed)
base_engine_t &	createEngine (seed_t seed, std::string const &kind_of_engine_to_make)
base_engine_t &	createEngine (seed_t seed, std::string const &kind_of_engine_to_make, label_t const &engine_label)
seed_t	get_seed_value (fhicl::ParameterSet const &pset, char const key[]="seed", seed_t const implicit_seed=-1)

Static Public Member Functions
static cet::exempt_ptr< Consumer >	non_module_context ()

Protected Member Functions
CurrentProcessingContext const *	currentContext () const
void	validateConsumedProduct (BranchType const bt, ProductInfo const &pi)
void	prepareForJob (fhicl::ParameterSet const &pset)
void	showMissingConsumes () const

Private Member Functions
std::string	ParticleStatus (int StatusCode)
std::string	ReactionChannel (int ccnc, int mode)
void	FillHistograms (simb::MCTruth mc)

Private Attributes
std::string	fNdkFile
std::ifstream *	fEventFile
TStopwatch	fStopwatch
std::string	fNDKModuleLabel
	keep track of how long it takes to run the job More...
TH1F *	fGenerated [6]
	Spectra as generated. More...
TH1F *	fVertexX
	vertex location of generated events in x More...
TH1F *	fVertexY
	vertex location of generated events in y More...
TH1F *	fVertexZ
	vertex location of generated events in z More...
TH2F *	fVertexXY
	vertex location in xy More...
TH2F *	fVertexXZ
	vertex location in xz More...
TH2F *	fVertexYZ
	vertex location in yz More...
TH1F *	fDCosX
	direction cosine in x More...
TH1F *	fDCosY
	direction cosine in y More...
TH1F *	fDCosZ
	direction cosine in z More...
TH1F *	fMuMomentum
	momentum of outgoing muons More...
TH1F *	fMuDCosX
	direction cosine of outgoing mu in x More...
TH1F *	fMuDCosY
	direction cosine of outgoing mu in y More...
TH1F *	fMuDCosZ
	direction cosine of outgoing mu in z More...
TH1F *	fEMomentum
	momentum of outgoing electrons More...
TH1F *	fEDCosX
	direction cosine of outgoing e in x More...
TH1F *	fEDCosY
	direction cosine of outgoing e in y More...
TH1F *	fEDCosZ
	direction cosine of outgoing e in z More...
TH1F *	fCCMode
	CC interaction mode. More...
TH1F *	fNCMode
	CC interaction mode. More...
TH1F *	fECons
	histogram to determine if energy is conserved in the event More...

Detailed Description

A module to check the results from the Monte Carlo generator.

Definition at line 64 of file NDKGen_module.cc.

Member Typedef Documentation

using art::EDProducer::ModuleType = EDProducer

inherited

Definition at line 34 of file EDProducer.h.

template<typename UserConfig , typename KeysToIgnore = void>

using art::EDProducer::Table = ProducerBase::Table<UserConfig, KeysToIgnore>

inherited

Definition at line 43 of file EDProducer.h.

using art::EDProducer::WorkerType = WorkerT<EDProducer>

inherited

Definition at line 35 of file EDProducer.h.

Constructor & Destructor Documentation

evgen::NDKGen::NDKGen ( fhicl::ParameterSet const &  pset )

explicit

Definition at line 125 of file NDKGen_module.cc.

  {
    this->reconfigure(pset); 
    fStopwatch.Start();

    produces< std::vector<simb::MCTruth> >();
    produces< sumdata::RunData, art::InRun >();

    fEventFile = new std::ifstream(fNdkFile.c_str());
    if(!fEventFile->good())
      exit(0);
    // create a default random engine; obtain the random seed from NuRandomService,
    // unless overridden in configuration with key "Seed"
    art::ServiceHandle<rndm::NuRandomService>()
      ->createEngine(*this, pset, "Seed");

  }

evgen::NDKGen::~NDKGen ( )

virtual

Definition at line 144 of file NDKGen_module.cc.

  {  
    fStopwatch.Stop();
  }

Member Function Documentation

void evgen::NDKGen::beginJob ( )

virtual

Reimplemented from art::EDProducer.

Definition at line 157 of file NDKGen_module.cc.

References geo::GeometryCore::DetHalfHeight(), geo::GeometryCore::DetHalfWidth(), geo::GeometryCore::DetLength(), art::TFileDirectory::make(), x, y, and z.

                       {
   
    // Get access to the TFile service.
    art::ServiceHandle<art::TFileService> tfs;

    fGenerated[0] = tfs->make<TH1F>("fGenerated_necc","",  100, 0.0, 20.0);
    fGenerated[1] = tfs->make<TH1F>("fGenerated_nebcc","", 100, 0.0, 20.0);
    fGenerated[2] = tfs->make<TH1F>("fGenerated_nmcc","",  100, 0.0, 20.0);
    fGenerated[3] = tfs->make<TH1F>("fGenerated_nmbcc","", 100, 0.0, 20.0);
    fGenerated[4] = tfs->make<TH1F>("fGenerated_nnc","",   100, 0.0, 20.0);
    fGenerated[5] = tfs->make<TH1F>("fGenerated_nbnc","",  100, 0.0, 20.0);
    
    fDCosX = tfs->make<TH1F>("fDCosX", ";dx/ds", 200, -1., 1.);
    fDCosY = tfs->make<TH1F>("fDCosY", ";dy/ds", 200, -1., 1.);
    fDCosZ = tfs->make<TH1F>("fDCosZ", ";dz/ds", 200, -1., 1.);

    fMuMomentum = tfs->make<TH1F>("fMuMomentum", ";p_{#mu} (GeV/c)", 500, 0., 50.);
    fMuDCosX    = tfs->make<TH1F>("fMuDCosX", ";dx/ds;", 200, -1., 1.);
    fMuDCosY    = tfs->make<TH1F>("fMuDCosY", ";dy/ds;", 200, -1., 1.);
    fMuDCosZ    = tfs->make<TH1F>("fMuDCosZ", ";dz/ds;", 200, -1., 1.);

    fEMomentum  = tfs->make<TH1F>("fEMomentum", ";p_{e} (GeV/c)", 500, 0., 50.);
    fEDCosX     = tfs->make<TH1F>("fEDCosX", ";dx/ds;", 200, -1., 1.);
    fEDCosY     = tfs->make<TH1F>("fEDCosY", ";dy/ds;", 200, -1., 1.);
    fEDCosZ     = tfs->make<TH1F>("fEDCosZ", ";dz/ds;", 200, -1., 1.);

    fCCMode = tfs->make<TH1F>("fCCMode", ";CC Interaction Mode;", 5, 0., 5.);
    fCCMode->GetXaxis()->SetBinLabel(1, "QE");
    fCCMode->GetXaxis()->SetBinLabel(2, "Res");
    fCCMode->GetXaxis()->SetBinLabel(3, "DIS");
    fCCMode->GetXaxis()->SetBinLabel(4, "Coh");
    fCCMode->GetXaxis()->SetBinLabel(5, "kInverseMuDecay");
    fCCMode->GetXaxis()->CenterLabels();

    fNCMode = tfs->make<TH1F>("fNCMode", ";NC Interaction Mode;", 5, 0., 5.);
    fNCMode->GetXaxis()->SetBinLabel(1, "QE");
    fNCMode->GetXaxis()->SetBinLabel(2, "Res");
    fNCMode->GetXaxis()->SetBinLabel(3, "DIS");
    fNCMode->GetXaxis()->SetBinLabel(4, "Coh");
    fNCMode->GetXaxis()->SetBinLabel(5, "kNuElectronElastic");
    fNCMode->GetXaxis()->CenterLabels();

    //fDeltaE = tfs->make<TH1F>("fDeltaE", ";#Delta E_{#nu} (GeV);", 200, -1., 1.); 
    fECons  = tfs->make<TH1F>("fECons", ";#Delta E(#nu,lepton);", 500, -5., 5.);

    art::ServiceHandle<geo::Geometry> geo;
    double x = 2.1*geo->DetHalfWidth();
    double y = 2.1*geo->DetHalfHeight();
    double z = 2.*geo->DetLength();
    int xdiv = TMath::Nint(2*x/5.);
    int ydiv = TMath::Nint(2*y/5.);
    int zdiv = TMath::Nint(2*z/5.);
    
    fVertexX = tfs->make<TH1F>("fVertexX", ";x (cm)", xdiv,  -x, x);
    fVertexY = tfs->make<TH1F>("fVertexY", ";y (cm)", ydiv,  -y, y);
    fVertexZ = tfs->make<TH1F>("fVertexZ", ";z (cm)", zdiv, -0.2*z, z);
    
    fVertexXY = tfs->make<TH2F>("fVertexXY", ";x (cm);y (cm)", xdiv,     -x, x, ydiv, -y, y);
    fVertexXZ = tfs->make<TH2F>("fVertexXZ", ";z (cm);x (cm)", zdiv, -0.2*z, z, xdiv, -x, x);
    fVertexYZ = tfs->make<TH2F>("fVertexYZ", ";z (cm);y (cm)", zdiv, -0.2*z, z, ydiv, -y, y);

  }

void evgen::NDKGen::beginRun ( art::Run &  run )

virtual

Reimplemented from art::EDProducer.

Definition at line 221 of file NDKGen_module.cc.

References geo::GeometryCore::DetectorName(), and art::Run::put().

  {

    // grab the geometry object to see what geometry we are using
    art::ServiceHandle<geo::Geometry> geo;
    std::unique_ptr<sumdata::RunData> runcol(new sumdata::RunData(geo->DetectorName()));

    run.put(std::move(runcol));

    return;
  }

template<typename T , BranchType = InEvent>

ProductToken<T> art::Consumer::consumes ( InputTag const &  )

inherited

template<typename T , art::BranchType BT>

art::ProductToken<T> art::Consumer::consumes ( InputTag const &  it )

inherited

Definition at line 147 of file Consumer.h.

References art::InputTag::instance(), art::InputTag::label(), and art::InputTag::process().

{
  if (!moduleContext_)
    return ProductToken<T>::invalid();

  consumables_[BT].emplace_back(ConsumableType::Product,
                                TypeID{typeid(T)},
                                it.label(),
                                it.instance(),
                                it.process());
  return ProductToken<T>{it};
}

template<typename T , art::BranchType BT>

void art::Consumer::consumesMany ( )

inherited

Definition at line 162 of file Consumer.h.

{
  if (!moduleContext_)
    return;

  consumables_[BT].emplace_back(ConsumableType::Many, TypeID{typeid(T)});
}

template<typename Element , BranchType = InEvent>

ViewToken<Element> art::Consumer::consumesView ( InputTag const &  )

inherited

template<typename T , art::BranchType BT>

art::ViewToken<T> art::Consumer::consumesView ( InputTag const &  it )

inherited

Definition at line 172 of file Consumer.h.

References art::InputTag::instance(), art::InputTag::label(), and art::InputTag::process().

{
  if (!moduleContext_)
    return ViewToken<T>::invalid();

  consumables_[BT].emplace_back(ConsumableType::ViewElement,
                                TypeID{typeid(T)},
                                it.label(),
                                it.instance(),
                                it.process());
  return ViewToken<T>{it};
}

EngineCreator::base_engine_t & EngineCreator::createEngine ( seed_t  seed )

inherited

Definition at line 26 of file EngineCreator.cc.

References art::EngineCreator::rng().

Referenced by evgen::CosmicsGen::CosmicsGen(), rndm::NuRandomService::createEngine(), cluster::fuzzyCluster::fuzzyCluster(), cluster::HoughLineFinder::HoughLineFinder(), art::MixFilter< T >::initEngine_(), larg4::LArG4::LArG4(), evgen::LightSource::LightSource(), evgen::NeutronOsc::NeutronOsc(), evgen::NucleonDecay::NucleonDecay(), opdet::OpMCDigi::OpMCDigi(), opdet::OptDetDigitizer::OptDetDigitizer(), phot::PhotonLibraryPropagation::PhotonLibraryPropagation(), detsim::SimDriftElectrons::SimDriftElectrons(), evgen::SingleGen::SingleGen(), evgen::SNNueAr40CCGen::SNNueAr40CCGen(), ToyOneShowerGen::ToyOneShowerGen(), and trkf::Track3DKalman::Track3DKalman().

{
  return rng()->createEngine(placeholder_schedule_id(), seed);
}

EngineCreator::base_engine_t & EngineCreator::createEngine ( seed_t  seed, std::string const &  kind_of_engine_to_make  )

inherited

Definition at line 32 of file EngineCreator.cc.

References art::EngineCreator::rng().

{
  return rng()->createEngine(
    placeholder_schedule_id(), seed, kind_of_engine_to_make);
}

EngineCreator::base_engine_t & EngineCreator::createEngine ( seed_t  seed, std::string const &  kind_of_engine_to_make, label_t const &  engine_label  )

inherited

Definition at line 40 of file EngineCreator.cc.

References art::EngineCreator::rng().

{
  return rng()->createEngine(
    placeholder_schedule_id(), seed, kind_of_engine_to_make, engine_label);
}

CurrentProcessingContext const * art::EDProducer::currentContext ( ) const

protectedinherited

Definition at line 120 of file EDProducer.cc.

References art::EDProducer::current_context_.

  {
    return current_context_.get();
  }

void evgen::NDKGen::endJob ( )

virtual

Reimplemented from art::EDProducer.

Definition at line 234 of file NDKGen_module.cc.

  {
    fEventFile->close();
  }

void evgen::NDKGen::FillHistograms ( simb::MCTruth  mc )

private

look for the outgoing lepton in the particle stack just interested in the first one

Definition at line 497 of file NDKGen_module.cc.

References simb::MCNeutrino::CCNC(), DEFINE_ART_MODULE, simb::MCParticle::E(), energy, simb::MCTruth::GetNeutrino(), simb::MCTruth::GetParticle(), simb::kCC, art::left(), simb::MCNeutrino::Mode(), simb::MCTruth::NParticles(), simb::MCNeutrino::Nu(), simb::MCParticle::P(), part, simb::MCParticle::PdgCode(), simb::MCParticle::Px(), simb::MCParticle::Py(), simb::MCParticle::Pz(), simb::MCParticle::Vx(), simb::MCParticle::Vy(), and simb::MCParticle::Vz().

  {
    // Decide which histograms to put the spectrum in
    int id = -1;
    if (mc.GetNeutrino().CCNC()==simb::kCC) {
      fCCMode->Fill(mc.GetNeutrino().Mode());
      if      (mc.GetNeutrino().Nu().PdgCode() ==  12) id = 0;
      else if (mc.GetNeutrino().Nu().PdgCode() == -12) id = 1;
      else if (mc.GetNeutrino().Nu().PdgCode() ==  14) id = 2;
      else if (mc.GetNeutrino().Nu().PdgCode() == -14) id = 3;
      else return;
    }
    else {
      fNCMode->Fill(mc.GetNeutrino().Mode());
      if (mc.GetNeutrino().Nu().PdgCode() > 0) id = 4;
      else                                     id = 5;
    }
    if (id==-1) abort();
    
    // Fill the specta histograms
    fGenerated[id]->Fill(mc.GetNeutrino().Nu().E() );
      
    //< fill the vertex histograms from the neutrino - that is always 
    //< particle 0 in the list
    simb::MCNeutrino       mcnu = mc.GetNeutrino();
    const simb::MCParticle nu   = mcnu.Nu();
    
    fVertexX->Fill(nu.Vx());
    fVertexY->Fill(nu.Vy());
    fVertexZ->Fill(nu.Vz());
    
    fVertexXY->Fill(nu.Vx(), nu.Vy());
    fVertexXZ->Fill(nu.Vz(), nu.Vx());
    fVertexYZ->Fill(nu.Vz(), nu.Vy());
    
    double mom = nu.P();
    if(std::abs(mom) > 0.){
      fDCosX->Fill(nu.Px()/mom);
      fDCosY->Fill(nu.Py()/mom);
      fDCosZ->Fill(nu.Pz()/mom);
    }


//     LOG_DEBUG("GENIEInteractionInformation") 
//       << std::endl
//       << "REACTION:  " << ReactionChannel(mc.GetNeutrino().CCNC(),mc.GetNeutrino().Mode()) 
//       << std::endl
//       << "-----------> Particles in the Stack = " << mc.NParticles() << std::endl
//       << std::setiosflags(std::ios::left) 
//       << std::setw(20) << "PARTICLE"
//       << std::setiosflags(std::ios::left) 
//       << std::setw(32) << "STATUS"
//       << std::setw(18) << "E (GeV)"
//       << std::setw(18) << "m (GeV/c2)"
//       << std::setw(18) << "Ek (GeV)"
//       << std::endl << std::endl;

//     const TDatabasePDG* databasePDG = TDatabasePDG::Instance();

//     // Loop over the particle stack for this event 
//     for(int i = 0; i < mc.NParticles(); ++i){
//       simb::MCParticle part(mc.GetParticle(i));
//       std::string name = databasePDG->GetParticle(part.PdgCode())->GetName();
//       int code = part.StatusCode();
//       std::string status = ParticleStatus(code);
//       double mass = part.Mass();
//       double energy = part.E(); 
//       double Ek = (energy-mass); // Kinetic Energy (GeV)
//       if(status=="kIStFinalStB4Interactions")
//      LOG_DEBUG("GENIEFinalState")
//        << std::setiosflags(std::ios::left) << std::setw(20) << name
//        << std::setiosflags(std::ios::left) << std::setw(32) <<status
//        << std::setw(18)<< energy
//        << std::setw(18)<< mass
//        << std::setw(18)<< Ek <<std::endl;
//       else 
//      LOG_DEBUG("GENIEFinalState") 
//        << std::setiosflags(std::ios::left) << std::setw(20) << name
//        << std::setiosflags(std::ios::left) << std::setw(32) << status
//        << std::setw(18) << energy
//        << std::setw(18) << mass <<std::endl; 

    std::cout << "REACTION:  " << ReactionChannel(mc.GetNeutrino().CCNC(),mc.GetNeutrino().Mode()) << std::endl;
    std::cout << "-----------> Particles in the Stack = " << mc.NParticles() << std::endl;
    std::cout << std::setiosflags(std::ios::left) 
              << std::setw(20) << "PARTICLE"
              << std::setiosflags(std::ios::left) 
              << std::setw(32) << "STATUS"
              << std::setw(18) << "E (GeV)"
              << std::setw(18) << "m (GeV/c2)"
              << std::setw(18) << "Ek (GeV)"
              << std::endl << std::endl;
    
    const TDatabasePDG* databasePDG = TDatabasePDG::Instance();
    
    // Loop over the particle stack for this event 
    for(int i = 0; i < mc.NParticles(); ++i){
      simb::MCParticle part(mc.GetParticle(i));
      std::string name;
      if (part.PdgCode() == 18040)
        name = "Ar40 18040";
      else if (part.PdgCode() != -99999 )
        {
          name = databasePDG->GetParticle(part.PdgCode())->GetName(); 
        }
      
      int code = part.StatusCode();
      std::string status = ParticleStatus(code);
      double mass = part.Mass();
      double energy = part.E(); 
      double Ek = (energy-mass); // Kinetic Energy (GeV)
      
      std::cout << std::setiosflags(std::ios::left) << std::setw(20) << name
                << std::setiosflags(std::ios::left) << std::setw(32) <<status
                << std::setw(18)<< energy
                << std::setw(18)<< mass
                << std::setw(18)<< Ek <<std::endl;  
    }

    if(mc.GetNeutrino().CCNC() == simb::kCC){
  
      for(int i = 0; i < mc.NParticles(); ++i){
        simb::MCParticle part(mc.GetParticle(i));
        if(abs(part.PdgCode()) == 11){
          fEMomentum->Fill(part.P());
          fEDCosX->Fill(part.Px()/part.P());
          fEDCosY->Fill(part.Py()/part.P());
          fEDCosZ->Fill(part.Pz()/part.P());
          fECons->Fill(nu.E() - part.E());
          break;
        }
        else if(abs(part.PdgCode()) == 13){
          fMuMomentum->Fill(part.P());
          fMuDCosX->Fill(part.Px()/part.P());
          fMuDCosY->Fill(part.Py()/part.P());
          fMuDCosZ->Fill(part.Pz()/part.P());
          fECons->Fill(nu.E() - part.E());
          break;
        }
      }// end loop over particles
    }//end if CC interaction

    return;
  }

EngineCreator::seed_t EngineCreator::get_seed_value ( fhicl::ParameterSet const &  pset, char const  key[] = "seed", seed_t const  implicit_seed = -1  )

inherited

Definition at line 49 of file EngineCreator.cc.

References fhicl::ParameterSet::get().

Referenced by art::MixFilter< T >::initEngine_().

{
  auto const& explicit_seeds = pset.get<std::vector<int>>(key, {});
  return explicit_seeds.empty() ? implicit_seed : explicit_seeds.front();
}

template<typename PROD , BranchType B>

ProductID art::EDProducer::getProductID ( std::string const &  instanceName = {} ) const

inlineinherited

Definition at line 123 of file EDProducer.h.

References art::EDProducer::moduleDescription_.

  {
    return ProducerBase::getProductID<PROD, B>(moduleDescription_,
                                               instanceName);
  }

template<typename PROD , BranchType B>

ProductID art::ProducerBase::getProductID ( ModuleDescription const &  moduleDescription, std::string const &  instanceName  ) const

inherited

Definition at line 56 of file ProducerBase.h.

References B, and art::ModuleDescription::moduleLabel().

Referenced by art::ProducerBase::modifiesEvent().

  {
    auto const& pd =
      get_ProductDescription<PROD>(B, md.moduleLabel(), instanceName);
    return pd.productID();
  }

template<typename T , BranchType = InEvent>

ProductToken<T> art::Consumer::mayConsume ( InputTag const &  )

inherited

template<typename T , art::BranchType BT>

art::ProductToken<T> art::Consumer::mayConsume ( InputTag const &  it )

inherited

Definition at line 190 of file Consumer.h.

References art::InputTag::instance(), art::InputTag::label(), and art::InputTag::process().

{
  if (!moduleContext_)
    return ProductToken<T>::invalid();

  consumables_[BT].emplace_back(ConsumableType::Product,
                                TypeID{typeid(T)},
                                it.label(),
                                it.instance(),
                                it.process());
  return ProductToken<T>{it};
}

template<typename T , art::BranchType BT>

void art::Consumer::mayConsumeMany ( )

inherited

Definition at line 205 of file Consumer.h.

{
  if (!moduleContext_)
    return;

  consumables_[BT].emplace_back(ConsumableType::Many, TypeID{typeid(T)});
}

template<typename Element , BranchType = InEvent>

ViewToken<Element> art::Consumer::mayConsumeView ( InputTag const &  )

inherited

template<typename T , art::BranchType BT>

art::ViewToken<T> art::Consumer::mayConsumeView ( InputTag const &  it )

inherited

Definition at line 215 of file Consumer.h.

References art::InputTag::instance(), art::InputTag::label(), and art::InputTag::process().

{
  if (!moduleContext_)
    return ViewToken<T>::invalid();

  consumables_[BT].emplace_back(ConsumableType::ViewElement,
                                TypeID{typeid(T)},
                                it.label(),
                                it.instance(),
                                it.process());
  return ViewToken<T>{it};
}

bool art::ProducerBase::modifiesEvent ( ) const

inlineinherited

Definition at line 40 of file ProducerBase.h.

References art::ProducerBase::getProductID().

    {
      return true;
    }

cet::exempt_ptr< art::Consumer > art::Consumer::non_module_context ( )

staticinherited

Definition at line 76 of file Consumer.cc.

Referenced by art::RootOutput::beginSubRun(), art::OutputModule::doBeginRun(), art::OutputModule::doBeginSubRun(), art::OutputModule::doEndRun(), art::OutputModule::doEndSubRun(), art::ProducingService::doPostReadEvent(), art::ProducingService::doPostReadRun(), art::ProducingService::doPostReadSubRun(), art::OutputModule::doWriteEvent(), art::ProcessPackage< L >::postScheduleSignal(), art::BeginEndPackage< Level::Run >::Begin::postScheduleSignal(), art::BeginEndPackage< Level::Run >::End::postScheduleSignal(), art::BeginEndPackage< Level::SubRun >::Begin::postScheduleSignal(), art::BeginEndPackage< Level::SubRun >::End::postScheduleSignal(), art::ProcessPackage< L >::preScheduleSignal(), art::BeginEndPackage< Level::Run >::Begin::preScheduleSignal(), art::BeginEndPackage< Level::SubRun >::Begin::preScheduleSignal(), art::EventProcessor::readEvent(), art::EventProcessor::readRun(), art::EmptyEvent::readRun_(), art::EventProcessor::readSubRun(), and art::EmptyEvent::readSubRun_().

{
  static Consumer invalid{InvalidTag{}};
  return &invalid;
}

std::string evgen::NDKGen::ParticleStatus ( int  StatusCode )

private

Definition at line 445 of file NDKGen_module.cc.

References sim::ParticleStatusName().

  {
    int code = StatusCode;
    std::string ParticleStatusName;

    switch(code)
      {
      case 0:
        ParticleStatusName = "kIStInitialState";
        break;
      case 1:
        ParticleStatusName = "kIStFinalState";
        break;
      case 11:
        ParticleStatusName = "kIStNucleonTarget";
        break;
      default:
        ParticleStatusName = "Status Unknown";
      }
    return ParticleStatusName;
  }

void art::Consumer::prepareForJob ( fhicl::ParameterSet const &  pset )

protectedinherited

Definition at line 89 of file Consumer.cc.

References fhicl::ParameterSet::get_if_present().

Referenced by art::EDProducer::doBeginJob(), art::EDFilter::doBeginJob(), and art::EDAnalyzer::doBeginJob().

{
  if (!moduleContext_)
    return;

  pset.get_if_present("errorOnMissingConsumes", requireConsumes_);
  for (auto& consumablesPerBranch : consumables_) {
    cet::sort_all(consumablesPerBranch);
  }
}

void evgen::NDKGen::produce ( art::Event &  evt )

virtual

Implements art::EDProducer.

Definition at line 239 of file NDKGen_module.cc.

References simb::MCTruth::Add(), simb::MCParticle::AddTrajectoryPoint(), geo::GeometryCore::DetHalfHeight(), geo::GeometryCore::DetHalfWidth(), geo::GeometryCore::DetLength(), E, art::EventID::event(), art::RandomNumberGenerator::getEngine(), art::Event::id(), in, simb::kUnknown, geo::TPCGeo::LocalToWorld(), simb::MCTruth::NParticles(), geo::GeometryCore::NTPC(), art::Event::put(), simb::MCTruth::SetOrigin(), and geo::GeometryCore::TPC().

  {

    std::cout << std::endl;
    std::cout<<"------------------------------------------------------------------------------"<<std::endl;
    //std::cout << "run    : " << evt.Header().Run() << std::endl;
    //std::cout << "subrun : " << evt.Header().Subrun() << std::endl;
    //std::cout << "event  : " << evt.Header().Event() << std::endl;
    std::cout << "event  : " << evt.id().event() << std::endl;  

    // TODO: fill more quantities out, as below.
    /*
    double X; // vertex position from Ndk
    double Y; // vertex position from Ndk
    double Z; // vertex position from Ndk
    double PDGCODE = -9999.; 
    double CHANNEL = -9999.;
    int channel = -9999;
    double energy = 0.; // in MeV from Ndk
    double cosx = 0.;
    double cosy = 0.;
    double cosz = 0.;

    int partnumber = 0;
    */

    std::string name, k, dollar;
     

    // event dump format on file output by the two commands ....
    // gevgen_ndcy  -g 1000180400  -m 8 -n 400 -o ndk
    // gevdump -f ndk.1000.ghep.root > ndk.out
    std::string Name;
    int Idx, Ist, PDG, Mother1, Mother2, Daughter1 ,Daughter2;
    double Px, Py, Pz, E, m ;
    std::string p1,p2,p3,p4,p5,p6,p7,p8,p9,p10,p11,p12,p13,p14;


    int trackid = -1; // set track id to -i as these are all primary particles and have id <= 0
    std::string primary("primary");
    int FirstMother = -1;
    double Mass = -9999;
    int Status = -9999;
    
    double P; // momentum of MCParticle IN GeV/c

    // TODO: Could perhaps imagine using these in NDk. 
    /*
    int targetnucleusPdg = -9999;
    int hitquarkPdg = -9999;
    double Q2 = -9999;
    */
    TLorentzVector Neutrino;
    TLorentzVector Lepton;
    TLorentzVector Target;
    TLorentzVector q;
    TLorentzVector Hadron4mom;

    // TODO: Could perhaps imagine using these in NDk. 
    /*
    int Tpdg = 0;  // for target 
    double Tmass = 0;
    int Tstatus = 11;
    double Tcosx, Tcosy, Tcosz, Tenergy;
    */

    TLorentzVector Tpos;
    
        
    std::unique_ptr< std::vector<simb::MCTruth> > truthcol(new std::vector<simb::MCTruth>);
    simb::MCTruth truth;

    art::ServiceHandle<geo::Geometry> geo;
    art::ServiceHandle<art::RandomNumberGenerator> rng;
    CLHEP::HepRandomEngine &engine = rng->getEngine();
    CLHEP::RandFlat flat(engine);

    double fvCut (5.0); // force vtx to be this far from any wall.

    // Find boundary of active volume
    double minx = 1e9;
    double maxx = -1e9;
    double miny = 1e9;
    double maxy = -1e9;
    double minz = 1e9;
    double maxz = -1e9;
    for (size_t i = 0; i<geo->NTPC(); ++i)
    {
      double local[3] = {0.,0.,0.};
      double world[3] = {0.,0.,0.};
      const geo::TPCGeo &tpc = geo->TPC(i);
      tpc.LocalToWorld(local,world);
      if (minx>world[0]-geo->DetHalfWidth(i)) 
        minx = world[0]-geo->DetHalfWidth(i);
      if (maxx<world[0]+geo->DetHalfWidth(i))
        maxx = world[0]+geo->DetHalfWidth(i);
      if (miny>world[1]-geo->DetHalfHeight(i))
        miny = world[1]-geo->DetHalfHeight(i);
      if (maxy<world[1]+geo->DetHalfHeight(i))
        maxy = world[1]+geo->DetHalfHeight(i);
      if (minz>world[2]-geo->DetLength(i)/2.)
        minz = world[2]-geo->DetLength(i)/2.;
      if (maxz<world[2]+geo->DetLength(i)/2.)
        maxz = world[2]+geo->DetLength(i)/2.;
    }

    // Assign vertice position
    double X0 = 0.0 + flat.fire( minx+fvCut , maxx-fvCut );
    double Y0 = 0.0 + flat.fire( miny+fvCut , maxy-fvCut );
    double Z0 = 0.0 + flat.fire( minz+fvCut , maxz-fvCut );

    std::cout << "NDKGen_module: X, Y, Z of vtx: " << X0 << ", "<< Y0 << ", "<< Z0 << std::endl;
    
    int GenieEvt = -999;
    
    if(!fEventFile->good())
      std::cout << "NdkFile: Problem reading Ndk file" << std::endl; 
    
    while(getline(*fEventFile,k)){

      if (k.find("** Event:")!= std::string::npos) {
        std::istringstream in;
        in.clear();
        in.str(k);
        std::string dummy;   
        in>> dummy>> dummy>> dummy >> dummy>> dummy>> dummy>> dummy >> dummy>> dummy>> dummy >> GenieEvt;
        std::cout<<"Genie Evt "<< GenieEvt <<" art evt "<<evt.id().event()<<"\n";
      }

      if (GenieEvt+1 != static_cast<int>(evt.id().event()))
        continue;
      else {

        if (!k.compare(0,25,"GENIE Interaction Summary")) // testing for new event.
          break;
        if (k.compare(0,1,"|") || k.compare(1,2,"  ")) continue; // uninteresting line if it doesn't start with "|" and if second and third characters aren't spaces.
        if (k.find("Fin-Init") != std::string::npos) continue; // Meh.
        if (k.find("Ar") != std::string::npos) continue; // Meh.
        if (k.find("Cl") != std::string::npos) continue; // ignore chlorine nucleus in nnbar events
        if (k.find("HadrBlob") != std::string::npos) continue; // Meh.
        if (k.find("NucBindE") != std::string::npos) continue; // Meh. atmo
        if (k.find("FLAGS") != std::string::npos) break; // Event end. gevgen_ndcy
        if (k.find("Vertex") != std::string::npos) break; // Event end. atmo

        //      if (!k.compare(26,1,"3") || !k.compare(26,1,"1")) ; // New event or stable particles.
        if (!k.compare(26,1,"1"))  // New event or stable particles.
        {

          std::istringstream in;
          in.clear();
          in.str(k);

          in>>p1>> Idx >>p2>> Name >>p3>> Ist >>p4>> PDG >>p5>>Mother1 >> p6 >> Mother2 >>p7>> Daughter1 >>p8>> Daughter2 >>p9>>Px>>p10>>Py>>p11>>Pz>>p12>>E>>p13>> m>>p14;
          //std::cout<<std::setprecision(9)<<dollar<<"  "<<name<<"  "<<PDGCODE<<"  "<<energy<<"  "<<cosx<<" "<<cosy<<"  "<<cosz<<"  "<<partnumber<<std::endl;
          if (Ist!=1) continue;

          std::cout << "PDG = " << PDG << std::endl;    
        
          const TDatabasePDG* databasePDG = TDatabasePDG::Instance();
          const TParticlePDG* definition = databasePDG->GetParticle(PDG);
          Mass = definition->Mass(); // GeV
          if (E-Mass < 0.001) continue; // KE is too low.
        
          //      if(partnumber == -1)
          Status = 1;
          
          simb::MCParticle mcpart(trackid,
                                PDG,
                                primary,                    
                                FirstMother,
                                Mass,
                                Status
                                );
                
          P = std::sqrt(pow(E,2.) - pow(Mass,2.)); // GeV/c
          std::cout << "Momentum = " << P << std::endl;
        
          TLorentzVector pos(X0, Y0, Z0, 0);

          Tpos = pos; // for target
        
          TLorentzVector mom(Px,Py,Pz, E);
        
          mcpart.AddTrajectoryPoint(pos,mom);
          truth.Add(mcpart);
        
        
        }// loop over particles in an event
        truth.SetOrigin(simb::kUnknown);
            
        //if (!k.compare(1,1,"FLAGS")) // end of event
        //  break;  

      }
    } // end while loop
    
    std::cout << "NDKGen.cxx: Putting " << truth.NParticles() << " tracks on stack." << std::endl; 
    truthcol->push_back(truth);
    //FillHistograms(truth);  
    evt.put(std::move(truthcol));
      
    return;
  }

std::string evgen::NDKGen::ReactionChannel ( int  ccnc, int  mode  )

private

Definition at line 469 of file NDKGen_module.cc.

  {
    std::string ReactionChannelName=" ";

    if(ccnc==0)
      ReactionChannelName = "kCC";
    else if(ccnc==1)
      ReactionChannelName = "kNC";
    else std::cout<<"Current mode unknown!! "<<std::endl;

    if(mode==0)
      ReactionChannelName += "_kQE";
    else if(mode==1)
      ReactionChannelName += "_kRes";
    else if(mode==2)
      ReactionChannelName += "_kDIS";
    else if(mode==3)
      ReactionChannelName += "_kCoh";
    else if(mode==4)
      ReactionChannelName += "_kNuElectronElastic";
    else if(mode==5)
      ReactionChannelName += "_kInverseMuDecay";
    else std::cout<<"interaction mode unknown!! "<<std::endl;

    return ReactionChannelName;
  }

void evgen::NDKGen::reconfigure

(

fhicl::ParameterSet const &

p

)

Definition at line 150 of file NDKGen_module.cc.

References fhicl::ParameterSet::get().

  {
    fNdkFile   =(p.get< std::string         >("NdkFile"));
    return;
  }

void art::Consumer::showMissingConsumes ( ) const

protectedinherited

Definition at line 125 of file Consumer.cc.

Referenced by art::EDProducer::doEndJob(), art::EDFilter::doEndJob(), art::EDAnalyzer::doEndJob(), and art::RootOutput::endJob().

{
  if (!moduleContext_)
    return;

  // If none of the branches have missing consumes statements, exit early.
  if (std::all_of(cbegin(missingConsumes_),
                  cend(missingConsumes_),
                  [](auto const& perBranch) { return perBranch.empty(); }))
    return;

  constexpr cet::HorizontalRule rule{60};
  mf::LogPrint log{"MTdiagnostics"};
  log << '\n'
      << rule('=') << '\n'
      << "The following consumes (or mayConsume) statements are missing from\n"
      << module_context(moduleDescription_) << '\n'
      << rule('-') << '\n';

  cet::for_all_with_index(
    missingConsumes_, [&log](std::size_t const i, auto const& perBranch) {
      for (auto const& pi : perBranch) {
        log << "  "
            << assemble_consumes_statement(static_cast<BranchType>(i), pi)
            << '\n';
      }
    });
  log << rule('=');
}

void art::Consumer::validateConsumedProduct ( BranchType const  bt, ProductInfo const &  pi  )

protectedinherited

Definition at line 101 of file Consumer.cc.

References art::errors::ProductRegistrationFailure.

{
  // Early exits if consumes tracking has been disabled or if the
  // consumed product is an allowed consumable.
  if (!moduleContext_)
    return;

  if (cet::binary_search_all(consumables_[bt], pi))
    return;

  if (requireConsumes_) {
    throw Exception(errors::ProductRegistrationFailure,
                    "Consumer: an error occurred during validation of a "
                    "retrieved product\n\n")
      << "The following consumes (or mayConsume) statement is missing from\n"
      << module_context(moduleDescription_) << ":\n\n"
      << "  " << assemble_consumes_statement(bt, pi) << "\n\n";
  }

  missingConsumes_[bt].insert(pi);
}

Member Data Documentation

TH1F* evgen::NDKGen::fCCMode

private

CC interaction mode.

Definition at line 112 of file NDKGen_module.cc.

TH1F* evgen::NDKGen::fDCosX

private

direction cosine in x

Definition at line 98 of file NDKGen_module.cc.

TH1F* evgen::NDKGen::fDCosY

private

direction cosine in y

Definition at line 99 of file NDKGen_module.cc.

TH1F* evgen::NDKGen::fDCosZ

private

direction cosine in z

Definition at line 100 of file NDKGen_module.cc.

TH1F* evgen::NDKGen::fECons

private

histogram to determine if energy is conserved in the event

Definition at line 117 of file NDKGen_module.cc.

TH1F* evgen::NDKGen::fEDCosX

private

direction cosine of outgoing e in x

Definition at line 108 of file NDKGen_module.cc.

TH1F* evgen::NDKGen::fEDCosY

private

direction cosine of outgoing e in y

Definition at line 109 of file NDKGen_module.cc.

TH1F* evgen::NDKGen::fEDCosZ

private

direction cosine of outgoing e in z

Definition at line 110 of file NDKGen_module.cc.

TH1F* evgen::NDKGen::fEMomentum

private

momentum of outgoing electrons

Definition at line 107 of file NDKGen_module.cc.

std::ifstream* evgen::NDKGen::fEventFile

private

Definition at line 83 of file NDKGen_module.cc.

TH1F* evgen::NDKGen::fGenerated[6]

private

Spectra as generated.

Definition at line 88 of file NDKGen_module.cc.

TH1F* evgen::NDKGen::fMuDCosX

private

direction cosine of outgoing mu in x

Definition at line 103 of file NDKGen_module.cc.

TH1F* evgen::NDKGen::fMuDCosY

private

direction cosine of outgoing mu in y

Definition at line 104 of file NDKGen_module.cc.

TH1F* evgen::NDKGen::fMuDCosZ

private

direction cosine of outgoing mu in z

Definition at line 105 of file NDKGen_module.cc.

TH1F* evgen::NDKGen::fMuMomentum

private

momentum of outgoing muons

Definition at line 102 of file NDKGen_module.cc.

TH1F* evgen::NDKGen::fNCMode

private

CC interaction mode.

Definition at line 113 of file NDKGen_module.cc.

std::string evgen::NDKGen::fNdkFile

private

Definition at line 82 of file NDKGen_module.cc.

std::string evgen::NDKGen::fNDKModuleLabel

private

keep track of how long it takes to run the job

Definition at line 86 of file NDKGen_module.cc.

TStopwatch evgen::NDKGen::fStopwatch

private

Definition at line 84 of file NDKGen_module.cc.

TH1F* evgen::NDKGen::fVertexX

private

vertex location of generated events in x

Definition at line 90 of file NDKGen_module.cc.

TH2F* evgen::NDKGen::fVertexXY

private

vertex location in xy

Definition at line 94 of file NDKGen_module.cc.

TH2F* evgen::NDKGen::fVertexXZ

private

vertex location in xz

Definition at line 95 of file NDKGen_module.cc.

TH1F* evgen::NDKGen::fVertexY

private

vertex location of generated events in y

Definition at line 91 of file NDKGen_module.cc.

TH2F* evgen::NDKGen::fVertexYZ

private

vertex location in yz

Definition at line 96 of file NDKGen_module.cc.

TH1F* evgen::NDKGen::fVertexZ

private

vertex location of generated events in z

Definition at line 92 of file NDKGen_module.cc.

---

The documentation for this class was generated from the following file:

• larsim/v06_53_00/source/larsim/EventGenerator/NDKGen_module.cc