23 #include "Framework/Algorithm/AlgFactory.h" 24 #include "Framework/EventGen/EventRecord.h" 25 #include "Framework/EventGen/EventRecordVisitorI.h" 26 #include "Framework/GHEP/GHepParticle.h" 27 #include "Framework/ParticleData/PDGLibrary.h" 28 #include "Framework/Utils/AppInit.h" 29 #include "Physics/NucleonDecay/NucleonDecayMode.h" 30 #include "Physics/NucleonDecay/NucleonDecayUtils.h" 46 #include "CLHEP/Random/RandFlat.h" 72 const genie::EventRecordVisitorI*
mcgen;
85 string sname =
"genie::EventGenerator";
86 string sconfig =
"NucleonDecay";
90 genie::PDGLibrary::Instance();
92 genie::AlgFactory* algf = genie::AlgFactory::Instance();
93 mcgen =
dynamic_cast<const genie::EventRecordVisitorI*
>(algf->GetAlgorithm(sname, sconfig));
95 throw cet::exception(
"NucleonDecay") <<
"Couldn't instantiate the nucleon decay generator";
97 int fDecayMode = p.get<
int>(
"DecayMode");
100 if (p.get<
int>(
"DecayedNucleon") > 0) {
dpdg = p.get<
int>(
"DecayedNucleon"); }
102 dpdg = genie::utils::nucleon_decay::DecayedNucleonPdgCode(gOptDecayMode);
105 produces<std::vector<simb::MCTruth>>();
106 produces<sumdata::RunData, art::InRun>();
109 genie::utils::app_init::RandGen(seed);
115 genie::EventRecord*
event =
new genie::EventRecord;
118 genie::Interaction*
interaction = genie::Interaction::NDecay(target, decay,
dpdg);
119 event->AttachSummary(interaction);
131 std::unique_ptr<std::vector<simb::MCTruth>> truthcol(
new std::vector<simb::MCTruth>);
144 if (minx > tpc.MinX()) minx = tpc.MinX();
145 if (maxx < tpc.MaxX()) maxx = tpc.MaxX();
146 if (miny > tpc.MinY()) miny = tpc.MinY();
147 if (maxy < tpc.MaxY()) maxy = tpc.MaxY();
148 if (minz > tpc.MinZ()) minz = tpc.MinZ();
149 if (maxz < tpc.MaxZ()) maxz = tpc.MaxZ();
153 double X0 =
flatDist.fire(minx, maxx);
154 double Y0 =
flatDist.fire(miny, maxy);
155 double Z0 =
flatDist.fire(minz, maxz);
157 TIter partitr(
event);
158 genie::GHepParticle*
part = 0;
163 std::string primary(
"primary");
165 while ((part = dynamic_cast<genie::GHepParticle*>(partitr.Next()))) {
168 trackid, part->Pdg(), primary, part->FirstMother(), part->Mass(), part->Status());
170 TLorentzVector pos(X0, Y0, Z0, 0);
171 TLorentzVector mom(part->Px(), part->Py(), part->Pz(), part->E());
172 tpart.AddTrajectoryPoint(pos, mom);
173 if (part->PolzIsSet()) {
175 part->GetPolarization(polz);
176 tpart.SetPolarization(polz);
178 tpart.SetRescatter(part->RescatterCode());
184 truthcol->push_back(truth);
186 e.
put(std::move(truthcol));
details::range_type< T > Iterate() const
Initializes the specified ID with the ID of the first cryostat.
base_engine_t & createEngine(seed_t seed)
genie::NucleonDecayMode_t gOptDecayMode
NucleonDecay(fhicl::ParameterSet const &p)
void SetOrigin(simb::Origin_t origin)
const genie::EventRecordVisitorI * mcgen
EDProducer(fhicl::ParameterSet const &pset)
Geometry information for a single TPC.
PutHandle< PROD > put(std::unique_ptr< PROD > &&edp, std::string const &instance={})
Functions for transforming GENIE objects into ART objects (and back)
void beginRun(art::Run &run) override
PutHandle< PROD > put(std::unique_ptr< PROD > &&edp, std::string const &instance={})
#define DEFINE_ART_MODULE(klass)
void produce(art::Event &e) override
An art service to assist in the distribution of guaranteed unique seeds to all engines within an art ...
NucleonDecay & operator=(NucleonDecay const &)=delete
void Add(simb::MCParticle const &part)
void SetEventGeneratorListAndTune(const std::string &evtlistname="", const std::string &tunename="${GENIE_XSEC_TUNE}")
cout<< "-> Edep in the target
std::string const & DetectorName() const
Returns a string with the name of the detector, as configured.
Event generator information.
Namespace collecting geometry-related classes utilities.
Event Generation using GENIE, cosmics or single particles.
art framework interface to geometry description
cet::coded_exception< error, detail::translate > exception
Event finding and building.
The data type to uniquely identify a cryostat.