Jpp  18.2.0-rc.1
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JCalibrateK40.cc
Go to the documentation of this file.
1 #include <string>
2 #include <iostream>
3 #include <iomanip>
4 #include <limits>
5 
6 #include "TROOT.h"
7 #include "TFile.h"
8 #include "TH1D.h"
9 #include "TH2D.h"
10 
13 
14 #include "JDAQ/JDAQTimesliceIO.hh"
16 #include "JDAQ/JDAQEvaluator.hh"
19 #include "JTrigger/JHitR0.hh"
20 #include "JTrigger/JMatchL0.hh"
21 #include "JTrigger/JHitToolkit.hh"
26 #include "JDetector/JDetector.hh"
33 #include "JSupport/JSupport.hh"
34 #include "JSupport/JMeta.hh"
36 #include "JTools/JRange.hh"
38 
40 #include "JROOT/JRootToolkit.hh"
41 #include "JROOT/JRootFileWriter.hh"
44 
45 #include "Jeep/JPrint.hh"
46 #include "Jeep/JParser.hh"
47 #include "Jeep/JMessage.hh"
48 
49 namespace {
50 
53  using JDETECTOR::JModule;
54 
55  // Methods for background estimation
56 
57  static const char* rates_t = "rates"; //!< singles rates from summary data
58  static const char* count_t = "counts"; //!< singles rates from L0 data
59  static const char* tails_t = "tails"; //!< tails of time residual distribution
60  static const char* rndms_t = "randoms"; //!< random time offsets (only L0 data)
61 
62  /**
63  * Maximal time-over-threshold.
64  */
65  static const double ToTmax_ns = 250;
66 
67  /**
68  * Auxiliary data structure for histogram management.
69  */
70  struct JHistogram :
71  public JCombinatorics
72  {
73  /**
74  * Default constructor.
75  */
76  JHistogram() :
77  h2s(NULL),
78  h1b(NULL),
79  h1L(NULL)
80  {}
81 
82 
83  /**
84  * Constructor.
85  *
86  * \param module module
87  * \param axis axis
88  */
89  JHistogram(const JModule& module, const JAbstractHistogram<Double_t>& axis)
90  {
91  using namespace JPP;
92 
93  // sort the PMT pairs by opening angle in the order largest angle -> smallest angle
94 
95  this->configure(module.size());
96 
97  this->sort(JPairwiseComparator(module));
98 
99  const Int_t nx = this->getNumberOfPairs();
100  const Double_t xmin = -0.5;
101  const Double_t xmax = nx - 0.5;
102 
103  h2s = new TH2D(MAKE_CSTRING(module.getID() << _2S), NULL, nx, xmin, xmax, axis.getNumberOfBins(), axis.getLowerLimit(), axis.getUpperLimit());
104  h1b = new TH1D(MAKE_CSTRING(module.getID() << _1B), NULL, nx, xmin, xmax);
105  h1L = new TH1D(MAKE_CSTRING(module.getID() << _1L), NULL, nx, xmin, xmax);
106 
107  h2s->Sumw2();
108  h1b->Sumw2();
109  h1L->Sumw2();
110  }
111 
112  /**
113  * Check validity.
114  *
115  * \return true if valid; else false
116  */
117  bool is_valid() const
118  {
119  return (h2s != NULL &&
120  h1b != NULL &&
121  h1L != NULL);
122  }
123 
124  TH2D* h2s; //!< time distribution per PMT pair
125  TH1D* h1b; //!< background rates per PMT pair
126  TH1D* h1L; //!< livetimes per PMT pair
127  };
128 
129 
130  /**
131  * Get coincidence probability of two PMTs within one module due to random background.
132  *
133  * \param E1 expected counts of 1st PMT
134  * \param E2 expected counts of 2nd PMT
135  * \param ED expected counts of module, excluding PMT 1 and PMT 2
136  * \param M_min multiplicity lower limit
137  * \param M_max multiplicity upper limit
138  * \return probability
139  */
140  inline double getP(const double E1, const double E2, const double ED, const int M_min, const int M_max)
141  {
142  // 1st, always calculate the probability for twofold multiplicity
143 
144  double P = E1 * E2 * exp(-(E1 + E2 + ED));
145 
146  // 2nd, calculate the probability of multiplicity at lower limit
147 
148  int MD = 1;
149 
150  for ( ; MD + 2 <= M_min; ++MD) {
151  P *= ED/MD;
152  }
153 
154  // 3rd, add the probabilities of multiplicities up to upper limit
155 
156  double P_sum = P;
157 
158  for ( ; MD + 2 <= M_max; ++MD) {
159  P *= ED/MD;
160  P_sum += P;
161  }
162 
163  return P_sum;
164  }
165 }
166 
167 
168 /**
169  * \file
170  *
171  * Auxiliary program to determine PMT parameters from K40 data.
172  *
173  * By default, the combinatorial background is estimated from the singles rate.\n
174  * In case L0 data are taken (and no summary data are available),
175  * the singles rate can be determined from the amount of L0 data.\n
176  * One can also estimate the background from the tails in the time residual distribution.\n
177  * In that case, option -B can be used to specify the minimal and maximal time offset in ns.\n
178  * For example -B "10 20".\n
179  * Note that the minimal and maximal time should be larger that the width of
180  * the time residual distribution and less than the time window of the L1 coincidence, respectively.\n
181  * The time window is applied symmetrically around a time offset of zero.\n
182  * Finally, if L0 data are available, one can estimate the background using
183  * the same procedure but with a random (read wrong) time calibration.
184  *
185  * The option -M can be used to specify a range of multiplicities.\n
186  * For example -M "2 2" will strictly select two-fold coincidences.
187  * Note that such a selection can bias the result.
188  *
189  * Also consult <a href="https://common.pages.km3net.de/jpp/JMonitor.PDF">documentation</a>.
190  * \author mdejong
191  */
192 int main(int argc, char **argv)
193 {
194  using namespace std;
195  using namespace JPP;
196  using namespace KM3NETDAQ;
197 
198  JMultipleFileScanner_t inputFile;
199  counter_type numberOfEvents;
200  string outputFile;
201  string detectorFile;
202  Double_t TMax_ns;
203  JRange<double> rateRange_Hz;
204  JRange<double> totRange_ns;
205  JRange<double> Tail_ns;
206  JRange<int> multiplicity;
207  double deadTime_us;
208  JROOTClassSelector selector;
209  string background;
210  JPreprocessor option;
211  int debug;
212 
213  try {
214 
215  JParser<> zap("Auxiliary program to determine PMT parameters from K40 data.");
216 
217  zap['f'] = make_field(inputFile, "input file.");
218  zap['o'] = make_field(outputFile, "output file.") = "calibrate_k40.root";
219  zap['n'] = make_field(numberOfEvents) = JLimit::max();
220  zap['a'] = make_field(detectorFile, "detector file.");
221  zap['T'] = make_field(TMax_ns, "time window [ns].") = 20.0;
222  zap['V'] = make_field(rateRange_Hz, "PMT rate range [Hz].") = JRange<double>(0.0, 20.0e3);
223  zap['t'] = make_field(totRange_ns, "PMT time-over-threshold range [ns].") = JRange<double>(4.0, ToTmax_ns);
224  zap['b'] = make_field(background, "background estimation method.") = rates_t, count_t, tails_t, rndms_t;
225  zap['B'] = make_field(Tail_ns, "time window used for background estimation.") = JRange<double>(15.0, 20.0);
226  zap['M'] = make_field(multiplicity, "multiplicity range of hits on DOM.") = JRange<int>(2, 31);
227  zap['D'] = make_field(deadTime_us, "L1 dead time (us)") = 0.0;
228  zap['C'] = make_field(selector, "timeslice selector, e.g. JDAQTimesliceL1.") = getROOTClassSelection<JDAQTimesliceTypes_t>();
229  zap['O'] = make_field(option, "hit pre-processing option.") = JPreprocessor::getOptions();
230  zap['d'] = make_field(debug, "debug flag.") = 1;
231 
232  zap(argc, argv);
233  }
234  catch(const exception &error) {
235  FATAL(error.what() << endl);
236  }
237 
238  //-----------------------------------------------------------
239  // check the input parameters
240  //-----------------------------------------------------------
241 
242 
243  if (selector == JDAQTimesliceL2::Class_Name() ||
244  selector == JDAQTimesliceSN::Class_Name()) {
245  FATAL("Option -C <selector> " << selector << " not compatible with calibration method." << endl);
246  }
247 
248  if (selector == JDAQTimeslice ::Class_Name() ||
249  selector == JDAQTimesliceL1::Class_Name()) {
250 
252 
253  try {
254  parameters = getTriggerParameters(inputFile);
255  }
256  catch(const JException& error) {
257  FATAL("No trigger parameters from input:" << error.what() << endl);
258  }
259 
260  if ((selector == JDAQTimeslice ::Class_Name() && parameters.writeL1.prescale > 0) ||
261  (selector == JDAQTimesliceL1::Class_Name())) {
262 
263  if (parameters.TMaxLocal_ns < TMax_ns) {
264  FATAL("Option -T <TMax_ns> = " << TMax_ns << " is larger than in the trigger " << parameters.TMaxLocal_ns << endl);
265  }
266 
267  if (background == rndms_t ||
268  background == count_t) {
269  FATAL("Option -C <selector> " << selector << " incompatible with option -b <background> " << background << endl);
270  }
271  }
272  }
273 
274  if (!multiplicity.is_valid() || multiplicity.getLowerLimit() < 2 || multiplicity.getUpperLimit() > NUMBER_OF_PMTS) {
275  FATAL("Invalid option -M <multiplicity> " << multiplicity << endl);
276  }
277 
278  if (!totRange_ns.is_valid()) {
279  FATAL("Invalid option -t <totRange_ns> " << totRange_ns << endl);
280  }
281 
282  if (background == tails_t) {
283 
284  if (!Tail_ns.is_valid()) {
285  FATAL("Invalid option -B <Tail_ns> " << Tail_ns << endl);
286  }
287 
288  if (Tail_ns.getUpperLimit() > TMax_ns) {
289  FATAL("Invalid option -B <Tail_ns> " << Tail_ns << "; upper limit larger than option -T <TMax_ns> " << TMax_ns << endl);
290  }
291  }
292 
293  //-----------------------------------------------------------
294  // load detector file
295  //-----------------------------------------------------------
296 
298 
299  try {
300  load(detectorFile, detector);
301  }
302  catch(const JException& error) {
303  FATAL(error);
304  }
305 
306  if (detector.empty()) {
307  FATAL("Empty detector." << endl);
308  }
309 
310  const JModuleRouter router(detector);
311 
312  //-----------------------------------------------------------
313  // determine time offsets for background method rndms_t
314  //-----------------------------------------------------------
315 
316  vector<double> t0; // time offsets for random background evaluation [ns]
317 
318  for (int i = 0; i != NUMBER_OF_PMTS; ++i) {
319  t0.push_back(i * 2 * TMax_ns);
320  }
321 
322  //-----------------------------------------------------------
323  // correction factor for rate due to specified time-over-threshold range
324  //-----------------------------------------------------------
325 
326  const JPMTAnalogueSignalProcessor cpu;
327 
328  const int NPE = 1;
329  const double RTU = (cpu.getIntegralOfChargeProbability(cpu.getNPE(totRange_ns.getLowerLimit()), cpu.getNPE(totRange_ns.getUpperLimit()), NPE)
330  /
331  cpu.getIntegralOfChargeProbability(cpu.getNPE(0.0), cpu.getNPE(ToTmax_ns), NPE));
332 
333 
334  //-----------------------------------------------------------
335  // initialise histograms
336  //-----------------------------------------------------------
337 
338  vector<JHistogram> zmap(detector.size());
339 
340  const Double_t ymin = -floor(TMax_ns) + 0.5;
341  const Double_t ymax = +floor(TMax_ns) - 0.5;
342  const Int_t ny = (Int_t) ((ymax - ymin) / 1.0);
343 
344  for (JDetector::iterator module = detector.begin(); module != detector.end(); ++module) {
345 
346  const JModuleAddress& address = router.getAddress(module->getID());
347 
348  if (!module->empty()) {
349 
350  NOTICE("Booking histograms for module " << module->getID() << endl);
351 
352  zmap[address.first] = JHistogram(*module, JAbstractHistogram<Double_t>(ny, ymin, ymax));
353  }
354  }
355 
356 
357  typedef JHitR0 hit_type;
358  typedef JSuperFrame2D<hit_type> JSuperFrame2D_t;
359  typedef JSuperFrame1D<hit_type> JSuperFrame1D_t;
360 
361  const JMatchL0<hit_type> match(TMax_ns); // time window self-coincidences [ns]
362 
363  const double deadTime_ns = deadTime_us * 1e3;
364 
365 
366  TFile out(outputFile.c_str(), "recreate");
367 
368  putObject(out, JMeta(argc, argv));
369 
370  counter_type counter = 0;
371 
372  for (JMultipleFileScanner_t::const_iterator i = inputFile.begin(); i != inputFile.end(); ++i) {
373 
374  STATUS("Processing: " << *i << endl);
375 
376  JSummaryFileRouter summary(*i);
377 
380 
381  for (JDAQHeader header; in.hasNext() && counter != numberOfEvents; ++counter) {
382 
383  STATUS("event: " << setw(10) << counter << '\r'); DEBUG(endl);
384 
385  const JDAQTimeslice* timeslice = in.next();
386 
387  if (header.getDetectorID() != timeslice->getDetectorID() ||
388  header.getRunNumber () != timeslice->getRunNumber ()) {
389 
390  header = timeslice->getDAQHeader();
391 
392  putObject(out, header);
393  }
394 
395  if (background == rates_t) {
396 
397  summary.update(timeslice->getDAQHeader());
398 
399  if (timeslice->getFrameIndex() != summary.getSummaryslice()->getFrameIndex()) {
400 
401  ERROR("Frame indices do not match at "
402  << "[counter = " << counter << "]"
403  << "[timeslice = " << timeslice->getFrameIndex() << "]"
404  << "[summaryslice = " << summary.getSummaryslice()->getFrameIndex() << "]" << endl);
405 
406  continue;
407  }
408  }
409 
410  for (JDAQTimeslice::const_iterator frame = timeslice->begin(); frame != timeslice->end(); ++frame) {
411 
412  if (router.hasModule(frame->getModuleID())) {
413 
414  const JModule& module = router.getModule(frame->getModuleID());
415 
416  //-----------------------------------------------------------
417  // determine background rates and veto per PMT
418  //-----------------------------------------------------------
419 
420  vector<double> rate_Hz(NUMBER_OF_PMTS, 0.0);
421  vector<bool> veto (NUMBER_OF_PMTS, false);
422 
423  if (background == count_t) {
424 
425  for (JDAQSuperFrame::const_iterator i = frame->begin(); i != frame->end(); ++i) {
426  rate_Hz[i->getPMT()] += RTU * 1e9 / getFrameTime();
427  }
428 
429  } else if (background == tails_t) {
430 
431  // see below
432 
433  } else if (background == rndms_t) {
434 
435  // see below
436 
437  } else if (background == rates_t) {
438 
439  if (summary.hasSummaryFrame(frame->getModuleID())) {
440 
441  const JDAQSummaryFrame& sum = summary.getSummaryFrame(frame->getModuleID());
442 
443  for (int i = 0; i != NUMBER_OF_PMTS; ++i){
444  rate_Hz[i] = sum.getRate (i) *RTU;
445  veto [i] = sum.getValue(i) == 0;
446  }
447 
448  } else {
449 
450  FATAL("Summary frame of module " << frame->getModuleID() << " not found at frame index " << timeslice->getFrameIndex() << endl);
451  }
452  }
453 
454  const JDAQFrameStatus status = frame->getDAQFrameStatus();
455 
456  // Set veto according DAQ or PMT status and rate;
457  // Hits of PMTs with veto will not be counted in livetime nor coincidences nor background
458 
459  for (int i = 0; i != NUMBER_OF_PMTS; ++i) {
460  if (!getDAQStatus(status, module, i) ||
461  !getPMTStatus(status, module, i) ||
462  !rateRange_Hz(rate_Hz[i])) {
463  veto[i] = true;
464  }
465  }
466 
467  const JHistogram& histogram = zmap[router.getAddress(frame->getModuleID()).first];
468  const JCombinatorics& combinatorics = histogram.getCombinatorics();
469 
470  TH2D* h2s = histogram.h2s;
471  TH1D* h1b = histogram.h1b;
472  TH1D* h1L = histogram.h1L;
473 
474  //-----------------------------------------------------------
475  // process data
476  //-----------------------------------------------------------
477 
478  JSuperFrame2D_t& buffer = JSuperFrame2D_t::demultiplex(*frame, module);
479 
480  buffer.preprocess(option, match);
481 
482  for (JSuperFrame2D_t::iterator i = buffer.begin(); i != buffer.end(); ++i) {
483  if (veto[i->getPMTAddress()]) {
484  i->reset();
485  }
486  }
487 
488  JSuperFrame1D_t& data = JSuperFrame1D_t::multiplex(buffer);
489 
490  DEBUG("Number of hits " << timeslice->getFrameIndex() << ":" << frame->getModuleID() << ' ' << frame->size() << ' ' << data.size() << endl);
491 
492  // Signal;
493  // Hits from PMTs that do not comply with rate cuts have been exluded above
494 
495  size_t numberOfSignalEvents = 0;
496 
497  double t1 = numeric_limits<double>::lowest();
498 
499  for (vector<hit_type>::const_iterator p = data.begin(); p != data.end(); ) {
500 
502 
503  while (++q != data.end() && q->getT() - p->getT() <= TMax_ns) {}
504 
505  const int N = distance(p,q);
506 
507  if (multiplicity(N)) {
508 
509  numberOfSignalEvents += 1;
510 
511  if (p->getT() > t1 + deadTime_ns) {
512 
513  for (vector<hit_type>::const_iterator __p = p; __p != q; ++__p) {
514  for (vector<hit_type>::const_iterator __q = __p; ++__q != q; ) {
515 
516  if (totRange_ns( __p->getToT()) && totRange_ns(__q->getToT())) {
517  h2s->Fill((double) combinatorics.getIndex(__p->getPMT(),__q->getPMT()),
518  JCombinatorics::getSign(__p->getPMT(),__q->getPMT()) * (__q->getT() - __p->getT()));
519  }
520  }
521  }
522 
523  t1 = p->getT();
524  }
525  }
526 
527  p = q;
528  }
529 
530  // Background;
531  // Note that rate cuts and veto have already been accounted for when data buffer was filled
532 
533  if (background == rndms_t) {
534 
535  for (vector<hit_type>::iterator i = data.begin(); i != data.end(); ++i) {
536  *i = hit_type(i->getPMT(), JHit(i->getT() + t0[i->getPMT()], i->getToT()));
537  }
538 
539  sort(data.begin(), data.end());
540 
541  double t1 = numeric_limits<double>::lowest();
542 
543  for (vector<hit_type>::const_iterator p = data.begin(); p != data.end(); ) {
544 
546 
547  while (++q != data.end() && q->getT() - p->getT() <= TMax_ns) {}
548 
549  const int N = distance(p,q);
550 
551  if (multiplicity(N)) {
552 
553  if (p->getT() > t1 + deadTime_ns) {
554 
555  for (vector<hit_type>::const_iterator __p = p; __p != q; ++__p) {
556  for (vector<hit_type>::const_iterator __q = __p; ++__q != q; ) {
557 
558  if (totRange_ns(__p->getToT()) && totRange_ns(__q->getToT())) {
559  h1b->Fill((double) combinatorics.getIndex(p->getPMT(),q->getPMT()), 1.0);
560  }
561  }
562  }
563 
564  t1 = p->getT();
565  }
566  }
567 
568  p = q;
569  }
570 
571  } else if (background == rates_t ||
572  background == count_t) {
573 
574  double Rs = 0.0; // total rate [Hz]
575 
576  for (int i = 0; i != NUMBER_OF_PMTS; ++i) {
577  if (!veto[i]) {
578  Rs += rate_Hz[i]; // [Hz]
579  }
580  }
581 
582  for (int i = 0; i != NUMBER_OF_PMTS; ++i) {
583  for (int j = i; ++j != NUMBER_OF_PMTS; ) {
584 
585  if (!veto[i] && !veto[j]) {
586 
587  const double R1 = rate_Hz[i]; // [Hz]
588  const double R2 = rate_Hz[j]; // [Hz]
589 
590  // evaluate expected counts within a time window of 2 * TMax_ns for the two PMTs and the other PMTs inside the optical module;
591  // expected rate due to random coincidences then is product of the probability for the specific observation and the number of trials
592  // the observation is one hit in PMT 1, one hit in PMT 2 and a number of hits in the other PMTs inside the same module
593  // corresponding to the given multiplicity range
594 
595  const double N = getP((R1) * 2 * TMax_ns * 1e-9,
596  (R2) * 2 * TMax_ns * 1e-9,
597  (Rs - R1 - R2) * 2 * TMax_ns * 1e-9,
598  multiplicity.getLowerLimit(),
599  multiplicity.getUpperLimit()) * getFrameTime() / (2*TMax_ns);
600 
601  h1b->Fill((double) combinatorics.getIndex(i,j), N);
602  }
603  }
604  }
605  }
606 
607  //-----------------------------------------------------------
608  // fill the livetime by PMT pairs
609  //-----------------------------------------------------------
610 
611  const double livetime_s = getFrameTime()*1e-9 * exp(-deadTime_ns * numberOfSignalEvents/getFrameTime());
612 
613  for (size_t i = 0; i < combinatorics.getNumberOfPairs(); ++i) {
614 
615  const int pmt1 = combinatorics.getPair(i).first;
616  const int pmt2 = combinatorics.getPair(i).second;
617 
618  if (!veto[pmt1] && !veto[pmt2]) {
619  h1L->Fill(i, livetime_s);
620  }
621  }
622  }
623  }
624  }
625  }
626  STATUS(endl);
627 
628  if (background == tails_t ) {
629 
630  for (vector<JHistogram>::iterator hist = zmap.begin() ; hist != zmap.end() ; ++hist) {
631 
632  if (hist->is_valid()) {
633 
634  JCombinatorics& combinatorics = *hist;
635 
636  TH2D* h2s = hist->h2s;
637  TH1D* h1b = hist->h1b;
638 
639  for (int i = 0; i != NUMBER_OF_PMTS; ++i) {
640  for (int j = i; ++j != NUMBER_OF_PMTS; ) {
641 
642  double Y = 0.0; // sum of counts in tail regions
643  double W = 0.0; // square of error of Y
644  int N = 0; // number of bins in tail regions
645  int k = combinatorics.getIndex(i,j); // PMT pair index (x-axis of h2s)
646 
647  for (int l = 1; l <= h2s->GetYaxis()->GetNbins(); ++l) {
648 
649  const Double_t dt = h2s->GetYaxis()->GetBinCenter(l) ;
650 
651  if (Tail_ns(fabs(dt))) {
652  Y += h2s->GetBinContent(k+1,l);
653  W += h2s->GetBinError(k+1,l) * h2s->GetBinError(k+1,l);
654  N++;
655  }
656  }
657 
658  h1b->SetBinContent(k+1, Y / N * (2.0*TMax_ns) / ((ymax - ymin)/ny) );
659  h1b->SetBinError (k+1, sqrt(W/N) * (2.0*TMax_ns) / ((ymax - ymin)/ny) );
660  }
661  }
662  }
663  }
664  }
665 
666  //---------------------------------------------
667  // save normalisation constants and store histograms
668  //---------------------------------------------
669 
670  TH1D weights_hist(weights_hist_t, NULL, 3, 0.0, 3.0);
671 
672  weights_hist.GetXaxis()->SetBinLabel(1, W1_overall_t); // [s]
673  weights_hist.GetXaxis()->SetBinLabel(2, WS_t); // [ns]
674  weights_hist.GetXaxis()->SetBinLabel(3, WB_t); // [ns]
675 
676  weights_hist.Fill(W1_overall_t, counter*getFrameTime()*1e-9); // [s]
677  weights_hist.Fill(WS_t, (ymax - ymin)/ny); // [ns]
678  weights_hist.Fill(WB_t, 2.0*TMax_ns); // [ns]
679 
680  out << weights_hist;
681 
682  for (vector<JHistogram>::iterator i = zmap.begin(); i != zmap.end(); ++i) {
683  if (i->is_valid()) {
684  out << *(i->h2s);
685  out << *(i->h1b);
686  out << *(i->h1L);
687  }
688  }
689 
690  out.Write();
691  out.Close();
692 }
Auxiliary class for ROOT I/O of application specific meta data.
Definition: JMeta.hh:70
const double xmax
Definition: JQuadrature.cc:24
Utility class to parse command line options.
Definition: JParser.hh:1514
Auxiliaries for pre-processing of hits.
General exception.
Definition: JException.hh:24
then fatal No hydrophone data file $HYDROPHONE_TXT fi sort gr k
Auxiliary class to convert pair of indices to unique index and back.
JAbscissa_t getLowerLimit() const
Get lower limit.
Definition: JRange.hh:202
int main(int argc, char *argv[])
Definition: Main.cc:15
ROOT TTree parameter settings of various packages.
double getRate(const int tdc, const double factor=1.0) const
Get count rate.
Basic data structure for L0 hit.
const JModule & getModule(const JObjectID &id) const
Get module parameters.
Data structure for a composite optical module.
Definition: JModule.hh:68
std::vector< T >::difference_type distance(typename std::vector< T >::const_iterator first, typename PhysicsEvent::const_iterator< T > second)
Specialisation of STL distance.
int getDetectorID() const
Get detector identifier.
static const char *const WB_t
Named bin for value of TMax_ns in JCalibrateK40.
L0 match criterion.
Definition: JMatchL0.hh:27
#define STATUS(A)
Definition: JMessage.hh:63
bool is_valid(const json &js)
Check validity of JSon data.
void update(const JDAQHeader &header)
Update router.
Detector data structure.
Definition: JDetector.hh:89
Auxiliary class to select ROOT class based on class name.
Router for direct addressing of module data in detector data structure.
#define R1(x)
*fatal Wrong number of arguments esac JCookie sh typeset Z DETECTOR typeset Z SOURCE_RUN typeset Z TARGET_RUN set_variable PARAMETERS_FILE $WORKDIR parameters
Definition: diff-Tuna.sh:38
#define MAKE_CSTRING(A)
Make C-string.
Definition: JPrint.hh:136
Long64_t counter_type
Type definition for counter.
then warning Cannot perform comparison test for histogram
then fatal Wrong number of arguments fi set_variable STRING $argv[1] set_variable DETECTORXY_TXT $WORKDIR $DETECTORXY_TXT tail read X Y CHI2 RMS printf optimum n $X $Y $CHI2 $RMS awk v Y
const JDAQSummaryFrame & getSummaryFrame() const
Get default summary frame.
Auxiliary class for multiplexing object iterators.
Simple data structure for histogram binning.
string outputFile
int getNumberOfBins() const
Get number of bins.
int getRunNumber() const
Get run number.
Data structure for detector geometry and calibration.
Tools for handling different hit types.
Acoustics hit.
then echo The file $DIR KM3NeT_00000001_00000000 root already please rename or remove it first
int getFrameIndex() const
Get frame index.
1-dimensional frame with time calibrated data from one optical module.
Q JDAQEvent livetime_s
Definition: JDataQuality.sh:57
int first
index of module in detector data structure
Scanning of objects from a single file according a format that follows from the extension of each fil...
static const char *const WS_t
Named bin for time residual bin width.
I/O formatting auxiliaries.
Auxiliary class to sort pairs of PMT addresses within optical module.
Detector file.
Definition: JHead.hh:226
Reduced data structure for L0 hit.
Definition: JHitR0.hh:25
const JDAQSummaryslice * getSummaryslice() const
Get summary slice.
Hit data structure.
Definition: JDAQHit.hh:34
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1989
int getIndex(const int first, const int second) const
Get index of pair of indices.
static const char *const _1B
Name extension for 1D background.
double getFrameTime()
Get frame time duration.
Definition: JDAQClock.hh:162
Data storage class for rate measurements of all PMTs in one module.
int getID() const
Get identifier.
Definition: JObjectID.hh:50
ROOT I/O of application specific meta data.
JAbscissa_t getUpperLimit() const
Get upper limit.
Definition: JRange.hh:213
#define NOTICE(A)
Definition: JMessage.hh:64
#define ERROR(A)
Definition: JMessage.hh:66
Data time slice.
File router for fast addressing of summary data.
virtual const pointer_type & next() override
Get next element.
Address of module in detector data structure.
Match operator for consecutive hits.
virtual bool hasNext() override
Check availability of next element.
bool putObject(TDirectory &dir, const TObject &object)
Write object to ROOT directory.
General purpose messaging.
double getP(const double expval, bool hit)
Get Poisson probability to observe a hit or not for given expectation value for the number of hits...
Definition: JFitToolkit.hh:41
#define FATAL(A)
Definition: JMessage.hh:67
Scanning of objects from multiple files according a format that follows from the extension of each fi...
Direct access to module in detector data structure.
static const char *const weights_hist_t
Histogram naming.
static const char *const _1L
Name extension for 1D live time.
const double xmin
Definition: JQuadrature.cc:23
Auxiliary base class for list of file names.
virtual const char * what() const override
Get error message.
Definition: JException.hh:64
bool getPMTStatus(const JStatus &status)
Test status of PMT.
then usage $script< input file >[option[primary[working directory]]] nWhere option can be N
Definition: JMuonPostfit.sh:40
void configure(const T &value, const JAbstractCollection< JAbscissa_t > &bounds, JBool< false > option)
Configuration of value.
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
Auxiliary class to define a range between two values.
Utility class to parse command line options.
bool hasSummaryFrame(const JDAQModuleIdentifier &module) const
Has summary frame.
const JDAQHeader & getDAQHeader() const
Get DAQ header.
Definition: JDAQHeader.hh:49
static const char *const W1_overall_t
Named bin for duration of the run.
const JModuleAddress & getAddress(const JObjectID &id) const
Get address of module.
bool hasModule(const JObjectID &id) const
Has module.
2-dimensional frame with time calibrated data from one optical module.
PMT analogue signal processor.
Object reading from a list of files.
int j
Definition: JPolint.hh:703
do set_variable DETECTOR_TXT $WORKDIR detector
bool getDAQStatus(const JDAQFrameStatus &frame, const JStatus &status)
Test status of DAQ.
KM3NeT DAQ constants, bit handling, etc.
static const int NUMBER_OF_PMTS
Total number of PMTs in module.
Definition: JDAQ.hh:26
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 JAcoustics sh $DETECTOR_ID source JAcousticsToolkit sh CHECK_EXIT_CODE typeset A EMITTERS get_tripods $WORKDIR tripod txt EMITTERS get_transmitters $WORKDIR transmitter txt EMITTERS for EMITTER in
Definition: JCanberra.sh:46
Auxiliary class for specifying the way of pre-processing of hits.
const JDAQFrameStatus & getDAQFrameStatus() const
Get DAQ frame status.
then fatal Wrong number of arguments fi set_variable DETECTOR $argv[1] set_variable STRING $argv[2] set_array QUANTILES set_variable FORMULA *[0] exp(-0.5 *(x-[1])*(x-[1])/([2]*[2]))" set_variable MODULE `getModule -a $DETECTOR -L "$STRING 0"` source JAcousticsToolkit.sh typeset -A TRIPODS get_tripods $WORKDIR/tripod.txt TRIPODS XMEAN
JTriggerParameters getTriggerParameters(const JMultipleFileScanner_t &file_list)
Get trigger parameters.
int debug
debug level
double getValue(const int tdc) const
Get value.
then $DIR JPlotNPE PDG P
Definition: JPlotNPE-PDG.sh:62
static const char *const _2S
Name extension for 2D counts.
#define DEBUG(A)
Message macros.
Definition: JMessage.hh:62
Template definition of histogram object interface.
Definition: JHistogram.hh:27