Jpp - the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Functions | Variables
JSirene.cc File Reference

Main program to simulate detector response to muons and showers. More...

#include <string>
#include <iostream>
#include <iomanip>
#include <vector>
#include <limits>
#include "TROOT.h"
#include "TFile.h"
#include "TH1D.h"
#include "TProfile.h"
#include "TRandom3.h"
#include "km3net-dataformat/offline/Head.hh"
#include "km3net-dataformat/offline/Evt.hh"
#include "km3net-dataformat/definitions/applications.hh"
#include "km3net-dataformat/definitions/trkmembers.hh"
#include "JLang/JSharedPointer.hh"
#include "JLang/Jpp.hh"
#include "JLang/JPredicate.hh"
#include "JSystem/JDate.hh"
#include "JPhysics/JCDFTable.hh"
#include "JPhysics/JPDFTypes.hh"
#include "JPhysics/JPDFToolkit.hh"
#include "JPhysics/JGeane.hh"
#include "JPhysics/JGeanz.hh"
#include "JPhysics/JRadiationSource.hh"
#include "JTools/JFunction1D_t.hh"
#include "JTools/JFunctionalMap_t.hh"
#include "JAAnet/JHit_t.hh"
#include "JAAnet/JTrk_t.hh"
#include "JAAnet/JAAnetToolkit.hh"
#include "JAAnet/JPDB.hh"
#include "JSirene/km3-opa_efrac.hh"
#include "JSirene/JSeaWater.hh"
#include "JSirene/JCDFTable1D.hh"
#include "JSirene/JCDFTable2D.hh"
#include "JSirene/JSireneToolkit.hh"
#include "JDetector/JDetector.hh"
#include "JDetector/JDetectorSubset.hh"
#include "JDetector/JDetectorToolkit.hh"
#include "JGeometry3D/JCylinder3D.hh"
#include "JSupport/JMultipleFileScanner.hh"
#include "JSupport/JFileRecorder.hh"
#include "JSupport/JMonteCarloFileSupportkit.hh"
#include "JSupport/JSupport.hh"
#include "JSupport/JMeta.hh"
#include "Jeep/JProperties.hh"
#include "Jeep/JPrint.hh"
#include "Jeep/JTimer.hh"
#include "Jeep/JParser.hh"
#include "Jeep/JMessage.hh"

Go to the source code of this file.

Functions

int main (int argc, char **argv)
 

Variables

int debug
 debug level More...
 
int numberOfBins = 200
 number of bins for average CDF integral of optical module More...
 
double safetyFactor = 1.7
 safety factor for average CDF integral of optical module More...
 

Detailed Description

Main program to simulate detector response to muons and showers.

sirene.png
Picture by Claudine Colnard

Note that CDF file descriptor should contain the wild card character JPHYSICS::WILD_CARD;
The file names are obtained by replacing JPHYSICS::WILD_CARD; with

The CDF tables can be produced with the script JMakePDF.sh: 

JMakePDF.sh -P

(option -h will print all available options). Note that the script will launch a large number of processes (JMakePDF and JMakePDG)
which may take a considerable amount of time to completion.
On a standard desktop, all jobs should be finished within 1/2 a day or so.

The same script should then be run with option -M to merge the PDF files, i.e: 

JMakePDF.sh -M

CDF tables are obtained by running the same script with option -C, i.e: 

JMakePDF.sh -C

The various PDFs can be drawn using the JDrawPDF or JDrawPDG applications.
The tabulated PDFs can be plotted using the JPlotPDF or JPlotPDG applications.
The tabulated CDFs can be plotted using the JPlotCDF or JPlotCDG applications.

Author
mdejong

Definition in file JSirene.cc.

Function Documentation

int main ( int  argc,
char **  argv 
)

Definition at line 175 of file JSirene.cc.

176 {
177  using namespace std;
178  using namespace JPP;
179 
180  string fileDescriptor;
181  JMultipleFileScanner<JAAnetTypes_t> inputFile;
182  JFileRecorder <JTYPELIST<JAAnetTypes_t, JMeta, JRootTypes_t>::typelist> outputFile;
183  JLimit_t& numberOfEvents = inputFile.getLimit();
184  string detectorFile;
185  double Ecut_GeV = 0.1; // minimal energy for generation of light from shower [GeV]
186  double Emin_GeV = 1.0; // minimal energy of muon for shower generation [GeV]
187  double Dmin_m = 0.1; // minimal distance for positioning [m]
188  double Tmax_ns = 0.1; // maximal time between hits on same PMT to be merged
189  int Nmax_npe = numeric_limits<int>::max(); // maximal number of photo-electrons
190  bool writeEMShowers;
191  bool keep;
192  size_t numberOfHits;
193  double factor;
194  UInt_t seed;
195 
196  try {
197 
198  JProperties properties;
199 
200  properties.insert(gmake_property(Ecut_GeV));
201  properties.insert(gmake_property(Emin_GeV));
202  properties.insert(gmake_property(Dmin_m));
203  properties.insert(gmake_property(Tmax_ns));
204  properties.insert(gmake_property(Nmax_npe));
205  properties.insert(gmake_property(numberOfBins));
206  properties.insert(gmake_property(safetyFactor));
207 
208  JParser<> zap("Main program to simulate detector response to muons and showers.");
209 
210  zap['@'] = make_field(properties) = JPARSER::initialised();
211  zap['F'] = make_field(fileDescriptor, "file name descriptor for CDF tables");
212  zap['f'] = make_field(inputFile) = JPARSER::initialised();
213  zap['o'] = make_field(outputFile) = "sirene.root";
214  zap['n'] = make_field(numberOfEvents) = JLimit::max();
215  zap['a'] = make_field(detectorFile) = "";
216  zap['s'] = make_field(writeEMShowers, "store generated EM showers in event");
217  zap['k'] = make_field(keep, "keep position of tracks");
218  zap['N'] = make_field(numberOfHits, "minimum number of hits to output event") = 1;
219  zap['U'] = make_field(factor, "scaling factor applied to light yields") = 1.0;
220  zap['S'] = make_field(seed) = 0;
221  zap['d'] = make_field(debug) = 1;
222 
223  zap(argc, argv);
224  }
225  catch(const exception &error) {
226  FATAL(error.what() << endl);
227  }
228 
229 
230  gRandom->SetSeed(seed);
231 
232 
233  const JMeta meta(argc, argv);
234 
235  const double Zbed = 0.0; // level of seabed [m]
236 
237  vector<JCDF_t> CDF;
238  vector<JCDG_t> CDG;
239 
240  CDF.push_back(JCDF_t(fileDescriptor, DIRECT_LIGHT_FROM_MUON));
241  CDF.push_back(JCDF_t(fileDescriptor, SCATTERED_LIGHT_FROM_MUON));
242  CDF.push_back(JCDF_t(fileDescriptor, DIRECT_LIGHT_FROM_DELTARAYS));
243  CDF.push_back(JCDF_t(fileDescriptor, SCATTERED_LIGHT_FROM_DELTARAYS));
244 
245  CDG.push_back(JCDG_t(fileDescriptor, DIRECT_LIGHT_FROM_EMSHOWER));
246  CDG.push_back(JCDG_t(fileDescriptor, SCATTERED_LIGHT_FROM_EMSHOWER));
247 
248 
249  double maximal_road_width = 0.0; // road width [m]
250  double maximal_distance = 0.0; // road width [m]
251 
252  for (size_t i = 0; i != CDF.size(); ++i) {
253 
254  DEBUG("Range CDF["<< CDF[i].type << "] " << CDF[i].function.intensity.getXmax() << " m" << endl);
255 
256  maximal_road_width = max(maximal_road_width, CDF[i].function.intensity.getXmax());
257  }
258 
259  for (size_t i = 0; i != CDG.size(); ++i) {
260 
261  DEBUG("Range CDG["<< CDG[i].type << "] " << CDG[i].function.intensity.getXmax() << " m" << endl);
262 
263  if (!is_scattered(CDF[i].type)) {
264  maximal_road_width = max(maximal_road_width, CDG[i].function.intensity.getXmax());
265  }
266 
267  maximal_distance = max(maximal_distance, CDG[i].function.intensity.getXmax());
268  }
269 
270  NOTICE("Maximal road width [m] " << maximal_road_width << endl);
271  NOTICE("Maximal distance [m] " << maximal_distance << endl);
272 
273 
274  if (detectorFile == "" || inputFile.empty()) {
275  STATUS("Nothing to be done." << endl);
276  return 0;
277  }
278 
279  JDetector detector;
280 
281  try {
282 
283  STATUS("Load detector... " << flush);
284 
285  load(detectorFile, detector);
286 
287  STATUS("OK" << endl);
288  }
289  catch(const JException& error) {
290  FATAL(error);
291  }
292 
293  // remove empty modules
294 
295  for (JDetector::iterator module = detector.begin(); module != detector.end(); ) {
296  if (!module->empty())
297  ++module;
298  else
299  module = detector.erase(module);
300  }
301 
302  vector< JSharedPointer<JRadiationInterface> > radiation, ionization;
303 
304  if (true) {
305 
306  STATUS("Setting up radiation tables... " << flush);
307 
308  //More accuracy can be achieved uncommenting the commented elements and its initializations, remember to do the same in JSeawater.hh (The code will run slower).
309 
310  const JRadiation hydrogen ( 1.0, 1.0, 40, 0.01, 0.1, 0.1);
311  const JRadiation oxygen ( 8.0, 16.0, 40, 0.01, 0.1, 0.1);
312  const JRadiation chlorine (17.0, 35.5, 40, 0.01, 0.1, 0.1);
313  const JRadiation sodium (11.0, 23.0, 40, 0.01, 0.1, 0.1);
314 /* const JRadiation calcium (20.0, 40.0, 40, 0.01, 0.1, 0.1);
315  const JRadiation magnesium(12.0, 24.3, 40, 0.01, 0.1, 0.1);
316  const JRadiation potassium(19.0, 39.0, 40, 0.01, 0.1, 0.1);
317  const JRadiation sulphur (16.0, 32.0, 40, 0.01, 0.1, 0.1);
318 */
319 
320  JSharedPointer<JRadiation> Hydrogen (new JRadiationFunction(hydrogen, 300, 0.2, 1.0e11));
321  JSharedPointer<JRadiation> Oxygen (new JRadiationFunction(oxygen, 300, 0.2, 1.0e11));
322  JSharedPointer<JRadiation> Chlorine (new JRadiationFunction(chlorine, 300, 0.2, 1.0e11));
323  JSharedPointer<JRadiation> Sodium (new JRadiationFunction(sodium, 300, 0.2, 1.0e11));
324 /* JSharedPointer<JRadiation> Calcium (new JRadiationFunction(calcium, 300, 0.2, 1.0e11));
325  JSharedPointer<JRadiation> Magnesium(new JRadiationFunction(magnesium,300, 0.2, 1.0e11));
326  JSharedPointer<JRadiation> Potassium(new JRadiationFunction(potassium,300, 0.2, 1.0e11));
327  JSharedPointer<JRadiation> Sulphur (new JRadiationFunction(sulphur, 300, 0.2, 1.0e11));
328 */
329 
330  JRadiationSource::source_type EErad(&JRadiation::TotalCrossSectionEErad, &JRadiation::EfromEErad);
331  JRadiationSource::source_type Brems(&JRadiation::TotalCrossSectionBrems, &JRadiation::EfromBrems);
332  JRadiationSource::source_type GNrad(&JRadiation::TotalCrossSectionGNrad, &JRadiation::EfromGNrad);
333  JRadiationSource::source_type Ion (&JRadiation::CalculateACoeff , NULL);
334 
335  radiation.push_back(new JRadiationSource(Oxygen, DENSITY_SEA_WATER * JSeaWater::O(), EErad));
336  radiation.push_back(new JRadiationSource(Chlorine, DENSITY_SEA_WATER * JSeaWater::Cl(), EErad));
337  radiation.push_back(new JRadiationSource(Hydrogen, DENSITY_SEA_WATER * JSeaWater::H(), EErad));
338  radiation.push_back(new JRadiationSource(Sodium, DENSITY_SEA_WATER * JSeaWater::Na(), EErad));
339 /* radiation.push_back(new JRadiationSource(Calcium, DENSITY_SEA_WATER * JSeaWater::Ca(), EErad));
340  radiation.push_back(new JRadiationSource(Magnesium,DENSITY_SEA_WATER * JSeaWater::Mg(), EErad));
341  radiation.push_back(new JRadiationSource(Potassium,DENSITY_SEA_WATER * JSeaWater::K(), EErad));
342  radiation.push_back(new JRadiationSource(Sulphur, DENSITY_SEA_WATER * JSeaWater::S(), EErad));
343 */
344 
345  radiation.push_back(new JRadiationSource(Oxygen, DENSITY_SEA_WATER * JSeaWater::O(), Brems));
346  radiation.push_back(new JRadiationSource(Chlorine, DENSITY_SEA_WATER * JSeaWater::Cl(), Brems));
347  radiation.push_back(new JRadiationSource(Hydrogen, DENSITY_SEA_WATER * JSeaWater::H(), Brems));
348  radiation.push_back(new JRadiationSource(Sodium, DENSITY_SEA_WATER * JSeaWater::Na(), Brems));
349 /* radiation.push_back(new JRadiationSource(Calcium, DENSITY_SEA_WATER * JSeaWater::Ca(), Brems));
350  radiation.push_back(new JRadiationSource(Magnesium,DENSITY_SEA_WATER * JSeaWater::Mg(), Brems));
351  radiation.push_back(new JRadiationSource(Potassium,DENSITY_SEA_WATER * JSeaWater::K(), Brems));
352  radiation.push_back(new JRadiationSource(Sulphur, DENSITY_SEA_WATER * JSeaWater::S(), Brems));
353 */
354 
355  radiation.push_back(new JRadiationSource(Oxygen, DENSITY_SEA_WATER * JSeaWater::O(), GNrad));
356  radiation.push_back(new JRadiationSource(Chlorine, DENSITY_SEA_WATER * JSeaWater::Cl(), GNrad));
357  radiation.push_back(new JRadiationSource(Hydrogen, DENSITY_SEA_WATER * JSeaWater::H(), GNrad));
358  radiation.push_back(new JRadiationSource(Sodium, DENSITY_SEA_WATER * JSeaWater::Na(), GNrad));
359 /* radiation.push_back(new JRadiationSource(Calcium, DENSITY_SEA_WATER * JSeaWater::Ca(), GNrad));
360  radiation.push_back(new JRadiationSource(Magnesium,DENSITY_SEA_WATER * JSeaWater::Mg(), GNrad));
361  radiation.push_back(new JRadiationSource(Potassium,DENSITY_SEA_WATER * JSeaWater::K(), GNrad));
362  radiation.push_back(new JRadiationSource(Sulphur, DENSITY_SEA_WATER * JSeaWater::S(), GNrad));
363 */
364  ionization.push_back(new JRadiationSource(Oxygen, DENSITY_SEA_WATER * JSeaWater::O(), Ion));
365  ionization.push_back(new JRadiationSource(Chlorine, DENSITY_SEA_WATER * JSeaWater::Cl(),Ion));
366  ionization.push_back(new JRadiationSource(Hydrogen, DENSITY_SEA_WATER * JSeaWater::H(), Ion));
367  ionization.push_back(new JRadiationSource(Sodium, DENSITY_SEA_WATER * JSeaWater::Na(),Ion));
368 /* ionization.push_back(new JRadiationSource(Calcium, DENSITY_SEA_WATER * JSeaWater::Ca(),Ion));
369  ionization.push_back(new JRadiationSource(Magnesium,DENSITY_SEA_WATER * JSeaWater::Mg(),Ion));
370  ionization.push_back(new JRadiationSource(Potassium,DENSITY_SEA_WATER * JSeaWater::K(), Ion));
371  ionization.push_back(new JRadiationSource(Sulphur, DENSITY_SEA_WATER * JSeaWater::S(), Ion));
372 */
373 
374  STATUS("OK" << endl);
375  }
376 
377 
378  Vec center(0,0,0);
379  Head header;
380 
381  try {
382 
383  header = inputFile.getHeader();
384 
385  JHead buffer(header);
386 
387  center = get<Vec>(buffer);
388 
389  buffer.simul.push_back(JAANET::simul());
390 
391  buffer.simul.rbegin()->program = APPLICATION_JSIRENE;
392  buffer.simul.rbegin()->version = getGITVersion();
393  buffer.simul.rbegin()->date = getDate();
394  buffer.simul.rbegin()->time = getTime();
395 
396  buffer.detector.push_back(JAANET::detector());
397 
398  buffer.detector.rbegin()->program = APPLICATION_JSIRENE;
399  buffer.detector.rbegin()->filename = detectorFile;
400 
401  buffer.push(&JHead::simul);
402  buffer.push(&JHead::detector);
403 
404  if (!keep) {
405 
406  buffer.coord_origin = coord_origin(0.0, 0.0, 0.0);
407  buffer.can.zmin += center.z;
408  buffer.can.zmax += center.z;
409 
410  buffer.push(&JHead::coord_origin);
411  buffer.push(&JHead::can);
412  }
413 
414  const JCircle2D circle(detector.begin(), detector.end());
415 
416  center += Vec(circle.getX(), circle.getY(), 0.0);
417 
418  copy(buffer, header);
419  }
420  catch(const JException& error) {
421  FATAL(error);
422  }
423 
424  if (!keep)
425  NOTICE("Offset applied to true tracks is: " << center << endl);
426  else
427  NOTICE("No offset applied to true tracks." << endl);
428 
429  JCylinder3D cylinder(detector.begin(), detector.end());
430 
431  cylinder.addMargin(maximal_distance);
432 
433  if (cylinder.getZmin() < Zbed) {
434  cylinder.setZmin(Zbed);
435  }
436 
437  NOTICE("Light generation volume: " << cylinder << endl);
438 
439  TProfile cpu("cpu", NULL, 14, 1.0, 8.0);
440  TH1D job("job", NULL, 400, 0.5, 400.5);
441 
442 
443  outputFile.open();
444 
445  if (!outputFile.is_open()) {
446  FATAL("Error opening file " << outputFile << endl);
447  }
448 
449  outputFile.put(meta);
450  outputFile.put(header);
451  outputFile.put(*gRandom);
452 
453 
454  JTimer timer;
455 
456  for (JMultipleFileScanner<Evt>& in = inputFile; in.hasNext(); ) {
457 
458  STATUS("event: " << setw(10) << in.getCounter() << '\r'); DEBUG(endl);
459 
460  job.Fill(1.0);
461 
462  Evt* evt = in.next();
463 
464  if (!keep) {
465  for (vector<Trk>::iterator track = evt->mc_trks.begin(); track != evt->mc_trks.end(); ++track) {
466  track->pos += center;
467  }
468  }
469 
470  Evt event(*evt); // output
471 
472  event.mc_hits.clear();
473 
474  JHits_t mc_hits; // temporary buffer
475 
476  timer.reset();
477  timer.start();
478 
479  for (vector<Trk>::const_iterator track = evt->mc_trks.begin(); track != evt->mc_trks.end(); ++track) {
480 
481  if (!track->is_finalstate()) continue; //only final state particles produce light
482 
483  if (is_muon(*track)) {
484 
485  // -----------------------------------------------
486  // muon
487  // -----------------------------------------------
488 
489  job.Fill(2.0);
490 
491  typedef JDetectorSubset<JModuleComparator3Z_t> JDetectorSubset_t;
492 
493  const JCylinder3D::intersection_type intersection = cylinder.getIntersection(getAxis(*track));
494 
495  double Zmin = intersection.first;
496  double Zmax = intersection.second;
497 
498  if (Zmax - Zmin <= Dmin_m) {
499  continue;
500  }
501 
502  JVertex vertex(0.0, track->t, track->E); // start of muon
503 
504  if (track->pos.z < Zbed) { // propagate muon through rock
505 
506  if (track->dir.z > 0.0)
507  vertex.step(gRock, (Zbed - track->pos.z) / track->dir.z);
508  else
509  continue;
510  }
511 
512  if (vertex.getZ() < Zmin) { // propagate muon through water
513  vertex.step(gWater, Zmin - vertex.getZ());
514  }
515 
516  if (vertex.getRange() <= Dmin_m) {
517  continue;
518  }
519 
520  job.Fill(3.0);
521 
522  const JDetectorSubset_t subdetector(detector, getAxis(*track), maximal_road_width);
523 
524  if (subdetector.empty()) {
525  continue;
526  }
527 
528  job.Fill(4.0);
529 
530  JTrack muon(vertex); // propagate muon trough detector
531 
532  while (vertex.getE() >= Emin_GeV && vertex.getZ() < Zmax) {
533 
534  const int Nr = radiation.size();
535 
536  double li[Nr]; // inverse interaction lengths
537  double ls = 1.0e-5; // minimal total inverse interaction length (100 km)^-1
538 
539  for (int i = 0; i != Nr; ++i) {
540  ls += li[i] = radiation[i]->getInverseInteractionLength(vertex.getE());
541  }
542 
543  const int Ni = ionization.size();
544 
545  double Ai[Ni];
546  double As=0;
547 
548  for (int i = 0; i != Ni; ++i) {
549  As +=Ai[i]= ionization[i]->getA(vertex.getE());
550  }
551 
552  const double ds = min(gRandom->Exp(1.0) / ls, vertex.getRange(As));
553 
554  vertex.step(ds,As); //Ionization continuous losses
555 
556  if (vertex.getE() >= Emin_GeV) {
557 
558  double Es = 0.0; // shower energy [GeV]
559  double y = gRandom->Uniform(ls);
560 
561  for (int i = 0; i != Nr; ++i) {
562 
563  y -= li[i];
564 
565  if (y < 0.0) {
566  Es = radiation[i]->getEnergyOfShower(vertex.getE()); // shower energy [GeV]
567  break;
568  }
569  }
570 
571  vertex.applyEloss(Es);
572 
573  if (Es >= Ecut_GeV) {
574  muon.push_back(vertex);
575  }
576  }
577  }
578 
579  // add muon end point
580 
581  if (vertex.getE() < Emin_GeV && vertex.getRange() > 0.0) {
582 
583  vertex.step(vertex.getRange()); //For energies below 1 GeV, step(ds) is used, Geane.hh default "a" parameter will be used.
584 
585  muon.push_back(vertex);
586  }
587 
588  // add information to output muon
589 
590  vector<Trk>::iterator trk = find_if(event.mc_trks.begin(),
591  event.mc_trks.end(),
592  make_predicate(&Trk::id, track->id));
593 
594  if (trk != event.mc_trks.end()) {
595  trk->len = (muon.rbegin()->getZ() < Zmax ? +1 : -1) * (muon.rbegin()->getZ() - muon.begin()->getZ());
596  trk->setusr(mc_usr_keys::energy_lost_in_can, muon.begin()->getE() - muon.rbegin()->getE());
597  }
598 
599  for (JDetector::const_iterator module = subdetector.begin(); module != subdetector.end(); ++module) {
600 
601  const double z0 = muon.begin()->getZ();
602  const double t0 = muon.begin()->getT();
603  const double R = module->getX();
604  const double Z = module->getZ() - R / getTanThetaC();
605 
606  if (Z >= muon.begin()->getZ() && Z <= muon.rbegin()->getZ()) {
607 
608  for (size_t i = 0; i != CDF.size(); ++i) {
609 
610  if (R < CDF[i].integral.getXmax()) {
611 
612  try {
613 
614  double W = 1.0; // mip
615 
616  if (is_deltarays(CDF[i].type)) {
617  W = getDeltaRaysFromMuon(muon.getE(Z)); // delta-rays
618  }
619 
620  const double NPE = CDF[i].integral.getNPE(R) * module->size() * factor * W;
621  const int N = gRandom->Poisson(NPE);
622 
623  if (N != 0) {
624 
625  double ns = 0.0;
626 
627  for (JModule::const_iterator pmt = module->begin(); pmt != module->end(); ++pmt) {
628 
629  const double R = pmt->getX();
630  const double Z = pmt->getZ() - z0;
631  const double theta = pmt->getTheta();
632  const double phi = fabs(pmt->getPhi());
633 
634  const double npe = CDF[i].function.getNPE(R, theta, phi) * factor * W;
635 
636  ns += npe;
637 
638  int n0 = min(getNumberOfPhotoElectrons(NPE, N, npe), Nmax_npe);
639 
640  job.Fill((double) (100 + CDF[i].type), (double) n0);
641 
642  while (n0 != 0) {
643 
644  const double t1 = CDF[i].function.getTime(R, theta, phi, gRandom->Rndm());
645  const int n1 = getNumberOfPhotoElectrons(n0);
646 
647  mc_hits.push_back(JHit_t(mc_hits.size() + 1,
648  pmt->getID(),
649  getHitType(CDF[i].type),
650  track->id,
651  t0 + (R * getTanThetaC() + Z) / C + t1,
652  n1));
653 
654  n0 -= n1;
655  }
656  }
657 
658  if (ns > NPE) {
659  job.Fill((double) (300 + CDF[i].type));
660  }
661  }
662  }
663  catch(const exception& error) {
664  job.Fill((double) (200 + CDF[i].type));
665  }
666  }
667  }
668  }
669  }
670 
671  for (JTrack::const_iterator vertex = muon.begin(); vertex != muon.end(); ++vertex) {
672 
673  const double z0 = vertex->getZ();
674  const double t0 = vertex->getT();
675  const double Es = vertex->getEs();
676 
677  int origin = track->id;
678 
679  if (writeEMShowers) {
680  origin = event.mc_trks.size() + 1;
681  }
682 
683  int number_of_hits = 0;
684 
685  JDetectorSubset_t::range_type range = subdetector.getRange(z0 - maximal_distance,
686  z0 + maximal_distance);
687 
688  for (JDetector::const_iterator module = range.begin(); module != range.end(); ++module) {
689 
690  const double dz = geanz.getMaximum(Es);
691  const double R = module->getX();
692  const double Z = module->getZ() - z0 - dz;
693  const double D = sqrt(R*R + Z*Z);
694  const double cd = Z / D;
695 
696  for (size_t i = 0; i != CDG.size(); ++i) {
697 
698  if (D < CDG[i].integral.getXmax()) {
699 
700  try {
701 
702  const double NPE = CDG[i].integral.getNPE(D, cd) * Es * module->size() * factor;
703  const int N = gRandom->Poisson(NPE);
704 
705  if (N != 0) {
706 
707  double ns = 0.0;
708 
709  for (JModule::const_iterator pmt = module->begin(); pmt != module->end(); ++pmt) {
710 
711  const double dz = geanz.getMaximum(Es);
712  const double R = pmt->getX();
713  const double Z = pmt->getZ() - z0 - dz;
714  const double theta = pmt->getTheta();
715  const double phi = fabs(pmt->getPhi());
716  const double D = sqrt(R*R + Z*Z);
717  const double cd = Z / D;
718 
719  const double npe = CDG[i].function.getNPE(D, cd, theta, phi) * Es * factor;
720 
721  ns += npe;
722 
723  int n0 = min(getNumberOfPhotoElectrons(NPE, N, npe), Nmax_npe);
724 
725  job.Fill((double) (100 + CDG[i].type), (double) n0);
726 
727  while (n0 != 0) {
728 
729  const double dz = geanz.getLength(Es, gRandom->Rndm());
730  const double Z = pmt->getZ() - z0 - dz;
731  const double D = sqrt(R*R + Z*Z);
732  const double cd = Z / D;
733 
734  const double t1 = CDG[i].function.getTime(D, cd, theta, phi, gRandom->Rndm());
735  const int n1 = getNumberOfPhotoElectrons(n0);
736 
737  mc_hits.push_back(JHit_t(mc_hits.size() + 1,
738  pmt->getID(),
739  getHitType(CDG[i].type),
740  origin,
741  t0 + (dz + D * getIndexOfRefraction()) / C + t1,
742  n1));
743 
744  n0 -= n1;
745 
746  number_of_hits += n1;
747  }
748  }
749 
750  if (ns > NPE) {
751  job.Fill((double) (300 + CDG[i].type));
752  }
753  }
754  }
755  catch(const exception& error) {
756  job.Fill((double) (200 + CDG[i].type));
757  }
758  }
759  }
760  }
761 
762  if (writeEMShowers && number_of_hits != 0) {
763 
764  event.mc_trks.push_back(JTrk_t(origin,
765  TRACK_TYPE_PHOTON,
766  track->id,
767  track->pos + track->dir * vertex->getZ(),
768  track->dir,
769  vertex->getT(),
770  Es));
771  }
772  }
773 
774  } else if (track->len>0) {
775 
776  // -----------------------------------------------
777  // decayed particles treated as mip (tau includes mip+deltaray)
778  // -----------------------------------------------
779 
780  job.Fill(6.0);
781 
782  const double z0 = 0.0;
783  const double z1 = z0 + track->len;
784  const double t0 = track->t;
785  const double E = track->E;
786 
787  const JTransformation3D transformation = getTransformation(*track);
788 
789  JModule buffer;
790 
791  for (JDetector::const_iterator module = detector.begin(); module != detector.end(); ++module) {
792 
793  const JPosition3D pos = transformation.transform(module->getPosition());
794 
795  const double R = pos.getX();
796  const double Z = pos.getZ() - R / getTanThetaC();
797 
798  if (Z < z0 ||
799  Z > z1 ||
800  R > maximal_road_width) {
801  continue;
802  }
803 
804  for (size_t i = 0; i != CDF.size(); ++i) {
805 
806  double W = 1.0; // mip
807 
808  if (is_deltarays(CDF[i].type)) {
809  if (is_tau(*track)) W = getDeltaRaysFromTau(E); // delta-rays
810  else continue;
811  }
812 
813  if (R < CDF[i].integral.getXmax()) {
814 
815  try {
816 
817  const double NPE = CDF[i].integral.getNPE(R) * module->size() * factor * W;
818  const int N = gRandom->Poisson(NPE);
819 
820  if (N != 0) {
821 
822  buffer = *module;
823 
824  buffer.transform(transformation);
825 
826  double ns = 0.0;
827 
828  for (JModule::const_iterator pmt = buffer.begin(); pmt != buffer.end(); ++pmt) {
829 
830  const double R = pmt->getX();
831  const double Z = pmt->getZ() - z0;
832  const double theta = pmt->getTheta();
833  const double phi = fabs(pmt->getPhi());
834 
835  const double npe = CDF[i].function.getNPE(R, theta, phi) * factor * W;
836 
837  ns += npe;
838 
839  int n0 = min(getNumberOfPhotoElectrons(NPE, N, npe), Nmax_npe);
840 
841  job.Fill((double) (120 + CDF[i].type), (double) n0);
842 
843  while (n0 != 0) {
844 
845  const double t1 = CDF[i].function.getTime(R, theta, phi, gRandom->Rndm());
846  const int n1 = getNumberOfPhotoElectrons(n0);
847 
848  mc_hits.push_back(JHit_t(mc_hits.size() + 1,
849  pmt->getID(),
850  getHitType(CDF[i].type),
851  track->id,
852  t0 + (R * getTanThetaC() + Z) / C + t1,
853  n1));
854 
855  n0 -= n1;
856  }
857  }
858 
859  if (ns > NPE) {
860  job.Fill((double) (320 + CDF[i].type));
861  }
862  }
863  }
864  catch(const exception& error) {
865  job.Fill((double) (220 + CDF[i].type));
866  }
867 
868  }
869  }
870  }
871 
872  if (!buffer.empty()) {
873  job.Fill(7.0);
874  }
875 
876  } else if (!is_neutrino(*track)) {
877 
878  if (JPDB::getInstance().hasPDG(track->type)) {
879 
880  // -----------------------------------------------
881  // electron or hadron
882  // -----------------------------------------------
883 
884  job.Fill(8.0);
885 
886  double E = track->E;
887 
888  try {
889  E = getKineticEnergy(E, JPDB::getInstance().getPDG(track->type).mass);
890  }
891  catch(const exception& error) {
892  ERROR(error.what() << endl);
893  }
894 
895  E = pythia(track->type, E);
896 
897  if (E < Ecut_GeV || cylinder.getDistance(getPosition(*track)) > Dmin_m) {
898  continue;
899  }
900 
901  const double z0 = 0.0;
902  const double t0 = track->t;
903 
904  const JTransformation3D transformation = getTransformation(*track);
905 
906  JModule buffer;
907 
908  for (JDetector::const_iterator module = detector.begin(); module != detector.end(); ++module) {
909 
910  const JPosition3D pos = transformation.transform(module->getPosition());
911 
912  const double dz = geanz.getMaximum(E);
913  const double R = pos.getX();
914  const double Z = pos.getZ() - z0 - dz;
915  const double D = sqrt(R*R + Z*Z);
916  const double cd = Z / D;
917 
918  if (D > maximal_distance) {
919  continue;
920  }
921 
922  for (size_t i = 0; i != CDG.size(); ++i) {
923 
924  if (D < CDG[i].integral.getXmax()) {
925 
926  try {
927 
928  const double NPE = CDG[i].integral.getNPE(D, cd) * E * module->size() * factor;
929  const int N = gRandom->Poisson(NPE);
930 
931  if (N != 0) {
932 
933  buffer = *module;
934 
935  buffer.transform(transformation);
936 
937  double ns = 0.0;
938 
939  for (JModule::const_iterator pmt = buffer.begin(); pmt != buffer.end(); ++pmt) {
940 
941  const double dz = geanz.getMaximum(E);
942  const double R = pmt->getX();
943  const double Z = pmt->getZ() - z0 - dz;
944  const double theta = pmt->getTheta();
945  const double phi = fabs(pmt->getPhi());
946  const double D = sqrt(R*R + Z*Z);
947  const double cd = Z / D;
948 
949  const double npe = CDG[i].function.getNPE(D, cd, theta, phi) * E * factor;
950 
951  ns += npe;
952 
953  int n0 = min(getNumberOfPhotoElectrons(NPE, N, npe), Nmax_npe);
954 
955  job.Fill((double) (140 + CDG[i].type), (double) n0);
956 
957  while (n0 != 0) {
958 
959  const double dz = geanz.getLength(E, gRandom->Rndm());
960  const double Z = pmt->getZ() - z0 - dz;
961  const double D = sqrt(R*R + Z*Z);
962  const double cd = Z / D;
963 
964  const double t1 = CDG[i].function.getTime(D, cd, theta, phi, gRandom->Rndm());
965  const int n1 = getNumberOfPhotoElectrons(n0);
966 
967  mc_hits.push_back(JHit_t(mc_hits.size() + 1,
968  pmt->getID(),
969  getHitType(CDG[i].type, true),
970  track->id,
971  t0 + (dz + D * getIndexOfRefraction()) / C + t1,
972  n1));
973 
974  n0 -= n1;
975  }
976  }
977 
978  if (ns > NPE) {
979  job.Fill((double) (340 + CDG[i].type));
980  }
981  }
982  }
983  catch(const exception& error) {
984  job.Fill((double) (240 + CDG[i].type));
985  }
986  }
987  }
988  }
989 
990  if (!buffer.empty()) {
991  job.Fill(9.0);
992  }
993 
994  } else {
995  job.Fill(21.0);
996  }
997  }
998  }
999 
1000  if (!mc_hits.empty()) {
1001 
1002  mc_hits.merge(Tmax_ns);
1003 
1004  event.mc_hits.resize(mc_hits.size());
1005 
1006  copy(mc_hits.begin(), mc_hits.end(), event.mc_hits.begin());
1007  }
1008 
1009  timer.stop();
1010 
1011  if (has_neutrino(event)) {
1012  cpu.Fill(log10(get_neutrino(event).E), (double) timer.usec_ucpu * 1.0e-3);
1013  }
1014 
1015  if (event.mc_hits.size() >= numberOfHits) {
1016 
1017  outputFile.put(event);
1018 
1019  job.Fill(10.0);
1020  }
1021  }
1022  STATUS(endl);
1023 
1024  outputFile.put(cpu);
1025  outputFile.put(job);
1026  outputFile.put(*gRandom);
1027  outputFile.close();
1028 }
JAANET::JHit_t
Auxiliary class to set-up Hit.
Definition: JHit_t.hh:25
mc_usr_keys::energy_lost_in_can
const char *const energy_lost_in_can
Definition: io_ascii.hh:45
JSUPPORT::JMeta
Auxiliary class for ROOT I/O of application specific meta data.
Definition: JMeta.hh:70
JSUPPORT::JFileRecorder
Object writing to file.
Definition: JFileRecorder.hh:41
JPARSER::JParser
Utility class to parse command line options.
Definition: JParser.hh:1500
JLANG::JException
General exception.
Definition: JException.hh:23
D
do echo Generating $dir eval D
Definition: JDrawLED.sh:50
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
JDETECTOR::JModule
Data structure for a composite optical module.
Definition: JModule.hh:57
CDF
then usage $script[detector file(input file)+[output file[CDF file descriptor]]] nNote that if more than one input file is all other arguments must be provided fi case set_variable CDF
Definition: JSirene.sh:33
JAANET::getTransformation
JTransformation3D getTransformation(const Trk &track)
Get transformation.
Definition: JAAnetToolkit.hh:265
JGEOMETRY2D::JCircle2D
Data structure for circle in two dimensions.
Definition: JCircle2D.hh:32
gmake_property
#define gmake_property(A)
macro to convert (template) parameter to JPropertiesElement object
Definition: JProperties.hh:1291
safetyFactor
double safetyFactor
safety factor for average CDF integral of optical module
Definition: JSirene.cc:65
JAANET::simul
Generator for simulation.
Definition: JHead.hh:460
JPHYSICS::getIndexOfRefraction
double getIndexOfRefraction()
Get average index of refraction of water corresponding to group velocity.
Definition: JPhysics/JConstants.hh:127
STATUS
#define STATUS(A)
Definition: JMessage.hh:63
JPHYSICS::SCATTERED_LIGHT_FROM_EMSHOWER
scattered light from EM shower
Definition: JPDFTypes.hh:41
JAANET::has_neutrino
bool has_neutrino(const Evt &evt)
Test whether given event has an incoming neutrino.
Definition: JAAnetToolkit.hh:418
JDETECTOR::JDetector
Detector data structure.
Definition: JDetector.hh:80
JAANET::is_muon
bool is_muon(const Trk &track)
Test whether given track is a (anti-)muon.
Definition: JAAnetToolkit.hh:321
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
JPHYSICS::getDeltaRaysFromTau
double getDeltaRaysFromTau(const double E)
Equivalent EM-shower energy due to delta-rays per unit tau track length.
Definition: JPDFToolkit.hh:103
JEEP::JProperties
Utility class to parse parameter values.
Definition: JProperties.hh:496
JPHYSICS::H
static const double H
Planck constant [eV s].
Definition: JPhysics/JConstants.hh:28
JPARSER::initialised
Empty structure for specification of parser element that is initialised (i.e. does not require input)...
Definition: JParser.hh:66
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.
JPHYSICS::gWater
static const JGeaneWater gWater
Function object for energy loss of muon in sea water.
Definition: JGeane.hh:323
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:116
JPHYSICS::DENSITY_SEA_WATER
static const double DENSITY_SEA_WATER
Fixed environment values.
Definition: JPhysics/JConstants.hh:37
JSIRENE::pythia
static const JPythia pythia
Function object for relative light yield as a function of GEANT particle code.
Definition: JPythia.hh:96
JAANET::getTime
double getTime(const Hit &hit)
Get true time of hit.
Definition: JAAnetToolkit.hh:79
outputFile
string outputFile
Definition: JDAQTimesliceSelector.cc:37
APPLICATION_JSIRENE
static const char *const APPLICATION_JSIRENE
detector simulation
Definition: applications.hh:19
JPHYSICS::DIRECT_LIGHT_FROM_MUON
direct light from muon
Definition: JPDFTypes.hh:29
JAANET::coord_origin
Coordinate origin.
Definition: JHead.hh:663
JPHYSICS::JRadiationFunction
Fast implementation of class JRadiation.
Definition: JRadiationSource.hh:35
JPHYSICS::JRadiationSource
Source of radiation.
Definition: JRadiationSource.hh:275
JAANET::is_neutrino
bool is_neutrino(const Trk &track)
Test whether given track is a neutrino.
Definition: JAAnetToolkit.hh:303
JPHYSICS::C
static const double C
Physics constants.
Definition: JPhysics/JConstants.hh:25
JPHYSICS::SCATTERED_LIGHT_FROM_MUON
scattered light from muon
Definition: JPDFTypes.hh:30
JSUPPORT::JLimit
Auxiliary class for defining the range of iterations of objects.
Definition: JLimit.hh:41
Vec
The Vec class is a straightforward 3-d vector, which also works in pyroot.
Definition: Vec.hh:12
JGEOMETRY3D::JCylinder3D
Cylinder object.
Definition: JCylinder3D.hh:39
JAANET::detector
Detector file.
Definition: JHead.hh:196
JEEP::getGITVersion
std::string getGITVersion(const std::string &tag)
Get GIT version for given GIT tag.
Definition: Jeep/JGITTags.hh:34
JLANG::JSharedPointer
The template JSharedPointer class can be used to share a pointer to an object.
Definition: JSharedPointer.hh:28
JPHYSICS::JRadiation
Auxiliary class for the calculation of the muon radiative cross sections.
Definition: JRadiation.hh:31
JLANG::getInstance
T & getInstance(const T &object)
Get static instance from temporary object.
Definition: JObject.hh:75
make_field
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1961
JPHYSICS::geanz
static const JGeanz geanz(1.85, 0.62, 0.54)
Function object for longitudinal EM-shower profile.
JAANET::getAxis
JAxis3D getAxis(const Trk &track)
Get axis.
Definition: JAAnetToolkit.hh:229
JEEP::JTimer
Auxiliary class for CPU timing and usage.
Definition: JTimer.hh:32
JAANET::getPosition
JPosition3D getPosition(const Vec &pos)
Get position.
Definition: JAAnetToolkit.hh:157
JPHYSICS::SCATTERED_LIGHT_FROM_DELTARAYS
scattered light from delta-rays
Definition: JPDFTypes.hh:36
NOTICE
#define NOTICE(A)
Definition: JMessage.hh:64
ERROR
#define ERROR(A)
Definition: JMessage.hh:66
JPHYSICS::DIRECT_LIGHT_FROM_EMSHOWER
direct light from EM shower
Definition: JPDFTypes.hh:40
JSIRENE::JTrack
Muon trajectory.
Definition: JSireneToolkit.hh:326
JAANET::JHits_t
Auxiliary data structure for list of hits with hit merging capability.
Definition: JHit_t.hh:105
JGEOMETRY3D::JCylinder3D::addMargin
void addMargin(const double D)
Add (safety) margin.
Definition: JCylinder3D.hh:173
debug
int debug
debug level
Definition: JSirene.cc:63
JSIRENE::getKineticEnergy
double getKineticEnergy(const double E, const double m)
Get kinetic energy of particle with given mass.
Definition: JSireneToolkit.hh:371
JSYSTEM::getDate
const char * getDate(const JDateAndTimeFormat option=ISO8601)
Get ASCII formatted date.
Definition: JSystem/JDate.hh:35
R
then usage $script[distance] fi case set_variable R
Definition: JDrawLED.sh:40
JDETECTOR::JDetectorSubset
Detector subset with binary search functionality.
Definition: JDetectorSubset.hh:119
JDETECTOR::JDetectorSubset_t
Detector subset without binary search functionality.
Definition: JDetectorSubset.hh:38
JAANET::JHead
Monte Carlo run header.
Definition: JHead.hh:1113
JSIRENE::JVertex
Vertex of energy loss of muon.
Definition: JSireneToolkit.hh:146
Head
The Head class reflects the header of Monte-Carlo event files, which consists of keys (also referred ...
Definition: Head.hh:67
FATAL
#define FATAL(A)
Definition: JMessage.hh:67
Trk::id
int id
track identifier
Definition: Trk.hh:15
JPHYSICS::DIRECT_LIGHT_FROM_DELTARAYS
direct light from delta-rays
Definition: JPDFTypes.hh:35
range
z range($ZMAX-$ZMIN)< $MINIMAL_DZ." fi fi typeset -Z 4 STRING typeset -Z 2 FLOOR JPlot1D -f $
Definition: module-Z:fit.sh:116
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
JDETECTOR::load
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
Definition: JDetectorToolkit.hh:445
JPHYSICS::JGeanz::getMaximum
double getMaximum(const double E) const
Get depth of shower maximum.
Definition: JGeanz.hh:162
JSUPPORT::JMultipleFileScanner
General purpose class for object reading from a list of file names.
Definition: JMultipleFileScanner.hh:167
JPHYSICS::is_scattered
bool is_scattered(const int pdf)
Test if given PDF type corresponds to scattered light.
Definition: JPDFTypes.hh:191
JDETECTOR::JModule::transform
void transform(const JRotation3D &R, const JVector3D &pos)
Transformation of geometry (see method JGEOMETRY3D::JPosition3D::transform(const JRotation3D&amp;, const JVector3D&amp;)).
Definition: JModule.hh:382
JPHYSICS::is_deltarays
bool is_deltarays(const int pdf)
Test if given PDF type corresponds to Cherenkov light from delta-rays.
Definition: JPDFTypes.hh:177
JGEOMETRY3D::JVector3D::getX
double getX() const
Get x position.
Definition: JVector3D.hh:94
JPHYSICS::getDeltaRaysFromMuon
double getDeltaRaysFromMuon(const double E)
Equivalent EM-shower energy due to delta-rays per unit muon track length.
Definition: JPDFToolkit.hh:75
JAANET::copy
void copy(const Head &from, JHead &to)
Copy header from from to to.
Definition: JHead.cc:139
JSYSTEM::ls
Auxiliary data structure to list files in directory.
Definition: JSystemToolkit.hh:302
JPHYSICS::getTanThetaC
double getTanThetaC()
Get average tangent of Cherenkov angle of water corresponding to group velocity.
Definition: JPhysics/JConstants.hh:149
JGEOMETRY3D::JPosition3D
Data structure for position in three dimensions.
Definition: JPosition3D.hh:36
JSUPPORT::JLimit::getLimit
const JLimit & getLimit() const
Get limit.
Definition: JLimit.hh:73
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
JAANET::get_neutrino
const Trk & get_neutrino(const Evt &evt)
Get incoming neutrino.
Definition: JAAnetToolkit.hh:429
detector
do set_variable DETECTOR_TXT $WORKDIR detector
Definition: detector-XY:fit1d.sh:37
JAANET::is_tau
bool is_tau(const Trk &track)
Test whether given track is a (anti-)tau.
Definition: JAAnetToolkit.hh:329
in
then fatal Wrong number of arguments fi set_variable DETECTOR $argv[1] set_variable INPUT_FILE $argv[2] eval JPrintDetector a $DETECTOR O IDENTIFIER eval JPrintDetector a $DETECTOR O SUMMARY source JAcoustics sh $DETECTOR_ID CHECK_EXIT_CODE typeset A TRIPODS get_tripods $WORKDIR tripod txt TRIPODS for EMITTER in
Definition: JCanberra.sh:38
JGEOMETRY3D::JTransformation3D::transform
JPosition3D transform(const JPosition3D &pos) const
Transform position.
Definition: JTransformation3D.hh:49
N
then usage $script[input file[working directory[option]]] nWhere option can be N
Definition: JMuonPostfit.sh:37
numberOfBins
int numberOfBins
number of bins for average CDF integral of optical module
Definition: JSirene.cc:64
JGEOMETRY3D::JVector3D::getZ
double getZ() const
Get z position.
Definition: JVector3D.hh:115
Evt::mc_trks
std::vector< Trk > mc_trks
MC: list of MC truth tracks.
Definition: Evt.hh:46
JAANET::JTrk_t
Auxiliary class to set-up Trk.
Definition: JTrk_t.hh:19
Evt
The Evt class respresent a Monte Carlo (MC) event as well as an offline event.
Definition: Evt.hh:19
E
then usage $script[input file[working directory[option]]] nWhere option can be E
Definition: JMuonPostfit.sh:37
DEBUG
#define DEBUG(A)
Message macros.
Definition: JMessage.hh:62

Variable Documentation

int debug

debug level

Definition at line 63 of file JSirene.cc.

int numberOfBins = 200

number of bins for average CDF integral of optical module

Definition at line 64 of file JSirene.cc.

double safetyFactor = 1.7

safety factor for average CDF integral of optical module

Definition at line 65 of file JSirene.cc.

Generated by   doxygen 1.8.5