StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
EventT.cxx
1 #include "StEvent.h"
2 #include "StPrimaryVertex.h"
3 #include "StEventInfo.h"
4 #include "StEventSummary.h"
5 #include "StTrack.h"
6 #include "StTrackNode.h"
7 #include "StPrimaryTrack.h"
8 #include "StGlobalTrack.h"
9 #include "StTrackDetectorInfo.h"
10 #include "StTrackGeometry.h"
11 #include "StSvtHit.h"
12 #include "StSsdHit.h"
13 #include "TGeoMatrix.h"
14 #include "StarRoot/THelixTrack.h"
15 #include "EventT.h"
16 #include "TrackT.h"
17 #include "HitT.h"
18 #include "TKey.h"
19 #include "TDirectory.h"
20 #include "TClass.h"
21 //#include "StSvtPool/SvtMatchedTree/SvtMatchedTree.h"
22 #include "TRVector.h"
23 #include "TRSymMatrix.h"
24 #include "StSvtBarrelHitCollection.h"
25 #include "StSvtHitCollection.h"
26 #include "StSvtLadderHitCollection.h"
27 #include "StSsdLadderHitCollection.h"
28 #include "StSsdWaferHitCollection.h"
29 #include "StSsdHitCollection.h"
30 #include "StDbUtilities/St_svtRDOstrippedC.h"
31 #include "StDedxPidTraits.h"
32 // #include "StDbUtilities/StSvtCoordinateTransform.hh"
33 // #include "StDbUtilities/StSvtLocalCoordinate.hh"
34 // #include "StDbUtilities/StSvtWaferCoordinate.hh"
35 #include "StDbUtilities/St_svtHybridDriftVelocityC.h"
36 ClassImp(EventTHeader);
37 ClassImp(EventT);
38 ClassImp(TrackT);
39 ClassImp(HitT);
40 
41 TClonesArray *EventT::fgTracks = 0;
42 TClonesArray *EventT::fgHits = 0;
43 THashList *EventT::fRotList = 0;
44 
45 
46 static Int_t _debug = 0;
47 //______________________________________________________________________________
48 EventT::EventT() : fIsValid(kFALSE)
49 {
50  // Create an EventT object.
51  // When the constructor is invoked for the first time, the class static
52  // variable fgTracks is 0 and the TClonesArray fgTracks is created.
53 
54  if (!fgTracks) fgTracks = new TClonesArray("TrackT", 1000);
55  fTracks = fgTracks;
56  fNtrack = 0;
57  if (!fgHits) fgHits = new TClonesArray("HitT", 1000);
58  fHits = fgHits;
59  fNhit = 0;
60 }
61 
62 //______________________________________________________________________________
63 EventT::~EventT()
64 {
65  Clear();
66  SafeDelete(fRotList);
67 }
68 
69 //______________________________________________________________________________
70 Int_t EventT::Build(StEvent *pEventT, UInt_t MinNoHits, Double_t pCut) {
71  static const Int_t NoFitPointCutForGoodTrackT = 15;
72  Int_t iok = 1;
73  fIsValid = kFALSE;
74  if (! pEventT) return iok;
75  UInt_t NprimVtx = pEventT->numberOfPrimaryVertices();
76  if (! NprimVtx) return iok;
77  StPrimaryVertex *pVertex=0;
78  Int_t ibest = -1;
79  Int_t nBestTracks = -1;
80  Int_t nGoodTpcTracks;
81  for (UInt_t ipr=0; ipr < NprimVtx ; ipr++) {
82  pVertex = pEventT->primaryVertex(ipr);
83  if (! pVertex) continue;
84  UInt_t nDaughters = pVertex->numberOfDaughters();
85  nGoodTpcTracks = 0;
86  for (UInt_t i=0; i < nDaughters; i++) {
87  StTrack* pTrackT = pVertex->daughter(i);
88  if ( pTrackT->fitTraits().numberOfFitPoints(kTpcId) >= NoFitPointCutForGoodTrackT) nGoodTpcTracks++;
89  }
90  if (nBestTracks < nGoodTpcTracks) {nBestTracks = nGoodTpcTracks; ibest = ipr;}
91  }
92  if (ibest < 0) return iok;
93  pVertex = pEventT->primaryVertex(ibest);
94  StSvtHitCollection* SvtHitCollection = pEventT->svtHitCollection();
95  StSsdHitCollection* SsdHitCollection = pEventT->ssdHitCollection();
96  St_svtRDOstrippedC *svtRDOs = St_svtRDOstrippedC::instance();
97  if (! SvtHitCollection && ! SsdHitCollection) { cout << "No SVT & SSD HitT Collections" << endl; return iok;}
98  const StThreeVectorF& xyzP = pVertex->position();
99  fVertex[0] = xyzP.x();
100  fVertex[1] = xyzP.y();
101  fVertex[2] = xyzP.z();
102  StMatrixF vCM = pVertex->covariantMatrix();
103  fCovariantMatrix[0] = vCM(1,1); // left triangular
104  fCovariantMatrix[1] = vCM(1,2);
105  fCovariantMatrix[2] = vCM(2,2);
106  fCovariantMatrix[3] = vCM(1,3);
107  fCovariantMatrix[4] = vCM(2,3);
108  fCovariantMatrix[5] = vCM(3,3);
109  fNPTracks = pVertex->numberOfDaughters();
110  //Save current Object count
111  Clear();
112 
113  StEventInfo* info = pEventT->info();
114  Int_t ev = 0, run = 0, time = 0;
115  if (info) {
116  ev = info->id();
117  run = info->runId();
118  time = info->time();
119  }
120  StEventSummary* summary = pEventT->summary();
121  Double32_t field = 0;
122  if (summary) field = summary->magneticField();
123  SetHeader(ev,run,time,field);
124  SetFlag(1);
125  // Create and Fill the TrackT objects
126  for (UInt_t t = 0; t < fNPTracks; t++) {
127  StTrack *pTrackT = pVertex->daughter(t);
128  if (! pTrackT) continue;
129  StTrackNode *node = pTrackT->node();
130  if (! node) continue;
131 #ifdef __USE_GLOBAL__
132  StGlobalTrack *gTrackT = (StGlobalTrack *) node->track(global);
133  if (! gTrackT) continue;
134 #endif
135  StTrackDetectorInfo* dInfo = pTrackT->detectorInfo();
136  if (! dInfo) continue;
137  static StDetectorId ids[2] = {kSvtId, kSsdId};
138  UInt_t Nsp = dInfo->numberOfPoints(ids[0]) + dInfo->numberOfPoints(ids[1]);
139  if (MinNoHits > 0 && Nsp < MinNoHits) continue;
140  UInt_t npoints = dInfo->numberOfPoints() + 100*(dInfo->numberOfPoints(ids[0]) + 10*dInfo->numberOfPoints(ids[1]));
141  UInt_t nPpoints =
142  pTrackT->numberOfPossiblePoints() +
143  100*(pTrackT->numberOfPossiblePoints(ids[0]) + 10*pTrackT->numberOfPossiblePoints(ids[1]));
144  StThreeVectorD g3 = pTrackT->geometry()->momentum();
145 #ifdef __USE_GLOBAL__
146  StThreeVectorD g3Gl = gTrackT->geometry()->momentum();
147 #endif
148  Double_t pT = g3.perp();
149  Double_t pMom = g3.mag();
150  if (pMom < pCut) continue;
151  TrackT *track = AddTrackT();
152  Double_t InvpT = 0;
153  Double_t TanL = 999999;
154  if (TMath::Abs(pT) > 1.e-7) {
155  InvpT = pTrackT->geometry()->charge()/pT;
156  TanL = g3.z()/pT;
157  }
158  track->SetInvpT(InvpT);
159  track->SetPhi(TMath::ATan2(g3.y(),g3.x()));
160  track->SetTanL(TanL);
161  static const Double_t EC = 2.9979251E-4;
162  Double_t Rho = - EC*InvpT*field;
163  track->SetRho(Rho);
164  Double_t I70 = 0;
165  Double_t TrackLength70 = 0;
166  StSPtrVecTrackPidTraits &traits = pTrackT->pidTraits();
167  UInt_t size = traits.size();
168  StDedxPidTraits *pid;
169  for (UInt_t i = 0; i < size; i++) {
170  if (! traits[i]) continue;
171  if ( traits[i]->IsZombie()) continue;
172  pid = dynamic_cast<StDedxPidTraits*>(traits[i]);
173  if (! pid || pid->method() != kTruncatedMeanId) continue;
174  I70 = pid->mean();
175  TrackLength70 = pid->length();
176  }
177  track->SetdEdx(I70,TrackLength70);
178 #ifdef __USE_GLOBAL__
179  Double_t pTGl = g3Gl.perp();
180  // Double_t pMomGl = g3Gl.mag();
181  Double_t InvpTGl = 0;
182  Double_t TanLGl = 999999;
183  if (TMath::Abs(pTGl) > 1.e-7) {
184  InvpTGl = gTrackT->geometry()->charge()/pTGl;
185  TanLGl = g3Gl.z()/pTGl;
186  }
187  track->SetInvpTGl(InvpTGl);
188  track->SetPhiGl(TMath::ATan2(g3Gl.y(),g3Gl.x()));
189  track->SetTanLGl(TanLGl);
190  Double_t RhoGl = - EC*InvpT*field;
191  track->SetRhoGl(RhoGl);
192  track->SetNPpoint(nPpoints);
193 #endif
194 
195  track->SetN(0);
196  track->SetNpoint(npoints);
197  track->SetNPpoint(nPpoints);
198 #if 0
199  const Double_t XyzDirRho[6] = {fVertex[0], fVertex[1], fVertex[2], track->GetTanL(), track->GetPhi(), track->GetRho()};
200  const Double_t XyzDirRhoGl[6] = {fVertex[0], fVertex[1], fVertex[2], track->GetTanL(), track->GetPhi(), track->GetRho()};
201 #endif
202  THashList *fRotList = RotMatrices();
203  if (! fRotList) continue;
204  TIter next(fRotList);
205  TGeoHMatrix *comb = 0;
206  while ((comb = (TGeoHMatrix *) next())) {
207  TString combName(comb->GetName());
208  if (! combName.BeginsWith("R")) continue;
209  Int_t Id;
210  sscanf(comb->GetName()+1,"%04i",&Id);
211  UInt_t Ladder = Id%100;
212  UInt_t Layer = Id/1000; if (Layer > 7) Layer = 7;
213  UInt_t Wafer = (Id - 1000*Layer)/100;
214  UInt_t Barrel = (Layer - 1)/2 + 1;
215  if (_debug) {
216  cout << comb->GetName() << "\tLayer/Ladder/Wafer = " << Layer << "/" << Ladder << "/" << Wafer << endl;
217  comb->Print();
218  }
219  static Double_t dz[2] = {3.00, 2.10};
220  static Double_t dx[2] = {3.00, 3.65};
221  Int_t k = 0; // svt
222  if (Layer > 6) k = 1;
223 #if 0
224  // check straight line
225  const Double_t *X = &XyzDirRho[0]; // Global
226  Double_t x0[3]; // Local
227  comb->MasterToLocal(X,x0);
228  Double_t CosL = 1./TMath::Sqrt(1. + XyzDirRho[3]*XyzDirRho[3]);
229  Double_t SinL = XyzDirRho[3]*CosL;
230  Double_t CosP = TMath::Cos(XyzDirRho[4]);
231  Double_t SinP = TMath::Sin(XyzDirRho[4]);
232  Double_t dir[3], d[3];
233  dir[0] = CosL*CosP;
234  dir[1] = CosL*SinP;
235  dir[2] = SinL;
236  comb->MasterToLocalVect(dir,d);
237  if (TMath::Abs(d[2]) < 1.e-7) continue;
238  Double_t s = - x0[2]/d[2]; if (_debug) cout << "straight line s " << s << endl;;
239  if (s < 0) continue;
240  TRVector uvP(3,x0);
241  TRVector tuvP(2,d[0]/d[2],d[1]/d[2]);
242  TRVector DD(3,d);
243  DD *= s;
244  uvP += DD; if (_debug) cout << "straight line uvP\t" << uvP << "\ttuvP \t" << tuvP << endl;
245  // svt ssd
246  if (TMath::Abs(uvP[0]) > dx[k] + 1.0) continue;
247  if (TMath::Abs(uvP[1]) > dz[k] + 1.0) continue;
248  // if (_debug) {
249 #endif
250  Double_t *rot = comb->GetRotationMatrix();
251  Double_t *tra = comb->GetTranslation();
252  const StThreeVectorD normal(rot[2], rot[5], rot[8]);
253  const StThreeVectorD middle(tra);
254  StPhysicalHelixD helixI = pTrackT->geometry()->helix();
255  Double_t sh = helixI.pathLength(middle, normal); if (_debug) cout << "StHelix sh " << sh << endl;
256  if (sh <= 0 || sh > 1e3) continue;
257  StThreeVectorD xyzG = helixI.at(sh); if (_debug) cout << "StHelix xyzG\t" << xyzG << endl;
258  Double_t xyzGPred[3] = {xyzG.x(),xyzG.y(),xyzG.z()};
259  Double_t uvPred[3];
260  comb->MasterToLocal(xyzGPred,uvPred);
261  TRVector xyzL(3,uvPred); if (_debug) cout << "StHelix xyzL\t" << xyzL << endl;
262  Double_t dirGPred[3] = {helixI.cx(sh),helixI.cy(sh),helixI.cz(sh)};
263  Double_t dxyzL[3];
264  comb->MasterToLocalVect(dirGPred,dxyzL);
265  Double_t tuvPred[2] = {dxyzL[0]/dxyzL[2], dxyzL[1]/dxyzL[2]};
266  if (_debug) cout << "StHelix tU/tV = " << tuvPred[0] << "\t" << tuvPred[1] << endl;
267 #ifdef __USE_GLOBAL__
268 
269  StPhysicalHelixD helixG = gTrackT->geometry()->helix();
270  sh = helixG.pathLength(middle, normal); if (_debug) cout << "StHelix sh " << sh << endl;
271  if (sh < 0 || sh > 1e3) continue;
272  xyzG = helixG.at(sh); if (_debug) cout << "StHelix xyzG\t" << xyzG << endl;
273  Double_t xyzGPredGl[3] = {xyzG.x(),xyzG.y(),xyzG.z()};
274  Double_t uvPredGl[3];
275  comb->MasterToLocal(xyzGPredGl,uvPredGl);
276  TRVector xyzLGl(3,uvPredGl); if (_debug) cout << "StHelix xyzL\t" << xyzL << endl;
277  Double_t dirGPredGl[3] = {helixG.cx(sh),helixG.cy(sh),helixG.cz(sh)};
278  comb->MasterToLocalVect(dirGPredGl,dxyzL);
279  Double_t tuvPredGl[2] = {dxyzL[0]/dxyzL[2], dxyzL[1]/dxyzL[2]};
280  if (_debug) cout << "StHelix tU/tV = " << tuvPredGl[0] << "\t" << tuvPredGl[1] << endl;
281 #endif
282 
283  if (TMath::Abs(uvPred[0]) > dx[k] + 1.0) continue;
284  if (TMath::Abs(uvPred[1]) > dz[k] + 1.0) continue;
285  StPtrVecHit &hvec = pTrackT->detectorInfo()->hits();
286  UInt_t NoHits = hvec.size();
287  StHit *hit = 0;
288  for (UInt_t j = 0; j < NoHits; j++) {
289  StHit *hitc = hvec[j];
290  if (! hitc) continue;
291  if (hitc->detector() == kSvtId) {
292  StSvtHit *htSvt = (StSvtHit *) hitc;
293  if (htSvt->barrel() == Barrel &&
294  htSvt->ladder() == Ladder &&
295  htSvt->wafer() == Wafer) {hit = hitc; break;}
296  }
297  if (hitc->detector() == kSsdId) {
298  StSsdHit *htSsd = (StSsdHit *) hitc;
299  if (htSsd->ladder() == Ladder &&
300  htSsd->wafer() == Wafer) {hit = hitc; break;}
301  }
302  }
303  HitT *ht = AddHitT();
304  UInt_t Hybrid = 1;
305  if (uvPred[0] >= 0) Hybrid = 2;
306  ht->SetId(Barrel,Layer,Ladder,Wafer,Hybrid);
307  ht->SetisTrack(t+1);
308  ht->SetUVPred (uvPred[0],uvPred[1]);
309  ht->SettUVPred(tuvPred[0],tuvPred[1]);
310  ht->SetXyzG(xyzGPred);
311  ht->SetDirG(dirGPred);
312 #ifdef __USE_GLOBAL__
313  ht->SetUVPredGl (uvPredGl[0],uvPredGl[1]);
314  ht->SettUVPredGl(tuvPredGl[0],tuvPredGl[1]);
315  ht->SetXyzGl(xyzGPredGl);
316  ht->SetDirGl(dirGPredGl);
317 #endif
318  if (hit) {
319  SetHitT(ht, hit, comb, track);
320  ht->SetisFitted(t+1);
321  }
322  StSPtrVecSvtHit *hitsvt = 0;
323  StSPtrVecSsdHit *hitssd = 0;
324  Int_t NoHitPerTrack = 0;
325  if (Layer < 7) { // svt
326  if (! SvtHitCollection) continue;
327  StSvtBarrelHitCollection* barrelCollection = SvtHitCollection->barrel(Barrel-1);
328  if (! barrelCollection) continue;
329  StSvtLadderHitCollection *ladderCollection = barrelCollection->ladder(Ladder-1);
330  if (! ladderCollection) continue;
331  if (svtRDOs && svtRDOs->svtRDOstrippedStatus(Barrel,Ladder,Wafer)) continue;
332  StSvtWaferHitCollection* waferCollection = ladderCollection->wafer(Wafer-1);
333  if (! waferCollection) continue;
334  hitsvt = &waferCollection->hits();
335  } else { // ssd
336  if (! SsdHitCollection) continue;
337  StSsdLadderHitCollection* ladderCollection = SsdHitCollection->ladder(Ladder-1);
338  if (! ladderCollection) continue;
339 
340  StSsdWaferHitCollection* waferCollection = ladderCollection->wafer(Wafer-1);
341  if (! waferCollection) continue;
342  hitssd = &waferCollection->hits();
343  }
344  if (! hitsvt && ! hitssd) continue;
345  if (hitsvt) NoHits = hitsvt->size();
346  if (hitssd) NoHits = hitssd->size();
347  ht->SetNofHits(NoHits);
348  if (! NoHits) continue;
349 
350  for (UInt_t l = 0; l < NoHits; l++) {
351  hit = 0;
352  if (hitsvt) hit = (*hitsvt)[l];
353  if (hitssd) hit = (*hitssd)[l];
354  if (hit) {
355  //if (hit->flag()>=4) continue;
356  //if (hit->flag()< 0) continue;
357  // cout << "hitFlag=" << hit->flag() << endl;
358  HitT *h = AddHitT();
359  h->SetHitFlag(UInt_t(hit->flag()));
360  h->SetUVPred (uvPred[0],uvPred[1]);
361  h->SettUVPred(tuvPred[0],tuvPred[1]);
362  h->SetXyzG(xyzGPred);
363  h->SetDirG(dirGPred);
364 #ifdef __USE_GLOBAL__
365 
366  h->SetUVPredGl (uvPredGl[0],uvPredGl[1]);
367  h->SettUVPredGl(tuvPredGl[0],tuvPredGl[1]);
368  h->SetXyzGl(xyzGPredGl);
369  h->SetDirGl(dirGPredGl);
370 #endif
371  SetHitT(h, hit, comb, track);
372  NoHitPerTrack++;
373  h->SetHitPerTrack(NoHitPerTrack);
374  // SetHitT(h, hit, comb, track, &TPDeriv);
375  }
376  }
377  }
378  }
379  fIsValid = kTRUE;
380  iok = 0;
381  return iok;
382 }
383 //______________________________________________________________________________
384 TrackT *EventT::AddTrackT()
385 {
386  // Add a new track to the list of tracks for this event.
387  // To avoid calling the very time consuming operator new for each track,
388  // the standard but not well know C++ operator "new with placement"
389  // is called. If tracks[i] is 0, a new TrackT object will be created
390  // otherwise the previous TrackT[i] will be overwritten.
391 
392  TClonesArray &tracks = *fTracks;
393  TrackT *track = new(tracks[fNtrack++]) TrackT();
394  //Save reference to last TrackT in the collection of Tracks
395  return track;
396 }
397 //______________________________________________________________________________
398 HitT *EventT::AddHitT()
399 {
400  // Add a new hit to the list of hits for this event.
401  // To avoid calling the very time consuming operator new for each hit,
402  // the standard but not well know C++ operator "new with placement"
403  // is called. If hits[i] is 0, a new HitT object will be created
404  // otherwise the previous HitT[i] will be overwritten.
405 
406  TClonesArray &hits = *fHits;
407  HitT *hit = new(hits[fNhit++]) HitT();
408  //Save reference to last HitT in the collection of Hits
409  return hit;
410 }
411 
412 //______________________________________________________________________________
413 void EventT::Clear(Option_t * /*option*/)
414 {
415  fTracks->Clear("C"); //will also call TrackT::Clear
416  fHits->Clear("C"); //will also call HitT::Clear
417 }
418 
419 //______________________________________________________________________________
420 void EventT::Reset(Option_t * /*option*/)
421 {
422  // Static function to reset all static objects for this event
423  // fgTracks->Delete(option);
424 
425  delete fgTracks; fgTracks = 0;
426  delete fgHits; fgHits = 0;
427 }
428 
429 //______________________________________________________________________________
430 void EventT::SetHeader(Int_t i, Int_t run, Int_t date, Double32_t field)
431 {
432  fNtrack = 0;
433  fNhit = 0;
434  fEvtHdr.Set(i, run, date, field);
435 }
436 //________________________________________________________________________________
437 void EventT::Print(Option_t *opt) const {
438  cout << "Run/EventT\t" << fEvtHdr.GetRun() << "/" << fEvtHdr.GetEvtNum() << "\tDate " << fEvtHdr.GetDate()
439  << "\tField " << fEvtHdr.GetField() << endl;
440  cout << "Total no. tracks " << GetTotalNoTracks() << "\tRecorded tracks " << GetNtrack()
441  << "\tRecorded hits " << GetNhit() << endl;
442  TRVector vertex(3,GetVertex());
443  TRSymMatrix cov(3,GetCovMatrix());
444  cout << "Primary vertex " << vertex << endl;
445  cout << "Its cov. matrix " << cov << endl;
446  for (UInt_t i = 0; i < GetNtrack(); i++) {cout << i << "\t"; GetTrackT(i)->Print();}
447  for (UInt_t i = 0; i < GetNhit(); i++) {cout << i << "\t"; GetHitT(i)->Print();}
448 
449 }
450 //________________________________________________________________________________
451 HitT *EventT::SetHitT(HitT *h, StHit *hit, TGeoHMatrix *comb, TrackT *track) {
452  struct RDO_t {
453  const Char_t *name;
454  Int_t ladder, barrel;
455  Int_t rdo;
456  };
457  const Int_t NRDOS = 36;
458  const RDO_t RDOS[NRDOS] = {
459  {"L01B1", 1,1, 1},{"L02B1", 2,1, 2},{"L03B1", 3,1, 4},{"L04B1", 4,1, 5},{"L05B1", 5,1, 7},
460  {"L06B1", 6,1, 8},{"L07B1", 7,1,10},{"L08B1", 8,1,11},{"L01B2", 1,2, 3},{"L02B2", 2,2, 3},
461  {"L03B2", 3,2, 3},{"L04B2", 4,2, 6},{"L05B2", 5,2, 6},{"L06B2", 6,2, 6},{"L07B2", 7,2, 9},
462  {"L08B2", 8,2, 9},{"L09B2", 9,2, 9},{"L10B2",10,2,12},{"L11B2",11,2,12},{"L12B2",12,2,12},
463  {"L01B3", 1,3, 1},{"L02B3", 2,3, 1},{"L03B3", 3,3, 2},{"L04B3", 4,3, 2},{"L05B3", 5,3, 4},
464  {"L06B3", 6,3, 4},{"L07B3", 7,3, 5},{"L08B3", 8,3, 5},{"L09B3", 9,3, 7},{"L10B3",10,3, 7},
465  {"L11B3",11,3, 8},{"L12B3",12,3, 8},{"L13B3",13,3,10},{"L14B3",14,3,10},{"L15B3",15,3,11},
466  {"L16B3",16,3,11}
467  };
468  UInt_t B = 0, L = 0, l = 0, W = 0, H = 0;
469  Int_t rdo = 0;
470  h->SetRDO(rdo);
471  if (hit->detector() == kSvtId) {
472  StSvtHit *ht = (StSvtHit *) hit;
473  B = ht->barrel();
474  L = ht->layer();
475  l = ht->ladder();
476  W = ht->wafer();
477  H = ht->hybrid();
478  h->SetId(B,L,l,W,H);
479  h->SetuvD(ht->localPosition(0), ht->localPosition(1));
480  h->SetAnode(ht->anode());
481  h->SetTimeB(ht->timebucket());
482  for (Int_t r = 0; r < NRDOS; r++) {
483  if (h->Ladder() == RDOS[r].ladder && h->Barrel() == RDOS[r].barrel) {
484  rdo = RDOS[r].rdo;
485  break;
486  }
487  }
488  if (rdo) h->SetRDO(rdo);
489  St_svtHybridDriftVelocityC *d = St_svtHybridDriftVelocityC::instance();
490  if (d && d->p(B,l,W,H)) {
491  h->SetLM(d->CalcU(B,l,W,H,ht->timebucket(),ht->anode()),
492  d->CalcV(H,ht->anode()));
493  h->SetuHat(d->uHat(B,l,W,H,ht->timebucket()));
494  }
495  }
496  h->SetUsedInFit(hit->usedInFit());
497  if (hit->detector() == kSsdId) {
498  StSsdHit *ht = (StSsdHit *) hit;
499  B = 4;
500  L = 7;
501  l = ht->ladder();
502  W = ht->wafer();
503  h->SetId(B,L,l,W,H);
504  h->SetLM(-ht->localPosition(0), ht->localPosition(1));
505  h->SetuHat(-ht->localPosition(0));
506  h->SetAnode(0);
507  h->SetTimeB(0);
508  }
509  StThreeVectorF position = hit->position();
510  Double_t xyzG[3] = {position.x(),position.y(),position.z()};
511  h->SetGC(xyzG[0],xyzG[1],xyzG[2]);
512  Double_t xyzL[3] = {0,0,0};
513  comb->MasterToLocal(xyzG,xyzL);
514  // if (TMath::Abs(xyzL[2]) > 0.1) continue;
515  Double_t uvw[3] = {h->GetU(),h->GetV(),0};
516  comb->LocalToMaster(uvw,xyzG);
517  h->Set(xyzG,uvw);
518  Double_t *rot = comb->GetRotationMatrix();
519  h->SetWG(rot[2],rot[5],rot[8]);
520  // Int_t IdH = GetIndexOfHitT(h);
521  Int_t IdH = fNhit - 1;
522  track->SetHitTId(IdH);
523  Double_t invpT = track->GetInvpT();
524  if (TMath::Abs(invpT) < 1e-7) invpT = 1e-7;
525  h->SetpT(1./invpT);
526  h->SetMom(track->GetMomentum());
527  h->SetWG(rot[2],rot[5],rot[8]);
528  TGeoHMatrix *rotL = (TGeoHMatrix *) RotMatrices()->FindObject(Form("WL%s",comb->GetName()+1));
529  Double_t xyzLadder[3] = {0,0,0};
530  if (rotL) {
531 
532  rotL->LocalToMaster(uvw,xyzLadder);
533  h->SetL(xyzLadder[0],xyzLadder[1],xyzLadder[2]);
534  Double_t uvwP[3] = {h->GetPredU(),h->GetPredV(),0};
535  rotL->LocalToMaster(uvwP,xyzLadder);
536  h->SetXyzL(xyzLadder);
537 #ifdef __USE_GLOBAL__
538 
539  Double_t uvwPGl[3] = {h->GetPredGlU(),h->GetPredGlV(),0};
540  rotL->LocalToMaster(uvwPGl,xyzLadder);
541  h->SetXyzGlL(xyzLadder);
542 #endif
543  } else {
544 
545  cout << Form("WL%s",comb->GetName()+1) << " has not been found" << endl;
546  h->SetL(xyzLadder[0],xyzLadder[1],xyzLadder[2]);
547  h->SetXyzL(xyzLadder);
548 #ifdef __USE_GLOBAL__
549  h->SetXyzGlL(xyzLadder);
550 #endif
551  }
552  return h;
553 }
554 //________________________________________________________________________________
555 void TrackT::Print(Option_t *opt) const {
556  cout << "TrackT: InvpT " << fInvpT << "\tTanL " << fTanL
557  << "\tPhi " << fPhi << "\tRho " << fRho
558  << "\tNpoint " << fNpoint << "\tNsp " << fNsp << endl;
559  for (UInt_t i = 0; i < fNsp; i++) cout << "\t" << fIdHitT[i];
560  cout << endl;
561 }
562 //________________________________________________________________________________
563 void HitT::SetId(Int_t B, Int_t L, Int_t l, Int_t W, Int_t H) {
564  struct Geom_t {
565  Int_t Barrel;
566  Int_t Layer;
567  Int_t NoLadders;
568  Int_t NoWafers;
569  };
570  const Int_t NoLayers = 7;
571  //Barrel, Layer Nladder Nwafer
572  static const Geom_t SvtSsdConfig[NoLayers] =
573  { {1, 1, 8, 4}, // even
574  {1, 2, 8, 4}, // odd
575  {2, 3, 12, 6}, // event
576  {2, 4, 12, 6}, // odd
577  {3, 5, 16, 7}, // even
578  {3, 6, 16, 7}, // odd
579  {4, 7, 20, 16} // Ssd
580  };
581  static const Int_t ssdSector[20] = {// 100*sector + ladder
582  101, 102,
583  203, 204, 205, 206, 207, 208, 209,
584  310, 311, 312,
585  413, 414, 415, 416, 417, 418, 419,
586  120
587  };
588  barrel = B; layer = L; ladder = l; wafer = W; hybrid = H;
589  if (barrel == 0) Id = 7000 + 100*wafer + ladder;
590  else Id = 1000*layer + 100*wafer + ladder;
591  sector = -1;
592  if (layer < 7) {
593  sector = 0;
594  if (ladder > SvtSsdConfig[layer-1].NoLadders/2) sector = 1;
595  } else {sector = ssdSector[ladder-1]/100 + 1; barrel = 4;}
596 }
597 #if 0
598 //________________________________________________________________________________
599 void HitT::SetId(StHit *shit) {
600  if (shit) {
601  if (shit->detector() == kSvtId) {
602  StSvtHit *hit = (StSvtHit *) shit;
603  B = hit->barrel();
604  L = hit->layer();
605  l = hit->ladder();
606  W = hit->wafer();
607  H = hit->hybrid();
608  }
609  if (shit->detector() == kSsdId) {
610  StSsdHit *hit = (StSsdHit *) shit;
611  B = 4;
612  L = 7;
613  l = hit->ladder();
614  W = hit->wafer();
615  }
616  }
617  SetId(B,L,l,W,H);
618 }
619 #endif
620 //________________________________________________________________________________
621 void HitT::Print(Option_t *opt) const {
622  cout << "HitT: Id " << Id << "\tpT = " << pT << "\tmomentum " << pMom << endl;
623  TRVector glob(3,&xG); cout << "Global :" << glob << endl;
624  cout << "Local u/v/w " << u << "/ " << v << "/ " << w << endl;
625  cout << "Prediction uP/vP " << uP << "/ " << vP << "\ttuP/tvP " << tuP << "/ " << tvP << endl;
626 }
627 //________________________________________________________________________________
628 void EventT::RestoreListOfRotations() {
629  if (fRotList) return;
630  if (! gDirectory) return;
631  fRotList = new THashList(100,0);
632  fRotList->SetOwner();
633  TIter nextkey(gDirectory->GetListOfKeys() );
634  TKey *key;
635  while ((key = (TKey*) nextkey())) {
636  TObject *obj = key->ReadObj();
637  if ( obj->IsA()->InheritsFrom( "TGeoHMatrix" ) ) {
638  fRotList->Add(obj);
639  }
640  }
641 }
Definition: StHit.h:125
Definition: TrackT.h:8
Definition: HitT.h:7
Definition: EventT.h:32
double cx(double s) const
pointing vector of helix at point s
Definition: StHelix.hh:203
pair< double, double > pathLength(double r) const
path length at given r (cylindrical r)
Definition: StHelix.cc:351