CbmRoot
CbmMvdClusterAna.cxx
Go to the documentation of this file.
1 // -------------------------------------------------------------------------
2 // ----- CbmMvdClusterAna source file -----
3 // ----- Created 27/04/15 by S. Amar-Youcef -----
4 // ------------------------------------------------------------------------
5 
6 //-- Include from Cbm --//
7 #include "CbmMvdClusterAna.h"
8 // #include "CbmStsTrack.h"
9 #include "CbmMvdCluster.h"
10 #include "CbmMvdDigi.h"
11 #include "CbmMvdDigiMatch.h"
12 #include "CbmMvdHit.h"
13 #include "CbmMvdHitMatch.h"
14 #include "CbmMvdPoint.h"
15 
16 // #include "CbmVertex.h"
17 #include "CbmLink.h"
18 #include "CbmMCTrack.h"
19 #include "CbmMatch.h"
20 #include "CbmTrackMatchNew.h"
21 
22 // #include "base/CbmLitToolFactory.h"
23 // #include "data/CbmLitTrackParam.h"
24 // #include "utils/CbmLitConverter.h"
25 
26 
27 //-- Include from Fair --//
28 #include "FairLogger.h"
29 #include "FairTrackParam.h"
30 
31 
32 //-- Include from Root --//
33 #include "TCanvas.h"
34 #include "TF1.h"
35 #include "TGeoBBox.h"
36 #include "TGeoManager.h"
37 #include "TH1.h"
38 #include "TH2.h"
39 #include "TLegend.h"
40 #include "TMath.h"
41 #include "TProfile.h"
42 #include "TROOT.h"
43 
44 //-- Include from C++ --//
45 #include <iomanip>
46 #include <iostream>
47 using namespace std;
48 using std::cout;
49 using std::endl;
50 using std::flush;
51 
52 // ----- Default constructor -------------------------------------------
54  : FairTask("MvdClusterAna")
55  , fMcPoints(NULL)
56  , fMvdDigis(NULL)
57  , fMvdClusters(NULL)
58  , fMvdHits(NULL)
59  , fMvdDigisMatch(NULL)
60  , fMvdClustersMatch(NULL)
61  , fMvdHitsMatch(NULL)
62  , fListMCTracks(NULL)
63  , fStsTrackArray(NULL)
64  , fStsTrackMatches(NULL)
65  , fMvdHisto1()
66  , fMvdHisto2()
67  , fProf()
68  , fNrMcPointsAll(-1)
69  , fNrHitsAll(-1)
70  , fMcperDigi()
71  , fMcperHit()
72  , fPixelpitch() {
73  ;
74 }
75 // -------------------------------------------------------------------------
76 
77 
78 // ----- Standard constructor ------------------------------------------
79 CbmMvdClusterAna::CbmMvdClusterAna(const char* name, Int_t iVerbose)
80  : FairTask(name, iVerbose)
81  , fMcPoints(NULL)
82  , fMvdDigis(NULL)
83  , fMvdClusters(NULL)
84  , fMvdHits(NULL)
85  , fMvdDigisMatch(NULL)
86  , fMvdClustersMatch(NULL)
87  , fMvdHitsMatch(NULL)
88  , fListMCTracks(NULL)
89  , fStsTrackArray(NULL)
90  , fStsTrackMatches(NULL)
91  , fMvdHisto1()
92  , fMvdHisto2()
93  , fProf()
94  , fNrMcPointsAll(-1)
95  , fNrHitsAll(-1)
96  , fMcperDigi()
97  , fMcperHit()
98  , fPixelpitch() {
99  ;
100 }
101 // -------------------------------------------------------------------------
102 
103 
104 // ----- Destructor ----------------------------------------------------
106 // -------------------------------------------------------------------------
107 
108 
109 // -------------------------------------------------------------------------
111  cout << "--------------------------------------------------------------------"
112  "-----"
113  << endl
114  << "-I- " << GetName() << "::Init: "
115  << " Start Initilisation " << endl
116  << "--------------------------------------------------------------------"
117  "-----"
118  << endl;
119 
120  FairRootManager* ioman = FairRootManager::Instance();
121  if (!ioman) {
122  cout << "-E- " << GetName() << "::Init: "
123  << "RootManager not instantised!" << endl;
124  return kFATAL;
125  }
126 
127  gGeoManager = (TGeoManager*) gROOT->FindObject("FAIRGeom");
128 
129  fMcPoints = (TClonesArray*) ioman->GetObject("MvdPoint");
130  fMvdDigis = (TClonesArray*) ioman->GetObject("MvdDigi");
131  fMvdClusters = (TClonesArray*) ioman->GetObject("MvdCluster");
132  fMvdHits = (TClonesArray*) ioman->GetObject("MvdHit");
133 
134  fMvdDigisMatch = (TClonesArray*) ioman->GetObject("MvdDigiMatch");
135  fMvdClustersMatch = (TClonesArray*) ioman->GetObject("MvdClusterMatch");
136  fMvdHitsMatch = (TClonesArray*) ioman->GetObject("MvdHitMatch");
137 
138  fListMCTracks = (TClonesArray*) ioman->GetObject("MCTrack");
139  // fStsTrackArray = (TClonesArray*) ioman->GetObject("StsTrack");
140  // fStsTrackMatches = (TClonesArray*) ioman->GetObject("StsTrackMatch");
141 
142  // fPrimVtx = (CbmVertex*) ioman->GetObject("PrimaryVertex");
143  /*
144  // Get pointer to PrimaryVertex object from IOManager if it exists
145  // The old name for the object is "PrimaryVertex" the new one
146  // "PrimaryVertex." Check first for the new name
147  fPrimVtx = dynamic_cast<CbmVertex*>(ioman->GetObject("PrimaryVertex."));
148  if (nullptr == fPrimVtx) {
149  fPrimVtx = dynamic_cast<CbmVertex*>(ioman->GetObject("PrimaryVertex"));
150  }
151  if (nullptr == fPrimVtx) {
152  LOG(fatal) << "No PrimaryVertex array!";
153  }
154 */
155  // fListMCTracks = (TClonesArray*) ioman->GetObject("MCTrack");
156  // fExtrapolator = CbmLitToolFactory::Instance()->CreateTrackExtrapolator("rk4");
157 
158 
159  fMvdHisto1[0] = new TH1F("Momentum", "Momentum", 100, 0, 30);
160  fMvdHisto1[1] = new TH1F("Angle", "Angle", 200, 0, 2);
161  fMvdHisto1[2] = new TH1F("DigisPerHit", "DigisPerHit", 100, 0, 100);
162  fMvdHisto1[3] = new TH1F("xResidual", "xResidual", 500, -50, 50);
163  fMvdHisto1[4] = new TH1F("yResidual", "yResidual", 500, -50, 50);
164  fMvdHisto1[5] = new TH1F("zResidual", "zResidual", 2000, -20, 20);
165  fMvdHisto1[6] = new TH1F("DistancePixX", "DistancePixX", 100, -2, 2);
166  fMvdHisto1[7] = new TH1F("DistancePixY", "DistancePixY", 100, -2, 2);
167  fMvdHisto1[8] = new TH1F("ClusterShape", "ClusterShape", 10, 0, 10);
168  fMvdHisto1[9] = new TH1F("ChargeSpectrum", "ChargeSpectrum", 50000, 0, 10000);
169  fMvdHisto1[10] = new TH1F("res_x_shape_0", "res_x_shape_0", 500, -50, 50);
170  fMvdHisto1[11] = new TH1F("res_x_shape_1", "res_x_shape_1", 500, -50, 50);
171  fMvdHisto1[12] = new TH1F("res_x_shape_2", "res_x_shape_2", 500, -50, 50);
172  fMvdHisto1[13] = new TH1F("res_x_shape_3", "res_x_shape_3", 500, -50, 50);
173  fMvdHisto1[14] = new TH1F("res_x_shape_4", "res_x_shape_4", 500, -50, 50);
174  fMvdHisto1[15] = new TH1F("res_x_shape_5", "res_x_shape_5", 500, -50, 50);
175  fMvdHisto1[16] = new TH1F("res_x_shape_6", "res_x_shape_6", 500, -50, 50);
176  fMvdHisto1[17] = new TH1F("res_x_shape_7", "res_x_shape_7", 500, -50, 50);
177  fMvdHisto1[18] = new TH1F("res_x_shape_8", "res_x_shape_8", 500, -50, 50);
178  fMvdHisto1[19] = new TH1F("res_x_shape_9", "res_x_shape_9", 500, -50, 50);
179  fMvdHisto1[20] = new TH1F("res_y_shape_0", "res_y_shape_0", 500, -50, 50);
180  fMvdHisto1[21] = new TH1F("res_y_shape_1", "res_y_shape_1", 500, -50, 50);
181  fMvdHisto1[22] = new TH1F("res_y_shape_2", "res_y_shape_2", 500, -50, 50);
182  fMvdHisto1[23] = new TH1F("res_y_shape_3", "res_y_shape_3", 500, -50, 50);
183  fMvdHisto1[24] = new TH1F("res_y_shape_4", "res_y_shape_4", 500, -50, 50);
184  fMvdHisto1[25] = new TH1F("res_y_shape_5", "res_y_shape_5", 500, -50, 50);
185  fMvdHisto1[26] = new TH1F("res_y_shape_6", "res_y_shape_6", 500, -50, 50);
186  fMvdHisto1[27] = new TH1F("res_y_shape_7", "res_y_shape_7", 500, -50, 50);
187  fMvdHisto1[28] = new TH1F("res_y_shape_8", "res_y_shape_8", 500, -50, 50);
188  fMvdHisto1[29] = new TH1F("res_y_shape_9", "res_y_shape_9", 500, -50, 50);
189  fMvdHisto1[30] = new TH1F("DigisPerMC", "DigisPerMC", 100, 0, 100);
190  fMvdHisto1[31] = new TH1F("HitsPerMC", "HitsPerMC", 11, 0, 11);
191  fMvdHisto1[32] =
192  new TH1F("McPerHit (merged)", "McPerHit (merged)", 11, 0, 11);
193  fMvdHisto1[33] = new TH1F("xPull", "xPull", 500, -50, 50);
194  fMvdHisto1[34] = new TH1F("yPull", "yPull", 500, -50, 50);
195  fMvdHisto1[35] = new TH1F("McPerDigi", "McPerDigi", 6, 0, 6);
196  fMvdHisto1[36] = new TH1F("McPerHit", "McPerHit", 11, 0, 11);
197 
198  fMvdHisto2[0] = new TH2F("dxpdyp", "dxpdyp", 100, -1, 1, 100, -1, 1);
199  fMvdHisto2[1] = new TH2F("dpnorm", "dpnorm", 100, -1, 1, 100, -1, 1);
200  fMvdHisto2[2] = new TH2F("mom_dx", "mom_dx", 100, 0, 30, 100, -50, 50);
201  fMvdHisto2[3] = new TH2F("mom_dy", "mom_dy", 100, 0, 30, 100, -50, 50);
202  fMvdHisto2[4] = new TH2F("ang_dx", "ang_dx", 100, 0, 2, 100, -50, 50);
203  fMvdHisto2[5] = new TH2F("ang_dy", "ang_dy", 100, 0, 2, 100, -50, 50);
204  fMvdHisto2[6] = new TH2F("dnr_dx", "dnr_dx", 100, 0, 100, 100, -50, 50);
205  fMvdHisto2[7] = new TH2F("dnr_dy", "dnr_dy", 100, 0, 100, 100, -50, 50);
206  fMvdHisto2[8] = new TH2F("dxp_dx", "dxp_dx", 50, 0, 1, 100, -50, 50);
207  fMvdHisto2[9] = new TH2F("dxp_dy", "dxp_dy", 50, 0, 1, 100, -50, 50);
208  fMvdHisto2[10] = new TH2F("cha_dx", "cha_dx", 100, 0, 10000, 100, -50, 50);
209  fMvdHisto2[11] = new TH2F("cha_dy", "cha_dy", 100, 0, 10000, 100, -50, 50);
210 
211  fMvdHisto2[12] = new TH2F("mom_dz", "mom_dz", 100, 0, 30, 2000, -20, 20);
212  fMvdHisto2[13] = new TH2F("ang_dz", "ang_dz", 200, 0, 2, 2000, -20, 20);
213  fMvdHisto2[14] = new TH2F("dx_dy", "dx_dy", 100, -50, 50, 100, -50, 50);
214 
215  fMvdHisto2[15] = new TH2F("mom_cha", "mom_cha", 100, 0, 3, 1000, 0, 10000);
216  fMvdHisto2[16] = new TH2F("ang_cha", "ang_cha", 200, 0, 2, 1000, 0, 10000);
217  fMvdHisto2[17] =
218  new TH2F("mom_chacut", "mom_chacut", 100, 0, 3, 1000, 0, 10000);
219  cout << "--------------------------------------------------------------------"
220  "-----"
221  << endl
222  << "-I- " << GetName() << "::Init: "
223  << " Finished Initilisation " << endl
224  << "--------------------------------------------------------------------"
225  "-----"
226  << endl;
227 
228  fNrMcPointsAll = 0;
229  fNrHitsAll = 0;
230 
231  for (Int_t i = 0; i < 6; i++) {
232  fMcperDigi[i] = 0;
233  }
234  for (Int_t i = 0; i < 11; i++) {
235  fMcperHit[i] = 0;
236  }
237 
238  return kSUCCESS;
239 }
240 // -------------------------------------------------------------------------
241 
242 
243 // -------------------------------------------------------------------------
244 void CbmMvdClusterAna::Exec(Option_t* /*opt*/) {
245  Int_t nMcpoints =
246  fMcPoints->GetEntriesFast(); // Number of Monte Carlo Points
247  Int_t nDigis = fMvdDigis->GetEntriesFast(); // Number of Mvd Digis
248  // Int_t nClusters = fMvdClusters->GetEntriesFast(); // Number of reconstructed Mvd Clusters
249  Int_t nHits = fMvdHits->GetEntriesFast(); // Number of reconstructed Mvd Hits
250 
251  // Int_t nDigisMatch = fMvdDigisMatch ->GetEntriesFast(); // Number of Matches from Digis to Monte Carlo
252  // Int_t nClustersMatch = fMvdClustersMatch ->GetEntriesFast(); // Number of Matches from Reconstructed Mvd Clusters to Mvd Hits (?)
253  // Int_t nHitsMatch = fMvdHitsMatch ->GetEntriesFast(); // Number of Matches from Reconstructed Mvd Hits to Monte Carlo
254 
255  // Int_t nMcTracks = fListMCTracks ->GetEntriesFast();
256  // Int_t nTracks = fStsTrackArray ->GetEntriesFast(); // Number of Tracks
257 
258  // cout<<"MC Points : "<< nMcpoints <<endl;
259  // cout<<"Hits : "<< nHits <<endl;
260  // cout<<"HitMatchs : "<< nHitsMatch <<endl;
261  // cout<<"Clusters : "<< nClusters <<endl;
262  // cout<<"ClusterMatch : "<< nClustersMatch<<endl;
263  // cout<<"Digis : "<< nDigis <<endl;
264  // cout<<"DigiMatchs : "<< nDigisMatch <<endl;
265  // cout<<"MC Tracks : "<< nMcTracks <<endl;
266  // cout<<"Tracks : "<< nTracks <<endl;
267  // cout<<"------"<<endl;
268  // -------------------
269  CbmMvdPoint* mvdPoint; // Monte Carlo Point
270  CbmMvdDigi* mvdDigi; // Digi
271  CbmMvdCluster* mvdCluster;
272  CbmMvdHit* mvdHit; // Digitized Hit
273 
274  CbmMvdDigiMatch* mvdDigiMatch; // Digi to MC point
275  CbmMatch* mvdClusterMatch; // Cluster to (?)
276  CbmMvdHitMatch* mvdHitMatch; // Hit to MC point
277 
278  CbmMCTrack* mcTrack;
279  // CbmStsTrack * stsTrack;
280  // CbmTrackMatchNew* trackMatch;
281  // CbmMatch * mvdMatch; // Reco Track to Hit Match
282  // -------------------
283  typedef map<pair<Int_t, Int_t>, Int_t>::iterator it_type;
284  map<pair<Int_t, Int_t>, Int_t> digiMap;
285  pair<Int_t, Int_t> digiCoor;
286  Int_t digiCharge;
287 
288  TVector3 cOrth;
289  TVector3 cVect;
290 
291  Int_t* digiList;
292  Double_t gloMC[3];
293  Double_t locMC[3];
294  Double_t gloHit[3];
295  Double_t locHit[3];
296  Double_t locRef[3];
297 
298  Double_t pitchx, pitchy;
299 
300  Int_t MAXHITS = nHits;
301  Int_t count;
302  Int_t charge;
303  UInt_t shape;
304  // bool correctshape;
305  Double_t ARR_momentum[MAXHITS];
306  Double_t ARR_angle[MAXHITS];
307  Int_t ARR_digis[MAXHITS];
308  Double_t ARR_dx[MAXHITS];
309  Double_t ARR_dy[MAXHITS];
310  Double_t ARR_dz[MAXHITS];
311  Double_t ARR_dxp[MAXHITS];
312  Double_t ARR_dyp[MAXHITS];
313  UInt_t ARR_shape[MAXHITS];
314  Int_t ARR_charge[MAXHITS];
315  Int_t XAXIS[1000];
316  Int_t YAXIS[1000];
317  Int_t xaxis, yaxis, xaxismin, xaxismax, yaxismin, yaxismax, xrel, yrel;
318  Int_t POS;
319  // -------------------
320  fNrMcPointsAll += nMcpoints;
321  fNrHitsAll += nHits;
322  // -------------------
323  mvdDigi = (CbmMvdDigi*) fMvdDigis->At(0);
324  pitchx = mvdDigi->GetPixelSizeX();
325  pitchy = mvdDigi->GetPixelSizeY();
326  fPixelpitch[0] = pitchx;
327  fPixelpitch[1] = pitchy;
328 
329  TGeoVolume* CurrentVolume;
330  TGeoBBox* VolumeShape;
331  // -------------------
332  std::map<std::pair<std::pair<Int_t, Int_t>, TString>, std::vector<int>>
333  DigisMap;
334  std::map<std::pair<std::pair<Int_t, Int_t>, TString>,
335  std::vector<int>>::iterator it;
336  std::pair<std::pair<Int_t, Int_t>, TString> DigiStation;
337  std::pair<Int_t, Int_t> Digi;
338  std::vector<int> McContrToHitList;
339  std::vector<int>
340  McContrList; // vector of Mc Points which contribute to a Digi
341  std::map<Int_t, Int_t> McInHit;
342  std::map<Int_t, Int_t>::iterator it2;
343  std::vector<int> DigisInMc(nMcpoints, 0);
344  std::map<Int_t, std::vector<int>> McsInHit;
345  std::map<Int_t, std::vector<int>> HitsInMc;
346  std::map<Int_t, std::vector<int>>::iterator it3;
347 
348  DigisMap.clear();
349  McInHit.clear();
350  McsInHit.clear();
351  HitsInMc.clear();
352 
353  // -------------------
354  // Analyze Digis:
355 
356  for (int iDigi = 0; iDigi < nDigis; iDigi++) {
357 
358  mvdDigi = (CbmMvdDigi*) fMvdDigis->At(iDigi);
359  mvdDigiMatch = (CbmMvdDigiMatch*) fMvdDigisMatch->At(iDigi);
360 
361  Int_t nMatchedIndex = mvdDigiMatch->GetMatchedLink().GetIndex();
362  mvdPoint = (CbmMvdPoint*) fMcPoints->At(
363  nMatchedIndex); // Get matched MC Point from Digi
364 
365  // mcTrack = (CbmMCTrack*) fListMCTracks->At ( mvdpoint->GetTrackID() ); // Get matched MC Track from MC Point
366 
367  gloMC[0] = (mvdPoint->GetXOut() + mvdPoint->GetX()) / 2.;
368  gloMC[1] = (mvdPoint->GetYOut() + mvdPoint->GetY()) / 2.;
369  gloMC[2] = (mvdPoint->GetZOut() + mvdPoint->GetZ()) / 2.;
370 
371  gGeoManager->FindNode(gloMC[0], gloMC[1], gloMC[2]);
372 
373  CurrentVolume = gGeoManager->GetCurrentVolume();
374  // VolumeShape = (TGeoBBox*)CurrentVolume->GetShape();
375 
376  Digi.first = mvdDigi->GetPixelX();
377  Digi.second = mvdDigi->GetPixelY();
378 
379  DigiStation.first = Digi;
380  DigiStation.second = CurrentVolume->GetName();
381 
382  McContrList.clear();
383 
384  for (Int_t iLink = 0; iLink < mvdDigiMatch->GetNofLinks(); iLink++) {
385  McContrList.push_back(mvdDigiMatch->GetLink(iLink).GetIndex());
386  }
387 
388  // it = DigisMap.find(DigiStation);
389  // if (it == DigisMap.end())
390  // {
391  DigisMap[DigiStation] = McContrList;
392  // }
393  // else // Should not be possible
394  // {
395  // }
396  }
397 
398  // -------------------
399  // Analyze Hits and Clusters
400 
401  for (Int_t iHit = 0; iHit < nHits; iHit++) {
402  cout << endl
403  << "run next Hit " << iHit << " of " << nHits << " Hits" << endl;
404  mvdHit = (CbmMvdHit*) fMvdHits->At(iHit);
405  mvdCluster = (CbmMvdCluster*) fMvdClusters->At(iHit);
406  mvdHitMatch = (CbmMvdHitMatch*) fMvdHitsMatch->At(iHit);
407  mvdClusterMatch = (CbmMatch*) fMvdClustersMatch->At(iHit);
408 
409  Int_t nHitMatch = mvdHitMatch->GetMatchedLink().GetIndex();
410  mvdPoint = (CbmMvdPoint*) fMcPoints->At(nHitMatch);
411 
412  mcTrack = (CbmMCTrack*) fListMCTracks->At(mvdPoint->GetTrackID());
413 
414  gloMC[0] = (mvdPoint->GetXOut() + mvdPoint->GetX()) / 2.;
415  gloMC[1] = (mvdPoint->GetYOut() + mvdPoint->GetY()) / 2.;
416  gloMC[2] = (mvdPoint->GetZOut() + mvdPoint->GetZ()) / 2.;
417 
418  gloHit[0] = mvdHit->GetX();
419  gloHit[1] = mvdHit->GetY();
420  gloHit[2] = mvdHit->GetZ();
421 
422  gGeoManager->FindNode(gloMC[0], gloMC[1], gloMC[2]);
423  gGeoManager->MasterToLocal(gloMC, locMC);
424  gGeoManager->MasterToLocal(gloHit, locHit);
425 
426  CurrentVolume = gGeoManager->GetCurrentVolume();
427  VolumeShape = (TGeoBBox*) CurrentVolume->GetShape();
428 
429  locRef[0] =
430  (mvdHit->GetIndexCentralX() + 0.5) * pitchx - VolumeShape->GetDX();
431  locRef[1] =
432  (mvdHit->GetIndexCentralY() + 0.5) * pitchy - VolumeShape->GetDY();
433  locRef[2] = 0;
434 
435  cOrth.SetXYZ(0., 0., mvdPoint->GetZOut() - mvdPoint->GetZ());
436  cVect.SetXYZ(mvdPoint->GetXOut() - mvdPoint->GetX(),
437  mvdPoint->GetYOut() - mvdPoint->GetY(),
438  mvdPoint->GetZOut() - mvdPoint->GetZ());
439 
440  digiMap = mvdCluster->GetPixelMap();
441  digiList = mvdCluster->GetDigiList();
442 
443  count = 0;
444  charge = 0;
445  shape = 0;
446 
447  McInHit.clear();
448  McContrToHitList.clear();
449  cout << endl << "ping" << endl;
450  for (it_type iterator = digiMap.begin(); iterator != digiMap.end();
451  iterator++) {
452  digiCoor = iterator->first;
453  digiCharge = iterator->second;
454  xaxis = digiCoor.first;
455  yaxis = digiCoor.second;
456 
457  XAXIS[count] = xaxis;
458  YAXIS[count] = yaxis;
459 
460  if (count == 0) {
461  xaxismin = xaxis;
462  xaxismax = xaxis;
463  yaxismin = yaxis;
464  yaxismax = yaxis;
465  } else {
466  if (xaxismin > xaxis) { xaxismin = xaxis; }
467  if (xaxismax < xaxis) { xaxismax = xaxis; }
468  if (yaxismin > yaxis) { yaxismin = yaxis; }
469  if (yaxismax < yaxis) { yaxismax = yaxis; }
470  }
471 
472  charge += digiCharge;
473  count++;
474 
475  DigiStation.first = digiCoor;
476  DigiStation.second = CurrentVolume->GetName();
477 
478  it = DigisMap.find(DigiStation);
479 
480  McContrList = it->second; // Mcs contributing to a Digi
481 
482  // if (it != DigisMap.end())
483  // {
484  fMcperDigi[McContrList
485  .size()]++; // Number of MC Points contributing to a Digi
486  cout << endl << "ping 2" << endl;
487  for (Int_t iMc = 0; iMc < McContrList.size(); iMc++) {
488  if (std::find(McContrToHitList.begin(),
489  McContrToHitList.end(),
490  McContrList[iMc])
491  == McContrToHitList.end()) {
492  McContrToHitList.push_back(McContrList[iMc]);
493  }
494 
495  it2 = McInHit.find(McContrList[iMc]);
496 
497  if (it2 != McInHit.end()) {
498  McInHit[McContrList[iMc]]++;
499  } else {
500  McInHit[McContrList[iMc]] = 1;
501  }
502  }
503  // }
504  // else // should not be possible
505  // {
506  // }
507  }
508 
509  McsInHit[iHit] = McContrToHitList;
510 
511  if (McInHit.size() < 11) {
512  fMcperHit[McInHit.size()]++;
513  } else {
514  fMcperHit[10]++;
515  }
516  count = 0;
517  cout << endl << "ping 3" << endl;
518  for (std::map<Int_t, Int_t>::iterator iterator = McInHit.begin();
519  iterator != McInHit.end();
520  iterator++) {
521  DigisInMc[iterator->first] += iterator->second;
522  count++;
523 
524  it3 = HitsInMc.find(iterator->first);
525 
526  if (it3 != HitsInMc.end()) {
527  McContrList = it3->second;
528  McContrList.push_back(iHit);
529  } else {
530  McContrList.clear();
531  McContrList.push_back(iHit);
532  }
533 
534  HitsInMc[iterator->first] = McContrList;
535  }
536  cout << endl << "ping 4" << endl;
537  if ((xaxismax - xaxismin) < 4 && (yaxismax - yaxismin) < 4) {
538  shape = 0;
539 
540  for (int i = 0; i < mvdCluster->GetNofDigis(); i++) {
541  xrel = XAXIS[i] - xaxismin;
542  yrel = YAXIS[i] - yaxismin;
543  POS = xrel + 5 * yrel;
544  shape += TMath::Power(2, POS);
545  }
546 
547  if (shape == 3) {
548  shape = 0;
549  } else if (shape == 99) {
550  shape = 1;
551  } else if (shape == 33) {
552  shape = 2;
553  } else if (shape == 7) {
554  shape = 3;
555  } else if (shape == 67) {
556  shape = 4;
557  } else if (shape == 97) {
558  shape = 5;
559  } else if (shape == 35) {
560  shape = 6;
561  } else if (shape == 98) {
562  shape = 7;
563  } else if (shape == 1) {
564  shape = 8;
565  } else {
566  shape = 9;
567  }
568  } else {
569  shape = 9;
570  }
571  cout << endl << "ping 5" << endl;
572  ARR_momentum[iHit] = mcTrack->GetP();
573  ARR_angle[iHit] = cVect.Angle(cOrth);
574  ARR_digis[iHit] = mvdCluster->GetNofDigis();
575  ARR_dx[iHit] = 10000 * (locHit[0] - locMC[0]);
576  ARR_dy[iHit] = 10000 * (locHit[1] - locMC[1]);
577  ARR_dz[iHit] = 10000 * (gloHit[2] - gloMC[2]);
578  ARR_dxp[iHit] =
579  -(int) (((locMC[0] + VolumeShape->GetDX()) / (1 * pitchx)) - 0.5)
580  + (double) (((locMC[0] + VolumeShape->GetDX()) / (1 * pitchx)) - 0.5);
581  ARR_dyp[iHit] =
582  -(int) (((locMC[1] + VolumeShape->GetDY()) / (1 * pitchy)) - 0.5)
583  + (double) (((locMC[1] + VolumeShape->GetDY()) / (1 * pitchy)) - 0.5);
584  ARR_shape[iHit] = shape;
585  ARR_charge[iHit] = charge;
586  cout << endl << "finished this hit" << endl;
587  }
588 
589  // -------------------
590  for (Int_t iHit = 0; iHit < nHits; iHit++) {
591  cout << endl << "fill stuff" << endl;
592  fMvdHisto1[0]->Fill(ARR_momentum[iHit]);
593  fMvdHisto1[1]->Fill(ARR_angle[iHit]);
594  fMvdHisto1[2]->Fill(ARR_digis[iHit]);
595  fMvdHisto1[3]->Fill(ARR_dx[iHit]);
596  fMvdHisto1[4]->Fill(ARR_dy[iHit]);
597  fMvdHisto1[5]->Fill(ARR_dz[iHit]);
598  fMvdHisto1[6]->Fill(ARR_dxp[iHit]);
599  fMvdHisto1[7]->Fill(ARR_dyp[iHit]);
600  fMvdHisto1[8]->Fill(ARR_shape[iHit]);
601  fMvdHisto1[9]->Fill(ARR_charge[iHit]);
602  // fMvdHisto1[10+ARR_shape[iHit]]->Fill(ARR_dx[iHit]);
603  // fMvdHisto1[20+ARR_shape[iHit]]->Fill(ARR_dy[iHit]);
604 
605  fMvdHisto2[0]->Fill(ARR_dxp[iHit], ARR_dyp[iHit], ARR_digis[iHit]);
606  fMvdHisto2[0]->Fill(ARR_dxp[iHit], -ARR_dyp[iHit], ARR_digis[iHit]);
607  fMvdHisto2[0]->Fill(-ARR_dxp[iHit], ARR_dyp[iHit], ARR_digis[iHit]);
608  fMvdHisto2[0]->Fill(-ARR_dxp[iHit], -ARR_dyp[iHit], ARR_digis[iHit]);
609  fMvdHisto2[1]->Fill(ARR_dxp[iHit], ARR_dyp[iHit]);
610  fMvdHisto2[1]->Fill(ARR_dxp[iHit], -ARR_dyp[iHit]);
611  fMvdHisto2[1]->Fill(-ARR_dxp[iHit], ARR_dyp[iHit]);
612  fMvdHisto2[1]->Fill(-ARR_dxp[iHit], -ARR_dyp[iHit]);
613  fMvdHisto2[2]->Fill(ARR_momentum[iHit], ARR_dx[iHit]);
614  fMvdHisto2[3]->Fill(ARR_momentum[iHit], ARR_dy[iHit]);
615  fMvdHisto2[4]->Fill(ARR_angle[iHit], ARR_dx[iHit]);
616  fMvdHisto2[5]->Fill(ARR_angle[iHit], ARR_dy[iHit]);
617  fMvdHisto2[6]->Fill(ARR_digis[iHit], ARR_dx[iHit]);
618  fMvdHisto2[7]->Fill(ARR_digis[iHit], ARR_dy[iHit]);
619  // fMvdHisto2[8]->Fill(ARR_dxp[iHit], ARR_dx[iHit]);
620  // fMvdHisto2[9]->Fill(ARR_dyp[iHit], ARR_dy[iHit]);
621  fMvdHisto2[10]->Fill(ARR_charge[iHit], ARR_dx[iHit]);
622  fMvdHisto2[11]->Fill(ARR_charge[iHit], ARR_dy[iHit]);
623  fMvdHisto2[12]->Fill(ARR_momentum[iHit], ARR_dz[iHit]);
624  fMvdHisto2[13]->Fill(ARR_angle[iHit], ARR_dz[iHit]);
625  fMvdHisto2[14]->Fill(ARR_dx[iHit], ARR_dy[iHit]);
626 
627  fMvdHisto1[33]->Fill(ARR_dx[iHit] / 3.8);
628  fMvdHisto1[34]->Fill(ARR_dy[iHit] / 4.8);
629 
630  fMvdHisto2[15]->Fill(ARR_momentum[iHit], ARR_charge[iHit]);
631  fMvdHisto2[16]->Fill(ARR_angle[iHit], ARR_charge[iHit]);
632 
633  if (ARR_angle[iHit] < 0.3) {
634  fMvdHisto2[17]->Fill(ARR_momentum[iHit], ARR_charge[iHit]);
635  }
636  }
637 
638  for (int iMc = 0; iMc < nMcpoints; iMc++) {
639  fMvdHisto1[30]->Fill(DigisInMc[iMc]);
640  }
641 
642  for (std::map<Int_t, std::vector<int>>::iterator iterator = HitsInMc.begin();
643  iterator != HitsInMc.end();
644  iterator++) {
645  McContrList = iterator->second;
646  fMvdHisto1[31]->Fill(McContrList.size());
647  }
648 
649  bool criteria;
650 
651  for (std::map<Int_t, std::vector<int>>::iterator iterator = McsInHit.begin();
652  iterator != McsInHit.end();
653  iterator++) {
654  criteria = true;
655 
656  McContrToHitList = iterator->second;
657 
658  for (int i = 0; i < McContrToHitList.size(); i++) {
659  it3 = HitsInMc.find(McContrToHitList[i]);
660 
661  if (it3 != HitsInMc.end()) {
662  if ((it3->second).size() != 1) { criteria = false; }
663  }
664  }
665 
666  if (criteria) { fMvdHisto1[32]->Fill(McContrToHitList.size()); }
667  }
668  // -------------------
669 }
670 // -------------------------------------------------------------------------
671 
672 
673 // -------------------------------------------------------------------------
675  cout << "======================" << endl;
676  cout << "'Mvd QA Output' Start!" << endl;
677  cout << "======================" << endl;
678  // -------------------
679  cout << "----------------------" << endl;
680  cout << "Pixelpitch: " << endl;
681  cout << "x: " << fPixelpitch[0] << endl;
682  cout << "y: " << fPixelpitch[1] << endl;
683  // -------------------
684  int sum;
685  cout << "----------------------" << endl;
686  cout << "MC points per Digi:" << endl;
687  sum = 0;
688  for (Int_t i = 0; i < 6; i++) {
689  sum += fMcperDigi[i];
690  }
691  for (Int_t i = 0; i < 6; i++) {
692  fMvdHisto1[35]->Fill(i, fMcperDigi[i]);
693  printf(
694  " %2i %10i %6.5f \n", i, fMcperDigi[i], 1. * fMcperDigi[i] / sum);
695  }
696  // -------------------
697  cout << "----------------------" << endl;
698  cout << "Digis per MC Point:" << endl;
699  for (Int_t i = 0; i < 11; i++) {
700  if (fMvdHisto1[30]->GetBinCenter(i) >= 0)
701  printf(" %2i %10i %6.5f \n",
702  (int) (fMvdHisto1[30]->GetBinCenter(i)),
703  (int) (fMvdHisto1[30]->GetBinContent(i)),
704  1. * (int) (fMvdHisto1[30]->GetBinContent(i)) * 1.
705  / fMvdHisto1[30]->GetEntries());
706  }
707  cout << " .." << endl;
708  // -------------------
709  cout << "----------------------" << endl;
710  cout << "Digis per Hit:" << endl;
711  for (Int_t i = 0; i < 11; i++) {
712  if (fMvdHisto1[2]->GetBinCenter(i) >= 0)
713  printf(" %2i %10i %6.5f \n",
714  (int) (fMvdHisto1[2]->GetBinCenter(i)),
715  (int) (fMvdHisto1[2]->GetBinContent(i)),
716  1. * (int) (fMvdHisto1[2]->GetBinContent(i)) * 1.
717  / fMvdHisto1[2]->GetEntries());
718  }
719  cout << " .." << endl;
720  // -------------------
721  cout << "----------------------" << endl;
722  cout << "Hits per MC Point:" << endl;
723  for (Int_t i = 0; i < 11; i++) {
724  if (fMvdHisto1[31]->GetBinCenter(i) >= 0)
725  printf(" %2i %10i %6.5f \n",
726  (int) (fMvdHisto1[31]->GetBinCenter(i)),
727  (int) (fMvdHisto1[31]->GetBinContent(i)),
728  1. * (int) (fMvdHisto1[31]->GetBinContent(i)) * 1.
729  / fMvdHisto1[31]->GetEntries());
730  }
731  cout << " .." << endl;
732  // -------------------
733  cout << "----------------------" << endl;
734  cout << "MC Points per Hit:" << endl;
735  sum = 0;
736  for (Int_t i = 0; i < 11; i++) {
737  sum += fMcperHit[i];
738  }
739  for (Int_t i = 0; i < 11; i++) {
740  fMvdHisto1[36]->Fill(i, fMcperHit[i]);
741  printf(" %2i %10i %6.5f \n", i, fMcperHit[i], 1. * fMcperHit[i] / sum);
742  }
743  // -------------------
744  cout << "----------------------" << endl;
745  cout << "MC Points per Hit (merged Clusters):" << endl;
746  for (Int_t i = 0; i < 11; i++) {
747  if (fMvdHisto1[32]->GetBinCenter(i) >= 0)
748  printf(" %2i %10i %6.5f \n",
749  (int) (fMvdHisto1[32]->GetBinCenter(i)),
750  (int) (fMvdHisto1[32]->GetBinContent(i)),
751  1. * (int) (fMvdHisto1[32]->GetBinContent(i)) * 1.
752  / fMvdHisto1[32]->GetEntries());
753  }
754  cout << " .." << endl;
755  // -------------------
756  cout << "----------------------" << endl;
757  // -------------------
758  Double_t xpar[3];
759  Double_t ypar[3];
760 
761  TF1* gFitx = new TF1("gausx", "gaus", -50, 50);
762  TF1* gFity = new TF1("gausy", "gaus", -50, 50);
763 
764  fMvdHisto1[3]->Fit(gFitx, "QRN0");
765  fMvdHisto1[4]->Fit(gFity, "QRN0");
766 
767  gFitx->GetParameters(&xpar[0]);
768  gFity->GetParameters(&ypar[0]);
769 
770  cout
771  << setw(40)
772  << "---------------------------------------------------------------------"
773  << endl;
774  cout << "Resolution:" << endl;
775  cout
776  << setw(40)
777  << "---------------------------------------------------------------------"
778  << endl;
779  cout << setw(9) << " Shape";
780  cout << setw(10) << " Mean(x)";
781  cout << setw(10) << " RMS(x)";
782  cout << setw(10) << " Mean(y)";
783  cout << setw(10) << " RMS(y)";
784  cout << setw(10) << " Mean(z)";
785  cout << setw(10) << " RMS(z)";
786  cout << endl;
787  cout
788  << setw(40)
789  << "---------------------------------------------------------------------";
790  cout << endl;
791 
792  for (int i = 0; i < 10; i++) {
793  printf("%6i %+5.4f %+5.4f %+5.4f %+5.4f \n",
794  i,
795  fMvdHisto1[10 + i]->GetMean(),
796  fMvdHisto1[10 + i]->GetRMS(),
797  fMvdHisto1[20 + i]->GetMean(),
798  fMvdHisto1[20 + i]->GetRMS());
799  }
800  cout
801  << setw(40)
802  << "---------------------------------------------------------------------"
803  << endl;
804  printf("%10s%+5.4f %+5.4f %+5.4f %+5.4f %+5.4f %+5.4f \n",
805  "All ",
806  fMvdHisto1[3]->GetMean(),
807  fMvdHisto1[3]->GetRMS(),
808  fMvdHisto1[4]->GetMean(),
809  fMvdHisto1[4]->GetRMS(),
810  fMvdHisto1[5]->GetMean(),
811  fMvdHisto1[5]->GetRMS());
812  printf("%10s%+5.4f %+5.4f %+5.4f %+5.4f \n",
813  "All (Fit) ",
814  xpar[1],
815  xpar[2],
816  ypar[1],
817  ypar[2]);
818 
819  cout
820  << setw(40)
821  << "---------------------------------------------------------------------"
822  << endl;
823  fMvdHisto1[33]->Fit(gFitx, "QRN0");
824  fMvdHisto1[34]->Fit(gFity, "QRN0");
825 
826  gFitx->GetParameters(&xpar[0]);
827  gFity->GetParameters(&ypar[0]);
828  cout << "Pulls:" << endl;
829  printf("%10s%+5.4f %+5.4f %+5.4f %+5.4f \n",
830  "All ",
831  fMvdHisto1[33]->GetMean(),
832  fMvdHisto1[33]->GetRMS(),
833  fMvdHisto1[34]->GetMean(),
834  fMvdHisto1[34]->GetRMS());
835  printf("%10s%+5.4f %+5.4f %+5.4f %+5.4f \n",
836  "All (Fit) ",
837  xpar[1],
838  xpar[2],
839  ypar[1],
840  ypar[2]);
841  cout << "(Using 3.8(x) and 4.8(y) as nominal resolution!)" << endl;
842  cout
843  << setw(40)
844  << "---------------------------------------------------------------------"
845  << endl;
846  cout
847  << setw(40)
848  << "---------------------------------------------------------------------"
849  << endl;
850  // -------------------
851  cout << "Hit Reco Efficiency: " << float(1. * fNrHitsAll / fNrMcPointsAll)
852  << "\t( " << fNrHitsAll << "\t" << fNrMcPointsAll << " )" << endl;
853  cout
854  << setw(40)
855  << "---------------------------------------------------------------------"
856  << endl;
857  // -------------------
858  cout << "======================" << endl;
859  cout << "'Mvd QA Output' End!" << endl;
860  cout << "======================" << endl;
861 
862  // -------------------
863  TCanvas* a = new TCanvas("Residuals", "Residuals");
864  a->Divide(3, 3);
865 
866  a->cd(1);
867  fMvdHisto1[3]->Draw();
868  fMvdHisto1[3]->GetXaxis()->SetTitle("Residual in x [um]");
869  fMvdHisto1[3]->GetYaxis()->SetTitle("Entries");
870 
871  a->cd(2);
872  fMvdHisto1[4]->Draw();
873  fMvdHisto1[4]->GetXaxis()->SetTitle("Residual in y [um]");
874  fMvdHisto1[4]->GetYaxis()->SetTitle("Entries");
875 
876  a->cd(3);
877  fMvdHisto1[5]->Draw();
878  fMvdHisto1[5]->GetXaxis()->SetTitle("Residual in z [um]");
879  fMvdHisto1[5]->GetYaxis()->SetTitle("Entries");
880 
881  a->cd(4);
882  fMvdHisto1[33]->Draw();
883  fMvdHisto1[33]->GetXaxis()->SetTitle("Pull in x [um]");
884  fMvdHisto1[33]->GetYaxis()->SetTitle("Entries");
885 
886  a->cd(5);
887  fMvdHisto1[34]->Draw();
888  fMvdHisto1[34]->GetXaxis()->SetTitle("Pull in y [um]");
889  fMvdHisto1[34]->GetYaxis()->SetTitle("Entries");
890 
891  int col[10] = {1, 2, 3, 4, 5, 6, 7, 8, 9, 11};
892 
893  TLegend* leg = new TLegend(0.1, 0.7, 0.48, 0.9);
894 
895  a->cd(7);
896  for (int i = 0; i < 10; i++) {
897  if (i == 0) {
898  fMvdHisto1[10 + i]->Draw();
899  fMvdHisto1[10 + i]->GetXaxis()->SetTitle("Residual in x [um]");
900  fMvdHisto1[10 + i]->GetYaxis()->SetTitle("Entries");
901  } else {
902  fMvdHisto1[10 + i]->Draw("same");
903  }
904  fMvdHisto1[10 + i]->SetLineColor(col[i]);
905  leg->AddEntry(fMvdHisto1[10 + i], Form("Shape %i", i));
906  }
907  leg->Draw();
908 
909  a->cd(8);
910  for (int i = 0; i < 10; i++) {
911  if (i == 0) {
912  fMvdHisto1[20 + i]->Draw();
913  fMvdHisto1[20 + i]->GetXaxis()->SetTitle("Residual in y [um]");
914  fMvdHisto1[20 + i]->GetYaxis()->SetTitle("Entries");
915  } else {
916  fMvdHisto1[20 + i]->Draw("same");
917  }
918  fMvdHisto1[20 + i]->SetLineColor(col[i]);
919  }
920  leg->Draw();
921  // -------------------
922  TCanvas* b = new TCanvas("Matches", "Matches");
923  b->Divide(3, 2);
924 
925  b->cd(1);
926  fMvdHisto1[35]->Draw();
927  fMvdHisto1[35]->GetXaxis()->SetTitle("Mc Points per Digi");
928  fMvdHisto1[35]->GetYaxis()->SetTitle("Entries");
929 
930  b->cd(2);
931  fMvdHisto1[30]->Draw();
932  fMvdHisto1[30]->GetXaxis()->SetTitle("Digis per Mc Point");
933  fMvdHisto1[30]->GetYaxis()->SetTitle("Entries");
934 
935  b->cd(3);
936  fMvdHisto1[2]->Draw();
937  fMvdHisto1[2]->GetXaxis()->SetTitle("Digis per Hit");
938  fMvdHisto1[2]->GetYaxis()->SetTitle("Entries");
939 
940  b->cd(4);
941  fMvdHisto1[31]->Draw();
942  fMvdHisto1[31]->GetXaxis()->SetTitle("Hits per Mc Point");
943  fMvdHisto1[31]->GetYaxis()->SetTitle("Entries");
944 
945  b->cd(5);
946  fMvdHisto1[36]->Draw();
947  fMvdHisto1[36]->GetXaxis()->SetTitle("Mc Points per Hit");
948  fMvdHisto1[36]->GetYaxis()->SetTitle("Entries");
949 
950  b->cd(6);
951  fMvdHisto1[32]->Draw();
952  fMvdHisto1[32]->GetXaxis()->SetTitle("Mc Points per Hit (merged)");
953  fMvdHisto1[32]->GetYaxis()->SetTitle("Entries");
954  // -------------------
955  TCanvas* d = new TCanvas("ResAna", "ResAna");
956  d->Divide(4, 2);
957 
958  fProf[0] = fMvdHisto2[2]->ProfileX("0", 1, -1, "s");
959  fProf[1] = fMvdHisto2[3]->ProfileX("1", 1, -1, "s");
960  fProf[2] = fMvdHisto2[4]->ProfileX("2", 1, -1, "s");
961  fProf[3] = fMvdHisto2[5]->ProfileX("3", 1, -1, "s");
962  fProf[4] = fMvdHisto2[6]->ProfileX("4", 1, -1, "s");
963  fProf[5] = fMvdHisto2[7]->ProfileX("5", 1, -1, "s");
964  fProf[6] = fMvdHisto2[10]->ProfileX("6", 1, -1, "s");
965  fProf[7] = fMvdHisto2[11]->ProfileX("7", 1, -1, "s");
966 
967  d->cd(1);
968  fProf[0]->Draw();
969  fProf[0]->GetXaxis()->SetTitle("Momentum [GeV]");
970  fProf[0]->GetYaxis()->SetTitle("Residual in x [um]");
971 
972  d->cd(5);
973  fProf[1]->Draw();
974  fProf[1]->GetXaxis()->SetTitle("Momentum [GeV]");
975  fProf[1]->GetYaxis()->SetTitle("Residual in y [um]");
976 
977  d->cd(2);
978  fProf[2]->Draw();
979  fProf[2]->GetXaxis()->SetTitle("Angle [rad]");
980  fProf[2]->GetYaxis()->SetTitle("Residual in x [um]");
981 
982  d->cd(6);
983  fProf[3]->Draw();
984  fProf[3]->GetXaxis()->SetTitle("Angle [rad]");
985  fProf[3]->GetYaxis()->SetTitle("Residual in y [um]");
986 
987  d->cd(3);
988  fProf[4]->Draw();
989  fProf[4]->GetXaxis()->SetTitle("Digis per Hit");
990  fProf[4]->GetYaxis()->SetTitle("Residual in x [um]");
991 
992  d->cd(7);
993  fProf[5]->Draw();
994  fProf[5]->GetXaxis()->SetTitle("Digis per Hit");
995  fProf[5]->GetYaxis()->SetTitle("Residual in y [um]");
996 
997  d->cd(4);
998  fProf[6]->Draw();
999  fProf[6]->GetXaxis()->SetTitle("Charge");
1000  fProf[6]->GetYaxis()->SetTitle("Residual in x [um]");
1001 
1002  d->cd(8);
1003  fProf[7]->Draw();
1004  fProf[7]->GetXaxis()->SetTitle("Charge");
1005  fProf[7]->GetYaxis()->SetTitle("Residual in y [um]");
1006 
1007  TCanvas* c = new TCanvas("Rest", "Rest");
1008  c->Divide(6, 2);
1009 
1010  c->cd(1);
1011  fMvdHisto1[0]->Draw();
1012  fMvdHisto1[0]->GetXaxis()->SetTitle("Momentum [GeV]");
1013  fMvdHisto1[0]->GetYaxis()->SetTitle("Entries");
1014 
1015  c->cd(2);
1016  fMvdHisto1[1]->Draw();
1017  fMvdHisto1[1]->GetXaxis()->SetTitle("Angle [rad]");
1018  fMvdHisto1[1]->GetYaxis()->SetTitle("Entries");
1019 
1020  c->cd(3);
1021  fMvdHisto1[2]->Draw();
1022  fMvdHisto1[2]->GetXaxis()->SetTitle("Digis per Hit");
1023  fMvdHisto1[2]->GetYaxis()->SetTitle("Entries");
1024 
1025  c->cd(4);
1026  fMvdHisto2[14]->Draw("colz");
1027  fMvdHisto2[14]->GetXaxis()->SetTitle("Residual in x [um]");
1028  fMvdHisto2[14]->GetYaxis()->SetTitle("Residual in y [um]");
1029  fMvdHisto2[14]->GetZaxis()->SetTitle("Entries");
1030 
1031  c->cd(5);
1032  fMvdHisto2[15]->Draw("colz");
1033  fMvdHisto2[15]->GetXaxis()->SetTitle("Momentum [GeV]");
1034  fMvdHisto2[15]->GetYaxis()->SetTitle("Charge");
1035  fMvdHisto2[15]->GetZaxis()->SetTitle("Entries");
1036 
1037  c->cd(6);
1038  fMvdHisto2[16]->Draw("colz");
1039  fMvdHisto2[16]->GetXaxis()->SetTitle("Angle [rad]");
1040  fMvdHisto2[16]->GetYaxis()->SetTitle("Charge");
1041  fMvdHisto2[16]->GetZaxis()->SetTitle("Entries");
1042 
1043  c->cd(7);
1044  fMvdHisto1[6]->Draw();
1045  fMvdHisto1[6]->GetXaxis()->SetTitle("Distance to Pixel in x [Pixelpitch]");
1046  fMvdHisto1[6]->GetYaxis()->SetTitle("Entries");
1047 
1048  c->cd(8);
1049  fMvdHisto1[7]->Draw();
1050  fMvdHisto1[7]->GetXaxis()->SetTitle("Distance to Pixel in y [Pixelpitch]");
1051  fMvdHisto1[7]->GetYaxis()->SetTitle("Entries");
1052 
1053  c->cd(9);
1054  fMvdHisto1[8]->Draw();
1055  fMvdHisto1[8]->GetXaxis()->SetTitle("Clustershape");
1056  fMvdHisto1[8]->GetYaxis()->SetTitle("Entries");
1057 
1058  c->cd(10);
1059  fMvdHisto1[9]->Draw();
1060  fMvdHisto1[9]->GetXaxis()->SetTitle("Charge");
1061  fMvdHisto1[9]->GetYaxis()->SetTitle("Entries");
1062 
1063  fProf[8] = fMvdHisto2[15]->ProfileX("8", 1, -1, "");
1064  fProf[9] = fMvdHisto2[16]->ProfileX("9", 1, -1, "");
1065  fProf[10] = fMvdHisto2[17]->ProfileX("10", 1, -1, "");
1066 
1067  c->cd(11);
1068  fProf[8]->Draw();
1069  fProf[8]->GetXaxis()->SetTitle("Momentum [GeV]");
1070  fProf[8]->GetYaxis()->SetTitle("Charge");
1071  fProf[10]->Draw("same");
1072  fProf[10]->SetLineColor(2);
1073 
1074  c->cd(12);
1075  fProf[9]->Draw();
1076  fProf[9]->GetXaxis()->SetTitle("Angle [rad]");
1077  fProf[9]->GetYaxis()->SetTitle("Charge");
1078 
1079  // fMvdHisto2[0]->Divide(fMvdHisto2[1]);
1080  // TCanvas* e=new TCanvas();
1081  // e->cd(1);
1082  // fMvdHisto2[0]->Draw("colz");
1083  // fMvdHisto2[0]->GetXaxis()->SetTitle("Distance in x MC to Pix [Pixelpitch]");
1084  // fMvdHisto2[0]->GetYaxis()->SetTitle("Distance in y MC to Pix [Pixelpitch]");
1085  // fMvdHisto2[0]->GetZaxis()->SetTitle("Cluster Multiplicity");
1086 
1087 
1088  //
1089 }
1090 // -------------------------------------------------------------------------
1091 
1092 //-----------------------------------------------------------------------------------------
CbmMvdClusterAna::fMvdHitsMatch
TClonesArray * fMvdHitsMatch
Definition: CbmMvdClusterAna.h:43
CbmMvdClusterAna::Init
InitStatus Init()
Definition: CbmMvdClusterAna.cxx:110
CbmMvdClusterAna::Finish
void Finish()
Definition: CbmMvdClusterAna.cxx:674
CbmMvdDigi::GetPixelY
Int_t GetPixelY()
Definition: CbmMvdDigi.cxx:143
CbmHit::GetZ
Double_t GetZ() const
Definition: CbmHit.h:70
CbmMatch::GetMatchedLink
const CbmLink & GetMatchedLink() const
Definition: CbmMatch.h:37
CbmMatch
Definition: CbmMatch.h:22
CbmMvdClusterAna.h
CbmPixelHit::GetX
Double_t GetX() const
Definition: CbmPixelHit.h:83
CbmPixelHit::GetY
Double_t GetY() const
Definition: CbmPixelHit.h:84
CbmMvdDigi.h
CbmMvdClusterAna::CbmMvdClusterAna
CbmMvdClusterAna()
Definition: CbmMvdClusterAna.cxx:53
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
CbmMvdCluster::GetPixelMap
std::map< std::pair< Int_t, Int_t >, Int_t > GetPixelMap()
Definition: CbmMvdCluster.h:47
CbmMvdClusterAna::fMvdHits
TClonesArray * fMvdHits
Definition: CbmMvdClusterAna.h:39
CbmMvdClusterAna::fNrMcPointsAll
int fNrMcPointsAll
Definition: CbmMvdClusterAna.h:53
GetRMS
static Double_t GetRMS(const char *name)
Definition: GenNoiseElectrons.cxx:34
CbmMvdCluster
Definition: CbmMvdCluster.h:27
CbmMvdHit::GetIndexCentralX
Int_t GetIndexCentralX() const
Definition: CbmMvdHit.h:64
CbmMatch.h
CbmMvdHit
Definition: CbmMvdHit.h:29
CbmMvdClusterAna::fMvdHisto2
TH2F * fMvdHisto2[50]
Definition: CbmMvdClusterAna.h:50
CbmCluster::GetNofDigis
Int_t GetNofDigis() const
Number of digis in cluster.
Definition: CbmCluster.h:69
CbmMvdClusterAna::fMvdClusters
TClonesArray * fMvdClusters
Definition: CbmMvdClusterAna.h:38
CbmMvdClusterAna::fMvdClustersMatch
TClonesArray * fMvdClustersMatch
Definition: CbmMvdClusterAna.h:42
CbmMvdPoint.h
CbmMvdDigi::GetPixelSizeX
Double_t GetPixelSizeX()
Definition: CbmMvdDigi.h:52
CbmMvdClusterAna::fMvdDigis
TClonesArray * fMvdDigis
Definition: CbmMvdClusterAna.h:37
CbmMvdPoint
Definition: CbmMvdPoint.h:28
CbmMvdClusterAna::fMvdHisto1
TH1F * fMvdHisto1[50]
Definition: CbmMvdClusterAna.h:49
CbmMvdClusterAna::fMvdDigisMatch
TClonesArray * fMvdDigisMatch
Definition: CbmMvdClusterAna.h:41
d
double d
Definition: P4_F64vec2.h:24
CbmMvdClusterAna::fMcPoints
TClonesArray * fMcPoints
Definition: CbmMvdClusterAna.h:36
CbmTrackMatchNew.h
ClassImp
ClassImp(CbmMvdClusterAna)
CbmMvdClusterAna
Definition: CbmMvdClusterAna.h:24
CbmMvdHit::GetIndexCentralY
Int_t GetIndexCentralY() const
Definition: CbmMvdHit.h:67
CbmMvdClusterAna::~CbmMvdClusterAna
~CbmMvdClusterAna()
Definition: CbmMvdClusterAna.cxx:105
CbmMvdPoint::GetZOut
Double_t GetZOut() const
Definition: CbmMvdPoint.h:72
CbmMvdClusterAna::fNrHitsAll
int fNrHitsAll
Definition: CbmMvdClusterAna.h:54
CbmMvdHitMatch
Definition: CbmMvdHitMatch.h:17
CbmMCTrack.h
shape
UInt_t shape
Definition: CbmMvdSensorDigiToHitTask.cxx:73
CbmMvdHitMatch.h
CbmMCTrack
Definition: CbmMCTrack.h:34
CbmMvdPoint::GetXOut
Double_t GetXOut() const
Definition: CbmMvdPoint.h:70
CbmMvdHit.h
CbmMvdDigi::GetPixelSizeY
Double_t GetPixelSizeY()
Definition: CbmMvdDigi.h:53
CbmMvdPoint::GetYOut
Double_t GetYOut() const
Definition: CbmMvdPoint.h:71
CbmMvdDigi
Definition: CbmMvdDigi.h:21
CbmMvdClusterAna::fListMCTracks
TClonesArray * fListMCTracks
Definition: CbmMvdClusterAna.h:45
CbmMvdDigi::GetPixelX
Int_t GetPixelX()
Definition: CbmMvdDigi.cxx:141
CbmMvdClusterAna::fMcperHit
Int_t fMcperHit[11]
Definition: CbmMvdClusterAna.h:57
CbmMvdClusterAna::Exec
void Exec(Option_t *opt)
Definition: CbmMvdClusterAna.cxx:244
CbmMCTrack::GetP
Double_t GetP() const
Definition: CbmMCTrack.h:100
CbmMvdClusterAna::fMcperDigi
Int_t fMcperDigi[6]
Definition: CbmMvdClusterAna.h:56
CbmMvdClusterAna::fProf
TProfile * fProf[50]
Definition: CbmMvdClusterAna.h:51
CbmMvdClusterAna::fPixelpitch
Double_t fPixelpitch[2]
Definition: CbmMvdClusterAna.h:59
CbmMvdCluster.h