BmnRoot 26.10.0-18-g4e3414f9
Software framework of the BM@N experiment
Loading...
Searching...
No Matches
EmbEffView.C
Go to the documentation of this file.
1/**
2 * Macro for tracks versions comparison
3 */
4R__ADD_INCLUDE_PATH($VMCWORKDIR)
5
6#if !defined(__CLING__) || defined(__MAKECLING__)
7
8#include <BmnGlobalTrack.h>
9#include <BmnHit.h>
10#include <CbmMCTrack.h>
11#include <DstEventHeader.h>
12#include <FairMCEventHeader.h>
13#include <FairMCPoint.h>
14#include <FairTrackParam.h>
15#include <Rtypes.h>
16#include <TCanvas.h>
17#include <TChain.h>
18#include <TClass.h>
19#include <TClonesArray.h>
20#include <TDatime.h>
21#include <TF1.h>
22#include <TFile.h>
23#include <TH1.h>
24#include <TH2.h>
25#include <TMath.h>
26#include <TPad.h>
27#include <TROOT.h>
28#include <TString.h>
29#include <TStyle.h>
30#include <TTree.h>
31#include <TVirtualPad.h>
34#include <ncurses.h>
36// #include <BmnHist.h>
37#include <physics/bmnrecotools/BmnRecoTools.h>
38
39#else
40#include "bmnloadlibs.C"
41#include "gconfig/basiclibs.C"
42
43#endif
44
45#define PAD_WIDTH_EMB 1920
46#define PAD_HEIGHT_EMB 1200
47// #define PAD_WIDTH_EMB 960
48// #define PAD_HEIGHT_EMB 540
49#define ROWS 2
50#define COLS 2
51#define EV_LIMIT 4
52using namespace std;
53using namespace TMath;
54
56{
57#if !(!defined(__CLING__) || defined(__MAKECLING__))
58 bmnloadlibs(); // load BmnRoot libraries
59#endif
60 gStyle->SetOptStat(0);
61 const Int_t nColors = 16;
62 Int_t ColorMap[nColors] = {1, 2, 3, 4, 5, 6, 7, 8, 9, 49, 11, 12, 13, 14, 15, 16};
63 Int_t runID = 4649;
64 Int_t fPeriodID = 7;
65 const Int_t boxWidth = 40;
66 const Int_t boxHeigth = 200;
67 const Int_t Bins = 400;
68 // TString inFileNameEve = "eve-p-nofield.root";
69 TString inFileNameBase = Form("dst-%d.root", runID); // Form("cbmdst_%d.root", runID);
70 // TString inFileNameEmb = "dst-p-nofield.root"; //Form("dst_MK_%d.root", runID); //~/BMNDST_%d.root
71 TString inFileNameEve = Form("eve-%d-isAll.root", runID);
72 // TString inFileNameBase = Form("dst-%d.root", runID); //Form("cbmdst_%d.root", runID);
73 // TString inFileNameEmb = Form("dst-p-nofield.root", runID); //Form("dst_MK_%d.root", runID);
74 // //~/BMNDST_%d.root TString inFileNameEmb = Form("dst-%d-eve-qgsm-4h-unreal.root", runID);
75 // //Form("dst_MK_%d.root", runID); //~/BMNDST_%d.root TString inFileNameEmb =
76 // Form("dst-%d-emb-qgsm-4h-isAll.root", runID); //Form("dst_MK_%d.root", runID); //~/BMNDST_%d.root
77 TString inFileNameEmb =
78 Form("dst-%d-emb-qgsm-4h-isAll-BaseOff.root", runID); // Form("dst_MK_%d.root", runID); //~/BMNDST_%d.root
79
80 TString fnames[COLS] = {inFileNameEve, inFileNameEmb};
81 TString treeNames[COLS] = {"bmndata", "bmndata"}; //{"BMN_DIGIT", "cbmsim", "cbmsim"};
82 TString EHNames[COLS] = {"MCEventHeader.", "DstEventHeader."};
83 TString hitNames[COLS] = {"StsPoint", "BmnInnerHits"};
84 TString hitNamesGEM[COLS] = {"StsPoint", "BmnGemStripHit"};
85 TString hitNamesSil[COLS] = {"SiliconPoint", "BmnSiliconHit"};
86 TString hitNamesCSC[COLS] = {"CSCPoint", "BmnCSCHit"};
87 TString trackNames[COLS] = {"MCTrack", "BmnGlobalTrack"};
88 TFile* files[COLS] = {NULL};
89 TTree* trees[COLS] = {NULL};
90 FairMCEventHeader* MCEH = NULL;
91 DstEventHeader* EHs[COLS] = {NULL};
92 TClonesArray* hits[COLS] = {NULL};
93 TClonesArray* tracks[COLS] = {NULL};
94 TClonesArray* silTrack[COLS] = {NULL};
95 TClonesArray* gemTrack[COLS] = {NULL};
96 TClonesArray* silHits[COLS] = {NULL};
97 TClonesArray* cscHits[COLS] = {NULL};
98 TClonesArray* gemHits[COLS] = {NULL};
99 TH2D* hHitMatch[COLS][ROWS] = {NULL};
100 hHitMatch[0][0] = new TH2D("hHM-XZ-MC", "Hits XZ MC", Bins, -boxWidth, boxWidth, Bins, 0, boxHeigth);
101 hHitMatch[0][1] = new TH2D("hHM-YZ-MC", "Hits YZ MC", Bins, -boxWidth, boxWidth, Bins, 0, boxHeigth);
102 hHitMatch[1][0] = new TH2D("hHM-XZ-Emb", "Hits XZ Emb", Bins, -boxWidth, boxWidth, Bins, 0, boxHeigth);
103 hHitMatch[1][1] = new TH2D("hHM-YZ-Emb", "Hits YZ Emb", Bins, -boxWidth, boxWidth, Bins, 0, boxHeigth);
104 // TH2D * h[COLS][ROWS][nColors] = {NULL};
105 // for (Int_t ic = 0; ic < nColors; ic++) {
106 // h[0][0][ic] = new TH2D(Form("h0-0-%d", ic), "Hits XZ MC", Bins, -boxWidth, boxWidth, Bins, 0, boxHeigth);
107 // h[0][1][ic] = new TH2D(Form("h0-1-%d", ic), "Hits YZ MC", Bins, -boxWidth, boxWidth, Bins, 0, boxHeigth);
108 // h[1][0][ic] = new TH2D(Form("h1-0-%d", ic), "Hits XZ Base", Bins, -boxWidth, boxWidth, Bins, 0,
109 // boxHeigth); h[1][1][ic] = new TH2D(Form("h1-1-%d", ic), "Hits YZ Base", Bins, -boxWidth, boxWidth, Bins,
110 // 0, boxHeigth); h[2][0][ic] = new TH2D(Form("h2-0-%d", ic), "Hits XZ Embedded", Bins, -boxWidth, boxWidth,
111 // Bins, 0, boxHeigth); h[2][1][ic] = new TH2D(Form("h2-1-%d", ic), "Hits YZ Embedded", Bins, -boxWidth,
112 // boxWidth, Bins, 0, boxHeigth);
113 // }
114 // TCanvas* can = new TCanvas("canHits", "", PAD_WIDTH_EMB * COLS, PAD_HEIGHT_EMB * ROWS);
115 // can->Divide(COLS, ROWS);
116 TCanvas* can = new TCanvas("canHits", "", PAD_WIDTH_EMB * 1, PAD_HEIGHT_EMB * ROWS);
117 can->Divide(1, ROWS);
118
119 /** HitMaker efficiency test */
120 /** open from xml*/
121 TString xmlConfFileName = fPeriodID == 7 ? "SiliconRunSpring2018.xml" : "SiliconRunSpring2017.xml";
122 xmlConfFileName = TString(getenv("VMCWORKDIR")) + "/parameters/silicon/XMLConfigs/" + xmlConfFileName;
123 printf("xmlConfFileName %s\n", xmlConfFileName.Data());
124 BmnSiliconStationSet* silStationSet = new BmnSiliconStationSet(xmlConfFileName);
125
126 xmlConfFileName = fPeriodID == 7 ? "GemRunSpring2018.xml" : "GemRunSpring2017.xml";
127 xmlConfFileName = TString(getenv("VMCWORKDIR")) + "/parameters/gem/XMLConfigs/" + xmlConfFileName;
128 printf("xmlConfFileName %s\n", xmlConfFileName.Data());
129 BmnGemStripStationSet* gemStationSet = new BmnGemStripStationSet(xmlConfFileName);
130
131 xmlConfFileName = fPeriodID == 7 ? "CSCRunSpring2018.xml" : "CSCRunSpring2017.xml";
132 xmlConfFileName = TString(getenv("VMCWORKDIR")) + "/parameters/csc/XMLConfigs/" + xmlConfFileName;
133 printf("xmlConfFileName %s\n", xmlConfFileName.Data());
134 BmnCSCStationSet* cscStationSet = new BmnCSCStationSet(xmlConfFileName);
135
136 /** create hitCount vectors */
137 vector<vector<vector<UInt_t>>> hitVectorSil;
138 vector<vector<vector<UInt_t>>> ptVectorSil;
139 vector<vector<vector<Double_t>>> effVectorSil;
140 for (Int_t iStation = 0; iStation < silStationSet->GetNStations(); iStation++) {
141 BmnSiliconStation* st = silStationSet->GetSiliconStation(iStation);
142 vector<vector<UInt_t>> row;
143 row.resize(st->GetNModules(), vector<UInt_t>(2, 0));
144 hitVectorSil.push_back(row);
145 ptVectorSil.push_back(row);
146 vector<vector<Double_t>> rowDbl;
147 rowDbl.resize(st->GetNModules(), vector<Double_t>(2, 0.0));
148 effVectorSil.push_back(rowDbl);
149 }
150
151 vector<vector<vector<UInt_t>>> hitVectorGem;
152 vector<vector<vector<UInt_t>>> ptVectorGem;
153 vector<vector<vector<Double_t>>> effVectorGem;
154 UInt_t sumMods = 0;
155 UInt_t maxLayers = 0;
156 vector<PadInfo*> canStripPads;
157 vector<PadInfo*> canStripPads2D;
158 vector<TH2D*> ptsHistGemSt;
159 vector<TH2D*> effHistGemSt;
160 vector<TH1D*> effHistGemStResX;
161 vector<TH1D*> effHistGemStResY;
162 vector<vector<vector<TH1D*>>> ptsHistGem;
163 vector<vector<vector<TH1D*>>> effHistGem;
164 vector<vector<vector<TH2D*>>> ptsHistGem2D;
165 vector<vector<vector<TH2D*>>> effHistGem2D;
166 const Int_t StBinX = 200;
167 const Int_t StBinY = 100;
168 const Double_t TanMaxX = 1.5;
169 const Double_t TanMaxY = 0.8;
170 const Double_t ResMax = 1;
171 for (Int_t iStation = 0; iStation < gemStationSet->GetNStations(); iStation++) {
172 BmnGemStripStation* st = gemStationSet->GetGemStation(iStation);
173 sumMods += st->GetNModules();
174 vector<vector<UInt_t>> row;
175 row.resize(st->GetNModules(), vector<UInt_t>(2, 0));
176 hitVectorGem.push_back(row);
177 ptVectorGem.push_back(row);
178 vector<vector<Double_t>> rowDbl;
179 rowDbl.resize(st->GetNModules(), vector<Double_t>(2, 0.0));
180 effVectorGem.push_back(rowDbl);
181 TString n = Form("h-pts-%d", iStation);
182 TH2D* hst = new TH2D(n.Data(), n.Data(), StBinX,
183 -TanMaxX, //-st->GetXMaxStation(),
184 TanMaxX, //-st->GetXMinStation(),
185 StBinY,
186 -TanMaxY, //,st->GetYMinStation(),
187 TanMaxY // st->GetYMaxStation()
188 );
189 hst->GetZaxis()->SetRangeUser(0, 1.2);
190 ptsHistGemSt.push_back(hst);
191
192 n = Form("h-eff-%d", iStation);
193 TH2D* hstEff = new TH2D(n.Data(), n.Data(), StBinX,
194 -TanMaxX, //-st->GetXMaxStation(),
195 TanMaxX, //-st->GetXMinStation(),
196 StBinY,
197 -TanMaxY, //,st->GetYMinStation(),
198 TanMaxY // st->GetYMaxStation()
199 );
200 hstEff->GetZaxis()->SetRangeUser(0, 1.2);
201 effHistGemSt.push_back(hstEff);
202
203 n = Form("h-resX-%d", iStation);
204 TH1D* hrx = new TH1D(n.Data(), n.Data(), StBinX, -ResMax, ResMax);
205 effHistGemStResX.push_back(hrx);
206
207 n = Form("h-resY-%d", iStation);
208 TH1D* hry = new TH1D(n.Data(), n.Data(), StBinX, -ResMax, ResMax);
209 effHistGemStResY.push_back(hry);
210
211 vector<vector<TH1D*>> hmods;
212 vector<vector<TH1D*>> hmodsEff;
213 vector<vector<TH2D*>> hmods2D;
214 vector<vector<TH2D*>> hmodsEff2D;
215 for (Int_t iMod = 0; iMod < st->GetNModules(); iMod++) {
216 BmnGemStripModule* mod = st->GetModule(iMod);
217 if (maxLayers < mod->GetNStripLayers() / 2)
218 maxLayers = mod->GetNStripLayers() / 2;
219 vector<TH1D*> hzones;
220 vector<TH1D*> hzonesEff;
221 vector<TH2D*> hzones2D;
222 vector<TH2D*> hzonesEff2D;
223 for (Int_t iZone = 0; iZone < mod->GetNStripLayers() / 2; iZone++) {
224 BmnGemStripLayer l = mod->GetStripLayer(2 * iZone);
225
226 TString name = Form("h-pts-%d-%d-%d", iStation, iMod, iZone);
227 TH2D* h2d = new TH2D(name.Data(), name.Data(), l.GetNStrips() / 4, -l.GetXMaxLayer(), -l.GetXMinLayer(),
228 l.GetNStrips() / 16, l.GetYMinLayer(), l.GetYMaxLayer());
229 h2d->GetZaxis()->SetRangeUser(0, 1.2);
230
231 name = name + "-projection";
232 TH1D* h = new TH1D(name.Data(), name.Data(), l.GetNStrips() / 4, -l.GetXMaxLayer(), -l.GetXMinLayer());
233 h->GetYaxis()->SetRangeUser(0, 1.2);
234 hzones.push_back(h);
235 hzones2D.push_back(h2d);
236
237 name = Form("h-eff-%d-%d-%d", iStation, iMod, iZone);
238 TH2D* heff2d = new TH2D(name.Data(), name.Data(), l.GetNStrips() / 4, -l.GetXMaxLayer(),
239 -l.GetXMinLayer(), l.GetNStrips() / 16, l.GetYMinLayer(), l.GetYMaxLayer());
240 heff2d->GetZaxis()->SetRangeUser(0, 1.2);
241
242 name = name + "-projection";
243 TH1D* heff =
244 new TH1D(name.Data(), name.Data(), l.GetNStrips() / 4, -l.GetXMaxLayer(), -l.GetXMinLayer());
245 hzonesEff.push_back(heff);
246 hzonesEff2D.push_back(heff2d);
247 }
248 hmods2D.push_back(hzones2D);
249 hmodsEff2D.push_back(hzonesEff2D);
250 hmods.push_back(hzones);
251 hmodsEff.push_back(hzonesEff);
252 }
253 ptsHistGem2D.push_back(hmods2D);
254 effHistGem2D.push_back(hmodsEff2D);
255 ptsHistGem.push_back(hmods);
256 effHistGem.push_back(hmodsEff);
257 }
258
259 TString name = "Gem_Eff_Canvas";
260 TCanvas* canGemStrip = new TCanvas(name, name, PAD_WIDTH_EMB * maxLayers, PAD_HEIGHT_EMB * sumMods);
261 name = "Gem_Eff_Canvas2D";
262 TCanvas* canGemStrip2D = new TCanvas(name, name, PAD_WIDTH_EMB * maxLayers, PAD_HEIGHT_EMB * sumMods);
263 canGemStrip->Divide(maxLayers, sumMods);
264 canGemStrip2D->Divide(maxLayers, sumMods);
265 canStripPads.resize(sumMods * maxLayers);
266 canStripPads2D.resize(sumMods * maxLayers);
267 Int_t modCtr = 0;
268 for (Int_t iStation = 0; iStation < gemStationSet->GetNStations(); iStation++) {
269 BmnGemStripStation* st = gemStationSet->GetGemStation(iStation);
270 for (Int_t iModule = 0; iModule < st->GetNModules(); iModule++) {
271 BmnGemStripModule* mod = st->GetModule(iModule);
272 for (Int_t iZone = 0; iZone < mod->GetNStripLayers() / 2; iZone++) {
273 PadInfo* p = new PadInfo();
274 p->opt = "";
275 p->current = effHistGem[iStation][iModule][iZone];
276 PadInfo* p2D = new PadInfo();
277 p2D->opt = "colz";
278 p2D->current = effHistGem2D[iStation][iModule][iZone];
279 Int_t iPad = modCtr * maxLayers + iZone;
280 canStripPads[iPad] = p;
281 canGemStrip->GetPad(iPad + 1)->SetGrid();
282 canStripPads2D[iPad] = p2D;
283 canGemStrip2D->GetPad(iPad + 1)->SetGrid();
284 }
285 modCtr++;
286 }
287 }
288
289 vector<vector<vector<UInt_t>>> hitVectorCSC;
290 vector<vector<vector<UInt_t>>> ptVectorCSC;
291 vector<vector<vector<Double_t>>> effVectorCSC;
292 for (Int_t iStation = 0; iStation < cscStationSet->GetNStations(); iStation++) {
293 BmnCSCStation* st = cscStationSet->GetStation(iStation);
294 vector<vector<UInt_t>> row;
295 row.resize(st->GetNModules(), vector<UInt_t>(2, 0));
296 hitVectorCSC.push_back(row);
297 ptVectorCSC.push_back(row);
298 vector<vector<Double_t>> rowDbl;
299 rowDbl.resize(st->GetNModules(), vector<Double_t>(2, 0.0));
300 effVectorCSC.push_back(rowDbl);
301 }
302
303 /** Open exp files */
304 for (Int_t i = 0; i < COLS; i++) {
305 files[i] = new TFile(fnames[i], "READ");
306 if (files[i]->IsOpen())
307 printf("file %s opened\n", fnames[i].Data());
308 else {
309 fprintf(stderr, "file %s open error\n", fnames[i].Data());
310 continue;
311 // return;
312 }
313 trees[i] = (TTree*)files[i]->Get(treeNames[i]);
314 trees[i]->SetBranchAddress(hitNames[i], &hits[i]);
315 trees[i]->SetBranchAddress(trackNames[i], &tracks[i]);
316 for (Int_t j = 0; j < ROWS; j++) {
317 hHitMatch[i][j]->SetMarkerStyle(i > 0 ? 27 : 43);
318 hHitMatch[i][j]->SetMarkerSize(i > 0 ? 3 : 3);
319 hHitMatch[i][j]->SetMarkerColor(ColorMap[i]);
320 // for (Int_t ic = 0; ic < nColors; ic++) {
321 // h[i][j][ic]->SetMarkerStyle(ic > 0 ? 41 : 2);
322 // h[i][j][ic]->SetMarkerSize(ic > 0 ? 2 : 1);
323 // h[i][j][ic]->SetMarkerColor(ColorMap[ic]);
324 // }
325 }
326 if (i) {
327 trees[i]->SetBranchAddress(EHNames[i], &EHs[i]);
328 trees[i]->SetBranchAddress("BmnSiliconTrack", &silTrack[i]);
329 trees[i]->SetBranchAddress("BmnGemTrack", &gemTrack[i]);
330 trees[i]->SetBranchAddress(hitNamesSil[i], &silHits[i]);
331 trees[i]->SetBranchAddress(hitNamesCSC[i], &cscHits[i]);
332 trees[i]->SetBranchAddress(hitNamesGEM[i], &gemHits[i]);
333 } else {
334 // trees[i]->SetBranchAddress("BmnSiliconTrack", &silTrack[i]);
335 // trees[i]->SetBranchAddress("BmnGemTrack", &gemTrack[i]);
336 trees[i]->SetBranchAddress(EHNames[i], &MCEH);
337 trees[i]->SetBranchAddress(hitNamesSil[i], &silHits[i]);
338 trees[i]->SetBranchAddress(hitNamesCSC[i], &cscHits[i]);
339 trees[i]->SetBranchAddress(hitNamesGEM[i], &gemHits[i]);
340 }
341 }
342 const UInt_t nEvs = trees[0]->GetEntries();
343
344 /** iterate over events for efficiency test */
345 printf("Gathering hits<->points correspondence:\n");
346 Long64_t iExpEvent = -1;
347 for (Int_t iEv = 0; iEv < nEvs; iEv++) {
348 DrawBar(iEv, nEvs);
349 trees[0]->GetEntry(iEv);
350 UInt_t EventId = MCEH->GetEventID();
351 // printf("iEv = %d EventId = %u \n", iEv, EventId);
352 BmnRecoTools::FillSetStsPoints(gemHits[0], gemStationSet, (vector<vector<vector<TH2*>>>&)ptsHistGem2D);
353 // BmnRecoTools::FillSetStsPoints(gemHits[0], gemStationSet, (vector<TH2* > &)ptsHistGemSt);
354 // FillSetStsPoints(gemHits[0], gemStationSet, ptVectorGem);
355 UInt_t i = 1;
356 while (iExpEvent < 1000000) {
357 // printf("iExpEvent %lld\n", iExpEvent);
358 trees[i]->GetEntry(++iExpEvent);
359 // printf("EHs[%lld]->GetEventId() = %u \n", iExpEvent, EHs[i]->GetEventId());
360 if (EventId == EHs[i]->GetEventId())
361 break;
362 if (EventId < EHs[i]->GetEventId()) {
363 printf("\n%d tree %u event lost!\n", i, EventId);
364 break;
365 }
366 }
367 BmnRecoTools::FillSetStsHits(gemHits[0], gemHits[i], gemStationSet,
368 (vector<vector<vector<TH2*>>>&)effHistGem2D);
369 // BmnRecoTools::FillSetStsHits(gemHits[0], gemHits[i], gemStationSet, (vector<TH2* > &)effHistGemSt,
370 // effHistGemStResX,effHistGemStResY);
371 // FillSetStsHits(gemHits[0], gemHits[i], gemStationSet, hitVectorGem);
372 // FillSetHits(gemHits[i], hitVectorGem);
373 // FillSetHits(silHits[i], hitVectorSil);
374 // FillSetHits(cscHits[i], hitVectorCSC);
375 }
376
377 // printf("\nHit Maker Efficiency\n");
378 // for (Int_t iStation = 0; iStation < gemStationSet->GetNStations(); iStation++) {
379 // printf("GEM Station = %d\n", iStation);
380 // BmnGemStripStation* st = gemStationSet->GetGemStation(iStation);
381 // for (Int_t iMod = 0; iMod < st->GetNModules(); iMod++) {
382 // BmnGemStripModule *mod = st->GetModule(iMod);
383 // for (Int_t iZone = 0; iZone < mod->GetNStripLayers() / 2; iZone++) {
384 // if (ptVectorGem[iStation][iMod][iZone] > 0)
385 // effVectorGem[iStation][iMod][iZone] =
386 // hitVectorGem[iStation][iMod][iZone] / (Double_t) ptVectorGem[iStation][iMod][iZone];
387 // printf("\tiMod = %d, Zone = %d Eff = %f\n", iMod, iZone, effVectorGem[iStation][iMod][iZone]);
388 // }
389 // }
390 // }
391 printf("\nHit Maker Efficiency BigZone HotZone \n");
392 for (Int_t iStation = 0; iStation < gemStationSet->GetNStations(); iStation++) {
393 printf("GEM Station = %d\n", iStation);
394 BmnGemStripStation* st = gemStationSet->GetGemStation(iStation);
395 effHistGemSt[iStation]->Divide(ptsHistGemSt[iStation]);
396 TCanvas* cc = new TCanvas("temp", "temp", 1800, 900);
397 effHistGemSt[iStation]->Draw("colz");
398 cc->SaveAs(Form("eff-lam-corr-station_%d.png", iStation));
399 delete cc;
400 // TCanvas* cr = new TCanvas("temp", "temp", 1800, 900);
401 // cr->Divide(2,1);
402 // cr->cd(1);
403 // effHistGemStResX[iStation]->Draw("");
404 // cr->cd(2);
405 // effHistGemStResY[iStation]->Draw("");
406 // TString residName = "eff-lam-tan-resid.pdf";
407 // if (iStation< gemStationSet->GetNStations() - 1){
408 // if (iStation == 0)
409 // residName+="(";
410 // else
411 // ;
412 // }
413 // else
414 // residName+=")";
415 // cr->Print(residName.Data());
416 // delete cr;
417 for (Int_t iMod = 0; iMod < st->GetNModules(); iMod++) {
418 BmnGemStripModule* mod = st->GetModule(iMod);
419 printf("\tiMod = %d\n", iMod);
420 for (Int_t iZone = 0; iZone < mod->GetNStripLayers() / 2; iZone++) {
421 // printf("\tiZone = %d\n", iZone);
422 Double_t ptInt = ptsHistGem2D[iStation][iMod][iZone]->Integral();
423 Double_t eff = 0.0;
424 if (ptInt > 0)
425 eff = effHistGem2D[iStation][iMod][iZone]->Integral() / ptInt;
426 if (iZone != 0)
427 printf(" & ");
428 printf(" %f", eff);
429 // effHistGem[iStation][iMod][iZone] =
430 // effHistGem2D[iStation][iMod][iZone]->ProjectionX(); ptsHistGem[iStation][iMod][iZone]
431 // = ptsHistGem2D[iStation][iMod][iZone]->ProjectionX();
432 effHistGem[iStation][iMod][iZone]->Divide(ptsHistGem[iStation][iMod][iZone]);
433 }
434 printf("\n");
435 }
436 }
437 // BmnHist::DrawRef(canGemStrip, &canStripPads);
438 // BmnHist::DrawRef(canGemStrip2D, &canStripPads2D);
439 //// canGemStrip->Print("eff-p.pdf");
440 //// canGemStrip->Print("eff-p.eps");
441 // canGemStrip->Print("eff-p.png");
442 // canGemStrip2D->Print("eff-lam-2d.png");
443
444 // printf("MC Points\n");
445 // for (Int_t iStation = 0; iStation < gemStationSet->GetNStations(); iStation++) {
446 // printf("GEM Station = %d\n", iStation);
447 // BmnGemStripStation* st = gemStationSet->GetGemStation(iStation);
448 // for (Int_t iMod = 0; iMod < st->GetNModules(); iMod++) {
449 // BmnGemStripModule *mod = st->GetModule(iMod);
450 // for (Int_t iZone = 0; iZone < mod->GetNStripLayers() / 2; iZone++) {
451 // printf("\tiMod = %d, Zone = %d nPoints = %d\n", iMod, iZone, ptVectorGem[iStation][iMod][iZone]);
452 // }
453 // }
454 // }
455 // printf("Exp Hits\n");
456 // for (Int_t iStation = 0; iStation < gemStationSet->GetNStations(); iStation++) {
457 // BmnGemStripStation* st = gemStationSet->GetGemStation(iStation);
458 // printf("GEM Station = %d\n", iStation);
459 // for (Int_t iMod = 0; iMod < st->GetNModules(); iMod++) {
460 // BmnGemStripModule *mod = st->GetModule(iMod);
461 // for (Int_t iZone = 0; iZone < mod->GetNStripLayers() / 2; iZone++) {
462 // printf("\tiMod = %d, Zone = %d nHits = %u\n", iMod, iZone, hitVectorGem[iStation][iMod][iZone]);
463 // }
464 // }
465 // }
466
467 // for (Int_t iStation = 0; iStation < silStationSet->GetNStations(); iStation++) {
468 // BmnSiliconStation* st = silStationSet->GetStation(iStation);
469 // printf("Sil Station = %d\n", iStation);
470 // for (Int_t iMod = 0; iMod < st->GetNModules(); iMod++) {
471 // printf("\tiMod = %d, nHits = %u\n", iMod, hitVectorSil[iStation][iMod]);
472 // }
473 // }
474 // for (Int_t iStation = 0; iStation < cscStationSet->GetNStations(); iStation++) {
475 // BmnCSCStation* st = cscStationSet->GetStation(iStation);
476 // printf("CSC Station = %d\n", iStation);
477 // for (Int_t iMod = 0; iMod < st->GetNModules(); iMod++) {
478 // printf("\tiMod = %d, nHits = %u\n", iMod, hitVectorCSC[iStation][iMod]);
479 // }
480 // }
481 return;
482
483 /** iterate interactively over events in 3 files */
484 cout << "#recorded entries = " << nEvs << " in " << fnames[0] << endl;
485 Int_t ch = 0;
486 Int_t iEv = 0;
487 do {
488 ch = getchar();
489 switch (ch) {
490 case 'a':
491 iEv = (iEv == 0) ? iEv : iEv - 1;
492 break;
493 case 'd':
494 iEv = (iEv < nEvs) ? iEv + 1 : iEv;
495 break;
496 default:
497 // iEv++;
498 break;
499 }
500 if (ch == '\n')
501 continue;
502 Int_t hitCount[COLS] = {0};
503 trees[0]->GetEntry(iEv);
504 UInt_t EventId = MCEH->GetEventID();
505 printf("iEv = %d EventId = %u \n", iEv, EventId);
506
507 for (Int_t i = 1; i < COLS; i++) {
508 trees[i]->GetEntry(iEv);
509 // Long64_t iExpEvent = -1;
510 // while (iExpEvent < 1000000) {
511 // // printf("iExpEvent %lld\n", iExpEvent);
512 // trees[i]->GetEntry(++iExpEvent);
513 // // printf("trees[i]->GetEntry(iExpEvent) %d\n",
514 // trees[i]->GetEntry(iExpEvent));
515 // // printf("check %lld \n", iExpEvent);
516 // if (EventId == EHs[i]->GetEventId())
517 // break;
518 // if (EventId < EHs[i]->GetEventId()) {
519 // printf("\n%d tree %u event lost!\n", i, EventId);
520 // break;
521 // }
522 // }
523 }
524 for (Int_t i = 0; i < COLS; i++) {
525 for (Int_t row = 0; row < ROWS; row++) {
526 hHitMatch[i][row]->Reset();
527 // for (Int_t ic = 0; ic < nColors; ic++)
528 // h[i][row][ic]->Reset();
529 }
530 if (i) { // printf("%d th tree hits %d\n", i, gemHits[i]->GetEntriesFast());
531 for (Int_t iHit = 0; iHit < gemHits[i]->GetEntriesFast(); iHit++) {
532 BmnHit* hit = (BmnHit*)gemHits[i]->UncheckedAt(iHit);
533 // h[i][0][0]->Fill(hit->GetX(), hit->GetZ());
534 // h[i][1][0]->Fill(hit->GetY(), hit->GetZ());
535 hHitMatch[i][0]->Fill(hit->GetX(), hit->GetZ());
536 hHitMatch[i][1]->Fill(hit->GetY(), hit->GetZ());
537 }
538 } else { // printf("%d th tree pts %d\n", i, gemHits[i]->GetEntriesFast());
539 for (Int_t iHit = 0; iHit < gemHits[i]->GetEntriesFast(); iHit++) {
540 FairMCPoint* hit = (FairMCPoint*)gemHits[i]->UncheckedAt(iHit);
541 hHitMatch[i][0]->Fill(hit->GetX(), hit->GetZ());
542 hHitMatch[i][1]->Fill(hit->GetY(), hit->GetZ());
543 }
544 }
545
546 //
547 // // hitCount[i] += hits[i]->GetEntriesFast();
548 //
549 // // if (tracks[1]->GetEntriesFast() < 2) {
550 // // iEv++;
551 // // continue;
552 // // }
553 // printf("%d branch %d tracks\n", i, tracks[i]->GetEntriesFast());
554 // if (i > 0) {
555 // for (Int_t iTrack = 0; iTrack < tracks[i]->GetEntriesFast(); iTrack++) {
556 // if (iTrack > nColors - 2) {
557 // printf("tracks overfull %d\n", iTrack);
558 // continue;
559 // }
560 // BmnGlobalTrack* track = (BmnGlobalTrack*) tracks[i]->UncheckedAt(iTrack);
561 // if (track->GetSilTrackIndex() != -1) {
562 // BmnTrack * silTr = (BmnTrack*)
563 // silTrack[i]->UncheckedAt(track->GetSilTrackIndex());
564 // // if (silTr->GetNHits() < 4)
565 // // continue;
566 // hitCount[i] += silTr->GetNHits();
567 // for (Int_t iHit = 0; iHit < silTr->GetNHits(); iHit++) {
568 // BmnHit *hit = (BmnHit *) silHits[i]->UncheckedAt(silTr->GetHitIndex(iHit));
569 // h[i][0][iTrack + 1]->Fill(hit->GetX(), hit->GetZ());
570 // h[i][1][iTrack + 1]->Fill(hit->GetY(), hit->GetZ());
571 // }
572 // }
573 // if (track->GetGemTrackIndex() != -1) {
574 // BmnTrack * gemTr = (BmnTrack*)
575 // gemTrack[i]->UncheckedAt(track->GetGemTrackIndex());
576 // // if (gemTr->GetNHits() < 4)
577 // // continue;
578 // hitCount[i] += gemTr->GetNHits();
579 // for (Int_t iHit = 0; iHit < gemTr->GetNHits(); iHit++) {
580 // BmnHit *hit = (BmnHit *) gemHits[i]->UncheckedAt(gemTr->GetHitIndex(iHit));
581 // h[i][0][iTrack + 1]->Fill(hit->GetX(), hit->GetZ());
582 // h[i][1][iTrack + 1]->Fill(hit->GetY(), hit->GetZ());
583 // }
584 // }
585 // // FairTrackParam* first = track->GetParamFirst();
586 // // h[i][2][iTrack]->Fill(first->GetTx());
587 // // h[i][3][iTrack]->Fill(first->GetTy());
588 // }
589 // } else {
590 // for (Int_t iTrack = 0; iTrack < tracks[i]->GetEntriesFast(); iTrack++) {
591 // CbmMCTrack* track = (CbmMCTrack*) tracks[i]->UncheckedAt(iTrack);
592 // // hitCount[i] += track->GetNPoints(kGEM);
593 // for (Int_t iPoint = 0; iPoint < gemHits[i]->GetEntriesFast(); iPoint++) {
594 // FairMCPoint* hit = (FairMCPoint*) gemHits[i]->UncheckedAt(iPoint);
595 // Int_t TrackID = hit->GetTrackID();
596 // if (TrackID != iTrack)
597 // continue;
598 // h[i][0][iTrack + 1]->Fill(hit->GetX(), hit->GetZ());
599 // h[i][1][iTrack + 1]->Fill(hit->GetY(), hit->GetZ());
600 // }
601 // for (Int_t iPoint = 0; iPoint < silHits[i]->GetEntriesFast(); iPoint++) {
602 // FairMCPoint* hit = (FairMCPoint*) silHits[i]->UncheckedAt(iPoint);
603 // Int_t TrackID = hit->GetTrackID();
604 // if (TrackID != iTrack)
605 // continue;
606 // h[i][0][iTrack + 1]->Fill(hit->GetX(), hit->GetZ());
607 // h[i][1][iTrack + 1]->Fill(hit->GetY(), hit->GetZ());
608 // }
609 // // FairTrackParam* first = track->GetParamFirst();
610 // // h[i][2][iTrack]->Fill(first->GetTx());
611 // // h[i][3][iTrack]->Fill(first->GetTy());
612 //
613 // }
614 // }
615 }
616 // if (/*hitCount[0] < EV_LIMIT || */hitCount[1] < EV_LIMIT || hitCount[2] < EV_LIMIT) {
617 // iEv++;
618 // continue;
619 // }
620 // for (Int_t col = 0; col < COLS; col++)
621 // for (Int_t row = 0; row < ROWS; row++) {
622 // TVirtualPad *pad = can->cd(row * COLS + col + 1);
623 // pad->Clear();
624 // for (Int_t iTrack = 0; iTrack < nColors; iTrack++)
625 // if (iTrack)
626 // h[col][row][iTrack]->Draw("SAME");
627 // else
628 // h[col][row][iTrack]->Draw("");
629 // }
630 for (Int_t row = 0; row < ROWS; row++) {
631 TVirtualPad* pad = can->cd(row + 1);
632 pad->Clear();
633 for (Int_t col = 0; col < COLS; col++) {
634 if (col)
635 hHitMatch[col][row]->Draw("SAME");
636 else
637 hHitMatch[col][row]->Draw("");
638 }
639 }
640 can->Update();
641 can->Modified();
642 can->SaveAs(Form("can-%d.pdf", runID));
643 } while (ch != 'q');
644
645 for (Int_t i = 0; i < COLS; i++) {
646 for (Int_t j = 0; j < ROWS; j++) {
647 delete hHitMatch[i][j];
648 }
649 }
650
651 return;
652}
void DrawBar(UInt_t iEv, UInt_t nEv)
Definition BmnMath.cxx:940
TCanvas * can(nullptr)
int i
Definition P4_F32vec4.h:22
int Int_t
unsigned int UInt_t
double Double_t
long long Long64_t
#define R__ADD_INCLUDE_PATH(PATH)
#define COLS
Definition EmbEffView.C:50
char * Form(const char *fmt,...)
R__EXTERN TStyle * gStyle
BmnCSCStation * GetStation(Int_t station_num)
Int_t GetNModules()
Double_t GetYMinLayer()
Double_t GetYMaxLayer()
Double_t GetXMinLayer()
Double_t GetXMaxLayer()
BmnGemStripLayer & GetStripLayer(Int_t num)
BmnGemStripStation * GetGemStation(Int_t station_num)
BmnGemStripModule * GetModule(Int_t module_num)
static void FillSetStsPoints(TClonesArray *pts, BmnGemStripStationSet *set, vector< vector< vector< TH2 * > > > &hitVec)
static void FillSetStsHits(TClonesArray *pts, TClonesArray *hits, BmnGemStripStationSet *set, vector< vector< vector< TH2 * > > > &hitVec)
BmnSiliconStation * GetSiliconStation(Int_t station_num)
Storage for pad content and it's options.
Definition PadInfo.h:20
string opt
Definition PadInfo.h:88
TH1 * current
Definition PadInfo.h:78
virtual void SetMarkerColor(Color_t mcolor=1)
virtual void SetMarkerStyle(Style_t mstyle=1)
virtual void SetMarkerSize(Size_t msize=1)
virtual void SetRangeUser(Double_t ufirst, Double_t ulast)
TAxis * GetZaxis()
TAxis * GetYaxis()
void Draw(Option_t *option="") override
void Reset(Option_t *option="") override
virtual Int_t Fill(const char *namex, const char *namey, Double_t w)
Int_t GetEntriesFast() const
void Divide(Int_t nx=1, Int_t ny=1, Float_t xmargin=0.01, Float_t ymargin=0.01, Int_t color=0) override
void SaveAs(const char *filename="", Option_t *option="") const override
TVirtualPad * GetPad(Int_t subpadnumber) const override
const char * Data() const
virtual Int_t GetEntry(Long64_t entry, Int_t getall=0)
virtual Long64_t GetEntries() const
Int_t SetBranchAddress(const char *bname, T **add, TBranch **ptr=nullptr)
virtual void SetGrid(Int_t valuex=1, Int_t valuey=1)=0
void Clear(Option_t *option="") override=0
#define PAD_HEIGHT_EMB
Definition EmbEffView.C:46
#define ROWS
Definition EmbEffView.C:49
#define PAD_WIDTH_EMB
Definition EmbEffView.C:45
void EmbEffView()
Definition EmbEffView.C:55
*brief user defined string literal for JSON values *sa std::size_t n
Definition json.hpp:20362
Definition json.hpp:4446