196{
200
202 JLimit_t& numberOfEvents = inputFile.getLimit();
218 double T_ns = 0.0;
220 string option;
222
223
224 try {
225
227
234
236
238
240 zap[
'f'] =
make_field(inputFile,
"input file (output of JXXXMuonReconstruction.sh)");
251
253 }
254 catch(const exception& error) {
256 }
257
260 }
261
264 }
265
268 }
269
271
273
275
277
280
281
282
283
285
288
290 ((
TRootCanvas *) cv->GetCanvasImp())->
Connect(
"CloseWindow()",
"TApplication", tp,
"Terminate()");
291 }
292
294
295 gROOT->SetStyle(
"gplot");
297
299
300 cv->SetFillStyle(4000);
302
305
307
309 p2->Draw();
310
311 const double Dmax = 1000.0;
312 const double Rmin = 0.0;
314 const double Tmin = min(parameters.
TMin_ns, -10.0);
315 const double Tmax = max(parameters.
TMax_ns, +100.0);
316 const double Amin = 0.002 * (Tmax - Tmin);
317 const double Amax = 0.8 * (Tmax - Tmin);
320
324
326
328 Xs[i] = 0.003 * (Xmax[i] -
Xmin[i]) * (0.5 * NUMBER_OF_PMTS);
329 }
330
333
335
337
338 H2[i].GetXaxis()->SetTitle(
Xlabel[i].c_str());
339 H2[i].GetYaxis()->SetTitle(
"#Deltat [ns]");
340
341 H2[i].GetXaxis()->CenterTitle(
true);
342 H2[i].GetYaxis()->CenterTitle(
true);
343
345
347
350
352
355 }
356
357
358 while (inputFile.hasNext()) {
359
360 cout <<
"event: " <<
setw(8) << inputFile.getCounter() <<
endl;
361
362 const Evt* evt = inputFile.next();
363
365
367
369 continue;
370 }
371
373
374 for (
const Hit&
hit : evt->hits) {
378 }
379
384
386
388
390
391 for (
const auto&
trk : evt->mc_trks) {
393 if (
trk.E > muon.
E) {
394
397
399 }
400 }
401 }
402 }
403
404
406
407 for (bool next = false; !next; ) {
408
411 }
412
414
417 else
419
423
424
425
426
427
429
430
431
433
434 for (JDataL0_t::const_iterator i = dataL0.begin(); i != dataL0.end(); ++i) {
435
436 JHitW0 hit(*i, summary.getRate(i->getPMTIdentifier(), parameters.
R_Hz));
437
439
442 }
443 }
444
445
446
447 sort(
data.begin(),
data.end(), JHitW0::compare);
448
450
451 double E_GeV = parameters.
E_GeV;
452
453
454
455
456
457
458
459
461
462 for (JDataW0_t::iterator
hit =
data.begin();
hit != __end; ++
hit) {
463
464 const double x =
hit->getX() - tz.getX();
465 const double y =
hit->getY() - tz.getY();
466 const double z =
hit->getZ();
467 const double R = sqrt(x*x + y*y);
468
470 }
471
472 const double z0 = tz.getZ();
473 const double z1 = 0.5 * (
zs.getLowerLimit() +
zs.getUpperLimit());
474
476
477
478
479
480 ostringstream os;
483
485
488
491 }
492
494 <<
FIXED(7,2) << tz.getX() <<
' '
495 <<
FIXED(7,2) << tz.getY() <<
' '
496 <<
FIXED(7,2) << tz.getZ() <<
' '
498
499 double chi2 = 0;
500
501 for (JDataW0_t::const_iterator
hit =
data.begin();
hit != __end; ++
hit) {
502
503 const double x =
hit->getX() - tz.getX();
504 const double y =
hit->getY() - tz.getY();
505 const double z =
hit->getZ() - tz.getZ();
506 const double R = sqrt(x*x + y*y);
507
509
511
513
514 const double theta = dir.getTheta();
515 const double phi = fabs(dir.getPhi());
516
517
518 const double E = E_GeV;
519 const double dt = T_ns.constrain(
hit->getT() - t1);
520
523
526 }
527
529
530 chi2 +=
H1.getChi2() -
H0.getChi2();
531
533 <<
setw(8) <<
hit->getModuleID() <<
'.' <<
FILL(2,
'0') << (
int)
hit->getPMTAddress() <<
FILL() <<
' '
535 <<
FIXED(7,2) << R <<
' '
536 <<
FIXED(7,4) << theta <<
' '
537 <<
FIXED(7,4) << phi <<
' '
539 <<
FIXED(7,3) <<
H1.getChi2() <<
' '
541
542 const double derivative =
H1.getDerivativeOfChi2() -
H0.getDerivativeOfChi2();
543
545
546 if (fabs(size) <
Amin) {
547 size = (size > 0.0 ? +
Amin : -
Amin);
548 }
else if (fabs(size) >
Amax) {
549 size = (size > 0.0 ? +
Amax : -
Amax);
550 }
551
553
555
557
559
560 a1.SetLineWidth(
graphics.lineWidth);
561 a1.SetLineStyle(
graphics.lineStyle);
562
563 arrow[i].push_back(a1);
564
566 }
567 }
568
569 double zmax = 0.0;
570
573 zmax =
H2[i].GetMaximum();
574 }
575 }
576
577 zmax *= 1.2;
578
580 H2[i].SetMaximum(zmax);
581 }
582
584 os <<
" Q = " <<
FIXED(4,0) <<
trk.lik
585 <<
'/' <<
FIXED(4,0) << -chi2;
587 os <<
" cos(#theta) = " <<
FIXED(6,3) <<
trk.dir.z;
589
591 os << " Monte Carlo";
594
595
596
597
598 TLatex title(0.05, 0.5, os.str().c_str());
599
600 title.SetTextAlign(12);
601 title.SetTextFont(42);
602 title.SetTextSize(0.6);
603
604 p2->cd();
605
606 title.Draw();
607
609
611
613
614 for (
auto& a1 :
arrow[i]) {
615 a1.Draw();
616 }
617
618 for (
auto&
m1 : marker[i]) {
620 }
621 }
622
625 }
626 }
627
628 cv->Update();
629
630
631
632
634
636
637 next = true;
638
639 } else {
640
641 static int count = 0;
642
643 if (count++ == 0) {
644 cout <<
endl <<
"Type '?' for possible options." <<
endl;
645 }
646
648
649 cout <<
"\n> " << flush;
650
652
653 case '?':
655 cout <<
"possible options: " <<
endl;
656 cout <<
'q' <<
" -> " <<
"exit application" <<
endl;
657 cout <<
'u' <<
" -> " <<
"update canvas" <<
endl;
658 cout <<
's' <<
" -> " <<
"save graphics to file" <<
endl;
659 cout <<
'M' <<
" -> " <<
"Monte Carlo true muon information" <<
endl;
660 cout <<
'F' <<
" -> " <<
"fit information" <<
endl;
662 cout <<
'L' <<
" -> " <<
"reload event selector" <<
endl;
663 }
664 cout <<
'r' <<
" -> " <<
"rewind input file" <<
endl;
665 cout <<
'R' <<
" -> " <<
"switch to ROOT mode (quit ROOT to continue)" <<
endl;
666 cout <<
'p' <<
" -> " <<
"print event information" <<
endl;
667 cout <<
' ' <<
" -> " <<
"next event (as well as any other key)" <<
endl;
668 break;
669
670 case 'q':
672 return 0;
673
674 case 'u':
675 cv->Update();
676 break;
677
678 case 's':
680 break;
681
682 case 'M':
685 else
688 break;
689
690 case 'F':
693 break;
694
695 case 'L':
699 }
700 break;
701
702 case 'R':
704 break;
705
706 case 'p':
710 for (
const auto&
trk : evt->mc_trks) {
712 }
713 for (
const auto&
trk : evt->trks) {
715 }
716 for (
const auto&
hit : evt->hits) {
718 }
719 }
720 break;
721
722 case 'r':
723 inputFile.rewind();
724
725 default:
726 next = true;
728 break;
729 }
730 }
731 }
732 }
733 }
734 }
736}
#define DEBUG(A)
Message macros.
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
#define MAKE_CSTRING(A)
Make C-string.
#define MAKE_STRING(A)
Make string.
#define gmake_property(A)
macros to convert (template) parameter to JPropertiesElement object
Utility class to parse parameter values.
Data structure for fit of straight line paralel to z-axis.
Data structure for direction in three dimensions.
Template definition of a multi-dimensional oscillation probability interpolation table.
JProperties getProperties(const JEquationParameters &eqpars=JOscProbInterface::getEquationParameters()) override final
Get properties of this class.
Auxiliary class for a hit with background rate value.
Data structure for size of TCanvas.
Wrapper class around ROOT TStyle.
File router for fast addressing of summary data.
Enable unbuffered terminal input.
Data structure for L0 hit.
static void setSlewing(const bool slewing)
Set slewing option.
Data structure for UTC time.
Auxiliary class to convert DAQ hit time to/from Monte Carlo hit time.
static const int JSTART_LENGTH_METRES
distance between projected positions on the track of optical modules for which the response does not ...
JDirection3D getDirection(const Vec &dir)
Get direction.
bool hasW(const Trk &trk, const int i)
Check availability of value.
bool has_muon(const Evt &evt)
Test whether given event has a muon.
void setW(Trk &trk, const int i, const double value)
Set associated value.
double getW(const Trk &track, const size_t index, const double value)
Get track information.
JPosition3D getPosition(const Vec &pos)
Get position.
bool is_muon(const Trk &track)
Test whether given track is a (anti-)muon.
std::string getPath(const std::string &file_name)
Get path, i.e. part before last JEEP::PATHNAME_SEPARATOR if any.
double getAngle(const JQuaternion3D &first, const JQuaternion3D &second)
Get space angle between quanternions.
std::string replace(const std::string &input, const std::string &target, const std::string &replacement)
Replace tokens in string.
static const double PI
Mathematical constants.
const double getInverseSpeedOfLight()
Get inverse speed of light.
double getTanThetaC()
Get average tangent of Cherenkov angle of water corresponding to group velocity.
const double getSpeedOfLight()
Get speed of light.
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
KM3NeT DAQ data structures and auxiliaries.
static const char WILDCARD
static const int NUMBER_OF_PMTS
Total number of PMTs in module.
The Evt class respresent a Monte Carlo (MC) event as well as an offline event.
int frame_index
from the raw data
int run_id
DAQ run identifier.
int det_id
detector identifier from DAQ
ULong64_t trigger_counter
trigger counter
void print(std::ostream &out=std::cout) const
Print event.
TTimeStamp t
UTC time of the timeslice, or the event_time for MC. (default: 01 Jan 1970 00:00:00)
Auxiliary data structure for sequence of same character.
Auxiliary data structure for floating point format specification.
Auxiliary data structure for muon PDF.
JFunction1D_t::result_type result_type
Empty structure for specification of parser element that is initialised (i.e. does not require input)...
Data structure for fit parameters.
double TTS_ns
transition-time spread [ns]
double TMin_ns
minimal time w.r.t. Cherenkov hypothesis [ns]
double roadWidth_m
road width [m]
double TMax_ns
maximal time w.r.t. Cherenkov hypothesis [ns]
double VMax_npe
maximum number of of photo-electrons
double R_Hz
default rate [Hz]
size_t numberOfPrefits
number of prefits
Auxiliary class for defining the range of iterations of objects.
static counter_type max()
Get maximum counter value.
Auxiliary data structure for floating point format specification.
The Trk class represents a Monte Carlo (MC) particle as well as a reconstructed track/shower.
double E
Energy [GeV] (either MC truth or reconstructed)
double t
track time [ns] (when the particle is at pos )
double len
length, if applicable [m]
Range of reconstruction stages.