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
Functions
JReconstruction/JEvD.cc File Reference

Program to display hit probabilities. More...

#include <string>
#include <iostream>
#include <iomanip>
#include <vector>
#include <algorithm>
#include "TROOT.h"
#include "TApplication.h"
#include "TCanvas.h"
#include "TStyle.h"
#include "TH2D.h"
#include "TArrow.h"
#include "TLatex.h"
#include "TMarker.h"
#include "km3net-dataformat/offline/Head.hh"
#include "km3net-dataformat/offline/Evt.hh"
#include "km3net-dataformat/online/JDAQ.hh"
#include "km3net-dataformat/tools/time_converter.hh"
#include "JROOT/JStyle.hh"
#include "JROOT/JCanvas.hh"
#include "JROOT/JRootToolkit.hh"
#include "JDAQ/JDAQEventIO.hh"
#include "JDAQ/JDAQSummarysliceIO.hh"
#include "JDetector/JDetector.hh"
#include "JDetector/JDetectorToolkit.hh"
#include "JDetector/JModuleRouter.hh"
#include "JTrigger/JHitL0.hh"
#include "JTrigger/JBuildL0.hh"
#include "JSupport/JSupport.hh"
#include "JSupport/JParallelFileScanner.hh"
#include "JSupport/JSummaryFileRouter.hh"
#include "JSupport/JMonteCarloFileSupportkit.hh"
#include "JAAnet/JAAnetToolkit.hh"
#include "JPhysics/JPDF_t.hh"
#include "JFit/JLine1Z.hh"
#include "JFit/JModel.hh"
#include "JReconstruction/JHitW0.hh"
#include "JReconstruction/JEvt.hh"
#include "JReconstruction/JEvtToolkit.hh"
#include "JReconstruction/JMuonGandalfParameters_t.hh"
#include "JReconstruction/JEventSelector.hh"
#include "JLang/JPredicate.hh"
#include "JLang/JComparator.hh"
#include "JMath/JMathToolkit.hh"
#include "JSystem/JKeypress.hh"
#include "JSystem/JProcess.hh"
#include "Jeep/JFunctionAdaptor.hh"
#include "Jeep/JPrint.hh"
#include "Jeep/JParser.hh"
#include "Jeep/JMessage.hh"

Go to the source code of this file.

Functions

int main (int argc, char **argv)
 

Detailed Description

Program to display hit probabilities.

Author
mdejong

Definition in file JReconstruction/JEvD.cc.

Function Documentation

int main ( int  argc,
char **  argv 
)

Definition at line 99 of file JReconstruction/JEvD.cc.

100 {
101  using namespace std;
102  using namespace JPP;
103  using namespace KM3NETDAQ;
104 
105  typedef JParallelFileScanner< JTypeList<JDAQEvent, JEvt> > JParallelFileScanner_t;
106  typedef JParallelFileScanner_t::multi_pointer_type multi_pointer_type;
107 
108  JParallelFileScanner_t inputFile;
109  JLimit_t& numberOfEvents = inputFile.getLimit();
110  string detectorFile;
111  string pdfFile;
112  string outputFile;
114  int application;
115  JEventSelector event_selector;
116  JCanvas canvas;
117  bool batch;
118  double arrowSize;
119  string arrowType;
120  double arrowScale;
121  int debug;
122 
123 
124  try {
125 
126  parameters.numberOfPrefits = 1;
127 
128  JParser<> zap("Program to display hit probabilities.");
129 
130  zap['w'] = make_field(canvas, "size of canvas <nx>x<ny> [pixels]") = JCanvas(1200, 600);
131  zap['f'] = make_field(inputFile, "input file (output of JXXXMuonReconstruction.sh)");
132  zap['a'] = make_field(detectorFile);
133  zap['n'] = make_field(numberOfEvents) = JLimit::max();
134  zap['P'] = make_field(pdfFile);
135  zap['o'] = make_field(outputFile, "graphics output file name") = MAKE_STRING("display_" << WILDCARD << ".gif");
137  zap['A'] = make_field(application) = JMUONGANDALF, JMUONENERGY, JMUONSTART;
138  zap['L'] = make_field(event_selector) = JPARSER::initialised();
139  zap['S'] = make_field(arrowSize) = 0.003;
140  zap['T'] = make_field(arrowType) = "|->";
141  zap['F'] = make_field(arrowScale) = 250.0;
142  zap['B'] = make_field(batch, "batch processing");
143  zap['d'] = make_field(debug) = 1;
144 
145  zap(argc, argv);
146  }
147  catch(const exception& error) {
148  FATAL(error.what() << endl);
149  }
150 
151  if (batch && outputFile == "") {
152  FATAL("Missing output file name " << outputFile << " in batch mode." << endl);
153  }
154 
155  if (!batch && outputFile == "") {
156  outputFile = MAKE_STRING(WILDCARD << ".gif");
157  }
158 
159  if (outputFile.find(WILDCARD) == string::npos) {
160  FATAL("Output file name " << outputFile << " has no wild card '" << WILDCARD << "'" << endl);
161  }
162 
163 
165 
166  try {
167  load(detectorFile, detector);
168  }
169  catch(const JException& error) {
170  FATAL(error);
171  }
172 
173  const JModuleRouter router(detector);
174 
175  JSummaryFileRouter summary(inputFile, parameters.R_Hz);
176 
177  const JMuonPDF_t pdf(pdfFile, parameters.TTS_ns);
178 
179  const JTimeRange T_ns(parameters.TMin_ns, parameters.TMax_ns);
180 
181  typedef vector<JHitL0> JDataL0_t;
182  typedef vector<JHitW0> JDataW0_t;
183 
184  const JBuildL0<JHitL0> buildL0;
185 
186 
187  // Monte Carlo
188 
189  JVector3D center(0.0, 0.0, 0.0);
190 
191  try {
192  center = get<JPosition3D>(getHeader(inputFile));
193  } catch(const exception& error) {}
194 
195 
196  // ROOT
197 
198  gROOT->SetBatch(batch);
199 
200  TApplication* tp = new TApplication("user", NULL, NULL);
201  TCanvas* cv = new TCanvas("display", "", canvas.x, canvas.y);
202 
203  JSinglePointer<TStyle> gStyle(new JStyle("gplot", cv->GetWw(), cv->GetWh()));
204 
205  gROOT->SetStyle("gplot");
206  gROOT->ForceStyle();
207 
208  const size_t NUMBER_OF_PADS = 3;
209 
210  cv->SetFillStyle(4000);
211  cv->SetFillColor(kWhite);
212  cv->Divide(NUMBER_OF_PADS, 1);
213 
214  const double Dmax = getMaximalDistance(detector);
215  const double Rmin = 0.0;
216  const double Rmax = min(parameters.roadWidth_m, 0.4 * Dmax);
217  const double Tmin = min(parameters.TMin_ns, -10.0);
218  const double Tmax = max(parameters.TMax_ns, +100.0);
219  const double Amin = 0.002 * (Tmax - Tmin); // minimal arrow length [ns]
220  const double Amax = 0.8 * (Tmax - Tmin); // maximal arrow length [ns]
221  const double ymin = min(-Amax, Tmin - 0.3 * Amax);
222  const double ymax = max(+Amax, Tmax + 0.5 * Amax);
223 
224  const string Xlabel[NUMBER_OF_PADS] = { "R [m]", "#phi [rad]", "z [m]" };
225  const double Xmin [NUMBER_OF_PADS] = { Rmin, -PI, -0.5 * Dmax };
226  const double Xmax [NUMBER_OF_PADS] = { Rmax, +PI, +0.5 * Dmax };
227 
228  double Xs[NUMBER_OF_PADS];
229 
230  for (size_t i = 0; i != NUMBER_OF_PADS; ++i) {
231  Xs[i] = 0.003 * (Xmax[i] - Xmin[i]) * (0.5 * NUMBER_OF_PMTS); // x-offset arrow as function of PMT number
232  }
233 
234  TH2D H2[NUMBER_OF_PADS];
235  TGraph G2[NUMBER_OF_PADS];
236 
237  for (size_t i = 0; i != NUMBER_OF_PADS; ++i) {
238 
239  H2[i] = TH2D(MAKE_CSTRING("h" << i), NULL, 100, Xmin[i] - Xs[i], Xmax[i] + Xs[i], 100, ymin, ymax);
240 
241  H2[i].GetXaxis()->SetTitle(Xlabel[i].c_str());
242  H2[i].GetYaxis()->SetTitle("#Deltat [ns]");
243 
244  H2[i].GetXaxis()->CenterTitle(true);
245  H2[i].GetYaxis()->CenterTitle(true);
246 
247  H2[i].SetStats(kFALSE);
248 
249  G2[i].Set(2);
250 
251  G2[i].SetPoint(0, H2[i].GetXaxis()->GetXmin(), 0.0);
252  G2[i].SetPoint(1, H2[i].GetXaxis()->GetXmax(), 0.0);
253 
254  cv->cd(i+1);
255 
256  H2[i].Draw("AXIS");
257  G2[i].Draw("SAME");
258  }
259 
260 
261  for (JTreeScanner<Evt> mc(inputFile); inputFile.hasNext(); ) {
262 
263  cout << "\revent: " << setw(8) << inputFile.getCounter() << flush;
264 
265  multi_pointer_type ps = inputFile.next();
266 
267  JDAQEvent* tev = ps;
268  JEvt* in = ps;
269  Evt* event = NULL;
270 
271  if (mc.getEntries() != 0) {
272  event = mc.getEntry(tev->getCounter()); // Monte Carlo true information
273  }
274 
275 
276  in->select(JHistory::is_event(application));
277 
278  if (!in->empty()) {
279 
280  sort(in->begin(), in->end(), qualitySorter);
281 
282  if (!event_selector(*tev, *in, event)) {
283  continue;
284  }
285 
286 
287  JDataL0_t dataL0;
288 
289  buildL0(*tev, router, true, back_inserter(dataL0));
290 
291  summary.update(*tev);
292 
293 
294  JFit muon; // Monte Carlo true muon
295 
296  if (event != NULL) {
297 
298  const time_converter converter = time_converter(*event, *tev);
299 
300  for (const auto& t1 : event->mc_trks) {
301  if (is_muon(t1)) {
302  if (t1.E > muon.getE()) {
303 
304  JTrack3E ta = getTrack(t1);
305 
306  ta.add(center);
307  ta.add(converter.putTime());
308 
309  muon = getFit(0, ta, 0.0, 0, t1.E, 1);
310 
311  muon.setW(JSTART_LENGTH_METRES, fabs(t1.len));
312  }
313  }
314  }
315  }
316 
317 
318  bool next_event = false; // goto next event
319  bool monte_carlo = false; // show Monte Carlo true muon
320  size_t index = 0; // index of fit
321 
322  while (!next_event) {
323 
324  JFit fit;
325 
326  if (!monte_carlo)
327  fit = (*in)[index];
328  else
329  fit = muon;
330 
331 
332  JRotation3D R (getDirection(fit));
333  JLine1Z tz(getPosition (fit).rotate(R), fit.getT());
334  JRange<double> Z_m;
335  /*
336  if (fit.hasW(JSTART_LENGTH_METRES)) {
337  Z_m = JRange<double>(0.0, fit.getW(JSTART_LENGTH_METRES)) + JRange<double>(parameters.ZMin_m, parameters.ZMax_m);
338  }
339  */
340  const JModel<JLine1Z> match(tz, parameters.roadWidth_m, T_ns, Z_m);
341 
342  // hit selection based on fit result
343 
344  JDataW0_t data;
345 
346  for (JDataL0_t::const_iterator i = dataL0.begin(); i != dataL0.end(); ++i) {
347 
348  JHitW0 hit(*i, summary.getRate(i->getPMTIdentifier()));
349 
350  hit.rotate(R);
351 
352  if (match(hit)) {
353  data.push_back(hit);
354  }
355  }
356 
357  // select first hit in PMT
358 
359  sort(data.begin(), data.end(), JHitW0::compare);
360 
361  JDataW0_t::iterator __end = unique(data.begin(), data.end(), equal_to<JDAQPMTIdentifier>());
362 
363  double E_GeV = parameters.E_GeV;
364  /*
365  if (fit.getE() > 0.1) {
366  E_GeV = fit.getE();
367  }
368  */
369 
370  // move fit to geometrical center of hits
371 
372  JRange<double> zs(make_array(data.begin(), __end, &JHitW0::getZ));
373 
374  const double z0 = tz.getZ();
375  const double z1 = 0.5 * (zs.getLowerLimit() + zs.getUpperLimit());
376 
377  tz.setZ(z1, getSpeedOfLight());
378 
379 
380  // graphics
381 
382  TLatex latex [NUMBER_OF_PADS];
383  vector<TArrow> arrow [NUMBER_OF_PADS];
384  vector<TMarker> marker[NUMBER_OF_PADS];
385 
386  for (int i = 0; i != NUMBER_OF_PADS; ++i) {
387 
388  latex[i].SetTextAlign(12);
389  latex[i].SetTextFont(42);
390  latex[i].SetTextSize(0.06);
391 
392  latex[i].SetX(H2[i].GetXaxis()->GetXmin() + 0.05 * (H2[i].GetXaxis()->GetXmax() - H2[i].GetXaxis()->GetXmin()));
393  latex[i].SetY(H2[i].GetYaxis()->GetXmax() - 0.05 * (H2[i].GetYaxis()->GetXmax() - H2[i].GetYaxis()->GetXmin()));
394  }
395 
396  if (fit.hasW(JSTART_LENGTH_METRES) && fit.getW(JSTART_LENGTH_METRES) > 0.0) {
397 
398  marker[2].push_back(TMarker(z0 - tz.getZ(), 0.0, kFullCircle));
399  marker[2].push_back(TMarker(z0 - tz.getZ() + fit.getW(JSTART_LENGTH_METRES), 0.0, kFullCircle));
400 
401  static_cast<TAttMarker&>(marker[2][0]) = TAttMarker(kRed, kFullCircle, 0.7);
402  static_cast<TAttMarker&>(marker[2][1]) = TAttMarker(kRed, kFullCircle, 0.7);
403  }
404 
405  double chi2 = 0;
406 
407  for (JDataW0_t::const_iterator hit = data.begin(); hit != __end; ++hit) {
408 
409  const double x = hit->getX() - tz.getX();
410  const double y = hit->getY() - tz.getY();
411  const double z = hit->getZ() - tz.getZ();
412  const double R = sqrt(x*x + y*y);
413 
414  const double t1 = tz.getT() + (z + R * getTanThetaC()) * getInverseSpeedOfLight();
415 
416  JDirection3D dir(hit->getDX(), hit->getDY(), hit->getDZ()); // PMT orientation
417 
418  dir.rotate(JRotation3Z(-atan2(y,x))); // rotate PMT axis to x-z plane
419 
420  const double theta = dir.getTheta();
421  const double phi = fabs(dir.getPhi()); // rotational symmetry of Cherenkov cone
422 
423  //const double E = gWater.getE(E_GeV, z); // correct for energy loss
424  const double E = E_GeV;
425  const double dt = T_ns.constrain(hit->getT() - t1);
426 
427  JMuonPDF_t::result_type H1 = pdf.calculate(E, R, theta, phi, dt);
428  JMuonPDF_t::result_type H0(hit->getR() * 1e-9, 0.0, T_ns);
429 
430  if (H1.V >= parameters.VMax_npe) {
431  H1 *= parameters.VMax_npe / H1.V;
432  }
433 
434  H1 += H0; // signal + background
435 
436  chi2 += H1.getChi2() - H0.getChi2();
437 
438  const double derivative = H1.getDerivativeOfChi2() - H0.getDerivativeOfChi2();
439 
440  double size = derivative * arrowScale; // size of arrow
441 
442  if (fabs(size) < Amin) {
443  size = (size > 0.0 ? +Amin : -Amin);
444  } else if (fabs(size) > Amax) {
445  size = (size > 0.0 ? +Amax : -Amax);
446  }
447 
448  const double X[NUMBER_OF_PADS] = { R, atan2(y,x), z - R/getTanThetaC() };
449 
450  const double xs = (double) (NUMBER_OF_PMTS - 2 * hit->getPMTAddress()) / (double) NUMBER_OF_PMTS;
451 
452  for (size_t i = 0; i != NUMBER_OF_PADS; ++i) {
453  arrow[i].push_back(TArrow(X[i] + xs*Xs[i], dt, X[i] + xs*Xs[i], dt + size, arrowSize, arrowType.c_str()));
454  }
455  }
456 
457  latex[0].SetTitle(MAKE_CSTRING("" << FILL(6,'0') << tev->getRunNumber() << FILL() <<
458  ":" << tev->getFrameIndex() <<
459  " " << tev->getCounter()));
460 
461  if (!monte_carlo)
462  latex[1].SetTitle(MAKE_CSTRING("[" << index << "]" << " " <<
463  "Q = " << FIXED(4,0) << getQuality(chi2) << " " <<
464  "E = " << SCIENTIFIC(7,1) << fit.getE() << " [GeV]"));
465  else
466  latex[1].SetTitle(MAKE_CSTRING("Q = " << FIXED(4,0) << getQuality(chi2)));
467 
468  if (muon.getStatus() >= 0) {
469 
470  if (!monte_carlo)
471  latex[2].SetTitle(MAKE_CSTRING("#Delta#alpha = " << FIXED(6,2) << getAngle(getDirection(muon), getDirection(fit)) << " [deg]"));
472  else
473  latex[2].SetTitle("Monte Carlo");
474  }
475 
476 
477  // draw
478 
479  for (int i = 0; i != NUMBER_OF_PADS; ++i) {
480 
481  cv->cd(i+1);
482 
483  latex[i].Draw();
484 
485  for (auto& a1 : arrow[i]) {
486  a1.Draw();
487  }
488 
489  for (auto& m1 : marker[i]) {
490  m1.Draw();
491  }
492  }
493 
494  cv->Update();
495 
496 
497  // action
498 
499  if (batch) {
500 
501  cv->SaveAs(replace(outputFile, WILDCARD, MAKE_STRING(inputFile.getCounter())).c_str());
502 
503  next_event = true;
504 
505  } else {
506 
507  for (int c = 0; c == 0; ) {
508 
509  if (!batch) {
510 
511  static bool first = true;
512 
513  if (first) {
514  cout << endl << "Type '?' for possible options." << endl;
515  }
516 
517  first = false;
518  }
519 
520  cout << "\n> " << flush;
521 
522  switch((c = JKeypress(true).get())) {
523 
524  case '?':
525  cout << endl;
526  cout << "possible options: " << endl;
527  cout << 'q' << " -> " << "exit application" << endl;
528  cout << 'u' << " -> " << "update canvas" << endl;
529  cout << 's' << " -> " << "save graphics to file" << endl;
530  cout << '+' << " -> " << "next fit" << endl;
531  cout << '-' << " -> " << "previous fit" << endl;
532  cout << 'M' << " -> " << "Monte Carlo true muon information" << endl;
533  cout << 'F' << " -> " << "fit information" << endl;
534  if (event_selector.is_valid()) {
535  cout << 'L' << " -> " << "reload event selector" << endl;
536  }
537  cout << 'r' << " -> " << "rewind input file" << endl;
538  cout << 'R' << " -> " << "switch to ROOT mode (quit ROOT to continue)" << endl;
539  cout << ' ' << " -> " << "next event (as well as any other key)" << endl;
540  c = 0;
541  break;
542 
543  case 'q':
544  cout << endl;
545  return 0;
546 
547  case 'u':
548  cv->Update();
549  c = 0;
550  break;
551 
552  case 's':
553  cv->SaveAs(replace(outputFile, WILDCARD, MAKE_STRING(inputFile.getCounter())).c_str());
554  c = 0;
555  break;
556 
557  case '+':
558  monte_carlo = false;
559  index = (index != in->size() - 1 ? index + 1 : 0);
560  break;
561 
562  case '-':
563  monte_carlo = false;
564  index = (index != 0 ? index - 1 : in->size() - 1);
565  break;
566 
567  case 'M':
568  if (muon.getStatus() >= 0)
569  monte_carlo = true;
570  else
571  ERROR(endl << "No Monte Carlo muon available." << endl);
572  break;
573 
574  case 'F':
575  monte_carlo = false;
576  break;
577 
578  case 'r':
579  inputFile.rewind();
580  next_event = true;
581  break;
582 
583  case 'L':
584  if (event_selector.is_valid()) {
585  execute(MAKE_STRING("make -f " << getPath(argv[0]) << "/JMakeEventSelector libs"), 3);
586  event_selector.reload();
587  }
588  break;
589 
590  case 'R':
591  tp->Run(kTRUE);
592  c = 0;
593  break;
594 
595  default:
596  next_event = true;
597  break;
598  }
599  }
600  }
601  }
602  }
603  }
604  cout << endl;
605 }
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.
General exception.
Definition: JException.hh:24
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 getQuality(const double chi2, const int NDF)
Get quality of fit.
JTrack3E getTrack(const Trk &track)
Get track.
Template specialisation of L0 builder for JHitL0 data type.
Definition: JBuildL0.hh:102
std::string getPath(const std::string &file_name)
Get path, i.e. part before last JEEP::PATHNAME_SEPARATOR if any.
Definition: JeepToolkit.hh:148
Detector data structure.
Definition: JDetector.hh:89
bool is_muon(const Trk &track)
Test whether given track is a (anti-)muon.
Router for direct addressing of module data in detector data structure.
Rotation matrix.
Definition: JRotation3D.hh:111
*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
General purpose class for parallel reading of objects from a single file or multiple files...
#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
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
int getRunNumber() const
Get run number.
3D track with energy.
Definition: JTrack3E.hh:30
Acoustic single fit.
bool qualitySorter(const JFit &first, const JFit &second)
Comparison of fit results.
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
int getFrameIndex() const
Get frame index.
JTime & add(const JTime &value)
Addition operator.
Head getHeader(const JMultipleFileScanner_t &file_list)
Get Monte Carlo header.
Rotation around Z-axis.
Definition: JRotation3D.hh:85
The template JSinglePointer class can be used to hold a pointer to an object.
Auxiliary class for defining the range of iterations of objects.
Definition: JLimit.hh:41
Detector file.
Definition: JHead.hh:226
JDirection3D getDirection(const Vec &dir)
Get direction.
JFunction1D_t::result_type result_type
Definition: JPDF_t.hh:145
Data structure for vector in three dimensions.
Definition: JVector3D.hh:34
JDirection3D & rotate(const JRotation3D &R)
Rotate.
Acoustic event fit.
#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
Auxiliary class to test history.
Definition: JHistory.hh:105
static const int JMUONGANDALF
double getTheta() const
Get theta angle.
Definition: JVersor3D.hh:128
JPosition3D getPosition(const Vec &pos)
Get position.
#define ERROR(A)
Definition: JMessage.hh:66
double putTime() const
Get Monte Carlo time minus DAQ/trigger time.
static const double PI
Mathematical constants.
File router for fast addressing of summary data.
Auxiliary data structure for sequence of same character.
Definition: JManip.hh:328
Auxiliary data structure for muon PDF.
Definition: JPDF_t.hh:135
#define FATAL(A)
Definition: JMessage.hh:67
$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
const double getSpeedOfLight()
Get speed of light.
static const int JSTART_LENGTH_METRES
distance between first and last hits in metres from JStart.cc
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
JFit getFit(const int id, const JMODEL::JString &string)
Get fit parameters of string.
const double getInverseSpeedOfLight()
Get inverse speed of light.
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
double getMaximalDistance(const JDetector &detector, const bool option=false)
Get maximal distance between modules in detector.
no fit printf nominal n $STRING awk v X
double getTanThetaC()
Get average tangent of Cherenkov angle of water corresponding to group velocity.
const JLimit & getLimit() const
Get limit.
Definition: JLimit.hh:84
do set_variable DETECTOR_TXT $WORKDIR detector
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 data structure for floating point format specification.
Definition: JManip.hh:484
JTriggerCounter_t getCounter() const
Get trigger counter.
Wrapper class around ROOT TStyle.
Definition: JStyle.hh:20
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