StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
StTofpNtupleMaker.cxx
1 /***************************************************************************
2  *
3  * $Id: StTofpNtupleMaker.cxx,v 1.8 2018/02/26 23:26:51 smirnovd Exp $
4  *
5  * Author: Frank Geurts
6  ***************************************************************************
7  *
8  * Description: Maker fills TOFp Tree from StTofCollection
9  *
10  **************************************************************************/
11 #include <iostream>
12 #include "StEventTypes.h"
13 #include "Stypes.h"
14 #include "StThreeVectorF.hh"
15 #include "StHelix.hh"
16 #include "StTrackGeometry.h"
17 #include "StEventUtilities/StuRefMult.hh"
18 #include "PhysicalConstants.h"
19 #include "StPhysicalHelixD.hh"
20 #include "StTofUtil/StTofGeometry.h"
21 #include "TNtuple.h"
22 #include "TFile.h"
23 #include "StMemoryInfo.hh"
24 #include "StTimer.hh"
25 #include "StTofUtil/tofPathLength.hh"
26 #include "StTofpNtupleMaker.h"
27 
28 
29 
30 //---------------------------------------------------------------------------
33  mTupleFileName="tofntuple.root";
34  setValidAdcRange(30,1200);
35  setValidTdcRange(1,2047);
36  doPrintMemoryInfo = kFALSE;
37  doPrintCpuInfo = kFALSE;
38 }
39 
42 
43 
44 //---------------------------------------------------------------------------
47 
48  if (mTupleFileName!="") bookNtuples();
49 
50  mTofGeom = new StTofGeometry();
51  mTofGeom->initDaqMap();
52 
53  mAcceptedEvents = 0;
54  mPvpdEntries = 0;
55  mTofpEvents = 0;
56  mTofpEntries = 0;
57 
58  return kStOK;
59 }
60 
63  if (!(mTupleFileName=="")){
64  mTupleFile->Write();
65  mTupleFile->Close();
66  LOG_INFO << "StTofpNtupleMaker::Finish() ntuple file "
67  << mTupleFileName << " closed." << endm;
68  }
69 
70  LOG_INFO << "StTofpNtupleMaker -- statistics" << endm;
71  LOG_INFO << " accepted events : " << mAcceptedEvents << endm;
72  LOG_INFO << " pVPD entries : " << mPvpdEntries << endm;
73  LOG_INFO << " TOFp entries/events : " << mTofpEntries << "/" << mTofpEvents << endm;
74  return kStOK;
75 }
76 
77 
78 //---------------------------------------------------------------------------
81  LOG_INFO << "StTofpNtupleMaker -- welcome" << endm;
82 
83  StEvent *event = (StEvent *) GetInputDS("StEvent");
84 
85  //.........................................................................
86  // event selection ...
87  if (!event ||
88  !event->tofCollection() ||
89  !event->tofCollection()->dataPresent() ||
90  !event->primaryVertex()){
91  LOG_INFO << "StTofpNtupleMaker -- nothing to do ... bye-bye"<< endm;
92  return kStOK;
93  }
94 
95  mAcceptedEvents++;
96 
97  StTimer timer;
98  if (doPrintCpuInfo) timer.start();
99  if (doPrintMemoryInfo) StMemoryInfo::instance()->snapshot();
100 
101  // determine TOF configuration from run#
102  mYear2= (event->runId()<4000000);
103  mYear3= (event->runId()>4000000&&event->runId()<5000000);
104  mYear4= (event->runId()>5000000);
105 
106 
107  //.........................................................................
108  // Collect global data for both ntuples
109 
110  //-- Primary vertex & trigger information
111  float xvtx = event->primaryVertex()->position().x();
112  float yvtx = event->primaryVertex()->position().y();
113  float zvtx = event->primaryVertex()->position().z();
114  StL0Trigger* pTrigger = event->l0Trigger();
115  unsigned int triggerWord = 0;
116  if (pTrigger) triggerWord = pTrigger->triggerWord();
117 
118  //-- Get the ZDC and CTB data (for pVPD ntuple)
119  StTriggerDetectorCollection *theTriggers = event->triggerDetectorCollection();
120  float zdcSumEast(-1.), zdcSumWest(-1), ctbSum(-1);
121  if (theTriggers){
122  StCtbTriggerDetector &theCtb = theTriggers->ctb();
123  StZdcTriggerDetector &theZdc = theTriggers->zdc();
124  zdcSumEast = theZdc.adcSum(east);
125  zdcSumWest = theZdc.adcSum(west);
126  ctbSum = 0;
127  for (unsigned int islat=0; islat<theCtb.numberOfSlats(); islat++)
128  for (unsigned int itray=0; itray<theCtb.numberOfTrays(); itray++)
129  ctbSum += theCtb.mips(itray, islat, 0);
130  }
131 
132  //-- Number of primary tracks
133  // (note: different meaning in event.root and richtof.root)
134  int refmult(0);
135  bool richTofMuDST = (event->summary()->numberOfExoticTracks() == -999);
136  if (richTofMuDST)
137  refmult = event->summary()->numberOfGoodTracks();
138  else
139  refmult = uncorrectedNumberOfPrimaries(*event);
140  if (Debug()){
141  LOG_INFO << " #Tracks :" << event->summary()->numberOfTracks()
142  << "\n #goodPrimaryTracks:" << event->summary()->numberOfGoodPrimaryTracks()
143  << "\n #uncorr.prim.tracks :" << refmult << endm;
144  if (!richTofMuDST)
145  LOG_INFO << " #goodTracks (global):" << event->summary()->numberOfGoodTracks() << endm;
146  }
147 
148  //-- Check and fill local copy of TOFx+pVPD ADC and TDC data
149  StTofCollection *theTof = event->tofCollection();
150  getTofData(theTof);
151 
152  //-- Update TOFp slat counters
153  float sumAdcTofp(0.); int nAdcTofp(0), nTdcTofp(0), nAdcTdcTofp(0);
154  for (int i=0;i<NTOFP;i++){
155  sumAdcTofp+=mTofpAdc[i];
156  bool tdcValid=validTdc(mTofpTdc[i]);
157  bool adcValid=validAdc(mTofpAdc[i]);
158  if (tdcValid) nTdcTofp++;
159  if (adcValid) nAdcTofp++;
160  if (adcValid && tdcValid) nAdcTdcTofp++;
161  }
162  if (Debug())
163  LOG_INFO << " TOFp #Adc:" << nAdcTofp << " #Tdc:" << nTdcTofp << endm;
164 
165 
166 
167  //.........................................................................
168  // build pVPD ntuple
169  if (!(mTupleFileName=="")){
170  int k(0);
171  float tuple[31];
172  tuple[k++] = event->runId(); // the run number
173  tuple[k++] = event->id(); // the event number
174  tuple[k++] = triggerWord;
175  tuple[k++] = event->summary()->magneticField();
176  tuple[k++] = zvtx; // z-vertex
177  tuple[k++] = event->primaryVertex()->chiSquared(); // z-vertex chi2
178  tuple[k++] = ctbSum; // CTB sum
179  tuple[k++] = zdcSumEast; // ZDC sum east
180  tuple[k++] = zdcSumWest; // ZDC sum west
181  tuple[k++] = refmult; // STAR reference multiplicity
182  tuple[k++] = event->summary()->numberOfGoodPrimaryTracks(); // #primary tracks
183  tuple[k++] = sumAdcTofp; // TOFp sum
184  tuple[k++] = nTdcTofp; // TOFp hits
185  for (int i=0;i<NPVPD;i++) tuple[k++] = mPvpdTdc[i];
186  for (int i=0;i<NPVPD;i++) tuple[k++] = mPvpdAdc[i];
187  for (int i=0;i<NPVPD;i++) tuple[k++] = mPvpdAdcLoRes[i];
188 
189  LOG_INFO << " pVPD update ..." << endm;
190  mPvpdTuple->Fill(tuple);
191  mPvpdEntries++;
192  }
193 
194 
195 
196  //.........................................................................
197  // build TOFp ntuple
198 
199  //-- make sure tofSlats are available
200  if (event->tofCollection()->slatsPresent()){
201 
202  mTofpEvents++;
203  int entriesThisEvent(0);
204 
205  //-- Loop over the slat container and retrieve the relevant parameters
206  StSPtrVecTofSlat& slatTofVec = theTof->tofSlats();
207  for (size_t i = 0; i < slatTofVec.size(); i++) {
208  StTofSlat *thisSlat = slatTofVec[i];
209  StTrack *thisTrack = thisSlat->associatedTrack();
210  StTrackGeometry *theTrackGeometry =
211  (mOuterTrackGeometry)?thisTrack->outerGeometry():thisTrack->geometry();
212 
213  //- retrieve and recalculate parameters
214  double pathLength = tofPathLength(&event->primaryVertex()->position(),
215  &thisSlat->position(),
216  theTrackGeometry->helix().curvature());
217  const StThreeVectorF momentum = theTrackGeometry->momentum();
218 
219  //- dig out from the dedx and rich pid traits
220  float dedx(0.), cherang(0);
221  int dedx_np(0), cherang_nph(0);
222  StSPtrVecTrackPidTraits& traits = thisTrack->pidTraits();
223  for (unsigned int it=0;it<traits.size();it++){
224  if (traits[it]->detector() == kTpcId){
225  StDedxPidTraits *pid = dynamic_cast<StDedxPidTraits*>(traits[it]);
226  if (pid && pid->method() ==kTruncatedMeanId){
227  dedx = pid->mean();
228  dedx_np = pid->numberOfPoints();
229  }
230  } else if (traits[it]->detector() == kRichId){
231  StRichPidTraits *pid = dynamic_cast<StRichPidTraits*>(traits[it]);
232  if (pid){
233  StRichSpectra* pidinfo = pid->getRichSpectra();
234  if (pidinfo && pidinfo->getCherenkovPhotons()>2){
235  cherang = pidinfo->getCherenkovAngle();
236  cherang_nph = pidinfo->getCherenkovPhotons();
237  }
238  }
239  }
240  }
241 
242  //- TOFp Slat Ntuple entry for a single matched slat
243  mSlatData.run = event->runId();
244  mSlatData.evt = event->id();
245  mSlatData.trgword = triggerWord;
246  mSlatData.magfield = event->summary()->magneticField();
247  mSlatData.ctbsum = ctbSum;
248  mSlatData.zdcsum = zdcSumEast + zdcSumWest;
249  mSlatData.xvtx = xvtx;
250  mSlatData.yvtx = yvtx;
251  mSlatData.zvtx = zvtx;
252  mSlatData.zvtxchi2 = event->primaryVertex()->chiSquared();
253  mSlatData.refmult = refmult;
254  mSlatData.nprimary = event->summary()->numberOfGoodPrimaryTracks();
255  mSlatData.meanpt = event->summary()->meanPt();
256  mSlatData.te1 = (int)mPvpdTdc[0]; mSlatData.te2 = (int)mPvpdTdc[1];
257  mSlatData.te3 = (int)mPvpdTdc[2]; mSlatData.tw1 = (int)mPvpdTdc[3];
258  mSlatData.tw2 = (int)mPvpdTdc[4]; mSlatData.tw3 = (int)mPvpdTdc[5];
259  if (mYear3||mYear4) {
260  mSlatData.ae1 = (int)mPvpdAdcLoRes[0]; mSlatData.ae2 = (int)mPvpdAdcLoRes[1];
261  mSlatData.ae3 = (int)mPvpdAdcLoRes[2]; mSlatData.aw1 = (int)mPvpdAdcLoRes[3];
262  mSlatData.aw2 = (int)mPvpdAdcLoRes[4]; mSlatData.aw3 = (int)mPvpdAdcLoRes[5];
263  } else {
264  mSlatData.ae1 = (int)mPvpdAdc[0]; mSlatData.ae2 = (int)mPvpdAdc[1];
265  mSlatData.ae3 = (int)mPvpdAdc[2]; mSlatData.aw1 = (int)mPvpdAdc[3];
266  mSlatData.aw2 = (int)mPvpdAdc[4]; mSlatData.aw3 = (int)mPvpdAdc[5];
267  }
268  mSlatData.slat = thisSlat->slatIndex();//jj+1;
269  mSlatData.tdc = thisSlat->tdc();//mTofpTdc[jj];
270  mSlatData.adc = thisSlat->adc();//(int)mTofpAdc[jj];
271  mSlatData.hitprof = thisSlat->hitProf(); //newerSlatHitVec[ii].hitProfile;
272  mSlatData.matchflag = thisSlat->matchFlag(); //newerSlatHitVec[ii].matchFlag;
273  mSlatData.zhit = thisSlat->zHit(); //localHitPos;
274  mSlatData.trackId = (Int_t)thisTrack->key();
275  mSlatData.ntrackpoints= thisTrack->detectorInfo()->numberOfPoints(kTpcId);
276  mSlatData.nfitpoints = thisTrack->fitTraits().numberOfFitPoints(kTpcId);
277  mSlatData.r_last = thisTrack->detectorInfo()->lastPoint().perp();
278  mSlatData.chi2 = thisTrack->fitTraits().chi2(0);
279 
280  mSlatData.s = fabs(pathLength);
281  mSlatData.p = momentum.mag()* theTrackGeometry->charge();
282  mSlatData.pt = momentum.perp();
283  mSlatData.px = momentum.x();
284  mSlatData.py = momentum.y();
285  mSlatData.pz = momentum.z();
286  mSlatData.eta = momentum.pseudoRapidity();
287  mSlatData.dedx = dedx;
288  mSlatData.dedx_np = dedx_np;
289  mSlatData.cherang = cherang;
290  mSlatData.cherang_nph = cherang_nph;
291 
292  mSlatTuple->Fill();
293  mTofpEntries++;
294  entriesThisEvent++;
295  }
296 
297  LOG_INFO << " TOFp update: " << entriesThisEvent << " entries" <<endm;
298  }
299 
300 
301  //- debug info
302  if (doPrintMemoryInfo) {
303  StMemoryInfo::instance()->snapshot();
304  StMemoryInfo::instance()->print();
305  }
306  if (doPrintCpuInfo) {
307  timer.stop();
308  LOG_INFO << "CPU time for StEventMaker::Make(): "
309  << timer.elapsedTime() << " sec\n" << endm;
310  }
311 
312  LOG_INFO << "StTofpNtupleMaker -- bye-bye" << endm;
313  return kStOK;
314 }
315 
316 
317 
318 //---------------------------------------------------------------------------
321  if (!tofCollection) return kStERR;
322  StSPtrVecTofData &tofData = tofCollection->tofData();
323 
324  // perform consistency check
325  bool dataOK(true);
326  if (Debug()) LOG_INFO << "TOF raw data consistency test ..."<< endm;
327  for (int i=0;i<48;i++){
328  if (tofData[i]->dataIndex() != mTofGeom->daqToSlatId(i)) {
329  dataOK = false;
330  LOG_INFO << "getTofData===>WARNING: " << tofData[i]->dataIndex() << " " << mTofGeom->daqToSlatId(i) << endm;
331  }
332  //if (Debug()) LOG_INFO << *tofData[i];
333  }
334 
335  for (int i=0;i<NTOFP;i++){
336  mTofpAdc[i] = tofData[i]->adc();
337  mTofpTdc[i] = tofData[i]->tdc();
338  }
339 
340  if (mYear2){
341  // swap ADC channels 3 and 4 ... this should move to the DAQreader!
342  float tmp = mTofpAdc[3];
343  mTofpAdc[3] = mTofpAdc[2];
344  mTofpAdc[2] = tmp;
345  }
346 
347  for (int i=0;i<NPVPD;i++){
348  mPvpdAdc[i] = tofData[42+i]->adc();
349  mPvpdTdc[i] = tofData[42+i]->tdc();
350  if (mYear3||mYear4)
351  mPvpdAdcLoRes[i] = tofData[54+i]->adc();
352  }
353 
354  if (!dataOK) return kStWarn;
355 
356  return kStOK;
357 }
358 
359 //---------------------------------------------------------------------------
361 void StTofpNtupleMaker::bookNtuples(){
362  mTupleFile = new TFile(mTupleFileName.c_str(), "RECREATE");
363  LOG_INFO << "StTofpNtupleMaker::bookNtuples() file "
364  << mTupleFileName << " opened" << endm;
365 
366  // pVPD timing
367  string varList = "run:evt:trgwrd:magfield:zvtx:zvtxchi2:ctbsum"
368  ":zdceast:zdcwest:refmult:nprimary:tofpsum"
369  ":ntofp:te1:te2:te3:tw1:tw2:tw3:ael1:ael2:ael3"
370  ":awl1:awl1:awl3:ae1:ae2:ae3:aw1:aw2:aw3";
371  mPvpdTuple = new TNtuple("pvpd","tofp timing",varList.c_str());
372 
373  // TOFp calibration ntuple
374  mSlatTuple = new TTree("tofp","TOFp slat data");
375  mSlatTuple->Branch("run",&mSlatData.run,"run/I");
376  mSlatTuple->Branch("evt",&mSlatData.evt,"evt/I");
377  mSlatTuple->Branch("trgword",&mSlatData.trgword,"trgword/I");
378  mSlatTuple->Branch("magfield",&mSlatData.magfield,"magfield/F");
379  mSlatTuple->Branch("ctbsum",&mSlatData.ctbsum,"ctbsum/F");
380  mSlatTuple->Branch("zdcsum",&mSlatData.zdcsum,"zdcsum/F");
381  mSlatTuple->Branch("xvtx",&mSlatData.xvtx,"xvtx/F");
382  mSlatTuple->Branch("yvtx",&mSlatData.yvtx,"yvtx/F");
383  mSlatTuple->Branch("zvtx",&mSlatData.zvtx,"zvtx/F");
384  mSlatTuple->Branch("zvtxchi2",&mSlatData.zvtxchi2,"zvtx/F");
385  mSlatTuple->Branch("refmult",&mSlatData.refmult,"refmult/I");
386  mSlatTuple->Branch("nprimary",&mSlatData.nprimary,"nprimary/I");
387  mSlatTuple->Branch("meanpt",&mSlatData.meanpt,"meanpt/F");
388  mSlatTuple->Branch("tdcstart",&mSlatData.tdcstart,"tdcstart/F");
389  mSlatTuple->Branch("pvpd",&mSlatData.te1,"te1/I:te2:te3:tw1:tw2:tw3:ae1:ae2:ae3:aw1:aw2:aw3");
390  mSlatTuple->Branch("slat",&mSlatData.slat,"slat/I:tdc:adc:hitprof:matchflag:zhit/F:zhitinner/F:zhitouter/F:ss:theta_xy:theta_zr");
391  //mSlatTuple->Branch("track",&mSlatData.ntrackpoints,"ntrackpoints/I:nfitpoints:r_last/F:chi2:s:p:pt:px:py:pz:eta:dedx/F:dedx_np/I:cherangle/F:cherangle_nph/I");
392  mSlatTuple->Branch("track",&mSlatData.trackId,"trackId/I:ntrackpoints/I:nfitpoints:r_last/F:chi2:s:p:pt:px:py:pz:eta:dedx/F:dedx_np/I:cherangle/F:cherangle_nph/I");
393 
394 
395  // TOFp matching ntuple
396  mMatchTuple = new TTree("match","TOFp match data");
397  mMatchTuple->Branch("id",&mMatchData.daqid,"daqid/I");
398  mMatchTuple->Branch("slat",&mMatchData.nneighbors,"nneighbors/I:zlocal/F:philocal/F:hitprof/I");
399  mMatchTuple->Branch("track",&mMatchData.nfitpoints,"nfitpoints/I:ntrackpoints/I:maxpoints/I:r_last/F:p/F:pt/F");
400 
401  mNoMatchTuple = new TTree("nomatch","TOFp no-match data");
402  mNoMatchTuple->Branch("id",&mMatchData.daqid,"daqid/I");
403  mNoMatchTuple->Branch("slat",&mMatchData.nneighbors,"nneighbors/I:zlocal/F:philocal/F:hitprof/I");
404  mNoMatchTuple->Branch("track",&mMatchData.nfitpoints,"nfitpoints/I:ntrackpoints/I:maxpoints/I:r_last/F:p/F:pt/F");
405 
406  return;
407 }
408 
409 
410 //---------------------------------------------------------------------------
411 /***************************************************************************
412  *
413  * $Log: StTofpNtupleMaker.cxx,v $
414  * Revision 1.8 2018/02/26 23:26:51 smirnovd
415  * StTof: Remove outdated ClassImp macro
416  *
417  * Revision 1.7 2018/02/26 23:13:20 smirnovd
418  * Move embedded CVS log messages to the end of file
419  *
420  * Revision 1.6 2007/04/17 23:01:28 dongx
421  * replaced with standard STAR Loggers
422  *
423  * Revision 1.5 2004/04/10 04:36:25 dongx
424  * additional update for AdcLoRes in ntuple
425  *
426  * Revision 1.4 2004/04/09 19:26:25 dongx
427  * Add some missing updates for year4, add AdcLoRes in ntuple
428  *
429  * Revision 1.3 2004/04/01 19:19:00 dongx
430  * update for year4 run
431  *
432  * Revision 1.2 2003/12/04 06:54:08 geurts
433  * introduced variables relevant to TOFp flow analysis
434  * * trackId, px and py
435  * * xvtx and yvtx
436  *
437  * Revision 1.1 2003/08/07 23:55:47 geurts
438  * first release
439  *
440  */
Int_t Make()
get tofp slat, pvpd rawdata and global data from StEvent and store in flat TTrees (ntuples) ...
Bool_t doPrintCpuInfo
control debug memory data
void initDaqMap()
set-up the default Daq-to-SlatId and Slat-to-DaqId mappings
Time-of-Flight Geometry Utilities.
Int_t getTofData(StTofCollection *)
control debug timing data
Int_t Init()
initialize ntuple and daqmap, and reset counters
Definition: Stypes.h:42
Definition: Stypes.h:40
Int_t Finish()
write and close the ntuple file
StTofpNtupleMaker(const Char_t *name="tofpNtuple")
constructor sets default parameters
~StTofpNtupleMaker()
default empty destructor
Definition: Stypes.h:45