CbmRoot
CbmRichGeoTest.cxx
Go to the documentation of this file.
1 
8 #include "CbmRichGeoTest.h"
9 #include "CbmDigiManager.h"
10 #include "CbmDrawHist.h"
11 #include "CbmHistManager.h"
12 #include "CbmMCDataArray.h"
13 #include "CbmMCDataManager.h"
14 #include "CbmMCEventList.h"
15 #include "CbmMCTrack.h"
16 #include "CbmReport.h"
17 #include "CbmRichConverter.h"
18 #include "CbmRichDigi.h"
19 #include "CbmRichDigiMapManager.h"
20 #include "CbmRichDraw.h"
21 #include "CbmRichGeoManager.h"
22 #include "CbmRichHit.h"
23 #include "CbmRichPmt.h"
24 #include "CbmRichPoint.h"
25 #include "CbmRichRing.h"
26 #include "CbmRichRingFitterCOP.h"
28 #include "CbmStudyReport.h"
29 #include "CbmTrackMatchNew.h"
30 #include "FairGeoNode.h"
31 #include "FairGeoTransform.h"
32 #include "FairRootManager.h"
33 #include "FairRunAna.h"
34 #include "FairRuntimeDb.h"
35 #include "utils/CbmUtils.h"
36 
37 #include "TCanvas.h"
38 #include "TClonesArray.h"
39 #include "TEllipse.h"
40 #include "TF1.h"
41 #include "TH1D.h"
42 #include "TH2D.h"
43 #include "TH3D.h"
44 #include "TLatex.h"
45 #include "TLegend.h"
46 #include "TLine.h"
47 #include "TMCProcess.h"
48 #include "TMath.h"
49 #include "TPad.h"
50 #include "TStyle.h"
51 #include "TSystem.h"
52 #include <TFile.h>
53 
54 #include <iomanip>
55 #include <iostream>
56 #include <map>
57 #include <sstream>
58 #include <string>
59 #include <vector>
60 
61 using namespace std;
62 
64  : FairTask("RichGeoTestQa")
65  , fOutputDir("")
66  , fRichHits(NULL)
67  , fRichRings(NULL)
68  , fRichPoints(NULL)
69  , fMcTracks(NULL)
70  , fRichRingMatches(NULL)
71  , fEventList(NULL)
72  , fCopFit(NULL)
73  , fTauFit(NULL)
74  , fHM(NULL)
75  , fEventNum(0)
76  , fMinNofHits(0)
77  , fMinAaxis(0.)
78  , fMaxAaxis(0.)
79  , fMinBaxis(0.)
80  , fMaxBaxis(0.)
81  , fMinRadius(0.)
82  , fMaxRadius(0.)
83  , fNofDrawnRings(0)
84  , fDrawPmts(true)
85  , fDrawEventDisplay(true) {
86  fEventNum = 0;
87  fNofDrawnRings = 0;
88  fMinNofHits = 7;
89 
90  fMinAaxis = 3.;
91  fMaxAaxis = 7.;
92  fMinBaxis = 3.;
93  fMaxBaxis = 7.;
94  fMinRadius = 3.;
95  fMaxRadius = 7.;
96 }
97 
99 
100 InitStatus CbmRichGeoTest::Init() {
101  cout << "CbmRichGeoTest::Init" << endl;
102  FairRootManager* ioman = FairRootManager::Instance();
103  if (NULL == ioman) {
104  Fatal("CbmRichGeoTest::Init", "RootManager not instantised!");
105  }
106 
107  CbmMCDataManager* mcManager =
108  (CbmMCDataManager*) ioman->GetObject("MCDataManager");
109  if (mcManager == nullptr)
110  LOG(fatal) << "CbmRichGeoTest::Init() NULL MCDataManager.";
111 
112  fMcTracks = mcManager->InitBranch("MCTrack");
113  if (NULL == fMcTracks) { LOG(fatal) << "CbmRichGeoTest::Init No MCTrack!"; }
114 
115  fRichPoints = mcManager->InitBranch("RichPoint");
116  if (NULL == fRichPoints) {
117  LOG(fatal) << "CbmRichGeoTest::Init No RichPoint!";
118  }
119 
120  fRichRefPlanePoints = mcManager->InitBranch("RefPlanePoint");
121  if (NULL == fRichRefPlanePoints) {
122  LOG(fatal) << "CbmRichGeoTest::Init No RefPlanePoint!";
123  }
124 
125  fRichHits = (TClonesArray*) ioman->GetObject("RichHit");
126  if (NULL == fRichHits) { LOG(fatal) << "CbmRichGeoTest::Init No RichHit!"; }
127 
129  fDigiMan->Init();
131  LOG(fatal) << "CbmRichGeoTest::Init No RichDigi!";
132  }
133 
135  LOG(fatal) << "CbmRichGeoTest::Init No RichDigiMatch!";
136  }
137 
138  fRichRings = (TClonesArray*) ioman->GetObject("RichRing");
139  if (NULL == fRichRings) { LOG(fatal) << "CbmRichGeoTest::Init No RichRing!"; }
140 
141  fRichRingMatches = (TClonesArray*) ioman->GetObject("RichRingMatch");
142  if (NULL == fRichRingMatches) {
143  LOG(fatal) << "CbmRichGeoTest::Init No RichRingMatch!";
144  }
145 
146  fEventList = (CbmMCEventList*) ioman->GetObject("MCEventList.");
147  if (NULL == fEventList) {
148  LOG(fatal) << "CbmRichGeoTest::Init No MCEventList!";
149  }
150 
153 
154  InitHistograms();
155 
156  return kSUCCESS;
157 }
158 
159 void CbmRichGeoTest::Exec(Option_t* /*option*/) {
160  fEventNum++;
161  cout << "CbmRichGeoTest, event No. " << fEventNum << endl;
162 
163  ProcessMc();
164  RingParameters();
165  ProcessHits();
166 }
167 
169  // double xMin = -120.;
170  // double xMax = 120.;
171  // int nBinsX = 60;
172  // double yMin = -210;
173  // double yMax = 210.;
174  // int nBinsY = 105;
175 
176  double xMin = -120.;
177  double xMax = 120.;
178  int nBinsX = 300;
179  double yMin = -210;
180  double yMax = 210.;
181  int nBinsY = 105;
182 
183  fHM = new CbmHistManager();
184 
185  fHM->Create2<TH2D>("fhHitsXY",
186  "fhHitsXY;X [cm];Y [cm];Counter",
187  nBinsX,
188  xMin,
189  xMax,
190  nBinsY,
191  yMin,
192  yMax);
193  fHM->Create2<TH2D>("fhPointsXY",
194  "fhPointsXY;X [cm];Y [cm];Counter",
195  nBinsX,
196  xMin,
197  xMax,
198  nBinsY,
199  yMin,
200  yMax);
201  fHM->Create2<TH2D>("fhPointsXYNoRotation",
202  "fhPointsXYNoRotation;X [cm];Y [cm];Counter",
203  nBinsX,
204  xMin,
205  xMax,
206  nBinsY,
207  yMin,
208  yMax);
209  fHM->Create1<TH1D>("fhHitsZ", "fhHitsZ;Z [cm];Yield", 1000, 150, 250);
210  fHM->Create1<TH1D>("fhPointsZ", "fhPointsZ;Z [cm];Yield", 100, 190, 250);
211 
212  for (Int_t i = 0; i < 2; i++) {
213  stringstream ss;
214  if (i == 0) ss << "_hits";
215  if (i == 1) ss << "_points";
216  string t = ss.str();
217  if (i == 0)
218  fHM->Create1<TH1D>("fhNofHits" + t,
219  "fhNofHits" + t + ";Nof hits in ring;Yield",
220  50,
221  -.5,
222  49.5);
223  if (i == 1)
224  fHM->Create1<TH1D>("fhNofHits" + t,
225  "fhNofHits" + t + ";Nof points in ring;Yield",
226  400,
227  -.5,
228  399.5);
229  // ellipse fitting parameters
230  fHM->Create2<TH2D>("fhBoverAVsMom" + t,
231  "fhBoverAVsMom" + t + ";p [GeV/c];B/A;Yield",
232  40,
233  0.,
234  10,
235  100,
236  0,
237  1);
238  fHM->Create2<TH2D>("fhXcYcEllipse" + t,
239  "fhXcYcEllipse" + t + ";X [cm];Y [cm];Yield",
240  nBinsX,
241  xMin,
242  xMax,
243  nBinsY,
244  yMin,
245  yMax);
246  fHM->Create2<TH2D>("fhBaxisVsMom" + t,
247  "fhBaxisVsMom" + t + ";p [GeV/c];B axis [cm];Yield",
248  40,
249  0.,
250  10,
251  200,
252  0.,
253  10.);
254  fHM->Create2<TH2D>("fhAaxisVsMom" + t,
255  "fhAaxisVsMom" + t + ";p [GeV/c];A axis [cm];Yield",
256  40,
257  0.,
258  10,
259  200,
260  0.,
261  10.);
262  fHM->Create2<TH2D>("fhChi2EllipseVsMom" + t,
263  "fhChi2EllipseVsMom" + t + ";p [GeV/c];#Chi^{2};Yield",
264  40,
265  0.,
266  10.,
267  50,
268  0.,
269  0.5);
270  // circle fitting parameters
271  fHM->Create2<TH2D>("fhXcYcCircle" + t,
272  "fhXcYcCircle" + t + ";x [cm];y [cm];Yield",
273  nBinsX,
274  xMin,
275  xMax,
276  nBinsY,
277  yMin,
278  yMax);
279  fHM->Create2<TH2D>("fhRadiusVsMom" + t,
280  "fhRadiusVsMom" + t + ";p [GeV/c];Radius [cm];Yield",
281  40,
282  0.,
283  10,
284  200,
285  0.,
286  10.);
287  fHM->Create2<TH2D>("fhChi2CircleVsMom" + t,
288  "fhChi2CircleVsMom" + t + ";p [GeV/c];#Chi^{2};Yield",
289  40,
290  0.,
291  10.,
292  50,
293  0.,
294  .5);
295  fHM->Create2<TH2D>("fhDRVsMom" + t,
296  "fhDRVsMom" + t + ";p [GeV/c];dR [cm];Yield",
297  40,
298  0,
299  10,
300  200,
301  -2.,
302  2.);
303 
304  fHM->Create1<TH1D>("fhBaxisUpHalf" + t,
305  "fhBaxisUpHalf" + t + ";B axis [cm];Yield",
306  200,
307  0.,
308  10.);
309  fHM->Create1<TH1D>("fhBaxisDownHalf" + t,
310  "fhBaxisDownHalf" + t + ";B axis [cm];Yield",
311  200,
312  0.,
313  10.);
314  }
315 
316  fHM->Create1<TH1D>("fhNofPhotonsPerHit",
317  "fhNofPhotonsPerHit;Number of photons per hit;Yield",
318  10,
319  -0.5,
320  9.5);
321 
322  // Difference between Mc Points and Hits fitting.
323  fHM->Create2<TH2D>(
324  "fhDiffAaxis",
325  "fhDiffAaxis;Nof hits in ring;A_{point}-A_{hit} [cm];Yield",
326  40,
327  0.,
328  40.,
329  100,
330  -1.5,
331  1.5);
332  fHM->Create2<TH2D>(
333  "fhDiffBaxis",
334  "fhDiffBaxis;Nof hits in ring;B_{point}-B_{hit} [cm];Yield",
335  40,
336  0.,
337  40.,
338  100,
339  -1.5,
340  1.5);
341  fHM->Create2<TH2D>(
342  "fhDiffXcEllipse",
343  "fhDiffXcEllipse;Nof hits in ring;Xc_{point}-Xc_{hit} [cm];Yield",
344  40,
345  0.,
346  40.,
347  100,
348  -1.5,
349  1.5);
350  fHM->Create2<TH2D>(
351  "fhDiffYcEllipse",
352  "fhDiffYcEllipse;Nof hits in ring;Yc_{point}-Yc_{hit} [cm];Yield",
353  40,
354  0.,
355  40.,
356  100,
357  -1.5,
358  1.5);
359  fHM->Create2<TH2D>(
360  "fhDiffXcCircle",
361  "fhDiffXcCircle;Nof hits in ring;Xc_{point}-Xc_{hit} [cm];Yield",
362  40,
363  0.,
364  40.,
365  100,
366  -1.5,
367  1.5);
368  fHM->Create2<TH2D>(
369  "fhDiffYcCircle",
370  "fhDiffYcCircle;Nof hits in ring;Yc_{point}-Yc_{hit} [cm];Yield",
371  40,
372  0.,
373  40.,
374  100,
375  -1.5,
376  1.5);
377  fHM->Create2<TH2D>(
378  "fhDiffRadius",
379  "fhDiffRadius;Nof hits in ring;Radius_{point}-Radius_{hit} [cm];Yield",
380  40,
381  0.,
382  40.,
383  100,
384  -1.5,
385  1.5);
386 
387  // R, A, B distribution for different number of hits from 0 to 40.
388  fHM->Create2<TH2D>("fhRadiusVsNofHits",
389  "fhRadiusVsNofHits;Nof hits in ring;Radius [cm];Yield",
390  40,
391  0.,
392  40.,
393  100,
394  0.,
395  10.);
396  fHM->Create2<TH2D>("fhAaxisVsNofHits",
397  "fhAaxisVsNofHits;Nof hits in ring;A axis [cm];Yield",
398  40,
399  0.,
400  40.,
401  100,
402  0.,
403  10.);
404  fHM->Create2<TH2D>("fhBaxisVsNofHits",
405  "fhBaxisVsNofHits;Nof hits in ring;B axis [cm];Yield",
406  40,
407  0.,
408  40.,
409  100,
410  0.,
411  10.);
412  fHM->Create2<TH2D>("fhDRVsNofHits",
413  "fhDRVsNofHits;Nof hits in ring;dR [cm];Yield",
414  40,
415  0.,
416  40.,
417  200,
418  -2.,
419  2.);
420 
421  // Hits and points.
422  fHM->Create1<TH1D>(
423  "fhDiffXhit", "fhDiffXhit;X_{point}-X_{hit} [cm];Yield", 200, -.5, .5);
424  fHM->Create1<TH1D>(
425  "fhDiffYhit", "fhDiffYhit;Y_{point}-Y_{hit} [cm];Yield", 200, -.5, .5);
426 
427  // Fitting efficiency.
428  fHM->Create1<TH1D>(
429  "fhNofHitsAll", "fhNofHitsAll;Nof hits in ring;Yield", 50, 0., 50.);
430  fHM->Create1<TH1D>("fhNofHitsCircleFit",
431  "fhNofHitsCircleFit;Nof hits in ring;Yield",
432  50,
433  0.,
434  50.);
435  fHM->Create1<TH1D>("fhNofHitsEllipseFit",
436  "fhNofHitsEllipseFit;Nof hits in ring;Yield",
437  50,
438  0.,
439  50.);
440  fHM->Create1<TH1D>("fhNofHitsCircleFitEff",
441  "fhNofHitsCircleFitEff;Nof hits in ring;Efficiency [%]",
442  50,
443  0.,
444  50.);
445  fHM->Create1<TH1D>("fhNofHitsEllipseFitEff",
446  "fhNofHitsEllipseFitEff;Nof hits in ring;Efficiency [%]",
447  50,
448  0.,
449  50.);
450 
451  // Detector acceptance efficiency vs. (pt,y) and p
452  fHM->Create1<TH1D>("fhMcMomEl", "fhMcMomEl;p [GeV/c];Yield", 24, 0., 12.);
453  fHM->Create2<TH2D>("fhMcPtyEl",
454  "fhMcPtyEl;Rapidity;P_{t} [GeV/c];Yield",
455  25,
456  0.,
457  4.,
458  20,
459  0.,
460  3.);
461  fHM->Create1<TH1D>("fhAccMomEl", "fhAccMomEl;p [GeV/c];Yield", 24, 0., 12.);
462  fHM->Create2<TH2D>("fhAccPtyEl",
463  "fhAccPtyEl;Rapidity;P_{t} [GeV/c];Yield",
464  25,
465  0.,
466  4.,
467  20,
468  0.,
469  3.);
470 
471  fHM->Create1<TH1D>("fhMcMomPi", "fhMcMomPi;p [GeV/c];Yield", 24, 0., 12.);
472  fHM->Create2<TH2D>("fhMcPtyPi",
473  "fhMcPtyPi;Rapidity;P_{t} [GeV/c];Yield",
474  25,
475  0.,
476  4.,
477  20,
478  0.,
479  3.);
480  fHM->Create1<TH1D>("fhAccMomPi", "fhAccMomPi;p [GeV/c];Yield", 24, 0., 12.);
481  fHM->Create2<TH2D>("fhAccPtyPi",
482  "fhAccPtyPi;Rapidity;P_{t} [GeV/c];Yield",
483  25,
484  0.,
485  4.,
486  20,
487  0.,
488  3.);
489 
490  // Numbers in dependence on XY position onto the photodetector.
491  fHM->Create3<TH3D>("fhNofHitsXYZ",
492  "fhNofHitsXYZ;X [cm];Y [cm];Nof hits in ring",
493  nBinsX,
494  xMin,
495  xMax,
496  nBinsY,
497  yMin,
498  yMax,
499  50,
500  0.,
501  50);
502  fHM->Create3<TH3D>("fhNofPointsXYZ",
503  "fhNofPointsXYZ;X [cm];Y [cm];Nof points in ring",
504  nBinsX,
505  xMin,
506  xMax,
507  nBinsY,
508  yMin,
509  yMax,
510  50,
511  100.,
512  300.);
513  fHM->Create3<TH3D>("fhBoverAXYZ",
514  "fhBoverAXYZ;X [cm];Y [cm];B/A",
515  nBinsX,
516  xMin,
517  xMax,
518  nBinsY,
519  yMin,
520  yMax,
521  100,
522  0.,
523  1.);
524  fHM->Create3<TH3D>("fhBaxisXYZ",
525  "fhBaxisXYZ;X [cm];Y [cm];B axis [cm]",
526  nBinsX,
527  xMin,
528  xMax,
529  nBinsY,
530  yMin,
531  yMax,
532  80,
533  3.,
534  7.);
535  fHM->Create3<TH3D>("fhAaxisXYZ",
536  "fhAaxisXYZ;X [cm];Y [cm];A axis [cm]",
537  nBinsX,
538  xMin,
539  xMax,
540  nBinsY,
541  yMin,
542  yMax,
543  80,
544  3.,
545  7.);
546  fHM->Create3<TH3D>("fhRadiusXYZ",
547  "fhRadiusXYZ;X [cm];Y [cm];Radius [cm]",
548  nBinsX,
549  xMin,
550  xMax,
551  nBinsY,
552  yMin,
553  yMax,
554  80,
555  3.,
556  7.);
557  fHM->Create3<TH3D>("fhdRXYZ",
558  "fhdRXYZ;X [cm];Y [cm];dR [cm]",
559  nBinsX,
560  xMin,
561  xMax,
562  nBinsY,
563  yMin,
564  yMax,
565  100,
566  -1.,
567  1.);
568 
569  int nBinsX1 = 60;
570  int xMin1 = -120;
571  int xMax1 = 120;
572  int nBinsY1 = 25;
573  int yMin1 = 100;
574  int yMax1 = 200;
575  // Numbers in dependence X or Y position onto the photodetector plane
576  fHM->Create2<TH2D>("fhNofHitsVsX",
577  "fhNofHitsVsX;X [cm];Nof hits in ring",
578  nBinsX1,
579  xMin1,
580  xMax1,
581  50,
582  0.,
583  50);
584  fHM->Create2<TH2D>("fhNofHitsVsY",
585  "fhNofHitsVsY;Abs(Y) [cm];Nof hits in ring",
586  nBinsY1,
587  yMin1,
588  yMax1,
589  50,
590  0.,
591  50);
592 
593  fHM->Create2<TH2D>("fhNofPointsVsX",
594  "fhNofPointsVsX;X [cm];Nof points in ring",
595  nBinsX1,
596  xMin1,
597  xMax1,
598  50,
599  100.,
600  300.);
601  fHM->Create2<TH2D>("fhNofPointsVsY",
602  "fhNofPointsVsY;Abs(Y) [cm];Nof points in ring",
603  nBinsY1,
604  yMin1,
605  yMax1,
606  50,
607  100.,
608  300.);
609 
610  fHM->Create2<TH2D>("fhBoverAVsX",
611  "fhBoverAVsX;X [cm];B/A",
612  nBinsX1,
613  xMin1,
614  xMax1,
615  100,
616  0.,
617  1.);
618  fHM->Create2<TH2D>("fhBoverAVsY",
619  "fhBoverAVsY;Abs(Y) [cm];B/A",
620  nBinsY1,
621  yMin1,
622  yMax1,
623  100,
624  0.,
625  1.);
626 
627  fHM->Create2<TH2D>("fhBaxisVsX",
628  "fhBaxisVsX;X [cm];B axis [cm]",
629  nBinsX1,
630  xMin1,
631  xMax1,
632  80,
633  3.,
634  7.);
635  fHM->Create2<TH2D>("fhBaxisVsY",
636  "fhBaxisVsY;Abs(Y) [cm];B axis [cm]",
637  nBinsY1,
638  yMin1,
639  yMax1,
640  80,
641  3.,
642  7.);
643 
644  fHM->Create2<TH2D>("fhAaxisVsX",
645  "fhAaxisVsX;X [cm];A axis [cm]",
646  nBinsX1,
647  xMin1,
648  xMax1,
649  80,
650  3.,
651  7.);
652  fHM->Create2<TH2D>("fhAaxisVsY",
653  "fhAaxisVsY;Abs(Y) [cm];A axis [cm]",
654  nBinsY1,
655  yMin1,
656  yMax1,
657  80,
658  3.,
659  7.);
660 
661  fHM->Create2<TH2D>("fhRadiusVsX",
662  "fhRadiusVsX;X [cm];Radius [cm]",
663  nBinsX1,
664  xMin1,
665  xMax1,
666  80,
667  3.,
668  7.);
669  fHM->Create2<TH2D>("fhRadiusVsY",
670  "fhRadiusVsY;Abs(Y) [cm];Radius [cm]",
671  nBinsY1,
672  yMin1,
673  yMax1,
674  80,
675  3.,
676  7.);
677 
678  fHM->Create2<TH2D>(
679  "fhdRVsX", "fhdRVsX;X [cm];dR [cm]", nBinsX1, xMin1, xMax1, 100, -1., 1.);
680  fHM->Create2<TH2D>("fhdRVsY",
681  "fhdRVsY;Abs(Y) [cm];dR [cm]",
682  nBinsY1,
683  yMin1,
684  yMax1,
685  100,
686  -1.,
687  1.);
688 
689  fHM->Create1<TH1D>("fhPhotonEPlaneZ+",
690  "fhPhotonEPlaneZ+;Photon energy [eV];Counter",
691  100,
692  0.,
693  10.);
694  fHM->Create1<TH1D>("fhPhotonEPlaneZ-",
695  "fhPhotonEPlaneZ-;Photon energy [eV];Counter",
696  100,
697  0.,
698  10.);
699  fHM->Create1<TH1D>("fhPhotonEPmtPoint",
700  "fhPhotonEPmtPoint;Photon energy [eV];Counter",
701  100,
702  0.,
703  10.);
704  fHM->Create1<TH1D>("fhPhotonEPmtHit",
705  "fhPhotonEPmtHit;Photon energy [eV];Counter",
706  100,
707  0.,
708  10.);
709 
710  fHM->Create1<TH1D>("fhLambdaPlaneZ+",
711  "fhLambdaPlaneZ+;Photon wavelength [nm];Counter",
712  100,
713  0.,
714  700.);
715  fHM->Create1<TH1D>("fhLambdaPlaneZ-",
716  "fhLambdaPlaneZ-;Photon wavelength [nm];Counter",
717  100,
718  0.,
719  700.);
720  fHM->Create1<TH1D>("fhLambdaPmtPoint",
721  "fhLambdaPmtPoint;Photon wavelength [nm];Counter",
722  100,
723  0.,
724  700.);
725  fHM->Create1<TH1D>("fhLambdaPmtHit",
726  "fhLambdaPmtHit;Photon wavelength [nm];Counter",
727  100,
728  0.,
729  700.);
730 }
731 
733  Int_t nofEvents = fEventList->GetNofEvents();
734  for (Int_t iE = 0; iE < nofEvents; iE++) {
735  Int_t fileId = fEventList->GetFileIdByIndex(iE);
736  Int_t eventId = fEventList->GetEventIdByIndex(iE);
737 
738  Int_t nofMcTracks = fMcTracks->Size(fileId, eventId);
739  for (Int_t iT = 0; iT < nofMcTracks; iT++) {
740  const CbmMCTrack* mcTrack =
741  static_cast<CbmMCTrack*>(fMcTracks->Get(fileId, eventId, iT));
742  if (!mcTrack) continue;
743  Int_t motherId = mcTrack->GetMotherId();
744  Int_t pdg = TMath::Abs(mcTrack->GetPdgCode());
745  Bool_t isMcPrimaryElectron =
746  (pdg == 11 && motherId == -1)
747  || (mcTrack->GetGeantProcessId() == kPPrimary && pdg == 11);
748 
749  if (isMcPrimaryElectron) {
750  fHM->H1("fhMcMomEl")->Fill(mcTrack->GetP());
751  fHM->H2("fhMcPtyEl")->Fill(mcTrack->GetRapidity(), mcTrack->GetPt());
752  }
753 
754  if (pdg == 211 && motherId == -1) {
755  fHM->H1("fhMcMomPi")->Fill(mcTrack->GetP());
756  fHM->H2("fhMcPtyPi")->Fill(mcTrack->GetRapidity(), mcTrack->GetPt());
757  }
758  }
759 
760  Int_t nofPoints = fRichPoints->Size(fileId, eventId);
761  for (Int_t iP = 0; iP < nofPoints; iP++) {
762  const CbmRichPoint* point =
763  static_cast<CbmRichPoint*>(fRichPoints->Get(fileId, eventId, iP));
764  if (point == NULL) continue;
765  TVector3 inPos(point->GetX(), point->GetY(), point->GetZ());
766  TVector3 outPos;
767  CbmRichGeoManager::GetInstance().RotatePoint(&inPos, &outPos);
768  fHM->H2("fhPointsXYNoRotation")->Fill(point->GetX(), point->GetY());
769  fHM->H2("fhPointsXY")->Fill(outPos.X(), outPos.Y());
770  fHM->H1("fhPointsZ")->Fill(point->GetZ());
771 
772  TVector3 mom;
773  point->Momentum(mom);
774  Double_t momMag = mom.Mag();
775  fHM->H1("fhPhotonEPmtPoint")->Fill(1.e9 * momMag);
776  Double_t lambda = CbmRichPmt::getLambda(momMag);
777  fHM->H1("fhLambdaPmtPoint")->Fill(lambda);
778  }
779 
780  Int_t nofPlanePoints = fRichRefPlanePoints->Size(fileId, eventId);
781  for (Int_t iP = 0; iP < nofPlanePoints; iP++) {
782  const CbmRichPoint* point = static_cast<CbmRichPoint*>(
783  fRichRefPlanePoints->Get(fileId, eventId, iP));
784  if (point == NULL) continue;
785  TVector3 mom;
786  point->Momentum(mom);
787  Double_t momMag = mom.Mag();
788  Double_t lambda = CbmRichPmt::getLambda(momMag);
789  if (point->GetPz() > 0) {
790  fHM->H1("fhPhotonEPlaneZ+")->Fill(1.e9 * momMag);
791  fHM->H1("fhLambdaPlaneZ+")->Fill(lambda);
792  } else {
793  fHM->H1("fhPhotonEPlaneZ-")->Fill(1.e9 * momMag);
794  fHM->H1("fhLambdaPlaneZ-")->Fill(lambda);
795  }
796  }
797  }
798 }
799 
801  Int_t fileId = 0;
802  Int_t nofRings = fRichRings->GetEntriesFast();
803  for (Int_t iR = 0; iR < nofRings; iR++) {
804  CbmRichRing* ring = (CbmRichRing*) fRichRings->At(iR);
805  if (NULL == ring) continue;
806  const CbmTrackMatchNew* ringMatch =
807  static_cast<CbmTrackMatchNew*>(fRichRingMatches->At(iR));
808  if (NULL == ringMatch) continue;
809  Int_t mcEventId = ringMatch->GetMatchedLink().GetEntry();
810  Int_t mcTrackId = ringMatch->GetMatchedLink().GetIndex();
811  const CbmMCTrack* mcTrack =
812  static_cast<CbmMCTrack*>(fMcTracks->Get(fileId, mcEventId, mcTrackId));
813  if (NULL == mcTrack) continue;
814 
815  Int_t motherId = mcTrack->GetMotherId();
816  Int_t pdg = TMath::Abs(mcTrack->GetPdgCode());
817  Double_t momentum = mcTrack->GetP();
818  Double_t pt = mcTrack->GetPt();
819  Double_t rapidity = mcTrack->GetRapidity();
820  Bool_t isMcPrimaryElectron =
821  (pdg == 11 && motherId == -1)
822  || (mcTrack->GetGeantProcessId() == kPPrimary && pdg == 11);
823 
824  if (ring->GetNofHits() >= fMinNofHits) {
825  if (isMcPrimaryElectron) {
826  fHM->H1("fhAccMomEl")->Fill(momentum);
827  fHM->H2("fhAccPtyEl")->Fill(rapidity, pt);
828  }
829  if (pdg == 211 && motherId == -1) {
830  fHM->H1("fhAccMomPi")->Fill(momentum);
831  fHM->H2("fhAccPtyPi")->Fill(rapidity, pt);
832  }
833  }
834 
835  if (!isMcPrimaryElectron) continue; // only primary electrons
836 
837  CbmRichRingLight ringPoint;
838  int nofRichPoints = fRichPoints->Size(fileId, mcEventId);
839  for (int iPoint = 0; iPoint < nofRichPoints; iPoint++) {
840  const CbmRichPoint* richPoint =
841  static_cast<CbmRichPoint*>(fRichPoints->Get(fileId, mcEventId, iPoint));
842  if (NULL == richPoint) continue;
843  Int_t trackId = richPoint->GetTrackID();
844  if (trackId < 0) continue;
845  CbmMCTrack* mcTrackRich =
846  static_cast<CbmMCTrack*>(fMcTracks->Get(fileId, mcEventId, trackId));
847  if (NULL == mcTrackRich) continue;
848  int motherIdRich = mcTrackRich->GetMotherId();
849  if (motherIdRich == mcTrackId) {
850  TVector3 posPoint;
851  richPoint->Position(posPoint);
852  TVector3 detPoint;
853  CbmRichGeoManager::GetInstance().RotatePoint(&posPoint, &detPoint);
854  CbmRichHitLight hit(detPoint.X(), detPoint.Y());
855  ringPoint.AddHit(hit);
856  }
857  }
858  fHM->H1("fhNofHitsAll")->Fill(ring->GetNofHits());
859 
860  CbmRichRingLight ringHit;
862 
863  FitAndFillHistCircle(0, &ringHit, momentum); //hits
864  FitAndFillHistCircle(1, &ringPoint, momentum); // points
865  FillMcVsHitFitCircle(&ringHit, &ringPoint);
866 
867  double r = ringHit.GetRadius();
868  double xc = ringHit.GetCenterX();
869  double yc = ringHit.GetCenterY();
870 
871  if (ringHit.GetRadius() > fMinRadius && ringHit.GetRadius() < fMaxRadius) {
872  fHM->H1("fhNofHitsCircleFit")->Fill(ringHit.GetNofHits());
873  }
874  if (fDrawEventDisplay && fNofDrawnRings < 10) {
875  DrawRing(&ringHit, &ringPoint);
876  }
877 
878  FitAndFillHistEllipse(0, &ringHit, momentum); // hits
879  FitAndFillHistEllipse(1, &ringPoint, momentum); // points
880  FillMcVsHitFitEllipse(&ringHit, &ringPoint);
881 
882  if (ringHit.GetAaxis() > fMinAaxis && ringHit.GetAaxis() < fMaxAaxis
883  && ringHit.GetBaxis() > fMinBaxis && ringHit.GetAaxis() < fMaxBaxis) {
884  fHM->H1("fhNofHitsEllipseFit")->Fill(ringHit.GetNofHits());
885 
886  double np = ringPoint.GetNofHits();
887  double a = ringHit.GetAaxis();
888  double b = ringHit.GetBaxis();
889  double nh = ring->GetNofHits();
890 
891  fHM->H3("fhNofHitsXYZ")->Fill(xc, yc, nh);
892  fHM->H3("fhNofPointsXYZ")->Fill(xc, yc, np);
893  fHM->H3("fhBoverAXYZ")->Fill(xc, yc, b / a);
894  fHM->H3("fhBaxisXYZ")->Fill(xc, yc, b);
895  fHM->H3("fhAaxisXYZ")->Fill(xc, yc, a);
896  fHM->H3("fhRadiusXYZ")->Fill(xc, yc, r);
897 
898  fHM->H2("fhNofHitsVsX")->Fill(xc, nh);
899  fHM->H2("fhNofPointsVsX")->Fill(xc, np);
900  fHM->H2("fhBoverAVsX")->Fill(xc, b / a);
901  fHM->H2("fhBaxisVsX")->Fill(xc, b);
902  fHM->H2("fhAaxisVsX")->Fill(xc, a);
903  fHM->H2("fhRadiusVsX")->Fill(xc, r);
904 
905  fHM->H2("fhNofHitsVsY")->Fill(abs(yc), nh);
906  fHM->H2("fhNofPointsVsY")->Fill(abs(yc), np);
907  fHM->H2("fhBoverAVsY")->Fill(abs(yc), b / a);
908  fHM->H2("fhBaxisVsY")->Fill(abs(yc), b);
909  fHM->H2("fhAaxisVsY")->Fill(abs(yc), a);
910  fHM->H2("fhRadiusVsY")->Fill(abs(yc), r);
911 
912  for (int iH = 0; iH < ringHit.GetNofHits(); iH++) {
913  double xh = ringHit.GetHit(iH).fX;
914  double yh = ringHit.GetHit(iH).fY;
915  double dr = r - sqrt((xc - xh) * (xc - xh) + (yc - yh) * (yc - yh));
916  fHM->H3("fhdRXYZ")->Fill(xc, yc, dr);
917  fHM->H2("fhdRVsX")->Fill(xc, dr);
918  fHM->H2("fhdRVsY")->Fill(abs(yc), dr);
919  }
920  }
921  } // iR
922 }
923 
925  CbmRichRingLight* ring,
926  double momentum) {
927  fTauFit->DoFit(ring);
928  double axisA = ring->GetAaxis();
929  double axisB = ring->GetBaxis();
930  double xcEllipse = ring->GetCenterX();
931  double ycEllipse = ring->GetCenterY();
932  int nofHitsRing = ring->GetNofHits();
933  string t = "";
934  if (histIndex == 0) {
935  t = "_hits";
936  } else if (histIndex == 1) {
937  t = "_points";
938  }
939  if (axisA > fMinAaxis && axisA < fMaxAaxis && axisB > fMinBaxis
940  && axisB < fMaxBaxis) {
941  fHM->H1("fhBoverAVsMom" + t)->Fill(momentum, axisB / axisA);
942  fHM->H2("fhXcYcEllipse" + t)->Fill(xcEllipse, ycEllipse);
943  }
944  fHM->H1("fhNofHits" + t)->Fill(nofHitsRing);
945 
946  if (ycEllipse > 149 || ycEllipse < -149) {
947  fHM->H1("fhBaxisUpHalf" + t)->Fill(axisB);
948  } else {
949  fHM->H1("fhBaxisDownHalf" + t)->Fill(axisB);
950  }
951 
952  fHM->H2("fhBaxisVsMom" + t)->Fill(momentum, axisB);
953  fHM->H2("fhAaxisVsMom" + t)->Fill(momentum, axisA);
954  fHM->H2("fhChi2EllipseVsMom" + t)
955  ->Fill(momentum, ring->GetChi2() / ring->GetNofHits());
956  if (histIndex == 0) { // only hit fit
957  fHM->H2("fhAaxisVsNofHits")->Fill(nofHitsRing, axisA);
958  fHM->H2("fhBaxisVsNofHits")->Fill(nofHitsRing, axisB);
959  }
960 }
961 
963  CbmRichRingLight* ring,
964  double momentum) {
965  fCopFit->DoFit(ring);
966  double radius = ring->GetRadius();
967  double xcCircle = ring->GetCenterX();
968  double ycCircle = ring->GetCenterY();
969  int nofHitsRing = ring->GetNofHits();
970  string t = "";
971  if (histIndex == 0) {
972  t = "_hits";
973  } else if (histIndex == 1) {
974  t = "_points";
975  }
976  fHM->H1("fhXcYcCircle" + t)->Fill(xcCircle, ycCircle);
977  fHM->H1("fhRadiusVsMom" + t)->Fill(momentum, radius);
978  fHM->H1("fhChi2CircleVsMom" + t)
979  ->Fill(momentum, ring->GetChi2() / ring->GetNofHits());
980 
981  for (int iH = 0; iH < nofHitsRing; iH++) {
982  double xh = ring->GetHit(iH).fX;
983  double yh = ring->GetHit(iH).fY;
984  double dr = radius
985  - sqrt((xcCircle - xh) * (xcCircle - xh)
986  + (ycCircle - yh) * (ycCircle - yh));
987  fHM->H1("fhDRVsMom" + t)->Fill(momentum, dr);
988 
989  if (histIndex == 0) { // only hit fit
990  fHM->H2("fhDRVsNofHits")->Fill(nofHitsRing, dr);
991  }
992  }
993 
994  if (histIndex == 0) { // only hit fit
995  fHM->H2("fhRadiusVsNofHits")->Fill(nofHitsRing, radius);
996  }
997 }
998 
1000  CbmRichRingLight* ringMc) {
1001  fHM->H2("fhDiffAaxis")
1002  ->Fill(ring->GetNofHits(), ringMc->GetAaxis() - ring->GetAaxis());
1003  fHM->H2("fhDiffBaxis")
1004  ->Fill(ring->GetNofHits(), ringMc->GetBaxis() - ring->GetBaxis());
1005  fHM->H2("fhDiffXcEllipse")
1006  ->Fill(ring->GetNofHits(), ringMc->GetCenterX() - ring->GetCenterX());
1007  fHM->H2("fhDiffYcEllipse")
1008  ->Fill(ring->GetNofHits(), ringMc->GetCenterY() - ring->GetCenterY());
1009 }
1010 
1012  CbmRichRingLight* ringMc) {
1013  fHM->H2("fhDiffXcCircle")
1014  ->Fill(ring->GetNofHits(), ringMc->GetCenterX() - ring->GetCenterX());
1015  fHM->H2("fhDiffYcCircle")
1016  ->Fill(ring->GetNofHits(), ringMc->GetCenterY() - ring->GetCenterY());
1017  fHM->H2("fhDiffRadius")
1018  ->Fill(ring->GetNofHits(), ringMc->GetRadius() - ring->GetRadius());
1019 }
1020 
1022  int fileId = 0;
1023  Int_t nofHits = fRichHits->GetEntriesFast();
1024  for (Int_t iH = 0; iH < nofHits; iH++) {
1025  const CbmRichHit* hit = static_cast<CbmRichHit*>(fRichHits->At(iH));
1026  if (hit == NULL) continue;
1027 
1028  Int_t digiIndex = hit->GetRefId();
1029  if (digiIndex < 0) continue;
1030  const CbmRichDigi* digi = fDigiMan->Get<CbmRichDigi>(digiIndex);
1031  if (NULL == digi) continue;
1032  const CbmMatch* digiMatch =
1033  fDigiMan->GetMatch(ECbmModuleId::kRich, digiIndex);
1034  if (NULL == digiMatch) continue;
1035 
1036  vector<CbmLink> links = digiMatch->GetLinks();
1037  for (UInt_t i = 0; i < links.size(); i++) {
1038  Int_t pointId = links[i].GetIndex();
1039  Int_t eventId = links[i].GetEntry();
1040  if (pointId < 0) continue; // noise hit
1041 
1042  const CbmRichPoint* pMCpt =
1043  static_cast<CbmRichPoint*>(fRichPoints->Get(fileId, eventId, pointId));
1044  if (NULL == pMCpt) continue;
1045 
1046  TVector3 inPos(pMCpt->GetX(), pMCpt->GetY(), pMCpt->GetZ());
1047  TVector3 outPos;
1048  CbmRichGeoManager::GetInstance().RotatePoint(&inPos, &outPos);
1049  fHM->H1("fhDiffXhit")->Fill(hit->GetX() - outPos.X());
1050  fHM->H1("fhDiffYhit")->Fill(hit->GetY() - outPos.Y());
1051 
1052  TVector3 mom;
1053  pMCpt->Momentum(mom);
1054  Double_t momMag = mom.Mag();
1055  fHM->H1("fhPhotonEPmtHit")->Fill(1.e9 * momMag);
1056  Double_t lambda = CbmRichPmt::getLambda(momMag);
1057  fHM->H1("fhLambdaPmtHit")->Fill(lambda);
1058  }
1059 
1060  //fHM->H1("fhNofPhotonsPerHit")->Fill(hit->GetNPhotons());
1061  fHM->H2("fhHitsXY")->Fill(hit->GetX(), hit->GetY());
1062  fHM->H1("fhHitsZ")->Fill(hit->GetZ());
1063  }
1064 }
1065 
1067  CbmRichRingLight* ringPoint) {
1068  stringstream ss;
1069  ss << "Event" << fNofDrawnRings;
1070  fNofDrawnRings++;
1071  TCanvas* c = fHM->CreateCanvas(ss.str().c_str(), ss.str().c_str(), 500, 500);
1072  c->SetGrid(true, true);
1073  TH2D* pad = new TH2D(ss.str().c_str(),
1074  (ss.str() + ";X [cm];Y [cm]").c_str(),
1075  1,
1076  -15.,
1077  15.,
1078  1,
1079  -15.,
1080  15);
1081  pad->SetStats(false);
1082  pad->Draw();
1083 
1084  // find min and max x and y positions of the hits
1085  // in order to shift drawing
1086  double xmin = 99999., xmax = -99999., ymin = 99999., ymax = -99999.;
1087  for (int i = 0; i < ringHit->GetNofHits(); i++) {
1088  double hitX = ringHit->GetHit(i).fX;
1089  double hitY = ringHit->GetHit(i).fY;
1090  if (xmin > hitX) xmin = hitX;
1091  if (xmax < hitX) xmax = hitX;
1092  if (ymin > hitY) ymin = hitY;
1093  if (ymax < hitY) ymax = hitY;
1094  }
1095  double xCur = (xmin + xmax) / 2.;
1096  double yCur = (ymin + ymax) / 2.;
1097 
1098  //Draw circle and center
1099  TEllipse* circle = new TEllipse(ringHit->GetCenterX() - xCur,
1100  ringHit->GetCenterY() - yCur,
1101  ringHit->GetRadius());
1102  circle->SetFillStyle(0);
1103  circle->SetLineWidth(3);
1104  circle->Draw();
1105  TEllipse* center = new TEllipse(
1106  ringHit->GetCenterX() - xCur, ringHit->GetCenterY() - yCur, .5);
1107  center->Draw();
1108 
1109  // Draw hits
1110  for (int i = 0; i < ringHit->GetNofHits(); i++) {
1111  TEllipse* hitDr = new TEllipse(
1112  ringHit->GetHit(i).fX - xCur, ringHit->GetHit(i).fY - yCur, .5);
1113  hitDr->SetFillColor(kRed);
1114  hitDr->Draw();
1115  }
1116 
1117  // Draw MC Points
1118  for (int i = 0; i < ringPoint->GetNofHits(); i++) {
1119  TEllipse* pointDr = new TEllipse(
1120  ringPoint->GetHit(i).fX - xCur, ringPoint->GetHit(i).fY - yCur, 0.15);
1121  pointDr->SetFillColor(kBlue);
1122  pointDr->Draw();
1123  }
1124 
1125  //Draw information
1126  stringstream ss2;
1127  ss2 << "(r, n)=(" << setprecision(3) << ringHit->GetRadius() << ", "
1128  << ringHit->GetNofHits() << ")";
1129  TLatex* latex = new TLatex(-8., 8., ss2.str().c_str());
1130  latex->Draw();
1131 }
1132 
1134  Int_t nofMc = (Int_t) fHM->H1("fhMcMomEl")->GetEntries();
1135  TH1D* hist =
1136  (TH1D*) fHM->H1("fhNofHits_hits")->Clone("fhAccVsMinNofHitsHist");
1137  hist->GetXaxis()->SetTitle("Required min nof hits in ring");
1138  hist->GetYaxis()->SetTitle("Detector acceptance [%]");
1139  Double_t sum = 0.;
1140  for (int i = hist->GetNbinsX(); i > 1; i--) {
1141  sum += fHM->H1("fhNofHits_hits")->GetBinContent(i);
1142  hist->SetBinContent(i, 100. * sum / nofMc);
1143  }
1144  return hist;
1145 }
1146 
1147 void CbmRichGeoTest::DrawH2MeanRms(TH2* hist, const string& canvasName) {
1148  TCanvas* c =
1149  fHM->CreateCanvas(canvasName.c_str(), canvasName.c_str(), 1200, 600);
1150  c->Divide(2, 1);
1151  c->cd(1);
1152  DrawH2WithProfile(hist);
1153  c->cd(2);
1154  TH1D* py =
1155  (TH1D*) hist->ProjectionY((string(hist->GetName()) + "_py").c_str())
1156  ->Clone();
1157  DrawH1andFitGauss(py);
1158  py->Scale(1. / py->Integral());
1159  py->GetYaxis()->SetTitle("Yield");
1160 }
1161 
1164 
1165  {
1166  TCanvas* c =
1167  fHM->CreateCanvas("richgeo_hits_xy", "richgeo_hits_xy", 1200, 1200);
1168  CbmRichDraw::DrawPmtH2(fHM->H2("fhHitsXY"), c);
1169  }
1170 
1171  {
1172  TCanvas* c =
1173  fHM->CreateCanvas("richgeo_points_xy", "richgeo_points_xy", 1200, 1200);
1174  CbmRichDraw::DrawPmtH2(fHM->H2("fhPointsXY"), c);
1175  }
1176 
1177  {
1178  TCanvas* c = fHM->CreateCanvas("richgeo_points_xy_no_rotation",
1179  "richgeo_points_xy_no_rotation",
1180  1200,
1181  1200);
1182  CbmRichDraw::DrawPmtH2(fHM->H2("fhPointsXYNoRotation"), c);
1183  }
1184 
1185  {
1186  fHM->CreateCanvas("richgeo_hits_z", "richgeo_hits_z", 800, 800);
1187  fHM->H1("fhHitsZ")->Scale(1. / fHM->H1("fhHitsZ")->Integral());
1188  DrawH1(fHM->H1("fhHitsZ"), kLinear, kLinear, "hist");
1189  }
1190 
1191  {
1192  fHM->CreateCanvas("richgeo_points_z", "richgeo_points_z", 800, 800);
1193  fHM->H1("fhPointsZ")->Scale(1. / fHM->H1("fhPointsZ")->Integral());
1194  DrawH1(fHM->H1("fhPointsZ"), kLinear, kLinear, "hist");
1195  }
1196 
1197  for (int i = 0; i < 2; i++) {
1198  string t;
1199  if (i == 0) {
1200  t = "_hits";
1201  } else if (i == 1) {
1202  t = "_points";
1203  }
1204 
1205  DrawH2MeanRms((TH2D*) fHM->H2("fhBoverAVsMom" + t),
1206  "richgeo" + t + "_ellipse_boa_vs_mom");
1207  fHM->H2("fhBoverAVsMom" + t)->GetYaxis()->SetRangeUser(0.8, 1.0);
1208 
1209  {
1210  fHM->CreateCanvas(("richgeo" + t + "_ellipse_xc_yc").c_str(),
1211  ("richgeo" + t + "_ellipse_xc_yc").c_str(),
1212  800,
1213  800);
1214  DrawH2(fHM->H2("fhXcYcEllipse" + t));
1215  }
1216 
1217  DrawH2MeanRms((TH2D*) fHM->H2("fhChi2EllipseVsMom" + t),
1218  "richgeo" + t + "_chi2_ellipse_vs_mom");
1219  DrawH2MeanRms((TH2D*) fHM->H2("fhAaxisVsMom" + t),
1220  "richgeo" + t + "_a_vs_mom");
1221  DrawH2MeanRms((TH2D*) fHM->H2("fhBaxisVsMom" + t),
1222  "richgeo" + t + "_b_vs_mom");
1223 
1224  {
1225  TCanvas* c =
1226  fHM->CreateCanvas(("richgeo" + t + "_b_up_down_halves").c_str(),
1227  ("richgeo" + t + "_b_up_down_halves").c_str(),
1228  1200,
1229  600);
1230  c->Divide(2, 1);
1231  c->cd(1);
1233  (TH1D*) fHM->H1("fhBaxisUpHalf" + t)->Clone(), true, true, 3., 6.);
1234  c->cd(2);
1236  (TH1D*) fHM->H1("fhBaxisDownHalf" + t)->Clone(), true, true, 3., 6.);
1237  }
1238 
1239  {
1240  TCanvas* c = fHM->CreateCanvas(("richgeo" + t + "_circle").c_str(),
1241  ("richgeo" + t + "_circle").c_str(),
1242  1200,
1243  600);
1244  c->Divide(2, 1);
1245  c->cd(1);
1246  DrawH1andFitGauss((TH1D*) fHM->H1("fhNofHits" + t)->Clone());
1247  cout << "Number of hits/points = " << fHM->H1("fhNofHits" + t)->GetMean()
1248  << endl;
1249  //gPad->SetLogy(true);
1250  c->cd(2);
1251  DrawH2(fHM->H2("fhXcYcCircle" + t));
1252  }
1253 
1254  DrawH2MeanRms((TH2D*) fHM->H2("fhChi2CircleVsMom" + t),
1255  "richgeo" + t + "_chi2_circle_vs_mom");
1256  DrawH2MeanRms((TH2D*) fHM->H2("fhRadiusVsMom" + t),
1257  "richgeo" + t + "_r_vs_mom");
1258  DrawH2MeanRms((TH2D*) fHM->H2("fhDRVsMom" + t),
1259  "richgeo" + t + "_dr_vs_mom");
1260  fHM->H2("fhDRVsMom" + t)->GetYaxis()->SetRangeUser(-1.05, 1.05);
1261  } // for loop
1262 
1263  {
1264  fHM->CreateCanvas(
1265  "richgeo_nof_photons_per_hit", "richgeo_nof_photons_per_hit", 800, 800);
1266  fHM->H1("fhNofPhotonsPerHit")
1267  ->Scale(1. / fHM->H1("fhNofPhotonsPerHit")->Integral());
1268  DrawH1(fHM->H1("fhNofPhotonsPerHit"), kLinear, kLinear, "hist");
1269  }
1270 
1271  {
1272  TCanvas* c = fHM->CreateCanvas(
1273  "richgeo_diff_ellipse", "richgeo_diff_ellipse", 1200, 600);
1274  c->Divide(4, 2);
1275  c->cd(1);
1276  DrawH2(fHM->H2("fhDiffAaxis"));
1277  c->cd(2);
1278  DrawH2(fHM->H2("fhDiffBaxis"));
1279  c->cd(3);
1280  DrawH2(fHM->H2("fhDiffXcEllipse"));
1281  c->cd(4);
1282  DrawH2(fHM->H2("fhDiffYcEllipse"));
1283  c->cd(5);
1284  DrawH1(fHM->H2("fhDiffAaxis")->ProjectionY(), kLinear, kLog, "hist");
1285  c->cd(6);
1286  DrawH1(fHM->H2("fhDiffBaxis")->ProjectionY(), kLinear, kLog, "hist");
1287  c->cd(7);
1288  DrawH1(fHM->H2("fhDiffXcEllipse")->ProjectionY(), kLinear, kLog, "hist");
1289  c->cd(8);
1290  DrawH1(fHM->H2("fhDiffYcEllipse")->ProjectionY(), kLinear, kLog, "hist");
1291  }
1292 
1293  {
1294  TCanvas* c =
1295  fHM->CreateCanvas("richgeo_diff_circle", "richgeo_diff_circle", 900, 600);
1296  c->Divide(3, 2);
1297  c->cd(1);
1298  DrawH2(fHM->H2("fhDiffXcCircle"));
1299  c->cd(2);
1300  DrawH2(fHM->H2("fhDiffYcCircle"));
1301  c->cd(3);
1302  DrawH2(fHM->H2("fhDiffRadius"));
1303  c->cd(4);
1304  DrawH1(fHM->H2("fhDiffXcCircle")->ProjectionY(), kLinear, kLog, "hist");
1305  c->cd(5);
1306  DrawH1(fHM->H2("fhDiffYcCircle")->ProjectionY(), kLinear, kLog, "hist");
1307  c->cd(6);
1308  DrawH1(fHM->H2("fhDiffRadius")->ProjectionY(), kLinear, kLog, "hist");
1309  }
1310 
1311  {
1312  TCanvas* c =
1313  fHM->CreateCanvas("richgeo_hits_residual", "richgeo_hits", 1200, 600);
1314  c->Divide(2, 1);
1315  c->cd(1);
1316  fHM->H1("fhDiffXhit")->Scale(1. / fHM->H1("fhDiffXhit")->Integral());
1317  DrawH1(fHM->H1("fhDiffXhit"), kLinear, kLinear, "hist");
1318  c->cd(2);
1319  fHM->H1("fhDiffYhit")->Scale(1. / fHM->H1("fhDiffYhit")->Integral());
1320  DrawH1(fHM->H1("fhDiffYhit"), kLinear, kLinear, "hist");
1321  }
1322 
1323  {
1324  TCanvas* c =
1325  fHM->CreateCanvas("richgeo_fit_eff", "richgeo_fit_eff", 1800, 600);
1326  c->Divide(3, 1);
1327  c->cd(1);
1328  DrawH1({(TH1D*) fHM->H1("fhNofHitsAll")->Clone(),
1329  (TH1D*) fHM->H1("fhNofHitsCircleFit")->Clone(),
1330  (TH1D*) fHM->H1("fhNofHitsEllipseFit")->Clone()},
1331  {"All", "Circle fit", "Ellipse fit"},
1332  kLinear,
1333  kLog,
1334  true,
1335  0.7,
1336  0.7,
1337  0.99,
1338  0.99,
1339  "hist");
1340  TH1D* fhNofHitsCircleFitEff =
1341  Cbm::DivideH1(fHM->H1("fhNofHitsCircleFit"), fHM->H1("fhNofHitsAll"));
1342  TH1D* fhNofHitsEllipseFitEff =
1343  Cbm::DivideH1(fHM->H1("fhNofHitsEllipseFit"), fHM->H1("fhNofHitsAll"));
1344  c->cd(2);
1345  DrawH1(fhNofHitsCircleFitEff);
1346  TLatex* circleFitEffTxt = new TLatex(
1347  15,
1348  0.5,
1349  CalcEfficiency(fHM->H1("fhNofHitsCircleFit"), fHM->H1("fhNofHitsAll"))
1350  .c_str());
1351  cout << "Circle fit efficiency:" << circleFitEffTxt << "%" << endl;
1352  circleFitEffTxt->Draw();
1353  c->cd(3);
1354  DrawH1(fhNofHitsEllipseFitEff);
1355  TLatex* ellipseFitEffTxt = new TLatex(
1356  15,
1357  0.5,
1358  CalcEfficiency(fHM->H1("fhNofHitsEllipseFit"), fHM->H1("fhNofHitsAll"))
1359  .c_str());
1360  cout << "Ellipse fit efficiency:" << ellipseFitEffTxt << "%" << endl;
1361  ellipseFitEffTxt->Draw();
1362  }
1363 
1364  {
1365  TCanvas* c =
1366  fHM->CreateCanvas("richgeo_acc_el", "richgeo_acc_el", 1500, 500);
1367  c->Divide(3, 1);
1368  c->cd(1);
1369  DrawH1({(TH1D*) fHM->H1("fhMcMomEl"), (TH1D*) fHM->H1("fhAccMomEl")},
1370  {"MC", "ACC"},
1371  kLinear,
1372  kLog,
1373  true,
1374  0.8,
1375  0.8,
1376  0.99,
1377  0.99,
1378  "hist");
1379  c->cd(2);
1380  DrawH2(fHM->H2("fhMcPtyEl"));
1381  c->cd(3);
1382  DrawH2(fHM->H2("fhAccPtyEl"));
1383  }
1384 
1385  TH1D* pxEff = Cbm::DivideH1((TH1D*) fHM->H1("fhAccMomEl")->Clone(),
1386  (TH1D*) fHM->H1("fhMcMomEl")->Clone(),
1387  "",
1388  100.,
1389  "Geometrical acceptance [%]");
1390  TH2D* pyzEff = Cbm::DivideH2((TH2D*) fHM->H1("fhAccPtyEl")->Clone(),
1391  (TH2D*) fHM->H1("fhMcPtyEl")->Clone(),
1392  "",
1393  100.,
1394  "Geometrical acceptance [%]");
1395  {
1396  fHM->CreateCanvas(
1397  "richgeo_acc_eff_el_mom", "richgeo_acc_eff_el_mom", 800, 800);
1398  string effEl = CalcEfficiency((TH1D*) fHM->H1("fhAccMomEl")->Clone(),
1399  (TH1D*) fHM->H1("fhMcMomEl")->Clone());
1400  cout << "Geometrical acceptance electrons:" << effEl << "%" << endl;
1401  DrawH1({pxEff},
1402  {"e^{#pm} (" + effEl + "%)"},
1403  kLinear,
1404  kLinear,
1405  true,
1406  0.6,
1407  0.55,
1408  0.88,
1409  0.65);
1410  }
1411 
1412  {
1413  fHM->CreateCanvas(
1414  "richgeo_acc_eff_el_pty", "richgeo_acc_eff_el_pty", 800, 800);
1415  DrawH2(pyzEff);
1416  }
1417 
1418  {
1419  TCanvas* c =
1420  fHM->CreateCanvas("richgeo_acc_pi", "richgeo_acc_pi", 1500, 500);
1421  c->Divide(3, 1);
1422  c->cd(1);
1423  DrawH1({(TH1D*) fHM->H1("fhMcMomPi"), (TH1D*) fHM->H1("fhAccMomPi")},
1424  {"MC", "ACC"},
1425  kLinear,
1426  kLog,
1427  true,
1428  0.8,
1429  0.8,
1430  0.99,
1431  0.99,
1432  "hist");
1433  c->cd(2);
1434  DrawH2(fHM->H2("fhMcPtyPi"));
1435  c->cd(3);
1436  DrawH2(fHM->H2("fhAccPtyPi"));
1437  }
1438 
1439  TH1D* pxPiEff = Cbm::DivideH1((TH1D*) fHM->H1("fhAccMomPi")->Clone(),
1440  (TH1D*) fHM->H1("fhMcMomPi")->Clone(),
1441  "",
1442  100.,
1443  "Geometrical acceptance [%]");
1444  TH2D* pyzPiEff = Cbm::DivideH2((TH2D*) fHM->H1("fhAccPtyPi")->Clone(),
1445  (TH2D*) fHM->H1("fhMcPtyPi")->Clone(),
1446  "",
1447  100.,
1448  "Geometrical acceptance [%]");
1449  {
1450  fHM->CreateCanvas(
1451  "richgeo_acc_eff_pi_mom", "richgeo_acc_eff_pi_mom", 800, 800);
1452  string effPi = CalcEfficiency((TH1D*) fHM->H1("fhAccMomPi")->Clone(),
1453  (TH1D*) fHM->H1("fhMcMomPi")->Clone());
1454  cout << "Geometrical acceptance pions:" << effPi << "%" << endl;
1455  DrawH1({pxPiEff},
1456  {"#pi^{#pm} (" + effPi + "%)"},
1457  kLinear,
1458  kLinear,
1459  true,
1460  0.6,
1461  0.55,
1462  0.88,
1463  0.65);
1464  }
1465 
1466  {
1467  fHM->CreateCanvas(
1468  "richgeo_acc_eff_pi_pty", "richgeo_acc_eff_pi_pty", 800, 800);
1469  DrawH2(pyzPiEff);
1470  pyzPiEff->GetZaxis()->SetRangeUser(0, 100);
1471  }
1472 
1473  {
1474  TCanvas* c = fHM->CreateCanvas(
1475  "richgeo_acc_eff_el_zoom", "richgeo_acc_eff_el_zoom", 1000, 500);
1476  c->Divide(2, 1);
1477  c->cd(1);
1478  TH1D* fhMcMomElClone = (TH1D*) fHM->H1("fhMcMomEl")->Clone();
1479  TH1D* fhAccMomElClone = (TH1D*) fHM->H1("fhAccMomEl")->Clone();
1480  fhMcMomElClone->GetXaxis()->SetRangeUser(0., 3.);
1481  fhAccMomElClone->GetXaxis()->SetRangeUser(0., 3.);
1482  fhMcMomElClone->SetMinimum(0.);
1483  DrawH1({fhMcMomElClone, fhAccMomElClone},
1484  {"MC", "ACC"},
1485  kLinear,
1486  kLog,
1487  true,
1488  0.8,
1489  0.8,
1490  0.99,
1491  0.99);
1492  gPad->SetLogy(false);
1493  c->cd(2);
1494  TH1D* px_eff_clone = (TH1D*) pxEff->Clone();
1495  px_eff_clone->GetXaxis()->SetRangeUser(0., 3.);
1496  px_eff_clone->SetMinimum(0.);
1497  DrawH1(px_eff_clone);
1498  }
1499 
1500  // Draw number vs position onto the photodetector plane
1501  {
1502  TCanvas* c = fHM->CreateCanvas(
1503  "richgeo_numbers_vs_xy_hits", "richgeo_numbers_vs_xy_hits", 1800, 600);
1504  c->Divide(3, 1);
1505  c->cd(1);
1506  DrawH3Profile(fHM->H3("fhNofHitsXYZ"), true, false, 10, 30);
1507  c->cd(2);
1508  DrawH2WithProfile(fHM->H2("fhNofHitsVsX"), false, true);
1509  c->cd(3);
1510  DrawH2WithProfile(fHM->H2("fhNofHitsVsY"), false, true);
1511  }
1512 
1513  {
1514  TCanvas* c = fHM->CreateCanvas("richgeo_numbers_vs_xy_points",
1515  "richgeo_numbers_vs_xy_points",
1516  1800,
1517  600);
1518  c->Divide(3, 1);
1519  c->cd(1);
1520  DrawH3Profile(fHM->H3("fhNofPointsXYZ"), true, false, 100., 300.);
1521  c->cd(2);
1522  DrawH2WithProfile(fHM->H2("fhNofPointsVsX"), false, true);
1523  c->cd(3);
1524  DrawH2WithProfile(fHM->H2("fhNofPointsVsY"), false, true);
1525  }
1526 
1527  {
1528  TCanvas* c = fHM->CreateCanvas(
1529  "richgeo_numbers_vs_xy_boa", "richgeo_numbers_vs_xy_boa", 1800, 600);
1530  c->Divide(3, 1);
1531  c->cd(1);
1532  DrawH3Profile(fHM->H3("fhBoverAXYZ"), true, false, 0.75, 1.0);
1533  c->cd(2);
1534  DrawH2WithProfile(fHM->H2("fhBoverAVsX"), false, true);
1535  fHM->H2("fhBoverAVsX")->GetYaxis()->SetRangeUser(0.75, 1.0);
1536  c->cd(3);
1537  DrawH2WithProfile(fHM->H2("fhBoverAVsY"), false, true);
1538  fHM->H2("fhBoverAVsY")->GetYaxis()->SetRangeUser(0.75, 1.0);
1539  }
1540 
1541  {
1542  TCanvas* c = fHM->CreateCanvas(
1543  "richgeo_numbers_vs_xy_b", "richgeo_numbers_vs_xy_b", 1800, 600);
1544  c->Divide(3, 1);
1545  c->cd(1);
1546  DrawH3Profile(fHM->H3("fhBaxisXYZ"), true, false, 4., 5.);
1547  c->cd(2);
1548  DrawH2WithProfile(fHM->H2("fhBaxisVsX"), false, true);
1549  c->cd(3);
1550  DrawH2WithProfile(fHM->H2("fhBaxisVsY"), false, true);
1551  }
1552 
1553  {
1554  TCanvas* c = fHM->CreateCanvas(
1555  "richgeo_numbers_vs_xy_a", "richgeo_numbers_vs_xy_a", 1800, 600);
1556  c->Divide(3, 1);
1557  c->cd(1);
1558  DrawH3Profile(fHM->H3("fhAaxisXYZ"), true, false, 4.4, 5.7);
1559  c->cd(2);
1560  DrawH2WithProfile(fHM->H2("fhAaxisVsX"), false, true);
1561  c->cd(3);
1562  DrawH2WithProfile(fHM->H2("fhAaxisVsY"), false, true);
1563  }
1564 
1565  {
1566  TCanvas* c = fHM->CreateCanvas(
1567  "richgeo_numbers_vs_xy_r", "richgeo_numbers_vs_xy_r", 1800, 600);
1568  c->Divide(3, 1);
1569  c->cd(1);
1570  DrawH3Profile(fHM->H3("fhRadiusXYZ"), true, false, 4.2, 5.2);
1571  c->cd(2);
1572  DrawH2WithProfile(fHM->H2("fhRadiusVsX"), false, true);
1573  c->cd(3);
1574  DrawH2WithProfile(fHM->H2("fhRadiusVsY"), false, true);
1575  }
1576 
1577  {
1578  TCanvas* c = fHM->CreateCanvas(
1579  "richgeo_numbers_vs_xy_dr", "richgeo_numbers_vs_xy_dr", 1800, 600);
1580  c->Divide(3, 1);
1581  c->cd(1);
1582  DrawH3Profile(fHM->H3("fhdRXYZ"), false, false, 0., .5);
1583  c->cd(2);
1584  DrawH2WithProfile(fHM->H2("fhdRVsX"), false, false);
1585  c->cd(3);
1586  DrawH2WithProfile(fHM->H2("fhdRVsY"), false, false);
1587  }
1588 
1589  {
1590  fHM->CreateCanvas(
1591  "richgeo_acc_vs_min_nof_hits", "richgeo_acc_vs_min_nof_hits", 600, 600);
1592  TH1D* h = CreateAccVsMinNofHitsHist();
1593  h->GetXaxis()->SetRangeUser(0., 40.0);
1594  DrawH1(h);
1595  }
1596 
1597  DrawH2MeanRms(fHM->H2("fhRadiusVsNofHits"), "richgeo_hits_r_vs_nof_hits");
1598  DrawH2MeanRms(fHM->H2("fhAaxisVsNofHits"), "richgeo_hits_a_vs_nof_hits");
1599  DrawH2MeanRms(fHM->H2("fhBaxisVsNofHits"), "richgeo_hits_b_vs_nof_hits");
1600  DrawH2MeanRms(fHM->H2("fhDRVsNofHits"), "richgeo_hits_dr_vs_nof_hits");
1601  fHM->H2("fhDRVsNofHits")->GetYaxis()->SetRangeUser(-1.05, 1.05);
1602 
1603  {
1604  TCanvas* c =
1605  fHM->CreateCanvas("richgeo_hits_rab", "richgeo_hits_rab", 1500, 600);
1606  c->Divide(3, 1);
1607  c->cd(1);
1609  fHM->H2("fhRadiusVsNofHits")
1610  ->ProjectionY(
1611  (string(fHM->H2("fhRadiusVsNofHits")->GetName()) + "_py").c_str()),
1612  true,
1613  true,
1614  2.,
1615  8.);
1616  c->cd(2);
1618  fHM->H2("fhAaxisVsNofHits")
1619  ->ProjectionY(
1620  (string(fHM->H2("fhAaxisVsNofHits")->GetName()) + "_py").c_str()),
1621  true,
1622  true,
1623  2.,
1624  8.);
1625  c->cd(3);
1627  fHM->H2("fhBaxisVsNofHits")
1628  ->ProjectionY(
1629  (string(fHM->H2("fhBaxisVsNofHits")->GetName()) + "_py").c_str()),
1630  true,
1631  true,
1632  2.,
1633  8.);
1634  }
1635 
1636 
1637  {
1638  TCanvas* c = fHM->CreateCanvas(
1639  "richgeo_photon_energy", "richgeo_photon_energy", 1500, 750);
1640  c->Divide(2, 1);
1641  c->cd(1);
1642  DrawH1({fHM->H1("fhPhotonEPlaneZ+"),
1643  fHM->H1("fhPhotonEPlaneZ-"),
1644  fHM->H1("fhPhotonEPmtPoint"),
1645  fHM->H1("fhPhotonEPmtHit")},
1646  {"Sens plane Z+", "Sens plane Z-", "PMT point", "PMT hit"});
1647  c->cd(2);
1648  DrawH1({fHM->H1("fhLambdaPlaneZ+"),
1649  fHM->H1("fhLambdaPlaneZ-"),
1650  fHM->H1("fhLambdaPmtPoint"),
1651  fHM->H1("fhLambdaPmtHit")},
1652  {"Sens plane Z+", "Sens plane Z-", "PMT Point", "PMT hit"});
1653  }
1654 }
1655 
1656 
1658  fHM->Create3<TH3D>("fhPointsXYZ",
1659  "fhPointsXYZ;X [cm];Y [cm];Z [cm];Yield",
1660  100,
1661  -50,
1662  50,
1663  100,
1664  -300,
1665  300,
1666  100,
1667  100,
1668  300);
1669  fHM->Create3<TH3D>("fhHitsXYZ",
1670  "fhHitsXYZ;X [cm];Y [cm];Z [cm];Yield",
1671  100,
1672  -50,
1673  50,
1674  100,
1675  -300,
1676  300,
1677  100,
1678  100,
1679  300);
1680  vector<Int_t> pixels =
1682  vector<Int_t> pmts = CbmRichDigiMapManager::GetInstance().GetPmtIds();
1683 
1684  {
1685  TCanvas* c =
1686  fHM->CreateCanvas("richgeo_pixels_xy", "richgeo_pixels_xy", 1500, 1500);
1687  c->SetGrid(true, true);
1688  TH2D* pad = new TH2D(
1689  "richgeo_pixels_xy", ";X [cm];Y [cm]", 1, -120, 120, 1, -210, 210);
1690  //TH2D* pad = new TH2D("richgeo_pixels_xy", ";X [cm];Y [cm]", 1, -20, 20, 1, 140, 180);
1691  pad->SetStats(false);
1692  pad->Draw();
1693  DrawPmtPoint("xy", pixels, true);
1694  }
1695 
1696  {
1697  TCanvas* c =
1698  fHM->CreateCanvas("richgeo_pixels_xz", "richgeo_pixels_xz", 1500, 1500);
1699  c->SetGrid(true, true);
1700  TH2D* pad = new TH2D(
1701  "richgeo_pixels_xz", ";Z [cm];X [cm]", 1, 200, 250, 1, -120., 120.);
1702  pad->SetStats(false);
1703  pad->Draw();
1704  DrawPmtPoint("zx", pixels, true);
1705  }
1706 
1707  {
1708  fHM->CreateCanvas("richgeo_pixels_yz", "richgeo_pixels_yz", 1500, 1500);
1709  TH2D* pad = new TH2D(
1710  "richgeo_pixels_yz", ";Z [cm];Y [cm]", 1, 200, 250, 1, -220, 220);
1711  pad->SetStats(false);
1712  pad->Draw();
1713  DrawPmtPoint("zy", pixels, true);
1714  }
1715 
1716  {
1717  TCanvas* c =
1718  fHM->CreateCanvas("richgeo_pmts_xy", "richgeo_pmts_xy", 1500, 1500);
1719  c->SetGrid(true, true);
1720  TH2D* pad =
1721  new TH2D("richgeo_pmts_xy", ";X [cm];Y [cm]", 1, -120, 120, 1, -210, 210);
1722  pad->SetStats(false);
1723  pad->Draw();
1724  DrawPmtPoint("xy", pmts, false);
1725  }
1726 
1727 
1728  for (unsigned int i = 0; i < pixels.size(); i++) {
1729  CbmRichPixelData* pixelData =
1731  TVector3 inPos(pixelData->fX, pixelData->fY, pixelData->fZ);
1732  TVector3 outPos;
1733  CbmRichGeoManager::GetInstance().RotatePoint(&inPos, &outPos);
1734 
1735  fHM->H3("fhPointsXYZ")->Fill(inPos.X(), inPos.Y(), inPos.Z());
1736  fHM->H3("fhHitsXYZ")->Fill(outPos.X(), outPos.Y(), outPos.Z());
1737  }
1738 
1739  {
1740  fHM->CreateCanvas(
1741  "richgeo_pixels_points_xyz", "richgeo_pixels_points_xyz", 1500, 1500);
1742  fHM->H3("fhPointsXYZ")->Draw();
1743  }
1744 
1745  {
1746  fHM->CreateCanvas(
1747  "richgeo_pixels_hits_xyz", "richgeo_pixels_hits_xyz", 1500, 1500);
1748  fHM->H3("fhHitsXYZ")->Draw();
1749  }
1750 }
1751 
1752 void CbmRichGeoTest::DrawPmtPoint(const string& coordinates,
1753  const vector<Int_t>& ids,
1754  Bool_t isDrawPixel) {
1755  for (unsigned int i = 0; i < ids.size(); i++) {
1756  TVector3 inPos;
1757  Double_t boxHalfSize = 0.0;
1758  if (isDrawPixel) {
1759  CbmRichPixelData* pixelData =
1761  inPos.SetXYZ(pixelData->fX, pixelData->fY, pixelData->fZ);
1762  boxHalfSize = 0.15;
1763  } else {
1764  CbmRichPmtData* pmtData =
1766  inPos.SetXYZ(pmtData->fX, pmtData->fY, pmtData->fZ);
1767  boxHalfSize = 0.5 * pmtData->fHeight;
1768  }
1769  TVector3 outPos;
1770  CbmRichGeoManager::GetInstance().RotatePoint(&inPos, &outPos);
1771  TBox* boxOut = nullptr;
1772 
1773  if (coordinates == string("xy"))
1774  boxOut = new TBox(outPos.X() - boxHalfSize,
1775  outPos.Y() - boxHalfSize,
1776  outPos.X() + boxHalfSize,
1777  outPos.Y() + boxHalfSize);
1778  if (coordinates == string("zx"))
1779  boxOut = new TBox(outPos.Z() - boxHalfSize,
1780  outPos.X() - boxHalfSize,
1781  outPos.Z() + boxHalfSize,
1782  outPos.X() + boxHalfSize);
1783  if (coordinates == string("zy"))
1784  boxOut = new TBox(outPos.Z() - boxHalfSize,
1785  outPos.Y() - boxHalfSize,
1786  outPos.Z() + boxHalfSize,
1787  outPos.Y() + boxHalfSize);
1788  if (boxOut != NULL) {
1789  if (isDrawPixel) {
1790  boxOut->SetFillColor(kBlue);
1791  } else {
1792  boxOut->SetFillStyle(0);
1793  boxOut->SetLineColor(kBlue);
1794  }
1795  boxOut->Draw();
1796  }
1797 
1798  TBox* boxIn = nullptr;
1799  if (coordinates == string("xy"))
1800  boxIn = new TBox(inPos.X() - boxHalfSize,
1801  inPos.Y() - boxHalfSize,
1802  inPos.X() + boxHalfSize,
1803  inPos.Y() + boxHalfSize);
1804  if (coordinates == string("zx"))
1805  boxIn = new TBox(inPos.Z() - boxHalfSize,
1806  inPos.X() - boxHalfSize,
1807  inPos.Z() + boxHalfSize,
1808  inPos.X() + boxHalfSize);
1809  if (coordinates == string("zy"))
1810  boxIn = new TBox(inPos.Z() - boxHalfSize,
1811  inPos.Y() - boxHalfSize,
1812  inPos.Z() + boxHalfSize,
1813  inPos.Y() + boxHalfSize);
1814  if (boxIn != NULL) {
1815  if (isDrawPixel) {
1816  boxIn->SetFillColor(kRed);
1817  } else {
1818  boxIn->SetFillStyle(0);
1819  boxIn->SetLineColor(kRed);
1820  }
1821  boxIn->Draw();
1822  }
1823  }
1824 }
1825 
1826 
1828  DrawHist();
1829  if (fDrawPmts) { DrawPmts(); }
1830 
1831  TDirectory* oldir = gDirectory;
1832  TFile* outFile = FairRootManager::Instance()->GetOutFile();
1833  if (outFile != NULL) {
1834  outFile->cd();
1835  fHM->WriteToFile();
1836  }
1837  gDirectory->cd(oldir->GetPath());
1838  fHM->SaveCanvasToImage(fOutputDir, "png");
1839 }
1840 
1841 string CbmRichGeoTest::CalcEfficiency(TH1* histRec, TH1* histAcc) {
1842  if (histAcc->GetEntries() == 0) {
1843  return "0";
1844  } else {
1845  Double_t eff =
1846  100. * Double_t(histRec->GetEntries()) / Double_t(histAcc->GetEntries());
1847  return Cbm::NumberToString(eff, 2);
1848  }
1849 }
1850 
1851 void CbmRichGeoTest::DrawFromFile(const string& fileName,
1852  const string& outputDir) {
1853  fOutputDir = outputDir;
1854 
1855  if (fHM != NULL) delete fHM;
1856 
1857  fHM = new CbmHistManager();
1858  TFile* file = new TFile(fileName.c_str());
1859  fHM->ReadFromFile(file);
1860 
1861  DrawHist();
1862 
1863  fHM->SaveCanvasToImage(fOutputDir, "png,eps");
1864 }
1865 
CbmRichPoint.h
CbmHit::GetZ
Double_t GetZ() const
Definition: CbmHit.h:70
CbmRichGeoTest::fMaxAaxis
Double_t fMaxAaxis
Definition: CbmRichGeoTest.h:218
CbmMCTrack::GetMotherId
Int_t GetMotherId() const
Definition: CbmMCTrack.h:71
CbmRichPixelData
Definition: CbmRichDetectorData.h:17
CbmMatch::GetMatchedLink
const CbmLink & GetMatchedLink() const
Definition: CbmMatch.h:37
CbmMCDataManager::GetObject
CbmMCDataObject * GetObject(const char *name)
Definition: CbmMCDataManager.cxx:137
CbmMCEventList::GetFileIdByIndex
Int_t GetFileIdByIndex(UInt_t index)
File number by index @value File number for event at given index in list.
Definition: CbmMCEventList.cxx:106
CbmMatch
Definition: CbmMatch.h:22
CbmRichPmtData::fY
Double_t fY
Definition: CbmRichDetectorData.h:42
h
Generates beam ions for transport simulation.
Definition: CbmBeamGenerator.h:17
CbmMCDataManager.h
CbmRichDigiMapManager::GetPmtIds
std::vector< Int_t > GetPmtIds()
Definition: CbmRichDigiMapManager.cxx:285
CbmRichDigiMapManager::GetPixelAddresses
std::vector< Int_t > GetPixelAddresses()
Definition: CbmRichDigiMapManager.cxx:280
CbmHistManager::WriteToFile
void WriteToFile()
Write all histograms to current opened file.
Definition: core/base/CbmHistManager.cxx:103
CbmPixelHit::GetX
Double_t GetX() const
Definition: CbmPixelHit.h:83
CbmRichGeoTest::fHM
CbmHistManager * fHM
Definition: CbmRichGeoTest.h:210
CbmRichGeoTest::DrawPmts
void DrawPmts()
DrawPmts.
Definition: CbmRichGeoTest.cxx:1657
CbmRichGeoManager::GetInstance
static CbmRichGeoManager & GetInstance()
Definition: CbmRichGeoManager.h:29
CbmRichRingFitterEllipseTau
Here the ring is fitted with Taubin algorithm from A. Ayriyan, G. Ososkov, N. Chernov.
Definition: CbmRichRingFitterEllipseTau.h:35
CbmRichDigi
Definition: CbmRichDigi.h:25
CbmPixelHit::GetY
Double_t GetY() const
Definition: CbmPixelHit.h:84
CbmHistManager::Create3
void Create3(const std::string &name, const std::string &title, Int_t nofBinsX, Double_t minBinX, Double_t maxBinX, Int_t nofBinsY, Double_t minBinY, Double_t maxBinY, Int_t nofBinsZ, Double_t minBinZ, Double_t maxBinZ)
Helper function for creation of 3-dimensional histograms and profiles. Template argument is a real ob...
Definition: CbmHistManager.h:140
CbmRichRingFitterEllipseTau.h
Here the ring is fitted with Taubin algorithm from A. Ayriyan, G. Ososkov, N. Chernov.
sqrt
friend F32vec4 sqrt(const F32vec4 &a)
Definition: L1/vectors/P4_F32vec4.h:41
CbmRichGeoTest::Init
virtual InitStatus Init()
Inherited from FairTask.
Definition: CbmRichGeoTest.cxx:100
CbmMCDataArray::Size
Int_t Size(Int_t fileNumber, Int_t eventNumber)
Definition: CbmMCDataArray.cxx:133
CbmHistManager::ReadFromFile
void ReadFromFile(TFile *file)
Read histograms from file.
Definition: core/base/CbmHistManager.cxx:110
DrawH3Profile
TH2D * DrawH3Profile(TH3 *h, Bool_t drawMean, Bool_t doGaussFit, Double_t zMin, Double_t zMax)
Definition: CbmDrawHist.cxx:369
CbmDigiManager::Init
InitStatus Init()
Initialisation.
Definition: CbmDigiManager.cxx:71
CbmRichRingLight::GetChi2
float GetChi2() const
Definition: CbmRichRingLight.h:242
CbmHistManager::Create2
void Create2(const std::string &name, const std::string &title, Int_t nofBinsX, Double_t minBinX, Double_t maxBinX, Int_t nofBinsY, Double_t minBinY, Double_t maxBinY)
Helper function for creation of 2-dimensional histograms and profiles. Template argument is a real ob...
Definition: CbmHistManager.h:104
CbmRichGeoTest::fMcTracks
CbmMCDataArray * fMcTracks
Definition: CbmRichGeoTest.h:202
CbmRichGeoTest::Finish
virtual void Finish()
Inherited from FairTask.
Definition: CbmRichGeoTest.cxx:1827
CbmRichGeoTest::FillMcVsHitFitCircle
void FillMcVsHitFitCircle(CbmRichRingLight *ring, CbmRichRingLight *ringMc)
Calculate difference between circle parameters for two fittings using hits and MC points for fit and ...
Definition: CbmRichGeoTest.cxx:1011
CbmMCTrack::GetPdgCode
Int_t GetPdgCode() const
Definition: CbmMCTrack.h:70
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
CbmRichDigiMapManager::GetPmtDataById
CbmRichPmtData * GetPmtDataById(Int_t id)
Definition: CbmRichDigiMapManager.cxx:288
CbmRichDigiMapManager.h
CbmRichGeoTest::FitAndFillHistEllipse
void FitAndFillHistEllipse(Int_t histIndex, CbmRichRingLight *ring, Double_t momentum)
Fit ring using ellipse fitter and fill histograms.
Definition: CbmRichGeoTest.cxx:924
CbmRichGeoTest::fDrawEventDisplay
Bool_t fDrawEventDisplay
Definition: CbmRichGeoTest.h:226
CbmMCDataManager::InitBranch
CbmMCDataArray * InitBranch(const char *name)
Definition: CbmMCDataManager.cxx:106
Cbm::DivideH2
TH2D * DivideH2(TH2 *h1, TH2 *h2, const string &histName, double scale, const string &titleZaxis)
Definition: CbmUtils.cxx:92
CbmRichPmtData::fZ
Double_t fZ
Definition: CbmRichDetectorData.h:43
CbmHistManager::H2
TH2 * H2(const std::string &name) const
Return pointer to TH2 histogram.
Definition: CbmHistManager.h:190
CbmMCDataArray.h
CbmRichGeoTest::fMaxBaxis
Double_t fMaxBaxis
Definition: CbmRichGeoTest.h:220
CbmRichRing::GetNofHits
Int_t GetNofHits() const
Definition: CbmRichRing.h:40
CbmRichGeoTest
RICH geometry checking and testing.
Definition: CbmRichGeoTest.h:46
CbmRichRing
Definition: CbmRichRing.h:17
CbmHit::GetRefId
Int_t GetRefId() const
Definition: CbmHit.h:72
CbmRichPmtData
Definition: CbmRichDetectorData.h:26
CbmRichHitLight::fX
float fX
Definition: CbmRichRingLight.h:34
CbmRichGeoTest::ProcessHits
void ProcessHits()
Calculate residuals between hits and MC points and fill histograms.
Definition: CbmRichGeoTest.cxx:1021
CbmRichRing.h
CbmDrawHist.h
Helper functions for drawing 1D and 2D histograms and graphs.
CbmRichDraw.h
CbmRichGeoTest::fRichPoints
CbmMCDataArray * fRichPoints
Definition: CbmRichGeoTest.h:201
CbmRichGeoTest::fMinAaxis
Double_t fMinAaxis
Definition: CbmRichGeoTest.h:217
CbmHistManager.h
Histogram manager.
CbmDigiManager::IsPresent
static Bool_t IsPresent(ECbmModuleId systemId)
Presence of a digi branch.
Definition: CbmDigiManager.cxx:112
CbmRichGeoTest::FitAndFillHistCircle
void FitAndFillHistCircle(Int_t histIndex, CbmRichRingLight *ring, Double_t momentum)
Fit ring using circle fitter and fill histograms.
Definition: CbmRichGeoTest.cxx:962
CbmHistManager::H3
TH3 * H3(const std::string &name) const
Return pointer to TH3 histogram.
Definition: CbmHistManager.h:210
CbmDigiManager::Instance
static CbmDigiManager * Instance()
Static instance.
Definition: CbmDigiManager.h:93
CbmDigiManager::IsMatchPresent
static Bool_t IsMatchPresent(ECbmModuleId systemId)
Presence of a digi match branch.
Definition: CbmDigiManager.cxx:104
CbmRichRingFitterCOP
Here the ring is fitted with the COP algorithm from A. Ayriyan/G. Ososkov.
Definition: CbmRichRingFitterCOP.h:28
CbmRichGeoTest::fDigiMan
CbmDigiManager * fDigiMan
Definition: CbmRichGeoTest.h:200
CbmRichGeoManager.h
Cbm::DivideH1
TH1D * DivideH1(TH1 *h1, TH1 *h2, const string &histName, double scale, const string &titleYaxis)
Definition: CbmUtils.cxx:70
CbmRichDigi.h
CbmRichGeoTest::InitHistograms
void InitHistograms()
Initialize histograms.
Definition: CbmRichGeoTest.cxx:168
CbmRichGeoTest::fMaxRadius
Double_t fMaxRadius
Definition: CbmRichGeoTest.h:222
DrawH1
void DrawH1(TH1 *hist, HistScale logx, HistScale logy, const string &drawOpt, Int_t color, Int_t lineWidth, Int_t lineStyle, Int_t markerSize, Int_t markerStyle)
Definition: CbmDrawHist.cxx:49
CbmRichPmt::getLambda
static Double_t getLambda(Double_t momentum)
Definition: CbmRichPmt.cxx:36
CbmRichGeoTest::DrawHist
void DrawHist()
Draw histograms.
Definition: CbmRichGeoTest.cxx:1162
CbmRichRingLight::GetNofHits
int GetNofHits() const
Return number of hits in ring.
Definition: CbmRichRingLight.h:108
CbmMCEventList::GetNofEvents
std::size_t GetNofEvents() const
Number of events in the list @value Number of events.
Definition: CbmMCEventList.h:90
CbmRichDraw::DrawPmtH2
static void DrawPmtH2(TH2 *h, TCanvas *c, Bool_t usePmtBins=false)
Definition: CbmRichDraw.h:21
CbmRichGeoTest::fTauFit
CbmRichRingFitterEllipseTau * fTauFit
Definition: CbmRichGeoTest.h:208
CbmHistManager
Histogram manager.
Definition: CbmHistManager.h:41
CbmRichRingLight::GetAaxis
float GetAaxis() const
Definition: CbmRichRingLight.h:163
CbmRichGeoTest::fRichRings
TClonesArray * fRichRings
Definition: CbmRichGeoTest.h:198
CbmRichGeoTest::DrawFromFile
void DrawFromFile(const string &fileName, const string &outputDir)
Draw histogram from file.
Definition: CbmRichGeoTest.cxx:1851
CbmRichGeoTest::fRichHits
TClonesArray * fRichHits
Definition: CbmRichGeoTest.h:197
CbmRichGeoTest::CreateAccVsMinNofHitsHist
TH1D * CreateAccVsMinNofHitsHist()
Create histogram: RICH detector acceptance vs. minimum required number of hits in ring.
Definition: CbmRichGeoTest.cxx:1133
CbmDigiManager::Get
const Digi * Get(Int_t index) const
Get a digi object.
Definition: CbmDigiManager.h:52
CbmRichConverter.h
Convert internal data classes to cbmroot common data classes.
CbmRichGeoTest::fCopFit
CbmRichRingFitterCOP * fCopFit
Definition: CbmRichGeoTest.h:207
CbmRichPixelData::fY
Double_t fY
Definition: CbmRichDetectorData.h:21
CbmRichGeoTest::fMinRadius
Double_t fMinRadius
Definition: CbmRichGeoTest.h:221
CbmRichRingLight::GetCenterY
float GetCenterY() const
Definition: CbmRichRingLight.h:160
CbmMCDataArray::Get
TObject * Get(const CbmLink *lnk)
Definition: CbmMCDataArray.h:47
kLinear
@ kLinear
Definition: CbmDrawHist.h:79
CbmRichGeoTest::ProcessMc
void ProcessMc()
Fill MC histogram for detector acceptance calculation.
Definition: CbmRichGeoTest.cxx:732
CbmTrackMatchNew.h
CbmRichRingLight::AddHit
void AddHit(CbmRichHitLight hit)
Add new hit to the ring.
Definition: CbmRichRingLight.h:87
CbmDigiManager::GetMatch
const CbmMatch * GetMatch(ECbmModuleId systemId, UInt_t index) const
Get a match object.
Definition: CbmDigiManager.cxx:54
CbmRichGeoTest::fRichRingMatches
TClonesArray * fRichRingMatches
Definition: CbmRichGeoTest.h:203
CbmRichPixelData::fX
Double_t fX
Definition: CbmRichDetectorData.h:20
DrawH2WithProfile
void DrawH2WithProfile(TH2 *hist, Bool_t doGaussFit, Bool_t drawOnlyMean, const string &drawOpt2D, Int_t profileColor, Int_t profileLineWidth)
Definition: CbmDrawHist.cxx:296
CbmHistManager::Create1
void Create1(const std::string &name, const std::string &title, Int_t nofBins, Double_t minBin, Double_t maxBin)
Helper function for creation of 1-dimensional histograms and profiles. Template argument is a real ob...
Definition: CbmHistManager.h:81
CbmRichGeoManager::RotatePoint
void RotatePoint(TVector3 *inPos, TVector3 *outPos, Bool_t noTilting=false)
Definition: CbmRichGeoManager.cxx:407
CbmRichPmtData::fX
Double_t fX
Definition: CbmRichDetectorData.h:41
ECbmModuleId::kRich
@ kRich
Ring-Imaging Cherenkov Detector.
ClassImp
ClassImp(CbmConverterManager) InitStatus CbmConverterManager
Definition: CbmConverterManager.cxx:12
CbmMCTrack::GetGeantProcessId
UInt_t GetGeantProcessId() const
Definition: CbmMCTrack.h:69
CbmMCEventList
Container class for MC events with number, file and start time.
Definition: CbmMCEventList.h:38
CbmRichGeoTest::FillMcVsHitFitEllipse
void FillMcVsHitFitEllipse(CbmRichRingLight *ring, CbmRichRingLight *ringMc)
Calculate difference between ellipse parameters for two fitting using hits and MC points for fit and ...
Definition: CbmRichGeoTest.cxx:999
CbmHistManager::H1
TH1 * H1(const std::string &name) const
Return pointer to TH1 histogram.
Definition: CbmHistManager.h:170
CbmRichGeoTest::fEventList
CbmMCEventList * fEventList
Definition: CbmRichGeoTest.h:204
CbmRichRingLight::GetBaxis
float GetBaxis() const
Definition: CbmRichRingLight.h:164
CbmRichHitLight
Definition: CbmRichRingLight.h:14
CbmUtils.h
CbmRichRingFitterCOP.h
Here the ring is fitted with the COP algorithm from A. Ayriyan/G. Ososkov.
CbmRichGeoTest::~CbmRichGeoTest
virtual ~CbmRichGeoTest()
Standard destructor.
Definition: CbmRichGeoTest.cxx:98
CbmMCTrack::GetPt
Double_t GetPt() const
Definition: CbmMCTrack.h:99
CbmHistManager::CreateCanvas
TCanvas * CreateCanvas(const std::string &name, const std::string &title, Int_t width, Int_t height)
Create and draw TCanvas and store pointer to it.
Definition: core/base/CbmHistManager.cxx:267
CbmRichRingLight::GetRadius
float GetRadius() const
Definition: CbmRichRingLight.h:161
CbmMCDataManager
Task class creating and managing CbmMCDataArray objects.
Definition: CbmMCDataManager.h:27
CbmHistManager::SaveCanvasToImage
void SaveCanvasToImage(const std::string &outputDir, const std::string &options="png,eps")
Save all stored canvases to images.
Definition: core/base/CbmHistManager.cxx:276
CbmRichGeoTest::DrawH2MeanRms
void DrawH2MeanRms(TH2 *hist, const string &canvasName)
Definition: CbmRichGeoTest.cxx:1147
CbmMCEventList.h
CbmMCTrack::GetRapidity
Double_t GetRapidity() const
Definition: CbmMCTrack.cxx:177
CbmRichGeoTest::fRichRefPlanePoints
CbmMCDataArray * fRichRefPlanePoints
Definition: CbmRichGeoTest.h:199
CbmRichGeoTest::fMinBaxis
Double_t fMinBaxis
Definition: CbmRichGeoTest.h:219
CbmRichRingFitterCOP::DoFit
virtual void DoFit(CbmRichRingLight *ring)
Inherited from CbmRichRingFitterBase.
Definition: CbmRichRingFitterCOP.cxx:19
CbmMCTrack.h
CbmMCEventList::GetEventIdByIndex
Int_t GetEventIdByIndex(UInt_t index)
Event number by index @value Event number for event at given index in list.
Definition: CbmMCEventList.cxx:77
CbmRichHitLight::fY
float fY
Definition: CbmRichRingLight.h:35
CbmDigiManager.h
CbmRichGeoTest::RingParameters
void RingParameters()
Loop over all rings in array and fill ring parameters histograms.
Definition: CbmRichGeoTest.cxx:800
CbmMCTrack
Definition: CbmMCTrack.h:34
CbmRichGeoTest.h
RICH geometry checking and testing.
CbmRichConverter::CopyHitsToRingLight
static void CopyHitsToRingLight(const CbmRichRing *ring1, CbmRichRingLight *ring2)
Copy hits from CbmRichRing to CbmRichRingLight.
Definition: CbmRichConverter.h:41
CbmRichGeoTest::CbmRichGeoTest
CbmRichGeoTest()
Standard constructor.
Definition: CbmRichGeoTest.cxx:63
CbmRichRingLight::GetHit
CbmRichHitLight GetHit(int ind)
Return hit by the index.
Definition: CbmRichRingLight.h:114
CbmRichGeoTest::Exec
virtual void Exec(Option_t *option)
Inherited from FairTask.
Definition: CbmRichGeoTest.cxx:159
CbmRichGeoTest::fOutputDir
string fOutputDir
Definition: CbmRichGeoTest.h:195
CbmRichGeoTest::fMinNofHits
Int_t fMinNofHits
Definition: CbmRichGeoTest.h:214
CbmRichPmt.h
CbmTrackMatchNew
Definition: CbmTrackMatchNew.h:19
CbmMatch::GetLinks
const std::vector< CbmLink > & GetLinks() const
Definition: CbmMatch.h:36
DrawH2
void DrawH2(TH2 *hist, HistScale logx, HistScale logy, HistScale logz, const string &drawOpt)
Definition: CbmDrawHist.cxx:84
SetDefaultDrawStyle
void SetDefaultDrawStyle()
Definition: CbmDrawHist.cxx:33
CbmRichHit.h
CbmReport.h
CbmRichRingLight::GetCenterX
float GetCenterX() const
Definition: CbmRichRingLight.h:159
CbmRichGeoTest::DrawRing
void DrawRing(CbmRichRingLight *ringHit, CbmRichRingLight *ringPoint)
Draw ring in separate TCanvas.
Definition: CbmRichGeoTest.cxx:1066
CbmRichRingFitterEllipseTau::DoFit
virtual void DoFit(CbmRichRingLight *ring)
Inherited from CbmRichRingFitterBase.
Definition: CbmRichRingFitterEllipseTau.cxx:94
CbmRichGeoTest::CalcEfficiency
string CalcEfficiency(TH1 *histRec, TH1 *histAcc)
Calculate efficiency.
Definition: CbmRichGeoTest.cxx:1841
CbmRichPoint
Definition: CbmRichPoint.h:24
CbmRichDigiMapManager::GetInstance
static CbmRichDigiMapManager & GetInstance()
Definition: CbmRichDigiMapManager.h:29
CbmStudyReport.h
Base class for study reports.
CbmRichPmtData::fHeight
Double_t fHeight
Definition: CbmRichDetectorData.h:45
Cbm::NumberToString
std::string NumberToString(const T &value, int precision=1)
Definition: CbmUtils.h:23
CbmRichGeoTest::fNofDrawnRings
Int_t fNofDrawnRings
Definition: CbmRichGeoTest.h:224
CbmRichRingLight
Definition: CbmRichRingLight.h:39
CbmMCTrack::GetP
Double_t GetP() const
Definition: CbmMCTrack.h:100
kLog
@ kLog
Definition: CbmDrawHist.h:78
DrawH1andFitGauss
void DrawH1andFitGauss(TH1 *hist, Bool_t drawResults, Bool_t doScale, Double_t userRangeMin, Double_t userRangeMax)
Definition: CbmDrawHist.cxx:266
CbmRichGeoTest::DrawPmtPoint
void DrawPmtPoint(const string &coordinates, const vector< Int_t > &ids, Bool_t isDrawPixel)
Definition: CbmRichGeoTest.cxx:1752
CbmRichHit
Definition: CbmRichHit.h:19
CbmRichDigiMapManager::GetPixelDataByAddress
CbmRichPixelData * GetPixelDataByAddress(Int_t address)
Definition: CbmRichDigiMapManager.cxx:267
CbmRichGeoTest::fEventNum
Int_t fEventNum
Definition: CbmRichGeoTest.h:212
CbmRichPixelData::fZ
Double_t fZ
Definition: CbmRichDetectorData.h:22
CbmRichGeoTest::fDrawPmts
Bool_t fDrawPmts
Definition: CbmRichGeoTest.h:225