StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
StEbyET0.cxx
1 // //
3 // StEbyET0 redoes the space charge correction //
4 // //
6 
7 #include "StEvent/StEvent.h"
8 #include "StEvent/StTpcHit.h"
9 #include "StEvent/StTpcHitCollection.h"
10 #include "StEvent/StTriggerData.h"
11 #include "StDetectorDbMaker/St_EbyET0C.h"
12 /*
13 #include "StDetectorDbMaker/St_tpcPadConfigC.h"
14 #include "StDetectorDbMaker/St_tpcTimeBucketCorC.h"
15 #include "StDbUtilities/StTpcCoordinateTransform.hh"
16 #include "StTpcDb/StTpcDb.h"
17 */
18 #include "TFile.h"
19 #include "TNtupleD.h"
20 #include "StEbyET0.h"
21 
22 StEbyET0* StEbyET0::mInstance = 0;
23 
24 ClassImp(StEbyET0)
25 
26 //_____________________________________________________________________________
27 StEbyET0::StEbyET0() : mRunId(0), mEventId(0), mTime(0.), mTree(0) {}
28 //_____________________________________________________________________________
29 StEbyET0::~StEbyET0() {
30  finishCalib();
31  mInstance=0;
32 }
33 //_____________________________________________________________________________
34 void StEbyET0::finishCalib() {
35  if (!mTree) return;
36  mTree->Write();
37  delete (mTree->GetDirectory()->GetFile());
38 }
39 //_____________________________________________________________________________
40 double StEbyET0::getTime(StEvent* event, int mode) {
41 
42  // check for same event
43  if (event->runId() == mRunId && event->id() == mEventId) return mTime;
44 
45  mRunId = event->runId();
46  mEventId = event->id();
47  mTime = 0.; // default event time should always be zero, not correcting for some other global T0
48 
49  // check for calibration mode
50  // fill tree without mTime subtracted
51  if (mode == 1) fillTree(event);
52 
53  // determine event time
54  double info[12];
55  St_EbyET0C* ebyeTable = St_EbyET0C::instance();
56  for (int row = 0; row < ebyeTable->GetNRows(); row++) {
57  memset(info,0,12*sizeof(double));
58  int detector = ebyeTable->detector(row);
59  if (detector < 0) break; // end of active rows in the table
60 
61  double coordinate = -9e23;
62 
63  // determine which time coordinate to use
64  switch (detector) {
65  case 0 : // first 4 use average of trigger detectors
66  case 1 :
67  case 2 :
68  case 3 : getTriggerInfo(event,static_cast<trigDetType> (detector),info);
69  if (info[0] > 0 && info[1] > 0)
70  coordinate = 0.5*(info[0]+info[1]);
71  break;
72  case 4 : // next 4 use east trigger detectors
73  case 5 :
74  case 6 :
75  case 7 : getTriggerInfo(event,static_cast<trigDetType> (detector-4),info);
76  if (info[0] > 0)
77  coordinate = info[0];
78  break;
79  case 8 : // next 4 use west trigger detectors
80  case 9 :
81  case 10 :
82  case 11 : getTriggerInfo(event,static_cast<trigDetType> (detector-8),info);
83  if (info[1] > 0)
84  coordinate = info[1];
85  break;
86  }
87 
88  // determine if coordinate is in acceptable range
89  if (coordinate < -8e23 ||
90  (ebyeTable->min(row) < ebyeTable->max(row) &&
91  (ebyeTable->min(row) > coordinate ||
92  ebyeTable->max(row) < coordinate))) continue;
93 
94  // use only this coordinate to determine a time correction
95  mTime = ebyeTable->time(row,coordinate);
96  break;
97  }
98 
99  // check for QA mode
100  // fill tree with mTime added
101  if (mode == 2) fillTree(event);
102 
103  return mTime;
104 }
105 //_____________________________________________________________________________
106 void StEbyET0::fillTree(StEvent* event) {
107 
108  if (!mTree) {
109  new TFile(Form("EbyET0.%d.%d.root",mRunId,mEventId),"RECREATE");
110  mTree = new TNtupleD("ebyeT0tree","Event by Event T0 tree",
111  "run:event:time"
112  ":vpde:vpdw:epde:epdw:bbce:bbcw:zdce:zdcw"
113  ":tpctie:tpcnie:tpcrie"
114  ":tpctiw:tpcniw:tpcriw"
115  ":tpctoe:tpcnoe:tpcroe"
116  ":tpctow:tpcnow:tpcrow");
117  }
118 
119  // reduce AutoSave frequency for larger trees
120  int treeEntries = mTree->GetEntriesFast();
121  switch (treeEntries) {
122  case 0 : mTree->SetAutoSave( 5); break;
123  case 25 : mTree->SetAutoSave( 20); break;
124  case 250 : mTree->SetAutoSave(100); break;
125  case 1250 : mTree->SetAutoSave(250); break;
126  }
127 
128  double info[23];
129  memset(info,0,23*sizeof(double));
130  info[0] = (double) mRunId;
131  info[1] = (double) mEventId;
132  info[2] = (double) (event->time());
133  getTriggerInfo(event,kVPD,&(info[3]));
134  getTriggerInfo(event,kEPD,&(info[5]));
135  getTriggerInfo(event,kBBC,&(info[7]));
136  getTriggerInfo(event,kZDC,&(info[9]));
137  getTpcInfo(event,&(info[11]));
138  mTree->Fill(info);
139 
140  return;
141 }
142 //_____________________________________________________________________________
143 void StEbyET0::getTriggerInfo(StEvent* event, trigDetType trigDet, double* info) {
144 
145  StTriggerData* trigdata = event->triggerData();
146  if (!trigdata) return;
147  switch (trigDet) {
148  case kVPD : info[0] = (double) (trigdata->vpdEarliestTDC(StBeamDirection::east,0));
149  info[1] = (double) (trigdata->vpdEarliestTDC(StBeamDirection::west,0));
150  break;
151  case kEPD : info[0] = (double) (trigdata->epdEarliestTDC(StBeamDirection::east,0));
152  info[1] = (double) (trigdata->epdEarliestTDC(StBeamDirection::west,0));
153  break;
154  case kBBC : info[0] = (double) (trigdata->bbcEarliestTDC(StBeamDirection::east,0));
155  info[1] = (double) (trigdata->bbcEarliestTDC(StBeamDirection::west,0));
156  break;
157  case kZDC : info[0] = (double) (trigdata->zdcEarliestTDC(StBeamDirection::east,0));
158  info[1] = (double) (trigdata->zdcEarliestTDC(StBeamDirection::west,0));
159  break;
160  default : info[0] = -999.;
161  info[1] = -999.;
162  }
163 
164  return;
165 }
166 //_____________________________________________________________________________
167 void StEbyET0::getTpcInfo(StEvent* event, double* info) {
168 
169  // counters
170  double tie = 0.; // hit times sum
171  double tiw = 0.;
172  double toe = 0.;
173  double tow = 0.;
174  double nie = 0.; // hit counts
175  double niw = 0.;
176  double noe = 0.;
177  double now = 0.;
178  double rie = 0.; // hit rows sum
179  double riw = 0.;
180  double roe = 0.;
181  double row = 0.;
182 
183 /*
184  // parameters: prompt hit windows
185  static const double innerMinTB = 7.7;
186  static const double innerMaxTB = 10.5;
187  static const double outerMinTB = 2.0;
188  static const double outerMaxTB = 5.0;
189 
190  // coordinate transform to get time of hit
191  static StTpcCoordinateTransform* mTpcTransForm = 0;
192  if (! mTpcTransForm) mTpcTransForm = new StTpcCoordinateTransform(gStTpcDb);
193  StTpcCoordinateTransform &transform = *mTpcTransForm;
194 
195  StTpcHitCollection* tpcHits = event->tpcHitCollection();
196 
197  if (tpcHits) {
198  for (unsigned int i=0; i<tpcHits->numberOfSectors(); i++) {
199  StTpcSectorHitCollection* tpcHitsSector = tpcHits->sector(i);
200  for (unsigned int j=0; j<tpcHitsSector->numberOfPadrows(); j++) {
201  StSPtrVecTpcHit& tpcHitsVec = tpcHitsSector->padrow(j)->hits();
202  for (unsigned int k=0; k<tpcHitsVec.size(); k++) {
203  StTpcHit* hit = tpcHitsVec[k];
204  if (!hit || hit->flag() != 0) continue;
205  unsigned int sector = hit->sector();
206  unsigned int padrow = hit->padrow();
207  int pad = TMath::Nint(hit->pad());
208  bool innerRow = St_tpcPadConfigC::instance()->IsRowInner(sector,padrow);
209  double tb = hit->timeBucket();
210 
211  // copied from StTpcHitMoverMaker
212  if (St_tpcTimeBucketCorC::instance()->getNumRows()) {
213  Int_t io = (innerRow ? 1 : 0);
214  Double_t noTmbks = hit->maxTmbk() - hit->minTmbk() + 1;
215  tb += St_tpcTimeBucketCorC::instance()->CalcCorrection(io, noTmbks);// units are tb
216  }
217  StTpcPadCoordinate padcoord(sector, padrow, pad, tb);
218  StTpcLocalSectorCoordinate coorS;
219  transform(padcoord,coorS,kFALSE);
220  // transform converts tb => time => z, so we have to step back to time;
221  // ignoring adding a zoffset = StTpcDb::instance->Dimenstions()->zInner/OuterOffset()
222  // because constant shifts aren't relavant for the EbyET0 calibration
223  double time = coorS.position().z() / (StTpcDb::instance()->DriftVelocity()*1e-6);
224 
225  // apply the EbyE T0 time correction (if it is already determined)
226  time += mTime;
227 
228  if (innerRow) {
229  if (tb < innerMaxTB && tb > innerMinTB) { // use only prompt hit candidates
230  if (sector <= 12) {
231  niw += 1.0; tiw += time; riw += padrow;
232  } else {
233  nie += 1.0; tie += time; rie += padrow;
234  }
235  }
236  } else {
237  // exclude some known problematic rows
238  if ((sector == 5 || sector == 19) &&
239  abs((float) (St_tpcPadConfigC::instance()->padRows(sector) - padrow) - 7.5) < 8.0)
240  continue;
241  if (tb < outerMaxTB && tb > outerMinTB) { // use only prompt hit candidates
242  if (sector <= 12) {
243  now += 1.0; tow += time; row += padrow;
244  } else {
245  noe += 1.0; toe += time; roe += padrow;
246  }
247  }
248  } // innerRow
249  } // hits-in-padrow for-loop
250  } // padrow for-loop
251  } // sector for-loop
252  } // tpcHits
253 */
254 
255  // pass back the average time, number of hits, and average padrow
256  info[ 0] = (nie > 0 ? tie/nie : 0); info[ 1] = nie; info[ 2] = (nie > 0 ? rie/nie : 0);
257  info[ 3] = (niw > 0 ? tiw/niw : 0); info[ 4] = niw; info[ 5] = (niw > 0 ? riw/niw : 0);
258  info[ 6] = (noe > 0 ? toe/noe : 0); info[ 7] = noe; info[ 8] = (noe > 0 ? roe/noe : 0);
259  info[ 9] = (now > 0 ? tow/now : 0); info[10] = now; info[11] = (now > 0 ? row/now : 0);
260 
261  return;
262 }
263 //_____________________________________________________________________________
264 // $Id: StEbyET0.cxx,v 1.3 2021/04/27 20:58:52 genevb Exp $
265 // $Log: StEbyET0.cxx,v $
266 // Revision 1.3 2021/04/27 20:58:52 genevb
267 // Handle missing trigger data
268 //
269 // Revision 1.2 2021/03/20 02:38:13 genevb
270 // Remove StTpcDb dependence for StEventUtilities
271 //
272 // Revision 1.1 2021/03/19 01:44:47 genevb
273 // Introduce Event-by-Event T0 corrections
274 //
275 //
St_EbyET0C
Definition: St_EbyET0C.h:6
Generated by   doxygen 1.8.5