214 std::cout<<
"------------------------------------------------------------------------------"<<std::endl;
236 std::string
name,
k, dollar;
243 int Idx, Ist,
PDG, Mother1, Mother2, Daughter1 ,Daughter2;
244 double Px, Py, Pz,
E,
m ;
245 std::string
p1,p2,p3,p4,p5,p6,p7,p8,p9,p10,p11,p12,p13,p14;
249 std::string primary(
"primary");
250 int FirstMother = -1;
262 TLorentzVector Neutrino;
263 TLorentzVector Lepton;
264 TLorentzVector Target;
266 TLorentzVector Hadron4mom;
279 std::unique_ptr< std::vector<simb::MCTruth> > truthcol(
new std::vector<simb::MCTruth>);
282 art::ServiceHandle<geo::Geometry const> geo;
285 double const fvCut{5.0};
294 for (
size_t i = 0; i<geo->NTPC(); ++i)
296 double local[3] = {0.,0.,0.};
297 double world[3] = {0.,0.,0.};
300 if (minx>world[0]-geo->DetHalfWidth(i))
301 minx = world[0]-geo->DetHalfWidth(i);
302 if (maxx<world[0]+geo->DetHalfWidth(i))
303 maxx = world[0]+geo->DetHalfWidth(i);
304 if (miny>world[1]-geo->DetHalfHeight(i))
305 miny = world[1]-geo->DetHalfHeight(i);
306 if (maxy<world[1]+geo->DetHalfHeight(i))
307 maxy = world[1]+geo->DetHalfHeight(i);
308 if (minz>world[2]-geo->DetLength(i)/2.)
309 minz = world[2]-geo->DetLength(i)/2.;
310 if (maxz<world[2]+geo->DetLength(i)/2.)
311 maxz = world[2]+geo->DetLength(i)/2.;
315 double X0 = 0.0 + flat.fire( minx+fvCut , maxx-fvCut );
316 double Y0 = 0.0 + flat.fire( miny+fvCut , maxy-fvCut );
317 double Z0 = 0.0 + flat.fire( minz+fvCut , maxz-fvCut );
319 std::cout <<
"NDKGen_module: X, Y, Z of vtx: " << X0 <<
", "<< Y0 <<
", "<< Z0 << std::endl;
324 std::cout <<
"NdkFile: Problem reading Ndk file" << std::endl;
328 if (k.find(
"** Event:")!= std::string::npos) {
329 std::istringstream
in;
333 in>> dummy>> dummy>> dummy >> dummy>> dummy>> dummy>> dummy >> dummy>> dummy>> dummy >> GenieEvt;
337 if (GenieEvt+1 != static_cast<int>(
evt.id().
event()))
341 if (!k.compare(0,25,
"GENIE Interaction Summary"))
343 if (k.compare(0,1,
"|") || k.compare(1,2,
" "))
continue;
344 if (k.find(
"Fin-Init") != std::string::npos)
continue;
345 if (k.find(
"Ar") != std::string::npos)
continue;
346 if (k.find(
"Cl") != std::string::npos)
continue;
347 if (k.find(
"HadrBlob") != std::string::npos)
continue;
348 if (k.find(
"NucBindE") != std::string::npos)
continue;
349 if (k.find(
"FLAGS") != std::string::npos)
break;
350 if (k.find(
"Vertex") != std::string::npos)
break;
353 if (!k.compare(26,1,
"1"))
356 std::istringstream
in;
360 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;
362 if (Ist!=1)
continue;
364 std::cout <<
"PDG = " << PDG << std::endl;
366 const TDatabasePDG* databasePDG = TDatabasePDG::Instance();
367 const TParticlePDG* definition = databasePDG->GetParticle(PDG);
368 Mass = definition->Mass();
369 if (E-Mass < 0.001)
continue;
374 simb::MCParticle mcpart(trackid,
382 P = std::sqrt(pow(E,2.) - pow(Mass,2.));
383 std::cout <<
"Momentum = " << P << std::endl;
385 TLorentzVector pos(X0, Y0, Z0, 0);
389 TLorentzVector mom(Px,Py,Pz, E);
391 mcpart.AddTrajectoryPoint(pos,mom);
396 truth.SetOrigin(simb::kUnknown);
405 std::cout <<
"NDKGen.cxx: Putting " << truth.NParticles() <<
" tracks on stack." << std::endl;
406 truthcol->push_back(truth);
408 evt.put(std::move(truthcol));
process_name stream1 can override from command line with o or output services DetectorPropertiesService services DetectorPropertiesService services DetectorPropertiesService services DetectorPropertiesService physics analyzers pmtresponse NeutronTrackingCut services LArG4Parameters gaussian physics producers generator PDG
process_name opdaq physics producers generator physics producers generator physics producers generator Z0
Geometry information for a single TPC.
process_name opdaq physics producers generator physics producers generator Y0
tuple m
now if test mode generate materials, CRT shell, world, gdml header else just generate CRT shell for u...
CLHEP::HepRandomEngine & fEngine
art-managed random-number engine
standard_singlep gaussian distribution X0
BEGIN_PROLOG vertical distance to the surface Name
if &&[-z"$BASH_VERSION"] then echo Attempting to switch to bash bash shellSwitch exit fi &&["$1"= 'shellSwitch'] shift declare a IncludeDirectives for Dir in
void LocalToWorld(const double *tpc, double *world) const
Transform point from local TPC frame to world frame.
physics associatedGroupsWithLeft p1
BEGIN_PROLOG could also be cout