Jpp  18.2.1
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JAAnet/JEvD.cc
Go to the documentation of this file.
1 #include <string>
2 #include <iostream>
3 #include <iomanip>
4 #include <vector>
5 #include <algorithm>
6 
7 #include "TROOT.h"
8 #include "TApplication.h"
9 #include "TCanvas.h"
10 #include "TStyle.h"
11 #include "TH2D.h"
12 #include "TArrow.h"
13 #include "TLatex.h"
14 #include "TMarker.h"
15 
21 
22 #include "JROOT/JStyle.hh"
23 #include "JROOT/JCanvas.hh"
24 #include "JROOT/JRootToolkit.hh"
25 
27 
28 #include "JTrigger/JHitL0.hh"
29 
30 #include "JSupport/JSupport.hh"
34 #include "JAAnet/JAAnetToolkit.hh"
35 
36 #include "JPhysics/JPDF_t.hh"
37 
38 #include "JFit/JLine1Z.hh"
39 #include "JFit/JModel.hh"
42 
43 #include "JLang/JPredicate.hh"
44 #include "JLang/JComparator.hh"
45 
46 #include "JMath/JMathToolkit.hh"
47 #include "JSystem/JKeypress.hh"
48 #include "JSystem/JProcess.hh"
49 
50 #include "Jeep/JFunctionAdaptor.hh"
51 #include "Jeep/JPrint.hh"
52 #include "Jeep/JParser.hh"
53 #include "Jeep/JMessage.hh"
54 
55 
56 namespace JAANET {
57 
59 
60  /**
61  * Event selector.
62  *
63  * The default constructor will accept all events.\n
64  * A different method can dynamically be loaded from a shared library using class JEEP::JFunctionAdaptor.
65  */
66  struct JEventSelector :
67  public JFunctionAdaptor<bool, const Trk&, const Evt&>
68  {
69  /**
70  * Default event selection.
71  *
72  * \param trk track
73  * \param evt event
74  * \return true
75  */
76  static inline bool select(const Trk& trk, const Evt& evt)
77  {
78  return true;
79  }
80 
81 
82  /**
83  * Default constructor.
84  */
86  {
87  this->function = select;
88  this->symbol = "select";
89  }
90  };
91 
92 
93  /**
94  * Check availability of value.
95  *
96  * \param trk track
97  * \param i index
98  * \return true if available; else false
99  */
100  inline bool hasW(const Trk& trk, const int i)
101  {
102  return (i >= 0 && i < (int) trk.fitinf.size());
103  }
104 
105 
106  /**
107  * Get associated value.
108  *
109  * \param trk track
110  * \param i index
111  * \return value
112  */
113  inline double getW(const Trk& trk, const int i)
114  {
115  return trk.fitinf.at(i);
116  }
117 
118 
119  /**
120  * Get associated value.
121  *
122  * \param trk track
123  * \param i index
124  * \param value default value
125  * \return value
126  */
127  inline double getW(const Trk& trk, const int i, const double value)
128  {
129  if (hasW(trk,i))
130  return trk.fitinf.at(i);
131  else
132  return value;
133  }
134 
135 
136  /**
137  * Set associated value.
138  *
139  * \param trk track
140  * \param i index
141  * \param value value
142  */
143  void setW(Trk& trk, const int i, const double value)
144  {
145  if (i >= (int) trk.fitinf.size()) {
146  trk.fitinf.resize(i + 1, 0.0);
147  }
148 
149  trk.fitinf[i] = value;
150  }
151 }
152 
153 namespace {
154 
155  /**
156  * Wild card character for file name substition.
157  */
158  const char WILDCARD = '%';
159 
160 
161  /**
162  * Execute command in shell.
163  *
164  * \param command command
165  */
166  inline void execute(const std::string& command, int debug)
167  {
168  using namespace std;
169  using namespace JPP;
170 
171  JProcess process(command);
172 
173  istream in(process.getInputStreamBuffer());
174 
175  for (string buffer; getline(in, buffer); ) {
176  DEBUG(buffer << endl);
177  }
178  }
179 }
180 
181 
182 /**
183  * \file
184  *
185  * Program to display hit probabilities.
186  *
187  * \author mdejong
188  */
189 int main(int argc, char **argv)
190 {
191  using namespace std;
192  using namespace JPP;
193  using namespace KM3NETDAQ;
194 
195  JSingleFileScanner<Evt> inputFile;
196  JLimit_t& numberOfEvents = inputFile.getLimit();
197  string pdfFile;
198  string outputFile;
200  int application;
201  JEventSelector event_selector;
202  JCanvas canvas;
203  bool batch;
204  double arrowSize;
205  string arrowType;
206  double arrowScale;
207  int debug;
208 
209 
210  try {
211 
212  parameters.numberOfPrefits = 1;
213 
214  JParser<> zap("Program to display hit probabilities.");
215 
216  zap['w'] = make_field(canvas, "size of canvas <nx>x<ny> [pixels]") = JCanvas(1200, 600);
217  zap['f'] = make_field(inputFile, "input file (output of JXXXMuonReconstruction.sh)");
218  zap['n'] = make_field(numberOfEvents) = JLimit::max();
219  zap['P'] = make_field(pdfFile);
220  zap['o'] = make_field(outputFile, "graphics output file name") = MAKE_STRING("display_" << WILDCARD << ".gif");
223  zap['L'] = make_field(event_selector) = JPARSER::initialised();
224  zap['S'] = make_field(arrowSize) = 0.003;
225  zap['T'] = make_field(arrowType) = "|->";
226  zap['F'] = make_field(arrowScale) = 250.0;
227  zap['B'] = make_field(batch, "batch processing");
228  zap['d'] = make_field(debug) = 1;
229 
230  zap(argc, argv);
231  }
232  catch(const exception& error) {
233  FATAL(error.what() << endl);
234  }
235 
236  if (batch && outputFile == "") {
237  FATAL("Missing output file name " << outputFile << " in batch mode." << endl);
238  }
239 
240  if (!batch && outputFile == "") {
241  outputFile = MAKE_STRING(WILDCARD << ".gif");
242  }
243 
244  if (outputFile.find(WILDCARD) == string::npos) {
245  FATAL("Output file name " << outputFile << " has no wild card '" << WILDCARD << "'" << endl);
246  }
247 
248 
249  JSummaryFileRouter summary(inputFile, parameters.R_Hz);
250 
251  const JMuonPDF_t pdf(pdfFile, parameters.TTS_ns);
252 
253  const JTimeRange T_ns(parameters.TMin_ns, parameters.TMax_ns);
254 
255  typedef vector<JHitL0> JDataL0_t;
256  typedef vector<JHitW0> JDataW0_t;
257 
258 
259  // ROOT
260 
261  gROOT->SetBatch(batch);
262 
263  TApplication* tp = new TApplication("user", NULL, NULL);
264  TCanvas* cv = new TCanvas("display", "", canvas.x, canvas.y);
265 
266  JSinglePointer<TStyle> gStyle(new JStyle("gplot", cv->GetWw(), cv->GetWh()));
267 
268  gROOT->SetStyle("gplot");
269  gROOT->ForceStyle();
270 
271  const size_t NUMBER_OF_PADS = 3;
272 
273  cv->SetFillStyle(4000);
274  cv->SetFillColor(kWhite);
275  cv->Divide(NUMBER_OF_PADS, 1);
276 
277  const double Dmax = 1000.0;
278  const double Rmin = 0.0;
279  const double Rmax = min(parameters.roadWidth_m, 0.4 * Dmax);
280  const double Tmin = min(parameters.TMin_ns, -10.0);
281  const double Tmax = max(parameters.TMax_ns, +100.0);
282  const double Amin = 0.002 * (Tmax - Tmin); // minimal arrow length [ns]
283  const double Amax = 0.8 * (Tmax - Tmin); // maximal arrow length [ns]
284  const double ymin = min(-Amax, Tmin - 0.3 * Amax);
285  const double ymax = max(+Amax, Tmax + 0.5 * Amax);
286 
287  const string Xlabel[NUMBER_OF_PADS] = { "R [m]", "#phi [rad]", "z [m]" };
288  const double Xmin [NUMBER_OF_PADS] = { Rmin, -PI, -0.5 * Dmax };
289  const double Xmax [NUMBER_OF_PADS] = { Rmax, +PI, +0.5 * Dmax };
290 
291  double Xs[NUMBER_OF_PADS];
292 
293  for (size_t i = 0; i != NUMBER_OF_PADS; ++i) {
294  Xs[i] = 0.003 * (Xmax[i] - Xmin[i]) * (0.5 * NUMBER_OF_PMTS); // x-offset arrow as function of PMT number
295  }
296 
297  TH2D H2[NUMBER_OF_PADS];
298  TGraph G2[NUMBER_OF_PADS];
299 
300  for (size_t i = 0; i != NUMBER_OF_PADS; ++i) {
301 
302  H2[i] = TH2D(MAKE_CSTRING("h" << i), NULL, 100, Xmin[i] - Xs[i], Xmax[i] + Xs[i], 100, ymin, ymax);
303 
304  H2[i].GetXaxis()->SetTitle(Xlabel[i].c_str());
305  H2[i].GetYaxis()->SetTitle("#Deltat [ns]");
306 
307  H2[i].GetXaxis()->CenterTitle(true);
308  H2[i].GetYaxis()->CenterTitle(true);
309 
310  H2[i].SetStats(kFALSE);
311 
312  G2[i].Set(2);
313 
314  G2[i].SetPoint(0, H2[i].GetXaxis()->GetXmin(), 0.0);
315  G2[i].SetPoint(1, H2[i].GetXaxis()->GetXmax(), 0.0);
316 
317  cv->cd(i+1);
318 
319  H2[i].Draw("AXIS");
320  G2[i].Draw("SAME");
321  }
322 
323 
324  while (inputFile.hasNext()) {
325 
326  cout << "\revent: " << setw(8) << inputFile.getCounter() << flush;
327 
328  Evt* evt = inputFile.next();
329 
330  if (has_reconstructed_track<JPP_RECONSTRUCTION_TYPE>(*evt, rec_stages_range(application))) {
331 
332  Trk fit = get_best_reconstructed_track<JPP_RECONSTRUCTION_TYPE>(*evt, rec_stages_range(application));
333 
334  if (!event_selector(fit, *evt)) {
335  continue;
336  }
337 
338 
339  JDataL0_t dataL0;
340 
341  for (const Hit& hit : evt->hits) {
342  dataL0.push_back(JHitL0(JDAQPMTIdentifier(hit.dom_id, hit.channel_id),
343  JAxis3D(getPosition(hit.pos), getDirection(hit.dir)),
344  JHit(hit.t, hit.tot)));
345  }
346 
347  summary.update(JDAQChronometer(evt->det_id,
348  evt->run_id,
349  evt->frame_index,
350  JDAQUTCExtended(evt->t.GetSec(), evt->t.GetNanoSec() / 16)));
351 
352  const time_converter converter = time_converter(*evt);
353 
354  Trk muon; // Monte Carlo true muon
355 
356  if (has_muon(*evt)) {
357 
358  for (const auto& trk : evt->mc_trks) {
359  if (is_muon(trk)) {
360  if (trk.E > muon.E) {
361 
362  muon = trk;
363  muon.t += converter.putTime();
364 
365  setW(muon, JSTART_LENGTH_METRES, fabs(muon.len));
366  }
367  }
368  }
369  }
370 
371 
372  bool next_event = false; // goto next event
373  bool monte_carlo = false; // show Monte Carlo true muon
374 
375  while (!next_event) {
376 
377  Trk trk;
378 
379  if (!monte_carlo)
380  trk = fit;
381  else
382  trk = muon;
383 
384  JRotation3D R (getDirection(trk));
385  JLine1Z tz(getPosition (trk).rotate(R), trk.t);
386  JRange<double> Z_m;
387  /*
388  if (hasW(trk, JSTART_LENGTH_METRES)) {
389  Z_m = JRange<double>(0.0, getW(fit,JSTART_LENGTH_METRES)) + JRange<double>(parameters.ZMin_m, parameters.ZMax_m);
390  }
391  */
392  const JModel<JLine1Z> match(tz, parameters.roadWidth_m, T_ns, Z_m);
393 
394  // hit selection based on fit result
395 
396  JDataW0_t data;
397 
398  for (JDataL0_t::const_iterator i = dataL0.begin(); i != dataL0.end(); ++i) {
399 
400  JHitW0 hit(*i, summary.getRate(i->getPMTIdentifier()));
401 
402  hit.rotate(R);
403 
404  if (match(hit)) {
405  data.push_back(hit);
406  }
407  }
408 
409  // select first hit in PMT
410 
411  sort(data.begin(), data.end(), JHitW0::compare);
412 
413  JDataW0_t::iterator __end = unique(data.begin(), data.end(), equal_to<JDAQPMTIdentifier>());
414 
415  double E_GeV = parameters.E_GeV;
416  /*
417  if (trk.E > 0.1) {
418  E_GeV = trk.E;
419  }
420  */
421 
422  // move fit to geometrical center of hits
423 
424  JRange<double> zs(make_array(data.begin(), __end, &JHitW0::getZ));
425 
426  const double z0 = tz.getZ();
427  const double z1 = 0.5 * (zs.getLowerLimit() + zs.getUpperLimit());
428 
429  tz.setZ(z1, getSpeedOfLight());
430 
431 
432  // graphics
433 
434  TLatex latex [NUMBER_OF_PADS];
435  vector<TArrow> arrow [NUMBER_OF_PADS];
436  vector<TMarker> marker[NUMBER_OF_PADS];
437 
438  for (int i = 0; i != NUMBER_OF_PADS; ++i) {
439 
440  latex[i].SetTextAlign(12);
441  latex[i].SetTextFont(42);
442  latex[i].SetTextSize(0.06);
443 
444  latex[i].SetX(H2[i].GetXaxis()->GetXmin() + 0.05 * (H2[i].GetXaxis()->GetXmax() - H2[i].GetXaxis()->GetXmin()));
445  latex[i].SetY(H2[i].GetYaxis()->GetXmax() - 0.05 * (H2[i].GetYaxis()->GetXmax() - H2[i].GetYaxis()->GetXmin()));
446  }
447 
448  if (hasW(trk, JSTART_LENGTH_METRES) && getW(trk, JSTART_LENGTH_METRES) > 0.0) {
449 
450  marker[2].push_back(TMarker(z0 - tz.getZ(), 0.0, kFullCircle));
451  marker[2].push_back(TMarker(z0 - tz.getZ() + getW(trk, JSTART_LENGTH_METRES), 0.0, kFullCircle));
452 
453  static_cast<TAttMarker&>(marker[2][0]) = TAttMarker(kRed, kFullCircle, 0.7);
454  static_cast<TAttMarker&>(marker[2][1]) = TAttMarker(kRed, kFullCircle, 0.7);
455  }
456 
457  double chi2 = 0;
458 
459  for (JDataW0_t::const_iterator hit = data.begin(); hit != __end; ++hit) {
460 
461  const double x = hit->getX() - tz.getX();
462  const double y = hit->getY() - tz.getY();
463  const double z = hit->getZ() - tz.getZ();
464  const double R = sqrt(x*x + y*y);
465 
466  const double t1 = tz.getT() + (z + R * getTanThetaC()) * getInverseSpeedOfLight();
467 
468  JDirection3D dir(hit->getDX(), hit->getDY(), hit->getDZ()); // PMT orientation
469 
470  dir.rotate(JRotation3Z(-atan2(y,x))); // rotate PMT axis to x-z plane
471 
472  const double theta = dir.getTheta();
473  const double phi = fabs(dir.getPhi()); // rotational symmetry of Cherenkov cone
474 
475  //const double E = gWater.getE(E_GeV, z); // correct for energy loss
476  const double E = E_GeV;
477  const double dt = T_ns.constrain(hit->getT() - t1);
478 
479  JMuonPDF_t::result_type H1 = pdf.calculate(E, R, theta, phi, dt);
480  JMuonPDF_t::result_type H0(hit->getR() * 1e-9, 0.0, T_ns);
481 
482  if (H1.V >= parameters.VMax_npe) {
483  H1 *= parameters.VMax_npe / H1.V;
484  }
485 
486  H1 += H0; // signal + background
487 
488  chi2 += H1.getChi2() - H0.getChi2();
489 
490  const double derivative = H1.getDerivativeOfChi2() - H0.getDerivativeOfChi2();
491 
492  double size = derivative * arrowScale; // size of arrow
493 
494  if (fabs(size) < Amin) {
495  size = (size > 0.0 ? +Amin : -Amin);
496  } else if (fabs(size) > Amax) {
497  size = (size > 0.0 ? +Amax : -Amax);
498  }
499 
500  const double X[NUMBER_OF_PADS] = { R, atan2(y,x), z - R/getTanThetaC() };
501 
502  const double xs = (double) (NUMBER_OF_PMTS - 2 * hit->getPMTAddress()) / (double) NUMBER_OF_PMTS;
503 
504  for (size_t i = 0; i != NUMBER_OF_PADS; ++i) {
505  arrow[i].push_back(TArrow(X[i] + xs*Xs[i], dt, X[i] + xs*Xs[i], dt + size, arrowSize, arrowType.c_str()));
506  }
507  }
508 
509  latex[0].SetTitle(MAKE_CSTRING("" << FILL(6,'0') << evt->run_id << FILL() <<
510  ":" << evt->frame_index <<
511  " " << evt->trigger_counter));
512 
513  if (!monte_carlo)
514  latex[1].SetTitle(MAKE_CSTRING("Q = " << FIXED(4,0) << -chi2 << " " <<
515  "E = " << SCIENTIFIC(7,1) << trk.E << " [GeV]"));
516 
517  else
518  latex[1].SetTitle(MAKE_CSTRING("Q = " << FIXED(4,0) << -chi2));
519 
520  if (is_muon(muon)) {
521 
522  if (!monte_carlo)
523  latex[2].SetTitle(MAKE_CSTRING("#Delta#alpha = " << FIXED(6,2) << getAngle(getDirection(muon), getDirection(trk)) << " [deg]" <<
524  " z = " << FIXED(6,3) << trk.dir.z));
525  else
526  latex[2].SetTitle("Monte Carlo");
527  }
528 
529 
530  // draw
531 
532  for (int i = 0; i != NUMBER_OF_PADS; ++i) {
533 
534  cv->cd(i+1);
535 
536  latex[i].Draw();
537 
538  for (auto& a1 : arrow[i]) {
539  a1.Draw();
540  }
541 
542  for (auto& m1 : marker[i]) {
543  m1.Draw();
544  }
545  }
546 
547  cv->Update();
548 
549 
550  // action
551 
552  if (batch) {
553 
554  cv->SaveAs(replace(outputFile, WILDCARD, MAKE_STRING(inputFile.getCounter())).c_str());
555 
556  next_event = true;
557 
558  } else {
559 
560  for (int c = 0; c == 0; ) {
561 
562  if (!batch) {
563 
564  static bool first = true;
565 
566  if (first) {
567  cout << endl << "Type '?' for possible options." << endl;
568  }
569 
570  first = false;
571  }
572 
573  cout << "\n> " << flush;
574 
575  switch((c = JKeypress(true).get())) {
576 
577  case '?':
578  cout << endl;
579  cout << "possible options: " << endl;
580  cout << 'q' << " -> " << "exit application" << endl;
581  cout << 'u' << " -> " << "update canvas" << endl;
582  cout << 's' << " -> " << "save graphics to file" << endl;
583  cout << 'M' << " -> " << "Monte Carlo true muon information" << endl;
584  cout << 'F' << " -> " << "fit information" << endl;
585  if (event_selector.is_valid()) {
586  cout << 'L' << " -> " << "reload event selector" << endl;
587  }
588  cout << 'r' << " -> " << "rewind input file" << endl;
589  cout << 'R' << " -> " << "switch to ROOT mode (quit ROOT to continue)" << endl;
590  cout << 'p' << " -> " << "print event information" << endl;
591  cout << ' ' << " -> " << "next event (as well as any other key)" << endl;
592  c = 0;
593  break;
594 
595  case 'q':
596  cout << endl;
597  return 0;
598 
599  case 'u':
600  cv->Update();
601  c = 0;
602  break;
603 
604  case 's':
605  cv->SaveAs(replace(outputFile, WILDCARD, MAKE_STRING(inputFile.getCounter())).c_str());
606  c = 0;
607  break;
608 
609  case 'M':
610  if (is_muon(muon))
611  monte_carlo = true;
612  else
613  ERROR(endl << "No Monte Carlo muon available." << endl);
614  break;
615 
616  case 'F':
617  monte_carlo = false;
618  break;
619 
620  case 'r':
621  inputFile.rewind();
622  next_event = true;
623  break;
624 
625  case 'L':
626  if (event_selector.is_valid()) {
627  execute(MAKE_STRING("make -f " << getPath(argv[0]) << "/JMakeEventSelector libs"), 3);
628  event_selector.reload();
629  }
630  break;
631 
632  case 'R':
633  tp->Run(kTRUE);
634  c = 0;
635  break;
636 
637  case 'p':
638  cout << endl;
639  evt->print(cout);
640  if (debug >= debug_t) {
641  for (const auto& trk : evt->mc_trks) {
642  cout << "MC "; trk.print(cout); cout << endl;
643  }
644  for (const auto& trk : evt->trks) {
645  cout << "fit "; trk.print(cout); cout << endl;
646  }
647  for (const auto& hit : evt->hits) {
648  cout << "hit "; hit.print(cout); cout << endl;
649  }
650  }
651  c = 0;
652  break;
653 
654  default:
655  next_event = true;
656  break;
657  }
658  }
659  }
660  }
661  }
662  }
663  cout << endl;
664 }
static const int JMUONSTART
Utility class to parse command line options.
Definition: JParser.hh:1514
double getAngle(const JQuaternion3D &first, const JQuaternion3D &second)
Get space angle between quanternions.
debug
Definition: JMessage.hh:29
T getLowerLimit() const
Get lower limit.
Definition: JRange.hh:202
int main(int argc, char *argv[])
Definition: Main.cc:15
Vec pos
hit position
Definition: Hit.hh:25
ROOT TTree parameter settings of various packages.
Auxiliary methods for geometrical methods.
TString replace(const TString &target, const TRegexp &regexp, const T &replacement)
Replace regular expression in input by given replacement.
Definition: JPrintResult.cc:63
Data structure for direction in three dimensions.
Definition: JDirection3D.hh:33
then usage $script< input file >[option[primary[working directory]]] nWhere option can be E
Definition: JMuonPostfit.sh:40
double t
track time [ns] (when the particle is at pos )
Definition: Trk.hh:19
double z
Definition: Vec.hh:14
JEventSelector()
Default constructor.
Definition: JAAnet/JEvD.cc:85
std::string getPath(const std::string &file_name)
Get path, i.e. part before last JEEP::PATHNAME_SEPARATOR if any.
Definition: JeepToolkit.hh:148
void update(const JDAQHeader &header)
Update router.
bool is_muon(const Trk &track)
Test whether given track is a (anti-)muon.
void setW(Trk &trk, const int i, const double value)
Set associated value.
Definition: JAAnet/JEvD.cc:143
Rotation matrix.
Definition: JRotation3D.hh:111
Vec dir
track direction
Definition: Trk.hh:18
double getRate() const
Get default rate.
*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
Template specialisation of class JModel to match hit with muon trajectory along z-axis.
Definition: JFit/JModel.hh:34
double getZ(const JPosition3D &pos) const
Get point of emission of Cherenkov light along muon path.
Definition: JLine1Z.hh:134
Empty structure for specification of parser element that is initialised (i.e. does not require input)...
Definition: JParser.hh:83
Auxiliary class to convert DAQ hit time to/from Monte Carlo hit time.
Auxiliary data structure for floating point format specification.
Definition: JManip.hh:446
string outputFile
Data structure for UTC time.
double E
Energy [GeV] (either MC truth or reconstructed)
Definition: Trk.hh:20
Range of reconstruction stages.
Acoustics hit.
then echo The file $DIR KM3NeT_00000001_00000000 root already please rename or remove it first
#define MAKE_STRING(A)
Make string.
Definition: JPrint.hh:127
Event selector.
Definition: JAAnet/JEvD.cc:66
void print(std::ostream &out=std::cout) const
Print hit.
Definition: Hit.hh:60
Basic data structure for L0 hit.
Axis object.
Definition: JAxis3D.hh:38
Rotation around Z-axis.
Definition: JRotation3D.hh:85
result_type calculate(const double E, const double R, const double theta, const double phi, const double t1) const
Get PDF.
Definition: JPDF_t.hh:233
The template JSinglePointer class can be used to hold a pointer to an object.
Scanning of objects from a single file according a format that follows from the extension of each fil...
double getW(const Trk &track, const size_t index, const double value)
Get track information.
Auxiliary class for defining the range of iterations of objects.
Definition: JLimit.hh:41
static const int JMUONPREFIT
I/O formatting auxiliaries.
JAxis3D & rotate(const JRotation3D &R)
Rotate axis.
Definition: JAxis3D.hh:225
JDirection3D getDirection(const Vec &dir)
Get direction.
JFunction1D_t::result_type result_type
Definition: JPDF_t.hh:145
JDirection3D & rotate(const JRotation3D &R)
Rotate.
Keyboard settings for unbuffered input.
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1989
Enable unbuffered terminal input.
Definition: JKeypress.hh:32
const array_type< JValue_t > & make_array(const JValue_t(&array)[N])
Method to create array of values.
Definition: JVectorize.hh:54
static const int JMUONGANDALF
double getTheta() const
Get theta angle.
Definition: JVersor3D.hh:128
JPosition3D getPosition(const Vec &pos)
Get position.
T getUpperLimit() const
Get upper limit.
Definition: JRange.hh:213
Definition of hit and track types and auxiliary methods for handling Monte Carlo data.
std::vector< double > fitinf
place to store additional fit info, see km3net-dataformat/definitions/fitparameters.csv
Definition: Trk.hh:32
#define ERROR(A)
Definition: JMessage.hh:66
then awk string
std::istream & getline(std::istream &in, JString &object)
Read string from input stream until end of line.
Definition: JString.hh:478
double len
length, if applicable [m]
Definition: Trk.hh:22
double putTime() const
Get Monte Carlo time minus DAQ/trigger time.
static const double PI
Mathematical constants.
T constrain(argument_type x) const
Constrain value to range.
Definition: JRange.hh:350
File router for fast addressing of summary data.
double getY() const
Get y position.
Definition: JVector3D.hh:104
Auxiliary data structure for muon PDF.
Vec dir
hit direction; i.e. direction of the PMT
Definition: Hit.hh:26
General purpose messaging.
Auxiliary data structure for sequence of same character.
Definition: JManip.hh:328
Auxiliary include file for time conversion between DAQ/trigger hit and Monte Carlo hit...
Auxiliary data structure for muon PDF.
Definition: JPDF_t.hh:135
#define FATAL(A)
Definition: JMessage.hh:67
Streaming of input and output from Linux command.
Definition: JProcess.hh:29
$WORKDIR ev_configure_dqsimulator txt echo process $DQ_SIMULATOR $i $SOURCE_HOST[$index] csh c(setenv ROOTSYS $ROOTSYS &&source $JPP_DIR/setenv.csh $JPP_DIR &&($DQ_SIMULATOR\-u\$NAME\$\-H\$SERVER\$\-M\$LOGGER\$\-d $DEBUG</dev/null > &/dev/null &))'
then JCookie sh JDataQuality D $DETECTOR_ID R
Definition: JDataQuality.sh:41
Auxiliary class for a hit with background rate value.
Definition: JHitW0.hh:21
Definition: Hit.hh:8
const double getSpeedOfLight()
Get speed of light.
static const int JSTART_LENGTH_METRES
distance between first and last hits in metres from JStart.cc
bool has_muon(const Evt &evt)
Test whether given event has a muon.
double getT(const JVector3D &pos) const
Get arrival time of Cherenkov light at given position.
Definition: JLine1Z.hh:114
static const int JMUONSIMPLEX
Utility class to parse command line options.
double t
hit time (from tdc+calibration or MC truth)
Definition: Hit.hh:23
Data structure for L0 hit.
Definition: JHitL0.hh:27
const double getInverseSpeedOfLight()
Get inverse speed of light.
int dom_id
module identifier from the data (unique in the detector).
Definition: Hit.hh:14
void setZ(const double z, const double velocity)
Set z-position of vertex.
Definition: JLine1Z.hh:75
then if[[!-f $DETECTOR]] then JDetector sh $DETECTOR fi cat $WORKDIR trigger_parameters txt<< EOFtrigger3DMuon.enabled=1;trigger3DMuon.numberOfHits=5;trigger3DMuon.gridAngle_deg=1;ctMin=0.0;TMaxLocal_ns=15.0;EOF set_variable TRIGGEREFFICIENCY_TRIGGERED_EVENTS_ONLY INPUT_FILES=() for((i=1;$i<=$NUMBER_OF_RUNS;++i));do JSirene.sh $DETECTOR $JPP_DATA/genhen.km3net_wpd_V2_0.evt.gz $WORKDIR/sirene_ ${i}.root JTriggerEfficiency.sh $DETECTOR $DETECTOR $WORKDIR/sirene_ ${i}.root $WORKDIR/trigger_efficiency_ ${i}.root $WORKDIR/trigger_parameters.txt $JPP_DATA/PMT_parameters.txt INPUT_FILES+=($WORKDIR/trigger_efficiency_ ${i}.root) done for ANGLE_DEG in $ANGLES_DEG[*];do set_variable SIGMA_NS 3.0 set_variable OUTLIERS 3 set_variable OUTPUT_FILE $WORKDIR/matrix\[${ANGLE_DEG}\deg\].root $JPP_DIR/examples/JReconstruction-f"$INPUT_FILES[*]"-o $OUTPUT_FILE-S ${SIGMA_NS}-A ${ANGLE_DEG}-O ${OUTLIERS}-d ${DEBUG}--!fiif[[$OPTION=="plot"]];then if((0));then for H1 in h0 h1;do JPlot1D-f"$WORKDIR/matrix["${^ANGLES_DEG}" deg].root:${H1}"-y"1 2e3"-Y-L TR-T""-\^"number of events [a.u.]"-> o chi2
Definition: JMatrixNZ.sh:106
Data structure for fit of straight line paralel to z-axis.
Definition: JLine1Z.hh:27
unsigned int tot
tot value as stored in raw data (int for pyroot)
Definition: Hit.hh:17
double getX() const
Get x position.
Definition: JVector3D.hh:94
no fit printf nominal n $STRING awk v X
double getTanThetaC()
Get average tangent of Cherenkov angle of water corresponding to group velocity.
unsigned int channel_id
PMT channel id {0,1, .., 30} local to moduke.
Definition: Hit.hh:15
Object reading from a list of files.
const JLimit & getLimit() const
Get limit.
Definition: JLimit.hh:84
void print(std::ostream &out=std::cout) const
Print track.
Definition: Trk.hh:183
KM3NeT DAQ constants, bit handling, etc.
static const int NUMBER_OF_PMTS
Total number of PMTs in module.
Definition: JDAQ.hh:26
Function adaptor.
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
The Trk class represents a Monte Carlo (MC) particle as well as a reconstructed track/shower.
Definition: Trk.hh:14
Auxiliary data structure for floating point format specification.
Definition: JManip.hh:484
Wrapper class around ROOT TStyle.
Definition: JStyle.hh:20
static bool select(const Trk &trk, const Evt &evt)
Default event selection.
Definition: JAAnet/JEvD.cc:76
static const int JMUONENERGY
int debug
debug level
The Evt class respresent a Monte Carlo (MC) event as well as an offline event.
Definition: Evt.hh:20
Data structure for size of TCanvas.
Definition: JCanvas.hh:26
#define DEBUG(A)
Message macros.
Definition: JMessage.hh:62
bool hasW(const Trk &trk, const int i)
Check availability of value.
Definition: JAAnet/JEvD.cc:100