Jpp
JSirene.cc
Go to the documentation of this file.
1 #include <string>
2 #include <iostream>
3 #include <iomanip>
4 #include <vector>
5 #include <limits>
6 
7 #include "TROOT.h"
8 #include "TFile.h"
9 #include "TH1D.h"
10 #include "TProfile.h"
11 #include "TRandom3.h"
12 
13 #include "km3net-dataformat/offline/Head.hh"
14 #include "km3net-dataformat/offline/Evt.hh"
15 
16 #include "JLang/JSharedPointer.hh"
17 #include "JLang/Jpp.hh"
18 #include "JLang/JPredicate.hh"
19 #include "JSystem/JDate.hh"
20 
21 #include "JPhysics/JPDFToolkit.hh"
22 #include "JPhysics/JCDFTable.hh"
23 #include "JPhysics/JPDFTypes.hh"
24 #include "JPhysics/JRadiationSource.hh"
25 #include "JTools/JFunction1D_t.hh"
26 #include "JTools/JFunctionalMap_t.hh"
27 
28 #include "JAAnet/JHit_t.hh"
29 #include "JAAnet/JTrk_t.hh"
30 #include "JAAnet/JAAnetToolkit.hh"
31 #include "JAAnet/JPDB.hh"
32 
33 #include "JSirene/km3-opa_efrac.hh"
34 //#include "JSirene/pythia.hh"
35 //#include "JSirene/JPythia.hh"
36 #include "JSirene/JSeaWater.hh"
37 #include "JSirene/JCDFTable1D.hh"
38 #include "JSirene/JCDFTable2D.hh"
39 #include "JSirene/JSireneToolkit.hh"
40 
41 #include "JDetector/JDetector.hh"
42 #include "JDetector/JDetectorSubset.hh"
43 #include "JDetector/JDetectorToolkit.hh"
44 #include "JGeometry3D/JCylinder3D.hh"
45 
46 #include "JSupport/JMultipleFileScanner.hh"
47 #include "JSupport/JFileRecorder.hh"
48 #include "JSupport/JMonteCarloFileSupportkit.hh"
49 #include "JSupport/JSupport.hh"
50 #include "JSupport/JMeta.hh"
51 
52 #include "Jeep/JProperties.hh"
53 #include "Jeep/JPrint.hh"
54 #include "Jeep/JTimer.hh"
55 #include "Jeep/JParser.hh"
56 #include "Jeep/JMessage.hh"
57 
58 
59 int debug; //!< debug level
60 
61 namespace {
62 
63  using namespace JPP;
64 
65  typedef JHermiteSplineFunction1D_t JFunction1D_t;
66  typedef JMAPLIST<JPolint1FunctionalMap,
67  JPolint1FunctionalGridMap,
68  JPolint1FunctionalGridMap>::maplist J3DMap_t;
69  typedef JMAPLIST<JPolint1FunctionalMap,
70  JPolint1FunctionalMap,
71  JPolint1FunctionalGridMap,
72  JPolint1FunctionalGridMap>::maplist J4DMap_t;
73 
74  typedef JCDFTable<JFunction1D_t, J3DMap_t> JCDF4D_t; // muon
75  typedef JCDFTable<JFunction1D_t, J4DMap_t> JCDF5D_t; // shower
76 
77  typedef JCDFTable1D<JFunction1D_t::argument_type, JFunction1D_t::result_type> JCDF1D_t; // muon
78  typedef JCDFTable2D<JFunction1D_t::argument_type, JFunction1D_t::result_type> JCDF2D_t; // shower
79 
80 
81  const int NUMBER_OF_BINS = 200; //!< number of bins for CDF integral
82  const double SAFETY_FACTOR = 1.7; //!< safety factor for CDF integral
83 
84 
85  /**
86  * Auxiliary class for CDF tables.
87  */
88  template<class function_type, // CDF function
89  class integral_type> // CDF integral
90  struct JABC {
91  /**
92  * Constructor.
93  *
94  * \param file_descriptor file name descriptor
95  * \param type PDF type
96  */
97  JABC(const std::string& file_descriptor, const JPDFType_t type)
98  {
99  using namespace std;
100 
101  this->type = type;
102 
103  const string file_name = getFilename(file_descriptor, this->type);
104 
105  STATUS("loading input from file " << file_name << "... " << flush);
106 
107  try {
108  function.load(file_name.c_str());
109  }
110  catch(const JException& error) {
111  FATAL(error.what() << endl);
112  }
113 
114  integral = integral_type(function, NUMBER_OF_BINS, SAFETY_FACTOR);
115 
116  STATUS("OK" << endl);
117  }
118 
119  JPDFType_t type;
120  function_type function;
121  integral_type integral;
122  };
123 
124 
125  typedef JABC<JCDF4D_t, JCDF1D_t> JCDF_t; //!< muon light
126  typedef JABC<JCDF5D_t, JCDF2D_t> JCDG_t; //!< shower light
127 }
128 
129 
130 /**
131  * \file
132  *
133  * Main program to simulate detector response to muons and showers.
134  *
135  * \image html sirene.png "Picture by Claudine Colnard"
136  *
137  * Note that CDF file descriptor should contain the wild card character JPHYSICS::WILD_CARD;\n
138  * The file names are obtained by replacing JPHYSICS::WILD_CARD; with
139  * - JPHYSICS::DIRECT_LIGHT_FROM_MUON;
140  * - JPHYSICS::SCATTERED_LIGHT_FROM_MUON;
141  * - JPHYSICS::DIRECT_LIGHT_FROM_DELTARAYS;
142  * - JPHYSICS::SCATTERED_LIGHT_FROM_DELTARAYS;
143  * - JPHYSICS::DIRECT_LIGHT_FROM_EMSHOWER; and
144  * - JPHYSICS::SCATTERED_LIGHT_FROM_EMSHOWER,
145  * respectively.
146  *
147  * The CDF tables can be produced with the script <tt>JMakePDF.sh</tt>:
148  * <pre>
149  * JMakePDF.sh -P
150  * </pre>
151  * (option <tt>-h</tt> will print all available options).
152  * Note that the script will launch a large number of processes (<tt>JMakePDF</tt> and <tt>JMakePDG</tt>)\n
153  * which may take a considerable amount of time to completion.\n
154  * On a standard desktop, all jobs should be finished within 1/2 a day or so.
155  *
156  * The same script should then be run with option <tt>-M</tt> to merge the PDF files, i.e:
157  * <pre>
158  * JMakePDF.sh -M
159  * </pre>
160  *
161  * CDF tables are obtained by running the same script with option <tt>-C</tt>, i.e:
162  * <pre>
163  * JMakePDF.sh -C
164  * </pre>
165  *
166  * The various PDFs can be drawn using the <tt>JDrawPDF</tt> or <tt>JDrawPDG</tt> applications.\n
167  * The tabulated PDFs can be plotted using the <tt>JPlotPDF</tt> or <tt>JPlotPDG</tt> applications.\n
168  * The tabulated CDFs can be plotted using the <tt>JPlotCDF</tt> or <tt>JPlotCDG</tt> applications.
169  * \author mdejong
170  */
171 int main(int argc, char **argv)
172 {
173  using namespace std;
174  using namespace JPP;
175 
176  string fileDescriptor;
177  JMultipleFileScanner<JAAnetTypes_t> inputFile;
178  JFileRecorder <JTYPELIST<JAAnetTypes_t, JMeta, JRootTypes_t>::typelist> outputFile;
179  JLimit_t& numberOfEvents = inputFile.getLimit();
180  string detectorFile;
181  double Ecut_GeV = 0.1; // minimal energy for generation of light from shower [GeV]
182  double Emin_GeV = 1.0; // minimal energy of muon for shower generation [GeV]
183  double Dmin_m = 0.1; // minimal distance for positioning [m]
184  double Tmax_ns = 0.1; // maximal time between hits on same PMT to be merged
185  int Nmax_npe = numeric_limits<int>::max(); // maximal number of photo-electrons
186  bool writeEMShowers;
187  bool keep;
188  size_t numberOfHits;
189  double factor;
190  UInt_t seed;
191 
192  try {
193 
194  JProperties properties;
195 
196  properties.insert(gmake_property(Ecut_GeV));
197  properties.insert(gmake_property(Emin_GeV));
198  properties.insert(gmake_property(Dmin_m));
199  properties.insert(gmake_property(Tmax_ns));
200  properties.insert(gmake_property(Nmax_npe));
201 
202  JParser<> zap("Main program to simulate detector response to muons and showers.");
203 
204  zap['@'] = make_field(properties) = JPARSER::initialised();
205  zap['F'] = make_field(fileDescriptor, "file name descriptor for CDF tables");
206  zap['f'] = make_field(inputFile) = JPARSER::initialised();
207  zap['o'] = make_field(outputFile) = "sirene.root";
208  zap['n'] = make_field(numberOfEvents) = JLimit::max();
209  zap['a'] = make_field(detectorFile) = "";
210  zap['s'] = make_field(writeEMShowers, "store generated EM showers in event");
211  zap['k'] = make_field(keep, "keep position of tracks");
212  zap['N'] = make_field(numberOfHits, "minimum number of hits to output event") = 1;
213  zap['U'] = make_field(factor, "scaling factor applied to light yields") = 1.0;
214  zap['S'] = make_field(seed) = 0;
215  zap['d'] = make_field(debug) = 1;
216 
217  zap(argc, argv);
218  }
219  catch(const exception &error) {
220  FATAL(error.what() << endl);
221  }
222 
223 
224  gRandom->SetSeed(seed);
225 
226 
227  const JMeta meta(argc, argv);
228 
229  const double Zbed = 0.0; // level of seabed [m]
230 
231  vector<JCDF_t> CDF;
232  vector<JCDG_t> CDG;
233 
234  CDF.push_back(JCDF_t(fileDescriptor, DIRECT_LIGHT_FROM_MUON));
235  CDF.push_back(JCDF_t(fileDescriptor, SCATTERED_LIGHT_FROM_MUON));
236  CDF.push_back(JCDF_t(fileDescriptor, DIRECT_LIGHT_FROM_DELTARAYS));
237  CDF.push_back(JCDF_t(fileDescriptor, SCATTERED_LIGHT_FROM_DELTARAYS));
238 
239  CDG.push_back(JCDG_t(fileDescriptor, DIRECT_LIGHT_FROM_EMSHOWER));
240  CDG.push_back(JCDG_t(fileDescriptor, SCATTERED_LIGHT_FROM_EMSHOWER));
241 
242 
243  double maximal_road_width = 0.0; // road width [m]
244  double maximal_distance = 0.0; // road width [m]
245 
246  for (size_t i = 0; i != CDF.size(); ++i) {
247 
248  DEBUG("Range CDF["<< CDF[i].type << "] " << CDF[i].function.intensity.getXmax() << " m" << endl);
249 
250  maximal_road_width = max(maximal_road_width, CDF[i].function.intensity.getXmax());
251  }
252 
253  for (size_t i = 0; i != CDG.size(); ++i) {
254 
255  DEBUG("Range CDG["<< CDG[i].type << "] " << CDG[i].function.intensity.getXmax() << " m" << endl);
256 
257  if (!is_scattered(CDF[i].type)) {
258  maximal_road_width = max(maximal_road_width, CDG[i].function.intensity.getXmax());
259  }
260 
261  maximal_distance = max(maximal_distance, CDG[i].function.intensity.getXmax());
262  }
263 
264  NOTICE("Maximal road width [m] " << maximal_road_width << endl);
265  NOTICE("Maximal distance [m] " << maximal_distance << endl);
266 
267 
268  if (detectorFile == "" || inputFile.empty()) {
269  STATUS("Nothing to be done." << endl);
270  return 0;
271  }
272 
273  JDetector detector;
274 
275  try {
276 
277  STATUS("Load detector... " << flush);
278 
279  load(detectorFile, detector);
280 
281  STATUS("OK" << endl);
282  }
283  catch(const JException& error) {
284  FATAL(error);
285  }
286 
287 
288  vector< JSharedPointer<JRadiationInterface> > radiation;
289 
290  if (true) {
291 
292  STATUS("Setting up radiation tables... " << flush);
293 
294  const JRadiation hydrogen( 1.0, 1.0, 40, 0.01, 0.1, 0.1);
295  const JRadiation oxygen ( 8.0, 16.0, 40, 0.01, 0.1, 0.1);
296  const JRadiation chlorine(17.0, 35.0, 40, 0.01, 0.1, 0.1);
297 
298  JSharedPointer<JRadiation> Hydrogen(new JRadiationFunction(hydrogen, 300, 0.2, 1.0e11));
299  JSharedPointer<JRadiation> Oxygen (new JRadiationFunction(oxygen, 300, 0.2, 1.0e11));
300  JSharedPointer<JRadiation> Chlorine(new JRadiationFunction(chlorine, 300, 0.2, 1.0e11));
301 
302  JRadiationSource::source_type EErad(&JRadiation::TotalCrossSectionEErad, &JRadiation::EfromEErad);
303  JRadiationSource::source_type Brems(&JRadiation::TotalCrossSectionBrems, &JRadiation::EfromBrems);
304  JRadiationSource::source_type GNrad(&JRadiation::TotalCrossSectionGNrad, &JRadiation::EfromGNrad);
305 
306  radiation.push_back(new JRadiationSource(Oxygen, DENSITY_SEA_WATER * JSeaWater::O(), EErad));
307  radiation.push_back(new JRadiationSource(Chlorine, DENSITY_SEA_WATER * JSeaWater::Cl(), EErad));
308  radiation.push_back(new JRadiationSource(Hydrogen, DENSITY_SEA_WATER * JSeaWater::H(), EErad));
309 
310  radiation.push_back(new JRadiationSource(Oxygen, DENSITY_SEA_WATER * JSeaWater::O(), Brems));
311  radiation.push_back(new JRadiationSource(Chlorine, DENSITY_SEA_WATER * JSeaWater::Cl(), Brems));
312  radiation.push_back(new JRadiationSource(Hydrogen, DENSITY_SEA_WATER * JSeaWater::H(), Brems));
313 
314  radiation.push_back(new JRadiationSource(Oxygen, DENSITY_SEA_WATER * JSeaWater::O(), GNrad));
315  radiation.push_back(new JRadiationSource(Chlorine, DENSITY_SEA_WATER * JSeaWater::Cl(), GNrad));
316  radiation.push_back(new JRadiationSource(Hydrogen, DENSITY_SEA_WATER * JSeaWater::H(), GNrad));
317 
318  STATUS("OK" << endl);
319  }
320 
321 
322  Vec center(0,0,0);
323  Head header;
324 
325  try {
326 
327  header = inputFile.getHeader();
328 
329  JHead buffer(header);
330 
331  center = get<Vec>(buffer);
332 
333  buffer.simul.push_back(JAANET::simul());
334 
335  buffer.simul.rbegin()->program = "JSirene";
336  buffer.simul.rbegin()->version = getGITVersion();
337  buffer.simul.rbegin()->date = getDate();
338  buffer.simul.rbegin()->time = getTime();
339 
340  buffer.push(&JHead::simul);
341 
342  if (!keep) {
343 
344  buffer.coord_origin = coord_origin(0.0, 0.0, 0.0);
345  buffer.can.zmin += center.z;
346  buffer.can.zmax += center.z;
347 
348  buffer.push(&JHead::coord_origin);
349  buffer.push(&JHead::can);
350  }
351 
352  const JCircle2D circle(detector.begin(), detector.end());
353 
354  center += Vec(circle.getX(), circle.getY(), 0.0);
355 
356  copy(buffer, header);
357  }
358  catch(const JException& error) {
359  FATAL(error);
360  }
361 
362  if (!keep)
363  NOTICE("Offset applied to true tracks is: " << center << endl);
364  else
365  NOTICE("No offset applied to true tracks." << endl);
366 
367  JCylinder3D cylinder(detector.begin(), detector.end());
368 
369  cylinder.addMargin(maximal_distance);
370 
371  if (cylinder.getZmin() < Zbed) {
372  cylinder.setZmin(Zbed);
373  }
374 
375  NOTICE("Ligth generation volume: " << cylinder << endl);
376 
377  TProfile cpu("cpu", NULL, 14, 1.0, 8.0);
378  TH1D job("job", NULL, 400, 0.5, 400.5);
379 
380 
381  outputFile.open();
382 
383  if (!outputFile.is_open()) {
384  FATAL("Error opening file " << outputFile << endl);
385  }
386 
387  outputFile.put(meta);
388  outputFile.put(header);
389  outputFile.put(*gRandom);
390 
391 
392  JTimer timer;
393 
394  for (JMultipleFileScanner<Evt>& in = inputFile; in.hasNext(); ) {
395 
396  STATUS("event: " << setw(10) << in.getCounter() << '\r'); DEBUG(endl);
397 
398  job.Fill(1.0);
399 
400  Evt* evt = in.next();
401 
402  if (!keep) {
403  for (vector<Trk>::iterator track = evt->mc_trks.begin(); track != evt->mc_trks.end(); ++track) {
404  track->pos += center;
405  }
406  }
407 
408  Evt event(*evt); // output
409 
410  event.mc_hits.clear();
411 
412  JHits_t mc_hits; // temporary buffer
413 
414  timer.reset();
415  timer.start();
416 
417  for (vector<Trk>::const_iterator track = evt->mc_trks.begin(); track != evt->mc_trks.end(); ++track) {
418 
419  if (is_muon(*track)) {
420 
421  // -----------------------------------------------
422  // muon
423  // -----------------------------------------------
424 
425  job.Fill(2.0);
426 
427  typedef JDetectorSubset<JModuleComparator3Z_t> JDetectorSubset_t;
428 
429  const pair<double, double> intersection = cylinder.getIntersection(getAxis(*track));
430 
431  double Zmin = intersection.first;
432  double Zmax = intersection.second;
433 
434  if (Zmax - Zmin <= Dmin_m) {
435  continue;
436  }
437 
438  JVertex vertex(0.0, track->t, track->E); // start of muon
439 
440  if (track->pos.z < Zbed) { // propagate muon through rock
441 
442  if (track->dir.z > 0.0)
443  vertex.step(gRock, (Zbed - track->pos.z) / track->dir.z);
444  else
445  continue;
446  }
447 
448  if (vertex.getZ() < Zmin) { // propagate muon through water
449  vertex.step(gWater, Zmin - vertex.getZ());
450  }
451 
452  if (vertex.getRange() <= Dmin_m) {
453  continue;
454  }
455 
456  job.Fill(3.0);
457 
458  const JDetectorSubset_t subdetector(detector, getAxis(*track), maximal_road_width);
459 
460  if (subdetector.empty()) {
461  continue;
462  }
463 
464  job.Fill(4.0);
465 
466  JTrack muon(vertex); // propagate muon trough detector
467 
468  while (vertex.getE() >= Emin_GeV && vertex.getZ() < Zmax) {
469 
470  const int N = radiation.size();
471 
472  double li[N]; // inverse interaction lengths
473  double ls = 1.0e-5; // minimal total inverse interaction length (100 km)^-1
474 
475  for (int i = 0; i != N; ++i) {
476  ls += li[i] = radiation[i]->getInverseInteractionLength(vertex.getE());
477  }
478 
479  const double ds = min(gRandom->Exp(1.0) / ls, vertex.getRange());
480 
481  vertex.step(ds);
482 
483  if (vertex.getE() >= Emin_GeV) {
484 
485  double Es = 0.0; // shower energy [GeV]
486  double y = gRandom->Uniform(ls);
487 
488  for (int i = 0; i != N; ++i) {
489 
490  y -= li[i];
491 
492  if (y < 0.0) {
493  Es = radiation[i]->getEnergyOfShower(vertex.getE()); // shower energy [GeV]
494  break;
495  }
496  }
497 
498  vertex.applyEloss(Es);
499 
500  if (Es >= Ecut_GeV) {
501  muon.push_back(vertex);
502  }
503  }
504  }
505 
506  // add muon end point
507 
508  if (vertex.getE() < Emin_GeV && vertex.getRange() > 0.0) {
509 
510  vertex.step(vertex.getRange());
511 
512  muon.push_back(vertex);
513  }
514 
515  // add information to output muon
516 
517  vector<Trk>::iterator trk = find_if(event.mc_trks.begin(),
518  event.mc_trks.end(),
519  make_predicate(&Trk::id, track->id));
520 
521  if (trk != event.mc_trks.end()) {
522  trk->len = (muon.rbegin()->getZ() < Zmax ? +1 : -1) * (muon.rbegin()->getZ() - muon.begin()->getZ());
523  trk->setusr(mc_usr_keys::energy_lost_in_can, muon.begin()->getE() - muon.rbegin()->getE());
524  }
525 
526  for (JDetector::const_iterator module = subdetector.begin(); module != subdetector.end(); ++module) {
527 
528  const double z0 = muon.begin()->getZ();
529  const double t0 = muon.begin()->getT();
530  const double R = module->getX();
531  const double Z = module->getZ() - R / getTanThetaC();
532 
533  if (Z >= muon.begin()->getZ() && Z <= muon.rbegin()->getZ()) {
534 
535  for (size_t i = 0; i != CDF.size(); ++i) {
536 
537  if (R < CDF[i].integral.getXmax()) {
538 
539  try {
540 
541  double W = 1.0; // mip
542 
543  if (is_deltarays(CDF[i].type)) {
544  W = getDeltaRaysFromMuon(muon.getE(Z)); // delta-rays
545  }
546 
547  const double NPE = CDF[i].integral.getNPE(R) * module->size() * factor * W;
548  const int N = gRandom->Poisson(NPE);
549 
550  if (N != 0) {
551 
552  double ns = 0.0;
553 
554  for (JModule::const_iterator pmt = module->begin(); pmt != module->end(); ++pmt) {
555 
556  const double R = pmt->getX();
557  const double Z = pmt->getZ() - z0;
558  const double theta = pmt->getTheta();
559  const double phi = fabs(pmt->getPhi());
560 
561  const double npe = CDF[i].function.getNPE(R, theta, phi) * factor * W;
562 
563  ns += npe;
564 
565  int n0 = min(getNumberOfPhotoElectrons(NPE, N, npe), Nmax_npe);
566 
567  job.Fill((double) (100 + CDF[i].type), (double) n0);
568 
569  while (n0 != 0) {
570 
571  const double t1 = CDF[i].function.getTime(R, theta, phi, gRandom->Rndm());
572  const int n1 = getNumberOfPhotoElectrons(n0);
573 
574  mc_hits.push_back(JHit_t(mc_hits.size() + 1,
575  pmt->getID(),
576  getHitType(CDF[i].type),
577  track->id,
578  t0 + (R * getTanThetaC() + Z) / C + t1,
579  n1));
580 
581  n0 -= n1;
582  }
583  }
584 
585  if (ns > NPE) {
586  job.Fill((double) (300 + CDF[i].type));
587  }
588  }
589  }
590  catch(const exception& error) {
591  job.Fill((double) (200 + CDF[i].type));
592  }
593  }
594  }
595  }
596  }
597 
598  for (JTrack::const_iterator vertex = muon.begin(); vertex != muon.end(); ++vertex) {
599 
600  const double z0 = vertex->getZ();
601  const double t0 = vertex->getT();
602  const double Es = vertex->getEs();
603 
604  int origin = track->id;
605 
606  if (writeEMShowers) {
607  origin = event.mc_trks.size() + 1;
608  }
609 
610  int number_of_hits = 0;
611 
612  JDetectorSubset_t::range_type range = subdetector.getRange(z0 - maximal_distance,
613  z0 + maximal_distance);
614 
615  for (JDetector::const_iterator module = range.begin(); module != range.end(); ++module) {
616 
617  const double dz = geanz.getMaximum(Es);
618  const double R = module->getX();
619  const double Z = module->getZ() - z0 - dz;
620  const double D = sqrt(R*R + Z*Z);
621  const double cd = Z / D;
622 
623  for (size_t i = 0; i != CDG.size(); ++i) {
624 
625  if (D < CDG[i].integral.getXmax()) {
626 
627  try {
628 
629  const double NPE = CDG[i].integral.getNPE(D, cd) * Es * module->size() * factor;
630  const int N = gRandom->Poisson(NPE);
631 
632  if (N != 0) {
633 
634  double ns = 0.0;
635 
636  for (JModule::const_iterator pmt = module->begin(); pmt != module->end(); ++pmt) {
637 
638  const double dz = geanz.getMaximum(Es);
639  const double R = pmt->getX();
640  const double Z = pmt->getZ() - z0 - dz;
641  const double theta = pmt->getTheta();
642  const double phi = fabs(pmt->getPhi());
643  const double D = sqrt(R*R + Z*Z);
644  const double cd = Z / D;
645 
646  const double npe = CDG[i].function.getNPE(D, cd, theta, phi) * Es * factor;
647 
648  ns += npe;
649 
650  int n0 = min(getNumberOfPhotoElectrons(NPE, N, npe), Nmax_npe);
651 
652  job.Fill((double) (100 + CDG[i].type), (double) n0);
653 
654  while (n0 != 0) {
655 
656  const double dz = geanz.getLength(Es, gRandom->Rndm());
657  const double Z = pmt->getZ() - z0 - dz;
658  const double D = sqrt(R*R + Z*Z);
659  const double cd = Z / D;
660 
661  const double t1 = CDG[i].function.getTime(D, cd, theta, phi, gRandom->Rndm());
662  const int n1 = getNumberOfPhotoElectrons(n0);
663 
664  mc_hits.push_back(JHit_t(mc_hits.size() + 1,
665  pmt->getID(),
666  getHitType(CDG[i].type),
667  origin,
668  t0 + (dz + D * getIndexOfRefraction()) / C + t1,
669  n1));
670 
671  n0 -= n1;
672 
673  number_of_hits += n1;
674  }
675  }
676 
677  if (ns > NPE) {
678  job.Fill((double) (300 + CDG[i].type));
679  }
680  }
681  }
682  catch(const exception& error) {
683  job.Fill((double) (200 + CDG[i].type));
684  }
685  }
686  }
687  }
688 
689  if (writeEMShowers && number_of_hits != 0) {
690 
691  event.mc_trks.push_back(JTrk_t(origin,
692  TRACK_TYPE_PHOTON,
693  track->id,
694  track->pos + track->dir * vertex->getZ(),
695  track->dir,
696  vertex->getT(),
697  Es));
698  }
699  }
700 
701  } else if (is_tau(*track)) {
702 
703  // -----------------------------------------------
704  // tau
705  // -----------------------------------------------
706 
707  job.Fill(6.0);
708 
709  const double z0 = 0.0;
710  const double z1 = z0 + track->len;
711  const double t0 = track->t;
712  const double E = track->E;
713 
714  const JTransformation3D transformation = getTransformation(*track);
715 
716  JModule buffer;
717 
718  for (JDetector::const_iterator module = detector.begin(); module != detector.end(); ++module) {
719 
720  const JPosition3D pos = transformation.transform(module->getPosition());
721 
722  const double R = pos.getX();
723  const double Z = pos.getZ() - R / getTanThetaC();
724 
725  if (Z < z0 ||
726  Z > z1 ||
727  R > maximal_road_width) {
728  continue;
729  }
730 
731  for (size_t i = 0; i != CDF.size(); ++i) {
732 
733  if (R < CDF[i].integral.getXmax()) {
734 
735  try {
736 
737  double W = 1.0; // mip
738 
739  if (is_deltarays(CDF[i].type)) {
740  W = getDeltaRaysFromTau(E); // delta-rays
741  }
742 
743  const double NPE = CDF[i].integral.getNPE(R) * module->size() * factor * W;
744  const int N = gRandom->Poisson(NPE);
745 
746  if (N != 0) {
747 
748  buffer = *module;
749 
750  buffer.transform(transformation);
751 
752  double ns = 0.0;
753 
754  for (JModule::const_iterator pmt = buffer.begin(); pmt != buffer.end(); ++pmt) {
755 
756  const double R = pmt->getX();
757  const double Z = pmt->getZ() - z0;
758  const double theta = pmt->getTheta();
759  const double phi = fabs(pmt->getPhi());
760 
761  const double npe = CDF[i].function.getNPE(R, theta, phi) * factor * W;
762 
763  ns += npe;
764 
765  int n0 = min(getNumberOfPhotoElectrons(NPE, N, npe), Nmax_npe);
766 
767  job.Fill((double) (120 + CDF[i].type), (double) n0);
768 
769  while (n0 != 0) {
770 
771  const double t1 = CDF[i].function.getTime(R, theta, phi, gRandom->Rndm());
772  const int n1 = getNumberOfPhotoElectrons(n0);
773 
774  mc_hits.push_back(JHit_t(mc_hits.size() + 1,
775  pmt->getID(),
776  getHitType(CDF[i].type),
777  track->id,
778  t0 + (R * getTanThetaC() + Z) / C + t1,
779  n1));
780 
781  n0 -= n1;
782  }
783  }
784 
785  if (ns > NPE) {
786  job.Fill((double) (320 + CDF[i].type));
787  }
788  }
789  }
790  catch(const exception& error) {
791  job.Fill((double) (220 + CDF[i].type));
792  }
793 
794  }
795  }
796  }
797 
798  if (!buffer.empty()) {
799  job.Fill(7.0);
800  }
801 
802  } else if (!is_neutrino(*track)) {
803 
804  if (JPDB::getInstance().hasPDG(track->type)) {
805 
806  // -----------------------------------------------
807  // electron or hadron
808  // -----------------------------------------------
809 
810  job.Fill(8.0);
811 
812  double E = track->E;
813 
814  try {
815  E = getKineticEnergy(E, JPDB::getInstance().getPDG(track->type).mass);
816  }
817  catch(const exception& error) {
818  ERROR(error.what() << endl);
819  }
820 
821  E = pythia(track->type, E);
822 
823  if (E < Ecut_GeV || cylinder.getDistance(getPosition(*track)) > Dmin_m) {
824  continue;
825  }
826 
827  const double z0 = 0.0;
828  const double t0 = track->t;
829 
830  const JTransformation3D transformation = getTransformation(*track);
831 
832  JModule buffer;
833 
834  for (JDetector::const_iterator module = detector.begin(); module != detector.end(); ++module) {
835 
836  const JPosition3D pos = transformation.transform(module->getPosition());
837 
838  const double dz = geanz.getMaximum(E);
839  const double R = pos.getX();
840  const double Z = pos.getZ() - z0 - dz;
841  const double D = sqrt(R*R + Z*Z);
842  const double cd = Z / D;
843 
844  if (D > maximal_distance) {
845  continue;
846  }
847 
848  for (size_t i = 0; i != CDG.size(); ++i) {
849 
850  if (D < CDG[i].integral.getXmax()) {
851 
852  try {
853 
854  const double NPE = CDG[i].integral.getNPE(D, cd) * E * module->size() * factor;
855  const int N = gRandom->Poisson(NPE);
856 
857  if (N != 0) {
858 
859  buffer = *module;
860 
861  buffer.transform(transformation);
862 
863  double ns = 0.0;
864 
865  for (JModule::const_iterator pmt = buffer.begin(); pmt != buffer.end(); ++pmt) {
866 
867  const double dz = geanz.getMaximum(E);
868  const double R = pmt->getX();
869  const double Z = pmt->getZ() - z0 - dz;
870  const double theta = pmt->getTheta();
871  const double phi = fabs(pmt->getPhi());
872  const double D = sqrt(R*R + Z*Z);
873  const double cd = Z / D;
874 
875  const double npe = CDG[i].function.getNPE(D, cd, theta, phi) * E * factor;
876 
877  ns += npe;
878 
879  int n0 = min(getNumberOfPhotoElectrons(NPE, N, npe), Nmax_npe);
880 
881  job.Fill((double) (140 + CDG[i].type), (double) n0);
882 
883  while (n0 != 0) {
884 
885  const double dz = geanz.getLength(E, gRandom->Rndm());
886  const double Z = pmt->getZ() - z0 - dz;
887  const double D = sqrt(R*R + Z*Z);
888  const double cd = Z / D;
889 
890  const double t1 = CDG[i].function.getTime(D, cd, theta, phi, gRandom->Rndm());
891  const int n1 = getNumberOfPhotoElectrons(n0);
892 
893  mc_hits.push_back(JHit_t(mc_hits.size() + 1,
894  pmt->getID(),
895  getHitType(CDG[i].type, true),
896  track->id,
897  t0 + (dz + D * getIndexOfRefraction()) / C + t1,
898  n1));
899 
900  n0 -= n1;
901  }
902  }
903 
904  if (ns > NPE) {
905  job.Fill((double) (340 + CDG[i].type));
906  }
907  }
908  }
909  catch(const exception& error) {
910  job.Fill((double) (240 + CDG[i].type));
911  }
912  }
913  }
914  }
915 
916  if (!buffer.empty()) {
917  job.Fill(9.0);
918  }
919 
920  } else {
921 
922  job.Fill(21.0);
923  }
924  }
925  }
926 
927  if (!mc_hits.empty()) {
928 
929  mc_hits.merge(Tmax_ns);
930 
931  event.mc_hits.resize(mc_hits.size());
932 
933  copy(mc_hits.begin(), mc_hits.end(), event.mc_hits.begin());
934  }
935 
936  timer.stop();
937 
938  if (has_neutrino(event)) {
939  cpu.Fill(log10(get_neutrino(event).E), (double) timer.usec_ucpu * 1.0e-3);
940  }
941 
942  if (event.mc_hits.size() >= numberOfHits) {
943 
944  outputFile.put(event);
945 
946  job.Fill(10.0);
947  }
948  }
949  STATUS(endl);
950 
951  outputFile.put(cpu);
952  outputFile.put(job);
953  outputFile.put(*gRandom);
954  outputFile.close();
955 }
JAANET::JHead::simul
std::vector< JAANET::simul > simul
Definition: JHead.hh:1077
JPHYSICS::DIRECT_LIGHT_FROM_DELTARAYS
direct light from delta-rays
Definition: JPDFTypes.hh:35
JMeta.hh
JTOOLS::JPolint1FunctionalGridMap
Type definition of a 1st degree polynomial interpolation based on a JGridMap implementation.
Definition: JFunctionalMap_t.hh:91
JCDFTable1D.hh
Head.hh
JPredicate.hh
JSIRENE::JVertex::step
JVertex & step(const double ds)
Step.
Definition: JSireneToolkit.hh:266
Vec::z
double z
Definition: Vec.hh:14
Evt::mc_trks
std::vector< Trk > mc_trks
MC: list of MC truth tracks
Definition: Evt.hh:45
JAANET::JHit_t
Auxiliary class to set-up Hit.
Definition: JHit_t.hh:25
JDETECTOR::JModule::transform
void transform(const JRotation3D &R, const JVector3D &pos)
Transformation of geometry (see method JGEOMETRY3D::JPosition3D::transform(const JRotation3D&,...
Definition: JModule.hh:261
JEEP::JTimer::reset
void reset()
Reset timer.
Definition: JTimer.hh:76
JPHYSICS::JRadiationFunction
Fast implementation of class JRadiation.
Definition: JRadiationSource.hh:35
JFileRecorder.hh
JEEP::JTimer
Auxiliary class for CPU timing and usage.
Definition: JTimer.hh:32
JSUPPORT::JLimit
Auxiliary class for defining the range of iterations of objects.
Definition: JLimit.hh:41
JGEOMETRY3D::JTransformation3D
Transformation.
Definition: JTransformation3D.hh:24
JPHYSICS::JCDFTable
Multi-dimensional CDF table for arrival time of Cherenkov light.
Definition: JCDFTable.hh:51
JSIRENE::JVertex
Vertex of energy loss of muon.
Definition: JSireneToolkit.hh:146
JMessage.hh
JSIRENE::JVertex::getEs
double getEs() const
Get shower energy.
Definition: JSireneToolkit.hh:214
JPHYSICS::gWater
static const JGeaneWater gWater
Function object for energy loss of muon in sea water.
Definition: JGeane.hh:328
JPrint.hh
JAANET::getAxis
JAxis3D getAxis(const Trk &track)
Get axis.
Definition: JAAnetToolkit.hh:244
JTOOLS::DENSITY_SEA_WATER
static const double DENSITY_SEA_WATER
Fixed environment values.
Definition: JConstants.hh:34
JAANET::JHead::coord_origin
JAANET::coord_origin coord_origin
Definition: JHead.hh:1086
mc_usr_keys::energy_lost_in_can
const char *const energy_lost_in_can
Definition: io_ascii.hh:42
JSIRENE::getHitType
JHitType_t getHitType(const JPDFType_t pdf, const bool shower=false)
Get hit type corresponding to given PDF type.
Definition: JSireneToolkit.hh:112
JAANET::JPDB::getInstance
static const JPDB & getInstance()
Get particle data book.
Definition: JPDB.hh:121
JPHYSICS::is_deltarays
bool is_deltarays(const int pdf)
Test if given PDF type corresponds to Cherenkov light from delta-rays.
Definition: JPDFTypes.hh:172
JGEOMETRY2D::JCircle2D
Data structure for circle in two dimensions.
Definition: JCircle2D.hh:29
JGEOMETRY3D::JVector3D::getZ
double getZ() const
Get z position.
Definition: JVector3D.hh:114
JEEP::JTimer::start
void start()
Start timer.
Definition: JTimer.hh:89
JPARSER::initialised
Empty structure for specification of parser element that is initialised (i.e.
Definition: JParser.hh:63
JDETECTOR::load
void load(const JString &file_name, JDetector &detector)
Load detector from input file.
Definition: JDetectorToolkit.hh:476
JCDFTable.hh
JAANET::is_tau
bool is_tau(const Trk &track)
Test whether given track is a (anti-)tau.
Definition: JAAnetToolkit.hh:375
JSharedPointer.hh
JPHYSICS::DIRECT_LIGHT_FROM_EMSHOWER
direct light from EM shower
Definition: JPDFTypes.hh:40
JGEOMETRY3D::JTransformation3D::transform
JPosition3D transform(const JPosition3D &pos) const
Transform position.
Definition: JTransformation3D.hh:49
JPDFToolkit.hh
std::vector
Definition: JSTDTypes.hh:12
debug
int debug
debug level
Definition: JSirene.cc:59
JSIRENE::pythia
static const JPythia pythia
Function object for relative light yield as a function of GEANT particle code.
Definition: JPythia.hh:96
main
int main(int argc, char **argv)
Definition: JSirene.cc:171
Evt
The Evt class respresent a Monte Carlo (MC) event as well as an offline event.
Definition: Evt.hh:19
JSIRENE::JSeaWater::Cl
static double Cl()
Definition: JSeaWater.hh:28
Evt.hh
JPARSER::JParser
Utility class to parse command line options.
Definition: JParser.hh:1493
JTOOLS::JMAPLIST
Auxiliary class for recursive map list generation.
Definition: JMapList.hh:108
JPHYSICS::getDeltaRaysFromMuon
double getDeltaRaysFromMuon(const double E)
Equivalent EM-shower energy due to delta-rays per unit muon track length.
Definition: JPDFToolkit.hh:81
JAANET::getTransformation
JTransformation3D getTransformation(const Trk &track)
Get transformation.
Definition: JAAnetToolkit.hh:280
JDate.hh
JTOOLS::JHermiteSplineFunction1D_t
Type definition of a spline interpolation method based on a JCollection with double result type.
Definition: JFunction1D_t.hh:92
NOTICE
#define NOTICE(A)
Definition: JMessage.hh:64
JAANET::coord_origin
Coordinate origin.
Definition: JHead.hh:489
JPHYSICS::JCDFTable2D
Custom class for CDF table in 2 dimensions.
Definition: JCDFTable2D.hh:39
JTOOLS::C
static const double C
Speed of light in vacuum [m/ns].
Definition: JConstants.hh:22
JHit_t.hh
JGEOMETRY3D::JCylinder3D::addMargin
void addMargin(const double D)
Add (safety) margin.
Definition: JCylinder3D.hh:163
JAANET::JHits_t::merge
void merge(const double Tmax_ns)
Merge hits on same PMT that are within given time window.
Definition: JHit_t.hh:118
JAANET::copy
void copy(const Head &from, JHead &to)
Copy header from from to to.
Definition: JHead.cc:152
JAANET::is_muon
bool is_muon(const Trk &track)
Test whether given track is a (anti-)muon.
Definition: JAAnetToolkit.hh:367
JAANET::JHead
Monte Carlo run header.
Definition: JHead.hh:839
JCDFTable2D.hh
JPHYSICS::JRadiation::EfromBrems
double EfromBrems(const double E) const
Bremsstrahlung shower energy.
Definition: JRadiation.hh:222
JAANET::can::zmin
double zmin
Bottom [m].
Definition: JHead.hh:369
JSupport.hh
JPP
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
Definition: JAAnetToolkit.hh:37
ERROR
#define ERROR(A)
Definition: JMessage.hh:66
JFunction1D_t.hh
JPHYSICS::JGeanz::getMaximum
double getMaximum(const double E) const
Get depth of shower maximum.
Definition: JGeanz.hh:162
JGEOMETRY3D::JCylinder3D
Cylinder object.
Definition: JCylinder3D.hh:37
JPHYSICS::gRock
static const JGeane_t gRock(2.67e-1 *0.9 *DENSITY_ROCK, 3.40e-4 *1.2 *DENSITY_ROCK)
Function object for energy loss of muon in rock.
JSIRENE::getKineticEnergy
double getKineticEnergy(const double E, const double m)
Get kinetic energy of particle with given mass.
Definition: JSireneToolkit.hh:348
JDETECTOR::JDetectorSubset
Detector subset with binary search functionality.
Definition: JDetectorSubset.hh:118
JTOOLS::JPolint1FunctionalMap
Type definition of a 1st degree polynomial interpolation based on a JMap implementation.
Definition: JFunctionalMap_t.hh:55
JGEOMETRY3D::JPosition3D
Data structure for position in three dimensions.
Definition: JPosition3D.hh:35
JPHYSICS::JRadiation::TotalCrossSectionBrems
double TotalCrossSectionBrems(const double E) const
Bremsstrahlung cross section.
Definition: JRadiation.hh:145
Head
The Head class reflects the header of Monte-Carlo event files, which consists of keys (also referred ...
Definition: Head.hh:39
JTOOLS::getIndexOfRefraction
double getIndexOfRefraction()
Get average index of refraction of water.
Definition: JConstants.hh:111
JFunctionalMap_t.hh
JAANET::JTrk_t
Auxiliary class to set-up Trk.
Definition: JTrk_t.hh:19
JAANET::get_neutrino
const Trk & get_neutrino(const Evt &evt)
Get incoming neutrino.
Definition: JAAnetToolkit.hh:438
JAAnetToolkit.hh
JPHYSICS::SCATTERED_LIGHT_FROM_EMSHOWER
scattered light from EM shower
Definition: JPDFTypes.hh:41
JPHYSICS::getDeltaRaysFromTau
double getDeltaRaysFromTau(const double E)
Equivalent EM-shower energy due to delta-rays per unit tau track length.
Definition: JPDFToolkit.hh:109
JAANET::has_neutrino
bool has_neutrino(const Evt &evt)
Test whether given event has an incoming neutrino.
Definition: JAAnetToolkit.hh:427
JSIRENE::JVertex::applyEloss
void applyEloss(const double Es)
Apply energy loss energy.
Definition: JSireneToolkit.hh:252
JAANET::getTime
double getTime(const Hit &hit)
Get true time of hit.
Definition: JAAnetToolkit.hh:87
JSYSTEM::ls
Auxiliary class to list files.
Definition: JSystemToolkit.hh:300
JDETECTOR::JModule
Data structure for a composite optical module.
Definition: JModule.hh:49
JEEP::getFilename
std::string getFilename(const std::string &file_name)
Get file name part, i.e.
Definition: JeepToolkit.hh:86
JLANG::getGITVersion
const char * getGITVersion()
Get GIT version.
Definition: Jpp.hh:58
JMultipleFileScanner.hh
JEEP::JTimer::usec_ucpu
unsigned long long usec_ucpu
Definition: JTimer.hh:225
std::pair
Definition: JSTDTypes.hh:15
STATUS
#define STATUS(A)
Definition: JMessage.hh:63
JSUPPORT::JMultipleFileScanner::hasNext
virtual bool hasNext()
Check availability of next element.
Definition: JMultipleFileScanner.hh:350
JAANET::can::zmax
double zmax
Top [m].
Definition: JHead.hh:370
Jpp.hh
km3-opa_efrac.hh
JParser.hh
JDetectorToolkit.hh
JSIRENE::JVertex::getZ
double getZ() const
Get position.
Definition: JSireneToolkit.hh:181
JPHYSICS::JPDFType_t
JPDFType_t
PDF types.
Definition: JPDFTypes.hh:27
JTOOLS::getTanThetaC
double getTanThetaC()
Get average tangent of Cherenkov angle of water.
Definition: JConstants.hh:133
JTrk_t.hh
JPHYSICS::is_scattered
bool is_scattered(const int pdf)
Test if given PDF type corresponds to scattered light.
Definition: JPDFTypes.hh:186
JAANET::JHead::can
JAANET::can can
Definition: JHead.hh:1083
JAANET::JHead::push
void push(T JHead::*pd)
Push given data member to Head.
Definition: JHead.hh:944
JSireneToolkit.hh
make_field
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1954
gmake_property
#define gmake_property(A)
macro to convert (template) parameter to JPropertiesElement object
Definition: JProperties.hh:1319
JSUPPORT::JMultipleFileScanner
General purpose class for object reading from a list of file names.
Definition: JMultipleFileScanner.hh:167
JDETECTOR::JDetector
Detector data structure.
Definition: JDetector.hh:80
JAANET::is_neutrino
bool is_neutrino(const Trk &track)
Test whether given track is a neutrino.
Definition: JAAnetToolkit.hh:349
JEEP::JProperties
Utility class to parse parameter values.
Definition: JProperties.hh:496
JAANET::getPosition
JPosition3D getPosition(const Vec &v)
Get position.
Definition: JAAnetToolkit.hh:197
Trk::id
int id
track identifier
Definition: Trk.hh:14
JAANET::detector
Detector file.
Definition: JHead.hh:130
JPHYSICS::JGeanz::getLength
double getLength(const double E, const double P, const double eps=1.0e-3) const
Get shower length for a given integrated probability.
Definition: JGeanz.hh:122
JSUPPORT::JMeta
Auxiliary class for ROOT I/O of application specific meta data.
Definition: JMeta.hh:71
DEBUG
#define DEBUG(A)
Message macros.
Definition: JMessage.hh:62
JPHYSICS::JRadiation::EfromEErad
double EfromEErad(const double E) const
Pair production shower energy.
Definition: JRadiation.hh:242
JSIRENE::getNumberOfPhotoElectrons
int getNumberOfPhotoElectrons(const double NPE, const int N, const double npe)
Get number of photo-electrons of a hit given the expectation values of the number of photo-electrons ...
Definition: JSireneToolkit.hh:54
std
Definition: jaanetDictionary.h:36
JPHYSICS::JRadiation::TotalCrossSectionEErad
virtual double TotalCrossSectionEErad(const double E) const
Pair production cross section.
Definition: JRadiation.hh:118
JCylinder3D.hh
Head::coord_origin
Vec coord_origin() const
Get coordinate origin.
Definition: Head.hh:290
JLANG::make_predicate
JPredicate< JResult_t T::*, JComparison::eq > make_predicate(JResult_t T::*member, const JResult_t value)
Helper method to create predicate for data member.
Definition: JPredicate.hh:128
JPDB.hh
JSYSTEM::getDate
const char * getDate()
Get ASCII formatted date.
Definition: JDate.hh:23
JSIRENE::JVertex::getE
double getE() const
Get muon energy.
Definition: JSireneToolkit.hh:203
JEEP::JTimer::stop
void stop()
Stop timer.
Definition: JTimer.hh:113
JTimer.hh
JSIRENE::JTrack::getE
double getE(const double z) const
Get muon energy at given position.
Definition: JSireneToolkit.hh:322
Vec
The Vec class is a straightforward 3-d vector, which also works in pyroot.
Definition: Vec.hh:12
JPHYSICS::JRadiation::EfromGNrad
double EfromGNrad(const double E) const
Photo-nuclear shower energy.
Definition: JRadiation.hh:266
JProperties.hh
JSIRENE::JSeaWater::O
static double O()
Estimated mass fractions of chemical elements in sea water.
Definition: JSeaWater.hh:25
JMonteCarloFileSupportkit.hh
JRadiationSource.hh
JAANET::TRACK_TYPE_PHOTON
Definition: JParticleTypes.hh:70
JLANG::JException::what
virtual const char * what() const
Get error message.
Definition: JException.hh:48
JPHYSICS::DIRECT_LIGHT_FROM_MUON
direct light from muon
Definition: JPDFTypes.hh:29
JSIRENE::JVertex::getT
double getT() const
Get time.
Definition: JSireneToolkit.hh:192
JDetector.hh
JGEOMETRY3D::JVector3D::getX
double getX() const
Get x position.
Definition: JVector3D.hh:93
JDetectorSubset.hh
JSIRENE::JVertex::getRange
double getRange() const
Get range of muon.
Definition: JSireneToolkit.hh:226
Evt::mc_hits
std::vector< Hit > mc_hits
MC: list of MC truth hits.
Definition: Evt.hh:44
FATAL
#define FATAL(A)
Definition: JMessage.hh:67
outputFile
string outputFile
Definition: JDAQTimesliceSelector.cc:37
JSIRENE::JTrack
Muon trajectory.
Definition: JSireneToolkit.hh:303
JPHYSICS::SCATTERED_LIGHT_FROM_DELTARAYS
scattered light from delta-rays
Definition: JPDFTypes.hh:36
JSUPPORT::JFileRecorder
Object writing to file.
Definition: JFileRecorder.hh:41
JPHYSICS::JRadiationSource
Source of radiation.
Definition: JRadiationSource.hh:231
JSeaWater.hh
JPDFTypes.hh
JAANET::JHits_t
Auxiliary data structure for list of hits with hit merging capability.
Definition: JHit_t.hh:107
JAANET::simul
Generator for simulation.
Definition: JHead.hh:308
JLANG::JException
General exception.
Definition: JException.hh:23
JDETECTOR::JDetectorSubset_t
Detector subset without binary search functionality.
Definition: JDetectorSubset.hh:37
JPHYSICS::JRadiation
Auxiliary class for the calculation of the muon radiative cross sections.
Definition: JRadiation.hh:39
JSIRENE::JSeaWater::H
static double H()
Definition: JSeaWater.hh:26
JSUPPORT::JLimit::max
static counter_type max()
Get maximum counter value.
Definition: JLimit.hh:117
JLANG::JSharedPointer
The template JSharedPointer class can be used to share a pointer to an object.
Definition: JSharedPointer.hh:28
JPHYSICS::JCDFTable1D
Custom class for CDF table in 1 dimension.
Definition: JCDFTable1D.hh:38
JPHYSICS::JRadiation::TotalCrossSectionGNrad
virtual double TotalCrossSectionGNrad(const double E) const
Photo-nuclear cross section.
Definition: JRadiation.hh:197
JPHYSICS::geanz
static const JGeanz geanz(1.85, 0.62, 0.54)
Function object for longitudinal EM-shower profile.
JPHYSICS::SCATTERED_LIGHT_FROM_MUON
scattered light from muon
Definition: JPDFTypes.hh:30
Generated by   doxygen 1.8.16