StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
StFgtAVEfficiencyMaker.cxx
1 #include "StFgtAVEfficiencyMaker.h"
3 #include "StRoot/StEvent/StFgtCollection.h"
4 #include "StRoot/StEvent/StFgtHitCollection.h"
5 #include "StRoot/StEvent/StFgtHit.h"
6 #include "StRoot/StFgtUtil/geometry/StFgtGeom.h"
7 #include "StRoot/StEvent/StEvent.h"
8 #include "StRoot/StEvent/StEventInfo.h"
9 #include "StRoot/StFgtUtil/geometry/StFgtGeom.h"
10 #include <TH2D.h>
11 #include <TROOT.h>
12 #include <TStyle.h>
13 #include <TCanvas.h>
14 #include <utility>
15 #include <TArc.h>
16 #include <TLine.h>
17 #include <set>
18 //Double_t StFgtAVEfficiencyMaker::getRPhiRatio()
19 //{
20 
21 //};
22 Double_t getRPhiRatio(StSPtrVecFgtHitConstIterator hitIterBegin, StSPtrVecFgtHitConstIterator hitIterEnd)
23 {
24  Short_t quad, disc, strip;
25  Char_t layer;
26  Bool_t stripDead=false;
27  Int_t numR=0;
28  Int_t numPhi=0;
29 
30  StSPtrVecFgtHitConstIterator hitIter=hitIterBegin;
31  for(;hitIter!=hitIterEnd;hitIter++)
32  {
33  StFgtGeom::decodeGeoId((*hitIter)->getCentralStripGeoId(),disc, quad, layer, strip);
34  if(layer=='R')
35  numR++;
36  else
37  numPhi++;
38  }
39 
40  if(numR+numPhi>0)
41  return (numR-numPhi)/((Double_t)(numR+numPhi));
42  else
43  return -1;
44 };
45 
50 {
51  // cout <<" in eff make " <<endl;
52  StEvent* eventPtr = 0;
53  eventPtr = (StEvent*)GetInputDS("StEvent");
54 
55  Int_t ierr = StFgtQaMaker::Make();
56  Float_t x;
57  Float_t y;
58  Int_t prvGeoId=-1;
59 
60  StFgtHitCollection* clusterColD1=mFgtCollectionPtr->getHitCollection(0);
61  StFgtHitCollection* clusterColD6=mFgtCollectionPtr->getHitCollection(5);
62 
63  //look at the r phi ratio for each disk
64  for(int i=0;i<6;i++)
65  {
66  StFgtHitCollection* tmpClusterCol=mFgtCollectionPtr->getHitCollection(i);
67  Double_t ratio=getRPhiRatio(tmpClusterCol->getHitVec().begin(),tmpClusterCol->getHitVec().end());
68  rPhiRatioPlots[i]->Fill(ratio);
69  // cout << "ratio for disk: " << i << " is " << ratio <<" disk has: " << tmpClusterCol->getHitVec().size() << "hits" <<endl;
70  }
71 
72  const StSPtrVecFgtHit &hitVecD1=clusterColD1->getHitVec();
73  const StSPtrVecFgtHit &hitVecD6=clusterColD6->getHitVec();
74 
75 
76  set<Int_t> usedPoints;//saves the points that have been used for tracks (and shouldn't be reused)
77  Double_t D1Pos=StFgtGeom::getDiscZ(0);
78  Double_t D6Pos=StFgtGeom::getDiscZ(5);
79  Double_t zArm=D6Pos-D1Pos;
80  StSPtrVecFgtHitConstIterator hitIterD1,hitIterD6, hitIterD1R, hitIterD6R, hitIter, hitIter2;
81  cout <<" there are " << hitVecD1.size() << " hits in d1 "<<endl;
82  for(hitIterD1=hitVecD1.begin();hitIterD1 != hitVecD1.end();hitIterD1++)
83  {
84  Short_t quad, disc, strip;
85  Char_t layer;
86  Bool_t stripDead=false;
87  StFgtGeom::decodeGeoId((*hitIterD1)->getCentralStripGeoId(),disc, quad, layer, strip);
88  if(layer=='P')
89  cout <<"found phi d1 " <<endl;
90  else
91  cout <<"found r d1 " <<endl;
92  }
93  cout <<" there are " << hitVecD6.size() << " hits in d6 "<<endl;
94  for(hitIterD1=hitVecD6.begin();hitIterD1 != hitVecD6.end();hitIterD1++)
95  {
96  Short_t quad, disc, strip;
97  Char_t layer;
98  Bool_t stripDead=false;
99  StFgtGeom::decodeGeoId((*hitIterD1)->getCentralStripGeoId(),disc, quad, layer, strip);
100  if(layer=='P')
101  cout <<"found phi d6 " <<endl;
102  else
103  cout <<"found r d6 " <<endl;
104  }
105 
107  for(int iSeed1=0;iSeed1<5;iSeed1++)
108  {
109  for(int iSeed2=iSeed1+1;iSeed2<6;iSeed2++)
110  {
111  StFgtHitCollection* clusterColSeed1=mFgtCollectionPtr->getHitCollection(iSeed1);
112  StFgtHitCollection* clusterColSeed2=mFgtCollectionPtr->getHitCollection(iSeed2);
113  const StSPtrVecFgtHit &hitVecSeed1=clusterColSeed1->getHitVec();
114  const StSPtrVecFgtHit &hitVecSeed2=clusterColSeed2->getHitVec();
115  D1Pos=StFgtGeom::getDiscZ(iSeed1);
116  D6Pos=StFgtGeom::getDiscZ(iSeed2);
117  zArm=D6Pos-D1Pos;
118 
119  // cout <<"looking at " << hitVecD1.size() << " hits in D1 and " << hitVecD6.size() << " in D6 " <<endl;
120  // for(hitIterD1=hitVecD1.begin();hitIterD1 != hitVecD1.end();hitIterD1++)
121  for(hitIterD1=hitVecSeed1.begin();hitIterD1 != hitVecSeed1.end();hitIterD1++)
122  {
123  Short_t quad, disc, strip;
124  Char_t layer;
125  Bool_t stripDead=false;
126  StFgtGeom::decodeGeoId((*hitIterD1)->getCentralStripGeoId(),disc, quad, layer, strip);
128  if(quad>2)
129  continue;
130 
131  //do 1D 'fit' with r strips and the (x,y) thing
132  Int_t geoIdD1=(*hitIterD1)->getCentralStripGeoId();
133  StFgtGeom::decodeGeoId(geoIdD1,disc, quad, layer, strip);
134  if(layer!='P')
135  continue;
136 
137  Float_t phiD1=(*hitIterD1)->getPositionPhi();
138  for(hitIterD1R=hitVecSeed1.begin();hitIterD1R != hitVecSeed1.end();hitIterD1R++)
139  {
140  Int_t geoIdR=(*hitIterD1R)->getCentralStripGeoId();
141  StFgtGeom::decodeGeoId(geoIdR,disc, quad, layer, strip);
142  if(layer!='R')
143  continue;
144 
145  Float_t rD1=(*hitIterD1R)->getPositionR();
146  Float_t xD1=rD1*cos(phiD1);
147  Float_t yD1=rD1*sin(phiD1);
148  for(hitIterD6=hitVecSeed2.begin();hitIterD6 != hitVecSeed2.end();hitIterD6++)
149  {
150  Int_t geoIdD6=(*hitIterD6)->getCentralStripGeoId();
151  StFgtGeom::decodeGeoId(geoIdD6,disc, quad, layer, strip);
152  if(layer!='P')
153  continue;
154  Float_t phiD6=(*hitIterD6)->getPositionPhi();
155  for(hitIterD6R=hitVecSeed2.begin();hitIterD6R != hitVecSeed2.end();hitIterD6R++)
156  {
157  Int_t geoIdR=(*hitIterD6R)->getCentralStripGeoId();
158  StFgtGeom::decodeGeoId(geoIdR,disc, quad, layer, strip);
159  if(layer!='R')
160  continue;
161 
162  Float_t rD6=(*hitIterD6R)->getPositionR();
163  Double_t xD6=rD6*cos(phiD6);
164  Double_t yD6=rD6*sin(phiD6);
165 
167  vector< pair<Int_t,Double_t> > v_x;
168  vector< pair<Int_t,Double_t> > v_y;
169  vector< pair<Int_t,Double_t> > v_r;
170 
171  vector< pair<Int_t,Double_t> > v_xFail;
172  vector< pair<Int_t,Double_t> > v_yFail;
173  vector< pair<Int_t,Double_t> > v_rFail;
174 
175  vector< pair< Int_t, Double_t> > v_rCharge;
176  vector< pair< Int_t, Double_t> > v_phiCharge;
177 
178  vector<Int_t> v_clusterSizeR;
179  vector<Int_t> v_clusterSizePhi;
180 
181  vector<Int_t> v_geoIDsR;
182  vector<Int_t> v_geoIDsPhi;
183 
184  int iFound=0;
185  int iFoundR=0;
186 
187  //zarm is d6 position - D1
188  Double_t xIpExp=xD1+(xD6-xD1)*(-D1Pos)/zArm;
189  Double_t yIpExp=yD1+(yD6-yD1)*(-D1Pos)/zArm;
191  Double_t zIpExpX0=(D6Pos-(xD6/xD1)*D1Pos)*1/(1-xD6/xD1);
193  Double_t zIpExpY0=(D6Pos-yD6/yD1*D1Pos)*1/(1-yD6/yD1);
194  // cout <<"
195 
196  for(int iD=0;iD<6;iD++)
197  {
198  if(iD==iSeed1 || iD==iSeed2)
199  continue;
200 
201  Bool_t found=false;
202  Bool_t foundR=false;
203  //check for hit
204  Double_t diskZ=StFgtGeom::getDiscZ(iD);
205  //expected
206 
207  Double_t xPosExp=xD1+(xD6-xD1)*(diskZ-D1Pos)/zArm;
208  Double_t yPosExp=yD1+(yD6-yD1)*(diskZ-D1Pos)/zArm;
209  Double_t rPosExp=rD1+(rD6-rD1)*(diskZ-D1Pos)/zArm;
210 
211  //at z=0;
212 
213 
214  // cout <<"x1: " << xD1 << " y1: " << yD1 <<" x6: " << xD6 <<" y6: " << yD6 << " zpos: " << diskZ <<" arm: " << zArm<<endl;
215  // cout <<"expect hit at : " << xPosExp <<" / " << yPosExp <<" r: " << rPosExp <<endl;
216  StFgtHitCollection* clusterCol=mFgtCollectionPtr->getHitCollection(iD);
217  const StSPtrVecFgtHit &hitVec=clusterCol->getHitVec();
218  // cout <<" there are " << hitVec.size() << " hits in d" << iD <<endl;
219  for(hitIter=hitVec.begin();hitIter!=hitVec.end();hitIter++)
220  {
221  //do 1D 'fit' with r strips and the (x,y) thing
222  Int_t geoIdPhi=(*hitIter)->getCentralStripGeoId();
223  StFgtGeom::decodeGeoId(geoIdPhi,disc, quad, layer, strip);
224  if(usedPoints.find(geoIdPhi)!=usedPoints.end())
225  continue;
226  if(layer!='P')
227  continue;
228  Float_t phi=(*hitIter)->getPositionPhi();
229  Double_t phiCharge=(*hitIter)->charge();
230  Int_t clusterSizePhi=(*hitIter)->getStripWeightMap().size();
231  for(hitIter2=hitVec.begin();hitIter2!=hitVec.end();hitIter2++)
232  {
233  Int_t geoIdR=(*hitIter2)->getCentralStripGeoId();
234  StFgtGeom::decodeGeoId(geoIdR,disc, quad, layer, strip);
235  if(usedPoints.find(geoIdR)!=usedPoints.end())
236  continue;
237  if(layer!='R')
238  continue;
239  Float_t r=(*hitIter2)->getPositionR();
240  Double_t rCharge=(*hitIter2)->charge();
241  Int_t clusterSizeR=(*hitIter2)->getStripWeightMap().size();
242  // cout <<"checking with r:" << r <<endl;
243  if(fabs(r-rPosExp)<0.5)
244  {
245  foundR=true;
246  // cout <<"found r: " << r <<endl;
247  v_r.push_back(pair<Int_t,Double_t> (iD,r));
248  }
249 
250  x=r*cos(phi);
251  y=r*sin(phi);
252  // cout <<"checking with x: " << x << " y: " << y <<endl;
253  if(fabs(x-xPosExp) < 0.5 && fabs(y-yPosExp)<0.5) //found hit
254  {
255  found=true;
256  cout <<"found! " <<" pushing back: iD: " << iD << "x: " << x << " y "<< y <<endl;
257  v_x.push_back(pair<Int_t,Double_t>(iD,x));
258  v_y.push_back(pair<Int_t,Double_t>(iD,y));
259  v_rCharge.push_back(pair<Int_t, Double_t>(iD,rCharge));
260  v_phiCharge.push_back(pair<Int_t, Double_t>(iD,phiCharge));
261  //to avoid double counting
262  v_geoIDsPhi.push_back(geoIdPhi);
263  v_geoIDsR.push_back(geoIdR);
264  // cout <<" we have rcharge: " << rCharge<<" phi: " << phiCharge <<endl;
265  v_clusterSizeR.push_back(clusterSizeR);
266  v_clusterSizePhi.push_back(clusterSizePhi);
267  }
268  }
269  }
270  //only one per disk
271  if(found)
272  iFound++;
273  else
274  {
275  // cout <<"failed to find, pushing back " << xPosExp << " y: " << yPosExp <<endl;
276  v_xFail.push_back(pair<Int_t,Double_t>(iD,xPosExp));
277  v_yFail.push_back(pair<Int_t,Double_t>(iD,yPosExp));
278  }
279  if(foundR)
280  iFoundR++;
281  else
282  {
283  // cout <<"failed to find r " << rPosExp<<endl;
284  v_rFail.push_back(pair<Int_t, Double_t>(iD,rPosExp));
285  }
286 
287  }
288 
289  // cout << " Ifound: " << iFound <<endl;
290  if(iFound>=2)
291  {
292  hIp->Fill(xIpExp,yIpExp);
293  hIpZAtX0->Fill(zIpExpX0);
294  hIpZAtY0->Fill(zIpExpY0);
295  // cout <<" we have zIPX0: " << zIpExpX0 <<" and y: " << zIpExpY0 <<endl;
296  }
297 
298 
299  // if(iFound>=2 && fabs(xIpExp)<20 && fabs(yIpExp)<20 && fabs(zIpExpX0)<40 && fabs(zIpExpY0)<40) //at least one hit plus seed and pointing roughly to the interaction region
300  if(iFound>=2 && fabs(zIpExpX0)<100 && fabs(zIpExpY0)<100) //at least one hit plus seed and pointing roughly to the interaction region
301  {
302  if(v_x.size()>iFound)
303  {
304  cout<<"more hits than disks hit!!!" <<endl;
305  cout <<" vx_size: " << v_x.size() <<" ifound: " << iFound <<endl;
306  }
307  else
308  {
309 
310  // cout <<"filled hip with " << xIpExp << " / " << yIpExp <<endl;
311  for(int i=0;i<v_x.size();i++)
312  {
313  usedPoints.insert(v_geoIDsPhi[i]);
314  usedPoints.insert(v_geoIDsR[i]);
315  Int_t disk=v_x[i].first;
316  Double_t x=v_x[i].second;
317  Double_t y=v_y[i].second;
318  Double_t rCharge=v_rCharge[i].second;
319  Double_t phiCharge=v_phiCharge[i].second;
320  // cout <<"filling disk: " << disk <<" with: " << x <<" / " <<y <<endl;
321  radioPlotsEff[disk]->Fill(x,y);
322  chargeCorr[disk]->Fill(rCharge,phiCharge);
323  h_clusterSizeR[disk]->Fill(v_clusterSizeR[i]);
324  h_clusterSizePhi[disk]->Fill(v_clusterSizePhi[i]);
325  h_clusterChargeR[disk]->Fill(rCharge);
326  h_clusterChargePhi[disk]->Fill(phiCharge);
327  }
328  for(int i=0;i<v_xFail.size();i++)
329  {
330  Int_t disk=v_xFail[i].first;
331  Double_t x=v_xFail[i].second;
332  Double_t y=v_yFail[i].second;
333  // cout <<"filling disk fail: " << disk <<" with: " << x <<" / " <<y <<endl;
334  radioPlotsNonEff[disk]->Fill(x,y);
335  }
336  hitCounter++;
337  }
338  }
339 
340  if(iFoundR>=1) //at least one hit plus seed
341  {
342  if(v_r.size()>iFound)
343  {
344  // cout<<"more r hits than disks hit!!!" <<endl;
345  }
346  else
347  {
348  for(int i=0;i<v_r.size();i++)
349  {
350  Int_t disk=v_r[i].first;
351  Double_t r=v_r[i].second;
352  // cout <<"filling r disk: " << disk <<" with: " << r <<endl;
353  rEff[disk]->Fill(r);
354  }
355  for(int i=0;i<v_rFail.size();i++)
356  {
357  Int_t disk=v_rFail[i].first;
358  Double_t r=v_rFail[i].second;
359  // cout <<"filling r disk fail: " << disk <<" with: " << r <<endl;
360  rNonEff[disk]->Fill(r);
361  }
362  hitCounterR++;
363  }
364  }
365 
366 
367  //start over
368  iFound=0;
369  iFoundR=0;
370  v_x.clear();
371  v_y.clear();
372  v_r.clear();
373  v_xFail.clear();
374  v_yFail.clear();
375  v_rFail.clear();
376 
377  }
378  }
379  }
380  }
381 
382  }
383  }
384 
385  return ierr;
386 
387 };
388 
389 StFgtAVEfficiencyMaker::StFgtAVEfficiencyMaker( const Char_t* name):runningEvtNr(0),hitCounter(0),hitCounterR(0)
390 {
391  cout <<"AVE constructor!!" <<endl;
392  StFgtQaMaker( name, 0,0, "qName" );
393 
394 };
395 
396 StFgtAVEfficiencyMaker::~StFgtAVEfficiencyMaker()
397 {
398 
399  //delete histogram arrays
400 };
401 
402 
404  gStyle->SetPalette(1);
405  cout <<"cluster plotter finish funciton " <<endl;
406  Int_t ierr = kStOk;
407 
408  TCanvas* cRadio=new TCanvas("radioPlots","radioPlot",1000,1500);
409  TCanvas* cRadioHits=new TCanvas("radioPlotsHits","radioPlotHits",1000,1500);
410  TCanvas* cRadioNonHits=new TCanvas("radioPlotsNonHits","radioPlotNonHits",1000,1500);
411  cRadio->Divide(2,3); //6 discs
412  cRadioHits->Divide(2,3); //6 discs
413  cRadioNonHits->Divide(2,3); //6 discs
414  TCanvas* cRPRatio=new TCanvas("rPhiRatio","rPhiRatios",1000,1500);
415  cRPRatio->Divide(2,3); //6 discs
416 
417  TCanvas* cREff=new TCanvas("crEff","crEff",1000,1500);
418 
419  cREff->Divide(2,3); //6 discs
420 
421  for(Int_t iD=0;iD<kFgtNumDiscs;iD++)
422  {
423  // cRadio->cd(iD+1)->SetLogz();
424  cRadioHits->cd(iD+1);
425  radioPlotsEff[iD]->Draw("colz");
426  cRadioNonHits->cd(iD+1);
427  radioPlotsNonEff[iD]->Draw("colz");
428 
429  }
430  cRadioHits->SaveAs("radioPlotsHits.png");
431  cRadioNonHits->SaveAs("radioPlotsNonHits.png");
432 
433 
434  TCanvas* cClusterSizeR=new TCanvas("clusterSizeR","clusterSizeR",1000,1500);
435  cClusterSizeR->Divide(2,3);
436  TCanvas* cClusterSizePhi=new TCanvas("clusterSizePhi","clusterSizePhi",1000,1500);
437  cClusterSizePhi->Divide(2,3);
438  TCanvas* cChargeCorr=new TCanvas("chargeCorr","chargeCorr",1000,1500);
439  cChargeCorr->Divide(2,3);
440 
441  TCanvas* cClusterChargePhi=new TCanvas("clusterChargePhi","clusterChargePhi",1000,1500);
442  cClusterChargePhi->Divide(2,3);
443  TCanvas* cClusterChargeR=new TCanvas("clusterChargeR","clusterChargeR",1000,1500);
444  cClusterChargeR->Divide(2,3);
445 
446  TCanvas cIPProj;
447  hIp->Draw("colz");
448  cIPProj.SaveAs("ipProj.png");
449 
450  hIpZAtX0->Draw();
451  cIPProj.SaveAs("ipProjAtX0.png");
452 
453  hIpZAtY0->Draw();
454  cIPProj.SaveAs("ipProjAtY0.png");
455 
456 
457  for(Int_t iD=0;iD<kFgtNumDiscs;iD++)
458  {
459  cClusterSizeR->cd(iD+1);
460  h_clusterSizeR[iD]->Draw();
461  cClusterSizePhi->cd(iD+1);
462  h_clusterSizePhi[iD]->Draw();
463  cChargeCorr->cd(iD+1);
464  chargeCorr[iD]->Draw("colz");
465 
466  cClusterChargeR->cd(iD+1);
467  h_clusterChargeR[iD]->Draw();
468  cClusterChargePhi->cd(iD+1);
469  h_clusterChargePhi[iD]->Draw();
470 
471  }
472  cClusterSizeR->SaveAs("clusterSizeR.png");
473  cClusterSizePhi->SaveAs("clusterSizePhi.png");
474  cChargeCorr->SaveAs("chargeCorrelation.png");
475 
476  cClusterChargeR->SaveAs("clusterChargeR.png");
477  cClusterChargePhi->SaveAs("clusterChargePhi.png");
478 
479  cout <<"saving .." <<endl;
480 
481  for(Int_t iD=0;iD<kFgtNumDiscs;iD++)
482  {
483  cRadio->cd(iD+1);
484 
485  TH2D* tmpAllCounts=(TH2D*)radioPlotsEff[iD]->Clone("tmp");
486  radioPlotsEff[iD]->Add(radioPlotsNonEff[iD]);//all counts
487  // radioPlotsEff[iD]->Add(radioPlotsNonEff[iD],-1); //subtract non eff
488  for(int nx=0;nx<radioPlotsEff[iD]->GetNbinsX();nx++)
489  {
490  for(int ny=0;ny<radioPlotsEff[iD]->GetNbinsY();ny++)
491  {
492  Double_t denom=radioPlotsEff[iD]->GetBinContent(nx,ny);
493  if(denom>0)
494  radioPlotsEff[iD]->SetBinContent(nx,ny,tmpAllCounts->GetBinContent(nx,ny)/denom);
495  else
496  radioPlotsEff[iD]->SetBinContent(nx,ny,0.0);
497  }
498  }
499  // radioPlotsEff[iD]->Divide(tmpAllCounts);
500  radioPlotsEff[iD]->Draw("colz");
501  TArc *innerArc = new TArc(0,0,103.5+11.5,0.0,90.0);
502  TArc *outerArc = new TArc(0,0,394.0-11.5,0.0,90.0);
503  TLine *left = new TLine( 11.5, 103.5+11.5, 11.5, 394.0-11.5 );
504  TLine *right = new TLine( 103.5+11.5, 11.5, 394.0-11.5, 11.5 );
505  TLine *center = new TLine(
506  (103.5+11.5)*cos( 0.785398163 ),
507  (103.5+11.5)*sin( 0.785398163 ),
508  (394.0-11.5)*cos( 0.785398163 ),
509  (394.0-11.5)*sin( 0.785398163 )
510  );
511  TLine *weird = new TLine(
512  (394.0-11.5)*cos( 0.190240888 ),
513  (394.0-11.5)*sin( 0.190240888 ),
514  (394.0-11.5)*cos( 0.891863248 ),
515  (394.0-11.5)*sin( 0.891863248 )
516  );
517  innerArc->SetFillStyle(0);
518  outerArc->SetFillStyle(0);
519  innerArc->SetLineWidth(2);
520  outerArc->SetLineWidth(2);
521  left->SetLineWidth(2);
522  right->SetLineWidth(2);
523  center->SetLineWidth(2);
524  weird->SetLineWidth(2);
525  outerArc->Draw("only");
526  innerArc->Draw("only");
527  left->Draw();
528  right->Draw();
529  center->Draw();
530  weird->Draw();
531 
532 
533 
534 
535 
536  cRPRatio->cd(iD+1);
537  rPhiRatioPlots[iD]->Draw();
538  cREff->cd(iD+1);
539 
540  TH1D* tmpR=(TH1D*)rEff[iD]->Clone("tmpR");
541  rEff[iD]->Add(rNonEff[iD]);
542  for(int nx=0;nx<rEff[iD]->GetNbinsX();nx++)
543  {
544  Double_t denom=rEff[iD]->GetBinContent(nx);
545  if(denom>0)
546  rEff[iD]->SetBinContent(nx,tmpR->GetBinContent(nx)/denom);
547  else
548  rEff[iD]->SetBinContent(nx,0.0);
549  }
550  rEff[iD]->Draw();
551  }
552  cRadio->SaveAs("radioPlotsEff.png");
553  cRadio->SaveAs("radioPlotsEff.pdf");
554 
555 
556  cREff->SaveAs("rEff.png");
557  cREff->SaveAs("rEff.pdf");
558 
559  cRPRatio->SaveAs("rpRatio.png");
560  cRPRatio->SaveAs("rpRatio.pdf");
561 
562  cout <<"returning after finish" <<endl;
563  return ierr;
564 };
565 
566 
572  cout <<"AVE!!" <<endl;
573  myRootFile=new TFile("clusterEff.root","RECREATE");
574  // outTxtFile=new ofstream;
575  // outTxtFile->open("clusters.txt");
576 
577 
578  Int_t ierr = kStOk;
579 
580  char buffer[100];
581 
582  radioPlotsEff=new TH2D*[kFgtNumDiscs];
583  radioPlotsNonEff=new TH2D*[kFgtNumDiscs];
584  rPhiRatioPlots=new TH1D*[kFgtNumDiscs];
585  rEff=new TH1D*[kFgtNumDiscs];
586  rNonEff=new TH1D*[kFgtNumDiscs];
587 
588 
589  chargeCorr=new TH2D*[kFgtNumDiscs];
590  h_clusterSizeR=new TH1D*[kFgtNumDiscs];
591  h_clusterSizePhi=new TH1D*[kFgtNumDiscs];
592 
593  h_clusterChargeR=new TH1D*[kFgtNumDiscs];
594  h_clusterChargePhi=new TH1D*[kFgtNumDiscs];
595 
596 
597  hIp=new TH2D("Proj_to_IP","Proj_to_Ip",50,-100,100,50,-100,100);
598  hIpZAtX0=new TH1D("Proj_to_IPAtX0","Proj_toIPAtX0",50,-100,100);
599  hIpZAtY0=new TH1D("Proj_to_IPAtY0","Proj_toIPAtY0",50,-100,100);
600 
601  for(int iD=0;iD<kFgtNumDiscs;iD++)
602  {
603  // cout <<"id: " << iD <<endl;
604 
605  sprintf(buffer,"radioDiskEff_%d",iD);
606  radioPlotsEff[iD]=new TH2D(buffer,buffer,20,-50,50,20,-50,50);
607 
608  // cout <<"1" <<endl;
609  sprintf(buffer,"rEff_%d",iD);
610  rEff[iD]=new TH1D(buffer,buffer,100,0,50);
611  sprintf(buffer,"rNonEff_%d",iD);
612  rNonEff[iD]=new TH1D(buffer,buffer,100,0,50);
613  sprintf(buffer,"clusterSizeR_Disk_%d",iD);
614  h_clusterSizeR[iD]=new TH1D(buffer,buffer,20,0,20);
615  sprintf(buffer,"clusterSizePhi_Disk_%d",iD);
616  h_clusterSizePhi[iD]=new TH1D(buffer,buffer,20,0,20);
617  sprintf(buffer,"clusterChargeR_Disk_%d",iD);
618  h_clusterChargeR[iD]=new TH1D(buffer,buffer,100,0,50000);
619  sprintf(buffer,"clusterChargePhi_Disk_%d",iD);
620  h_clusterChargePhi[iD]=new TH1D(buffer,buffer,100,0,50000);
621  // cout <<"2" <<endl;
622  for(int nx=0;nx<radioPlotsEff[iD]->GetNbinsX();nx++)
623  {
624  for(int ny=0;ny<radioPlotsEff[iD]->GetNbinsY();ny++)
625  {
626  // radioPlotsEff[iD]->SetBinContent(nx,ny,0.1);//so that there is no divide by zero
627  }
628  }
629  // cout <<"3" <<endl;
630  sprintf(buffer,"r_phi_ChargeCorrelationInDisk_%d",iD);
631  chargeCorr[iD]=new TH2D(buffer,buffer,50,0,50000,50,0,50000);
632  sprintf(buffer,"radioDiskNonEff_%d",iD);
633  radioPlotsNonEff[iD]=new TH2D(buffer,buffer,20,-50,50,20,-50,50);
634  for(int nx=0;nx<rEff[iD]->GetNbinsX();nx++)
635  {
636  // rEff[iD]->SetBinContent(nx,0.1);
637  }
638  // cout <<"4" <<endl;
639  sprintf(buffer,"rPhiRatio_%d",iD);
640  rPhiRatioPlots[iD]=new TH1D(buffer,buffer,100,-2,10);
641  }
642  return ierr;
643 };
644 ClassImp(StFgtAVEfficiencyMaker);
virtual TObject * Clone(const char *newname="") const
the custom implementation fo the TObject::Clone
Definition: TDataSet.cxx:308
virtual Int_t Make()
Definition: Stypes.h:41