StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
StCentralityAnalyzer.cxx
1 
2 // Centrality headers
3 #include "StCentralityAnalyzer.h"
4 #include "Stypes.h"
5 #include "Param.h"
6 
7 // ROOT headers
8 #include "TString.h"
9 #include "TMath.h"
10 
11 // C++ headers
12 #include <iostream>
13 
14 ClassImp(StCentralityAnalyzer)
15 
16 
17 //________________
19  const Char_t* oFileName):
20 mOutFileName(oFileName), mOutFile(nullptr), mReader(reader), mDst(nullptr),
21  mTriggerIdCollection(), mEventsPassed(0), mDebug(kFALSE),
22  hRefMult(nullptr), hGRefMult(nullptr), hPrimVertNum(nullptr), hPrimVertZ(nullptr),
23  hZdcAdcEast(nullptr), hZdcAdcWest(nullptr), hZdcAdcSum(nullptr), hZdcCoincidenceRate(nullptr),
24  hEpdAdcEast(nullptr), hEpdAdcWest(nullptr), hEpdAdcSum(nullptr),
25  hPrimVertXvsY(nullptr), hPrimVertVpdVzDiff(nullptr), hRefMultVsTofTrayMult(nullptr),
26  hGRefMultVsTofTrayMult(nullptr), hRefMultVsTofMatched(nullptr), hRefMultVsBemcMatch(nullptr),
27  hRefMultVsRunNumber(nullptr), hGRefMultVsRunNumber(nullptr), hTofTrayMultVsRunNumber(nullptr),
28  hRefMultVsZdcCoincidenceRate(nullptr),hZdcAdcSumVsRunNumber(nullptr), hPrimTrackPtVsRunNumber(nullptr),
29  hPrimTrackNHitsVsRunNumber(nullptr), hPrimTrackDedxVsRunNumber(nullptr),
30  hGlobTrackPtVsRunNumber(nullptr), hGlobTrackNHitsVsRunNumber(nullptr),
31  hGlobTrackDedxVsRunNumber(nullptr), hAvgRefMultVsZdcCoincidenceRate(nullptr) {
32 
33  mVtxZ[0] = -70.; mVtxZ[1] = 70.;
34  mVtxR[0] = 0.; mVtxR[1] = 2.;
35  mVtxXShift = 0.; mVtxYShift = 0.;
36  mVpdVzDiff[0] = -5.; mVpdVzDiff[1] = 5.;
37  mRefMult[0] = 0.; mRefMult[1] = 1000.;
38  // Next numbers are for isobar runs
39  mRunIdBins = 57990;
40  mRunIdRange[0] = 19071030; mRunIdRange[1] = 19129020;
41  mMom[0] = 0.1; mMom[1] = 1e9;
42  mEta[0] = -1.; mEta[1] = 1.;
43  mNHits[0] = 15; mNHits[1] = 100;
44  mNHitsRatio[0] = 0.; mNHitsRatio[1] = 1.1;
45 }
46 
47 //________________
48 StCentralityAnalyzer::~StCentralityAnalyzer() {
49  /* empty */
50 }
51 
52 //________________
53 Int_t StCentralityAnalyzer::Init() {
54  if ( !mReader ) {
55  std::cout << "No StPicoDstReader has been provided" << std::endl;
56  return kStErr;
57  }
58 
59  mDst = mReader->picoDst();
60  if ( !mDst ) {
61  std::cout << "No StPicoDst has been retrieved from reader" << std::endl;
62  return kStErr;
63  }
64 
65  mOutFile = new TFile(mOutFileName, "recreate");
66  createHistograms();
67 
68  return StMaker::Init();
69 }
70 
71 //________________
73  if( mOutFile ) {
74  std::cout << "StCentralityAnalyzer::Finish - Writing histograms to the output file...";
75  mOutFile->Write();
76  mOutFile->Close();
77  std::cout << "\t [DONE]" << std::endl;
78  }
79 
80  return kStOK;
81 }
82 
83 //________________
85 
86  mEventsPassed++;
87  Bool_t readEvent = mReader->readPicoEvent(mEventsPassed);
88  //std::cout << "Working on event# " << mEventsPassed << std::endl;
89 
90  StPicoEvent* event = mDst->event();
91  TRandom3* rndm3 = new TRandom3(0);
92 
93  if ( !event ) {
94  std::cout << "No StPicoEvent has been retrieved from DST" << std::endl;
95  return kStSkip;
96  }
97 
98  // Skip bad runs
99  if ( isInBadRunList( event->runId() ) ) return kStOk;
100 
101  // Skip bad events
102  if ( !isGoodEvent(event) ) return kStOk;
103 
104  //std::cout << "Event vtxZ: " << event->primaryVertex().Z() << std::endl;
105  double VtxX = event->primaryVertex().X();
106  double VtxY = event->primaryVertex().Y();
107  double VtxZ = event->primaryVertex().Z();
108  double VzVpd = event->vzVpd();
109 
110  // Fill event information
111  Double_t refMult = event->refMult();
112 
113  // Luminosity correction
114  if(mUseLumCorr){
115  Double_t thislumcorr = calculateLumCorr(event->ZDCx());
116  double zeroToOne = rndm3->Rndm();
117  double nonIntegerRefMult = (Double_t)(refMult)-0.5+zeroToOne;
118  refMult = TMath::Nint(nonIntegerRefMult*thislumcorr);
119  }
120 
121  //fill tree
122  mTree_Vz = VtxZ;
123  mTree_ZDCx = event->ZDCx();
124  mTree_refMult = refMult;
125 
126  //correct refMult Vz dependence
127  if(mUseVzCorr){
128  refMult = getRefMultCorrVz(refMult, VtxZ);
129  mTree_refMultCor = refMult;
130  }
131  mTree_nBTOFMatched = event->nBTOFMatch();
132  mTree->Fill();
133 
134  //RefMult shape weight
135  //Double weight = getShapeWeight(VtxZ, refMult);
136  //it should be used like hRefMult->Fill(refMult, weight);
137 
138  Int_t grefMult = event->grefMult();
139 
140  double VtxZBinDouble = VtxZ/2. + 36.; //change this if Vz histogram labels change
141  if(VtxZ==73.0) hRefMultVtxZ[72]->Fill(refMult); //edge case
142  else if(VtxZ==-73.0) hRefMultVtxZ[0]->Fill(refMult); //edge case
143  else if(VtxZBinDouble - (int)VtxZBinDouble == 0.5){ //more likely edge case, filling both hists if Vz is on boundary
144  hRefMultVtxZ[(int)VtxZBinDouble]->Fill(refMult);
145  hRefMultVtxZ[(int)VtxZBinDouble+1]->Fill(refMult);
146  }
147  else hRefMultVtxZ[TMath::Nint(VtxZBinDouble)]->Fill(refMult); //most events by far
148 
149  hPrimVertXvsY->Fill(VtxX, VtxY);
150  hPrimVertZ->Fill(VtxZ);
151  hPrimVertVpdVzDiff->Fill(VtxZ-VzVpd);
152 
153  // Will be faster than direct call in case of simultanious usage
154  UInt_t runId = event->runId();
155 
156  hRefMult->Fill( refMult );
157  hGRefMult->Fill( grefMult );
158  hRefMultVsRunNumber->Fill( runId, refMult );
159  hRefMultVsTofTrayMult->Fill( refMult, event->btofTrayMultiplicity() );
160  hRefMultVsTofMatched->Fill( event->nBTOFMatch() , refMult );
161  hGRefMultVsTofTrayMult->Fill( grefMult, event->btofTrayMultiplicity() );
162  hGRefMultVsTofMatched->Fill( grefMult, event->nBTOFMatch() );
163  hRefMultVsRunNumber->Fill( runId, refMult );
164  hGRefMultVsRunNumber->Fill( runId, grefMult );
165  hRefMultVsZdcCoincidenceRate->Fill(event->ZDCx(),refMult);
166  hAvgRefMultVsZdcCoincidenceRate->Fill(event->ZDCx(), refMult);
167  hZdcCoincidenceRateVsRunNumber->Fill( runId, event->ZDCx() );
168  hTofTrayMultVsRunNumber->Fill( runId, event->btofTrayMultiplicity() );
169  hTofMatchedVsRunNumber->Fill( runId, event->nBTOFMatch() );
170 
171  std::vector<unsigned int> TheseTrigs = event->triggerIds();
172  for(int jTrigs=0; jTrigs<TheseTrigs.size(); jTrigs++){
173  int thisTrig = TheseTrigs[jTrigs];
174  if(thisTrig==600001) hAvgRefMultVsZdcCoincidenceRateForTrig[0]->Fill(event->ZDCx(), refMult);
175  if(thisTrig==600011) hAvgRefMultVsZdcCoincidenceRateForTrig[1]->Fill(event->ZDCx(), refMult);
176  if(thisTrig==600021) hAvgRefMultVsZdcCoincidenceRateForTrig[2]->Fill(event->ZDCx(), refMult);
177  if(thisTrig==600031) hAvgRefMultVsZdcCoincidenceRateForTrig[3]->Fill(event->ZDCx(), refMult);
178  }
179 
180  for ( UInt_t iTrk=0; iTrk<mDst->numberOfTracks(); iTrk++ ) {
181  StPicoTrack* track = mDst->track( iTrk );
182  if ( !track ) continue;
183 
184  if ( !isGoodTrack( track ) ) continue;
185 
186  // Fill track information
187  hGlobTrackNHitsVsRunNumber->Fill( runId, track->nHits() );
188  hGlobTrackPtVsRunNumber->Fill( runId, track->gPt() );
189  } // for ( Int_t iTrk=0; iTrk<mDst->numberOfTracks(); iTrk++ )
190 
191  return kStOK;
192 }
193 
194 //________________
195 void StCentralityAnalyzer::Clear(Option_t *opt) {
196  StMaker::Clear();
197 }
198 
199 //________________
200 //Informed by https://drupal.star.bnl.gov/STAR/system/files/Centrality_for_Run18_27GeVAuAu_ZaochenYe_20190827.pdf
201 Double_t StCentralityAnalyzer::calculateLumCorr(Double_t ZDCx) {
202  double f_ZDCx = m_LumCorr_a*ZDCx + m_LumCorr_b;
203  double LumCorr = m_LumCorr_bprime/(f_ZDCx);
204  return LumCorr;
205 }
206 
207 //________________
208 Bool_t StCentralityAnalyzer::isGoodEvent(StPicoEvent* ev) {
209  Bool_t goodEvent = ( isGoodVertex( ev->primaryVertex().X(),
210  ev->primaryVertex().Y(),
211  ev->primaryVertex().Z(),
212  ev->vzVpd() ) &&
213  isGoodTrigger( ev->triggerIds() ) );
214  Bool_t pileUpRejected = true;
215  if(mUsePileUp){
216  pileUpRejected = isNotPileUp( ev->refMult(), ev->nBTOFMatch() );
217  }
218  return ( goodEvent && pileUpRejected );
219 }
220 
221 //________________
222 void StCentralityAnalyzer::addTriggerId(const unsigned int& id) {
223  Bool_t isInList = ( std::find(mTriggerIdCollection.begin(), mTriggerIdCollection.end(), id) != mTriggerIdCollection.end() );
224  if ( !isInList ) {
225  mTriggerIdCollection.push_back( id );
226  }
227 }
228 
229 //________________
230 Bool_t StCentralityAnalyzer::isNotPileUp(UShort_t refMult, UShort_t btofMatched) {
231 
232  //double refmultcutmode=m_a0+m_a1*(btofMatched)+m_a2*pow(btofMatched,2)+m_a3*pow(btofMatched,3)+m_a4*pow(btofMatched,4);
233  double refmultcutmax = ( m_b0 + m_b1*(btofMatched) + m_b2*pow(btofMatched,2) +
234  m_b3*pow(btofMatched,3) + m_b4*pow(btofMatched,4) );
235  double refmultcutmin = ( m_c0 + m_c1*(btofMatched) + m_c2*pow(btofMatched,2) +
236  m_c3*pow(btofMatched,3) + m_c4*pow(btofMatched,4) );
237 
238  return ( refMult<refmultcutmax && refMult>refmultcutmin );
239 }
240 
241 //_________________
242 Bool_t StCentralityAnalyzer::isInBadRunList(unsigned int runId) {
243  Bool_t isInRuRu_200GeV = ( std::find( bad_run_list_ruru_200gev.begin(), bad_run_list_ruru_200gev.end(), runId) != bad_run_list_ruru_200gev.end() );
244  Bool_t isInZrZr_200GeV = ( std::find( bad_run_list_zrzr_200gev.begin(), bad_run_list_zrzr_200gev.end(), runId) != bad_run_list_zrzr_200gev.end() );
245 
246  return ( isInRuRu_200GeV || isInZrZr_200GeV );
247 }
248 
249 //________________
250 Bool_t StCentralityAnalyzer::isGoodTrigger(std::vector<unsigned int> triggers) {
251  Bool_t isInList = false;
252  for ( unsigned int i=0; i<triggers.size(); i++ ) {
253  if ( std::find(mTriggerIdCollection.begin(), mTriggerIdCollection.end(), triggers.at(i)) != mTriggerIdCollection.end() ) {
254  isInList = true;
255  break;
256  }
257  } //for ( unsigned int i=0; i<triggers.size(); i++ )
258  return isInList;
259 }
260 
261 //_________________
262 Bool_t StCentralityAnalyzer::isGoodVertex(Float_t x, Float_t y, Float_t z, Float_t vpdVz) {
263  Bool_t mIsGoodPositionZ = false;
264  Bool_t mIsGoodPositionR = false;
265  Bool_t mIsGoodVpdVzDiff = false;
266 
267  Float_t vpdVzDiff = z - vpdVz;
268  Float_t vtxPositionR = TMath::Sqrt( (x-mVtxXShift)*(x-mVtxXShift) +
269  (y-mVtxYShift)*(y-mVtxYShift) );
270 
271  mIsGoodPositionZ = ( (z > mVtxZ[0]) && (z < mVtxZ[1]) );
272  mIsGoodPositionR = ( vtxPositionR >= mVtxR[0] && vtxPositionR < mVtxR[1] );
273  mIsGoodVpdVzDiff = ( vpdVzDiff > mVpdVzDiff[0] && vpdVzDiff < mVpdVzDiff[1] );
274 
275  if(mDebug) {
276  Bool_t isGoodVertex = ( mIsGoodPositionZ && mIsGoodPositionR && mIsGoodVpdVzDiff );
277  std::cout << "IsGoodVertex: " << isGoodVertex << " : " << std::endl;
278  std::cout << "position z: " << z << " IsGood: " << mIsGoodPositionZ << std::endl
279  << "position r: " << vtxPositionR << " IsGood: " << mIsGoodPositionR << std::endl
280  << "vpdvz diff: " << vpdVzDiff << " IsGood: " << mIsGoodVpdVzDiff << std::endl;
281  } //if(mDebug)
282 
283  return ( mIsGoodPositionZ && mIsGoodPositionR && mIsGoodVpdVzDiff );
284 }
285 
286 //________________
287 Bool_t StCentralityAnalyzer::isGoodTrack(StPicoTrack *trk) {
288  Bool_t mGoodTrk = false;
289  Float_t hitRatio = ( (float)trk->nHitsFit() / (float)trk->nHitsMax() );
290  mGoodTrk = ( (trk->gPtot() >= mMom[0]) && (trk->gPtot() <= mMom[1]) &&
291  (trk->gPt() >= 0.1) &&
292  (trk->gMom().Eta() >= mEta[0]) &&
293  (trk->gMom().Eta() <= mEta[1]) &&
294  (trk->nHits() >= mNHits[0]) &&
295  (trk->nHits() <= mNHits[1]) &&
296  (hitRatio >= mNHitsRatio[0]) &&
297  (hitRatio <= mNHitsRatio[1]) );
298 
299  if(mDebug) {
300  std::cout << "IsGoodTrack: " << mGoodTrk << " : " << std::endl
301  << "primary : " << trk->isPrimary() << std::endl
302  << "gPtot : " << trk->gPtot() << std::endl
303  << "gPt : " << trk->gPt() << std::endl
304  << "gEta : " << trk->gMom().Eta() << std::endl;
305  if( trk->isPrimary() ) {
306  std::cout << "Ptot : " << trk->pMom().Mag() << std::endl
307  << "Pt : " << trk->pMom().Perp() << std::endl
308  << "Eta : " << trk->pMom().Eta() << std::endl;
309  } //if( trk->isPrimary() )
310 
311  std::cout << "nhits : " << trk->nHits() << std::endl
312  << "hitRatio: " << hitRatio << std::endl;
313  } //if(mDebug)
314  return mGoodTrk;
315 }
316 
317 //________________
318 void StCentralityAnalyzer::createHistograms() {
319  createEventHistograms();
320  createGlobalTrackHistograms();
321  createPrimaryTrackHistograms();
322 }
323 
324 //________________
325 void StCentralityAnalyzer::createEventHistograms() {
326  std::cout << "Creating event histograms...";
327 
328  mTree = new TTree("mTree","backup tree for refMult correction");
329  mTree->Branch("ZDCx", &mTree_ZDCx, "ZDCx/D");
330  mTree->Branch("Vz", &mTree_Vz, "Vz/D");
331  mTree->Branch("refMult", &mTree_refMult, "refMult/D");
332  mTree->Branch("refMultCor", &mTree_refMultCor, "refMultCor/D");
333  mTree->Branch("nBTOFMatched", &mTree_nBTOFMatched, "nBTOFMatched/D");
334 
335  Int_t refMultBins = 500;
336  Float_t refMult[2] = { 0.0, 500.0 }; //refMult hist should have bin edges on integer values in order to make centrality cuts exact
337  hRefMult = new TH1F( Form("hRefMult"),
338  Form("refMult;refMult;events"),
339  refMultBins, refMult[0], refMult[1]);
340  for ( Int_t iBin=0; iBin<73; iBin++ ) {
341  Float_t zmin = -73.;
342  Float_t zstep = 2.;
343  Float_t zrange[2] = { zmin + iBin*zstep, zmin + (iBin+1)*zstep};
344  hRefMultVtxZ[iBin] = new TH1F( Form("hRefMultVtxZ_%d",iBin),
345  Form("Reference multiplicity for %3.1f #leq z #leq %3.1f",zrange[0],zrange[1]),
346  refMultBins, refMult[0], refMult[1] );
347  } // for ( Int_t iBin=0; iBin<73; iBin++ )
348  hGRefMult = new TH1F( Form("hGRefMult"),
349  Form("gRefMult;gRefMult;events"),
350  2500, 0.0, 2500.0);
351  hPrimVertNum = new TH1F( Form("hPrimVertNum"),
352  Form("hPrimVertNum;Number of pVtx; events"),
353  15, -0.5, 14.5);
354  hPrimVertZ = new TH1F( Form("hPrimVertZ"), Form("hPrimVertZ;z (cm);events/4 cm"), 110, -220., 220.);
355  hPrimVertXvsY = new TH2F( Form("hPrimVertXvsY"),
356  Form("hPrimVertXvsY;x (cm);y (cm);events"),
357  40, -10., 10., 40, -10., 10.);
358  hPrimVertVpdVzDiff = new TH1F( Form("hPrimVertVpdVzDiff"),
359  Form("hPrimVertVpdVzDiff;Vz_{TPC}-Vz_{VPD} (cm);events"),
360  40, -20., 20.);
361  hRefMultVsTofTrayMult = new TH2F( Form("hRefMultVsTofTrayMult"),
362  Form("hRefMultVsTofTrayMult;refMult;btofTrayMultiplicity;events"),
363  650, -0.5, 649.5, 650, -0.5, 649.5 );
364  hRefMultVsTofMatched = new TH2F( Form("hRefMultVsTofMatched"),
365  Form("hRefMultVsTofTrayMult;btofMatched;refMult;events"),
366  650, -0.5, 649.5, 650, -0.5, 649.5 );
367  hGRefMultVsTofTrayMult = new TH2F( Form("hGRefMultVsTofTrayMult"),
368  Form("hGRefMultVsTofTrayMult;gRefMult;btofTrayMultiplicity;events"),
369  1950, -0.5, 1949.5,
370  650, -0.5, 649.5);
371  hGRefMultVsTofMatched = new TH2F( Form("hGRefMultVsTofMatched"),
372  Form("hGRefMultVsTofTrayMult;gRefMult;BTofMatched;events"),
373  1950, -0.5, 1949.5,
374  650, -0.5, 649.5);
375  hRefMultVsZdcCoincidenceRate = new TH2F( Form("hRefMultVsZdcCoincidenceRate"),
376  Form("hRefMultVsZdcCoincidenceRate;ZdcCoincidenceRate (Hz);refMult"),
377  700, 6000., 15000.,650,-0.5,649.5);
378  hRefMultVsRunNumber = new TProfile( Form("hRefMultVsRunNumber"),
379  Form("hRefMultVsRunNumber;runId;<refMult>"),
380  mRunIdBins, mRunIdRange[0], mRunIdRange[1] );
381  hAvgRefMultVsZdcCoincidenceRate = new TProfile( Form("hAvgRefMultVsZdcCoincidenceRate"),
382  Form("hAvgRefMultVsZdcCoincidenceRate;ZdcCoincidenceRate (Hz);<refMult>"),
383  700, 6000., 15000.);
384  int Trigs[4]={600001,600011,600021,600031};
385  for(int iTrig=0; iTrig<4; iTrig++){
386  hAvgRefMultVsZdcCoincidenceRateForTrig[iTrig] = new TProfile( Form("hAvgRefMultVsZdcCoincidenceRateForTrig_%d",Trigs[iTrig]),
387  Form("hAvgRefMultVsZdcCoincidenceRate_%d;ZdcCoincidenceRate (Hz);<refMult>",Trigs[iTrig]),
388  700, 6000., 15000.);
389  }
390 
391  hGRefMultVsRunNumber = new TProfile( Form("hGRefMultVsRunNumber"),
392  Form("hGRefMultVsRunNumber;runId;<gRefMult>"),
393  mRunIdBins, mRunIdRange[0], mRunIdRange[1] );
394  hZdcCoincidenceRateVsRunNumber = new TProfile( Form("hZdcCoincidenceRateVsRunNumber"),
395  Form("hZdcCoincidenceRateVsRunNumber;runId;<ZDC coincidence rate> (kHz)"),
396  mRunIdBins, mRunIdRange[0], mRunIdRange[1] );
397  hTofTrayMultVsRunNumber = new TProfile( Form("hTofTrayMultVsRunNumber"),
398  Form("hTofTrayMultVsRunNumber;runId;btofTrayMult"),
399  mRunIdBins, mRunIdRange[0], mRunIdRange[1] );
400  hTofMatchedVsRunNumber = new TProfile( Form("hTofMatchedVsRunNumber"),
401  Form("hTofTrayMultVsRunNumber;runId;BTofMatched"),
402  mRunIdBins, mRunIdRange[0], mRunIdRange[1] );
403 
404  std::cout << "\t[DONE]" << std::endl;
405 }
406 
407 //________________
408 void StCentralityAnalyzer::createGlobalTrackHistograms() {
409  std::cout << "Creating global track histograms...";
410  hGlobTrackNHitsVsRunNumber = new TProfile( Form("hGlobTrackNHitsVsRunNumber"),
411  Form("hGlobTrackNHitsVsRunNumber;runId;primary track <nHits>"),
412  mRunIdBins, mRunIdRange[0], mRunIdRange[1] );
413  hGlobTrackDedxVsRunNumber = new TProfile( Form("hGlobTrackDedxVsRunNumber"),
414  Form("hGlobTrackDedxVsRunNumber;runId;primary track <dE/dx> (keV/cm)"),
415  mRunIdBins, mRunIdRange[0], mRunIdRange[1] );
416  hGlobTrackPtVsRunNumber = new TProfile( Form("hGlobTrackPtVsRunNumber"),
417  Form("hGlobTrackPtVsRunNumber;runId;global track <p_T> (GeV/c)"),
418  mRunIdBins, mRunIdRange[0], mRunIdRange[1] );
419  std::cout << "\t[DONE]" << std::endl;
420 }
421 
422 //________________
423 void StCentralityAnalyzer::createPrimaryTrackHistograms() {
424  std::cout << "Creating primary track histograms...";
425  hPrimTrackPtVsRunNumber = new TProfile( Form("hPrimTrackPtVsRunNumber"),
426  Form("hPrimTrackPtVsRunNumber;runId;primary track <p_T> (GeV/c)"),
427  mRunIdBins, mRunIdRange[0], mRunIdRange[1] );
428  hPrimTrackNHitsVsRunNumber = new TProfile( Form("hPrimTrackNHitsVsRunNumber"),
429  Form("hPrimTrackNHitsVsRunNumber;runId;primary track <nHits>"),
430  mRunIdBins, mRunIdRange[0], mRunIdRange[1] );
431  hPrimTrackDedxVsRunNumber = new TProfile( Form("hPrimTrackDedxVsRunNumber"),
432  Form("hPrimTrackDedxVsRunNumber;runId;primary track <dE/dx> (keV/cm)"),
433  mRunIdBins, mRunIdRange[0], mRunIdRange[1] );
434  std::cout << "\t[DONE]" << std::endl;
435 }
436 
437 //_________________
438 const std::vector<unsigned int> StCentralityAnalyzer::bad_run_list_ruru_200gev = {
439  19120009,19102023,19102054,19103022,
440  19083049, 19083050, 19083051, 19083052, 19083053,
441  19083054, 19083055, 19083056, 19083057, 19083058,
442  19083059, 19083060, 19083061, 19083062, 19083063,
443  19083064, 19083065, 19083066, 19083067, 19084001,
444  19084002, 19084003, 19084004, 19084005, 19084006,
445  19084007, 19084008, 19084010, 19084011, 19084013,
446  19084022, 19084024, 19084025, 19084026, 19084027,
447  19084028, 19084029, 19084030, 19084031, 19084032,
448  19084033, 19084034, 19084035, 19084036, 19084037,
449  19084038, 19084039, 19084053, 19084055, 19084057,
450  19084059, 19084060, 19084061, 19084062, 19084063,
451  19084064, 19084065, 19084066, 19084067, 19084068,
452  19084070, 19084071, 19084072, 19085001, 19085002,
453  19085003, 19085004, 19085005, 19085006, 19085007,
454  19085008, 19085009, 19085010, 19085011, 19085012,
455  19085013, 19085014, 19085015, 19085016, 19085017,
456  19085018, 19085019, 19085020, 19085021, 19085023,
457  19085024, 19085025, 19085026, 19085058, 19086026,
458  19086060, 19086061, 19086062, 19086063, 19086064,
459  19086066, 19086067, 19086069, 19086070, 19086072,
460  19086073, 19086074, 19086076, 19086077, 19086080,
461  19087001, 19087012, 19087014, 19087015, 19087016,
462  19087017, 19087021, 19087022, 19087038, 19087042,
463  19088051, 19088052, 19088053, 19088055, 19090009,
464  19090010, 19090011, 19090012, 19090015, 19090016,
465  19090018, 19090019, 19090021, 19090022, 19090023,
466  19090024, 19090025, 19090032, 19092051, 19093042,
467  19093043, 19095061, 19096002, 19096005, 19096006,
468  19096057, 19097057, 19098017, 19098018, 19098020,
469  19100045, 19103007, 19103041, 19105024, 19105026,
470  19106023, 19106034, 19107045, 19110015, 19110039,
471  19112012, 19112029, 19115020, 19116035, 19120047,
472  19120048, 19122004, 19122005
473  /* Here will be some runIds in the format: 11111111, 22222222, etc */
474 };
475 
476 //_________________
477 const std::vector<unsigned int> StCentralityAnalyzer::bad_run_list_zrzr_200gev = {
478  19120009,19102023,19102054,19103022,
479  19083049, 19083050, 19083051, 19083052, 19083053,
480  19083054, 19083055, 19083056, 19083057, 19083058,
481  19083059, 19083060, 19083061, 19083062, 19083063,
482  19083064, 19083065, 19083066, 19083067, 19084001,
483  19084002, 19084003, 19084004, 19084005, 19084006,
484  19084007, 19084008, 19084010, 19084011, 19084013,
485  19084022, 19084024, 19084025, 19084026, 19084027,
486  19084028, 19084029, 19084030, 19084031, 19084032,
487  19084033, 19084034, 19084035, 19084036, 19084037,
488  19084038, 19084039, 19084053, 19084055, 19084057,
489  19084059, 19084060, 19084061, 19084062, 19084063,
490  19084064, 19084065, 19084066, 19084067, 19084068,
491  19084070, 19084071, 19084072, 19085001, 19085002,
492  19085003, 19085004, 19085005, 19085006, 19085007,
493  19085008, 19085009, 19085010, 19085011, 19085012,
494  19085013, 19085014, 19085015, 19085016, 19085017,
495  19085018, 19085019, 19085020, 19085021, 19085023,
496  19085024, 19085025, 19085026, 19085058, 19086026,
497  19086060, 19086061, 19086062, 19086063, 19086064,
498  19086066, 19086067, 19086069, 19086070, 19086072,
499  19086073, 19086074, 19086076, 19086077, 19086080,
500  19087001, 19087012, 19087014, 19087015, 19087016,
501  19087017, 19087021, 19087022, 19087038, 19087042,
502  19088051, 19088052, 19088053, 19088055, 19090009,
503  19090010, 19090011, 19090012, 19090015, 19090016,
504  19090018, 19090019, 19090021, 19090022, 19090023,
505  19090024, 19090025, 19090032, 19092051, 19093042,
506  19093043, 19095061, 19096002, 19096005, 19096006,
507  19096057, 19097057, 19098017, 19098018, 19098020,
508  19100045, 19103007, 19103041, 19105024, 19105026,
509  19106023, 19106034, 19107045, 19110015, 19110039,
510  19112012, 19112029, 19115020, 19116035, 19120047,
511  19120048, 19122004, 19122005
512  /* Here will be some runIds in the format: 11111111, 22222222, etc */
513 };
514 
515 //_________________
516 Double_t StCentralityAnalyzer::getRefMultCorrVz(Double_t RefMult, Double_t Vz){
517 
518  const Double_t RefMult_ref = m_vzCorr0; //RefMult at |z|<1.
519  const Double_t RefMult_z = m_vzCorr0
520  + m_vzCorr1*Vz + m_vzCorr2*pow(Vz, 2)
521  + m_vzCorr3*pow(Vz, 3) + m_vzCorr4*pow(Vz, 4)
522  + m_vzCorr5*pow(Vz, 5) + m_vzCorr6*pow(Vz, 6);
523 
524  Double_t ScaleFactor = 1.0;
525 
526  if(RefMult_z > 0.0){
527  ScaleFactor = RefMult_ref / RefMult_z;
528  }
529 
530  return RefMult * ScaleFactor;
531 
532 }
533 
534 //_________________
535 Double_t StCentralityAnalyzer::getShapeWeight(Double_t Vz, Double_t RefMult){
536 
537  // no shape correction for -9<=Vz<=9
538  if(Vz>=-9 && Vz<=9) return 1.;
539  //obtain index to load weight
540  Double_t VtxZBinDouble = Vz/2. + 17.;
541  Int_t VzIndex = 0;
542  if(Vz == 25.) VzIndex = 29;
543  else if(Vz == -35.) VzIndex = 0;
544  else VzIndex = TMath::Nint(VtxZBinDouble);
545  //handle VzIndex for Vz>9
546  if(VzIndex >= 22) VzIndex = VzIndex - 9;
547  //retrive shape weight
548  Double_t weight = ShapeWeightArray[mShapeIndex][VzIndex][TMath::Nint(RefMult)];
549  //handle bad weight
550  if(weight == 0 || TMath::IsNaN(weight)) weight = 1.;
551  return weight;
552 
553 }
Float_t gPtot() const
Return global total momentum (GeV/c)
Definition: StPicoTrack.h:70
Bool_t isPrimary() const
Return if track is primary.
Definition: StPicoTrack.h:178
Bool_t readPicoEvent(Long64_t iEvent)
Read next event in the chain.
Allows to read picoDst file(s)
UShort_t nBTOFMatch() const
Return number of tracks that matched TOF.
Definition: StPicoEvent.h:62
void Clear(Option_t *option="")
User defined functions.
StPicoDst * picoDst()
Return a pointer to picoDst (return NULL if no dst is found)
virtual void Clear(Option_t *option="")
User defined functions.
Definition: StMaker.cxx:634
TVector3 gMom() const
Return momentum (GeV/c) of the global tracks at the point of DCA to the primary vertex.
Definition: StPicoTrack.h:64
static StPicoEvent * event()
Return pointer to current StPicoEvent (class holding the event wise information)
Definition: StPicoDst.h:63
Holds information about track parameters.
Definition: StPicoTrack.h:35
Process and build distributions for centrality determination.
static UInt_t numberOfTracks()
Return number of tracks.
Definition: StPicoDst.h:104
Definition: Stypes.h:49
Float_t vzVpd() const
Return z position of the primary vertex estimated by VPD.
Definition: StPicoEvent.h:164
Float_t gPt() const
Return global tranverse momentum (GeV/c)
Definition: StPicoTrack.h:68
TVector3 pMom() const
Return momentum (GeV/c) of the primary track. Return (0,0,0) if not primary track.
Definition: StPicoTrack.h:62
Int_t refMult() const
Return RefMult (-0.5&lt;eta&lt;0.5)
Definition: StPicoEvent.h:83
static StPicoTrack * track(Int_t i)
Return pointer to i-th track.
Definition: StPicoDst.h:65
Definition: Stypes.h:40
Int_t nHitsFit() const
Return number of hits fit.
Definition: StPicoTrack.h:106
Stores global information about the event.
Definition: StPicoEvent.h:24
std::vector< unsigned int > triggerIds() const
Return trigger list of the current event.
Definition: StPicoEvent.h:65
Int_t nHitsMax() const
Return number of hits possible.
Definition: StPicoTrack.h:108
Int_t nHits() const
Return number of hits.
Definition: StPicoTrack.h:104
Definition: Stypes.h:44
TVector3 primaryVertex() const
Return primary vertex position.
Definition: StPicoEvent.h:52
Definition: Stypes.h:41