StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
LaserEvent.cxx
1 //$Id: LaserEvent.cxx,v 1.10 2018/04/11 02:43:21 smirnovd Exp $
2 // $Log: LaserEvent.cxx,v $
3 // Revision 1.10 2018/04/11 02:43:21 smirnovd
4 // Enable TPC/iTPC switch via St_tpcPadConfig
5 //
6 // This is accomplished by substituting St_tpcPadPlanes with St_tpcPadConfig.
7 // A sector ID is passed to St_tpcPadConfig in order to extract parameters for
8 // either TPC or iTPC
9 //
10 // Revision 1.9 2014/03/13 21:59:44 fisyak
11 // add cluster position in Local Sector Coordinate System
12 //
13 // Revision 1.8 2014/02/13 18:21:28 fisyak
14 // Add protection against cicling in fitting
15 //
16 // Revision 1.7 2011/01/10 20:36:12 fisyak
17 // Use sector/padrow in global => local transformation
18 //
19 // Revision 1.6 2008/06/02 13:48:02 fisyak
20 // Add t0 handlers for Tpx/Tpc time offsets
21 //
22 // Revision 1.5 2007/12/10 19:54:02 fisyak
23 // Add Id and Log, correct spelling error in README
24 //
25 #include <assert.h>
26 #include "TRandom.h"
27 #include "TDirectory.h"
28 #include "LaserEvent.h"
29 #if 1
30 #include "StProbPidTraits.h"
31 #include "StDedxPidTraits.h"
32 #endif
33 #include "StVertex.h"
34 #include "StPrimaryVertex.h"
35 #include "StTrack.h"
36 #include "StPrimaryTrack.h"
37 #include "StGlobalTrack.h"
38 #include "StTrackNode.h"
39 #include "StTpcDb/StTpcDb.h"
40 #include "StDbUtilities/StTpcCoordinateTransform.hh"
41 #include "StEventTypes.h"
42 #include "TGeoMatrix.h"
43 #include "StDetectorDbMaker/St_tpcPadConfigC.h"
44 ClassImp(LaserRaft);
45 ClassImp(EventHeader);
46 ClassImp(LaserEvent);
47 ClassImp(Vertex);
48 ClassImp(Track);
49 ClassImp(Hit);
50 TClonesArray *LaserEvent::fgVertices = 0;
51 TClonesArray *LaserEvent::fgTracks = 0;
52 TClonesArray *LaserEvent::fgHits = 0;
53 TClonesArray *LaserEvent::fgFit = 0;
54 //______________________________________________________________________________
55 LaserB::LaserB(const LaserRaft &laser) : Sector(laser.Sector), Raft(laser.Raft), Bundle(laser.Bundle), Mirror(laser.Mirror),
56  XyzL(laser.XyzL), XyzU(laser.XyzU), XyzB(laser.XyzB),
57  dirL(laser.dirL), dirU(laser.dirU),dirB(laser.dirB),
58  Theta(laser.Theta), Phi(laser.Phi){
59 }
60 //______________________________________________________________________________
61 LaserEvent::LaserEvent()
62 {
63  // Create an LaserEvent object.
64  // When the constructor is invoked for the first time, the class static
65  // variable fgTracks is 0 and the TClonesArray fgTracks is created.
66 
67  if (!fgVertices) fgVertices = new TClonesArray("Vertex", 1000);
68  fVertices = fgVertices;
69  fNvertex = 0;
70  if (!fgTracks) fgTracks = new TClonesArray("Track", 1000);
71  fTracks = fgTracks;
72  fNtrack = 0;
73 
74  if (!fgHits) fgHits = new TClonesArray("Hit", 1000);
75  fHits = fgHits;
76  fNhit = 0;
77 
78  if (!fgFit) fgFit = new TClonesArray("FitDV", 1000);
79  fgFit->ExpandCreate(12);
80  fFit = fgFit;
81 }
82 
83 //______________________________________________________________________________
84 LaserEvent::~LaserEvent()
85 {
86  Clear("C");
87 }
88 
89 //______________________________________________________________________________
90 Vertex *LaserEvent::AddVertex(StPrimaryVertex *vertex) {
91  if (! vertex) return 0;
92  TClonesArray &vertices = *fVertices;
93  Vertex *vx = new(vertices[fNvertex++]) Vertex(vertex);
94  const TGeoHMatrix &Tpc2Global = gStTpcDb->Tpc2GlobalMatrix();
95  Tpc2Global.MasterToLocal(vx->Xyz.xyz(),vx->XyzL.xyz());
96  return vx;
97 }
98 //______________________________________________________________________________
99 Track *LaserEvent::AddTrack(Int_t sector, StTrack *track, LaserB *laser, Double_t z) {
100  if (! track) return 0;
101  TClonesArray &tracks = *fTracks;
102  Track *t = new(tracks[fNtrack++]) Track(sector, track, laser, z);
103  return t;
104 }
105 //______________________________________________________________________________
106 Hit *LaserEvent::AddHit(StTpcHit *hit, Int_t trackKey) {
107  if (! hit) return 0;
108  TClonesArray &hits = *fHits;
109  Hit *t = new(hits[fNhit++]) Hit(hit,trackKey);
110  return t;
111 }
112 //________________________________________________________________________________
113 void LaserEvent::AddTrackFit(Track *t) {
114  if (! t) return;
115  Int_t sector = t->Laser.Sector;
116  Int_t bundle = t->Laser.Bundle;
117  Int_t mirror = t->Laser.Mirror;
118  Int_t s2 = (sector-1)/2;
119  if (s2 >= 0 && s2 < 12) {
120  TClonesArray &fits = *fFit;
121  FitDV *fit = (FitDV *) fits[s2];
122  fit->Sector = sector;
123  Int_t N = fit->N;
124  Double32_t x = t->Laser.XyzL.z();
125  Double32_t y = t->XyzPL.z() - x;
126  if (N < 42) {
127  fit->X[N] = x;
128  fit->Y[N] = y;
129  fit->Bundle[N] = bundle;
130  fit->Mirror[N] = mirror;
131  fit->N = N+1;
132  }
133  }
134 }
135 //______________________________________________________________________________
136 void LaserEvent::Clear(Option_t *option) {
137  fTracks->Clear(option);
138  fHits->Clear(option);
139  fVertices->Clear(option);
140  fgFit->Clear(option);
141 }
142 //______________________________________________________________________________
143 void LaserEvent::Reset() {
144  // Static function to reset all static objects for this event
145  SafeDelete(fgTracks);
146  SafeDelete(fgHits);
147  SafeDelete(fgVertices);
148  SafeDelete(fgFit);
149 }
150 
151 //______________________________________________________________________________
152 void LaserEvent::SetHeader(Int_t i, Int_t run, Int_t date, Int_t time)
153 {
154  fNvertex = 0;
155  fNtrack = 0;
156  fNhit = 0;
157  fEvtHdr.Set(i, run, date, time);
158 }
159 //______________________________________________________________________________
160 void LaserEvent::SetHeader(Int_t i, Int_t run, Int_t date, Int_t time,
161  Float_t tzero, Float_t drivel, Float_t clock)
162 {
163  SetHeader(i, run, date, time);
164  fEvtHdr.SetE(tzero, drivel, clock);
165 }
166 //______________________________________________________________________________
167 void LaserEvent::SetHeader(Int_t i, Int_t run, Int_t date, Int_t time,
168  Float_t tzero, Float_t drivel, Float_t clock, Float_t trigger)
169 {
170  SetHeader(i, run, date, time);
171  fEvtHdr.SetE(tzero, drivel, clock, trigger);
172 }
173 //______________________________________________________________________________
174 Track::Track(Int_t sector, StTrack *track, LaserB *theLaser, Double_t z) :
175  Flag(0),mType(kUndefinedVtxId), mSector(sector),
176  fpTInv(-999), thePath(0), dPhi(-999), dTheta(-999), zLastHit(z) {
177  StTrackNode* node = track->node();
178  StGlobalTrack* gTrack = static_cast<StGlobalTrack*>(node->track(global));
179  fgeoIn = *((StHelixModel *) gTrack->geometry());
180  fgeoOut = *((StHelixModel *) gTrack->outerGeometry());
181  fDca = *((StDcaGeometry *) gTrack->dcaGeometry());
182  StThreeVectorD g3 = fgeoOut.momentum();
183  fpTInv = fgeoOut.charge()/g3.perp();
184  fTheta = fgeoOut.momentum().theta();
185  fPhi = fgeoOut.momentum().phi();
186  StPhysicalHelixD helixO = fgeoOut.helix();
187  thePath = helixO.pathLength(Laser.XyzG.x(),Laser.XyzG.y());
188  if (theLaser) Laser = *theLaser;
189  if (gTrack) {
190  StPrimaryTrack *pTrack = static_cast<StPrimaryTrack*>(node->track(primary));
191  if (pTrack) {
192  StPrimaryVertex* vertex = (StPrimaryVertex*) pTrack->vertex();
193  if (vertex) {
194  mType = vertex->type();
195  Vertex = vertex->position();
196  Flag = 1;
197  }
198  }
199  mKey = gTrack->key();
200  mFlag = gTrack->flag();
201  mNumberOfPossiblePointsTpc = gTrack->numberOfPossiblePoints(kTpcId);
202  mImpactParameter = gTrack->impactParameter();
203  mLength = gTrack->length();
204 #if 1
205  StTrackFitTraits& fitTraits = gTrack->fitTraits();
206  mNumberOfFitPointsTpc = fitTraits.numberOfFitPoints(kTpcId);
207  mPrimaryVertexUsedInFit = fitTraits.primaryVertexUsedInFit();
208 
209  StSPtrVecTrackPidTraits &traits = gTrack->pidTraits();
210  unsigned int size = traits.size();
211  StDedxPidTraits *pid = 0;
212  for (unsigned int i = 0; i < size; i++) {
213  if (! traits[i]) continue;
214  if ( traits[i]->IsZombie()) continue;
215  pid = dynamic_cast<StDedxPidTraits*>(traits[i]);
216  if (pid && pid->detector() == kTpcId && pid->method() == kTruncatedMeanId) {
217  fNdEdx = pid->numberOfPoints(); // Number of points used in dE/dx calc
218  fdEdx = pid->mean();
219  break;
220  }
221  }
222 #endif
223  }
224 }
225 //______________________________________________________________________________
226 void Track::SetPredictions(TGeoHMatrix *Raft2Tpc, TGeoHMatrix *Bundle2Tpc, TGeoHMatrix *Mirror2Tpc) {
227  if ( Flag ) return;
228  if ( Laser.Sector < 1 || Laser.Sector > 24) return;
229  Flag = 2;
230  const TGeoHMatrix &Tpc2Global = gStTpcDb->Tpc2GlobalMatrix();
231  StPhysicalHelixD helixO = fgeoOut.helix();
232  thePath = helixO.pathLength(Laser.XyzG.x(),Laser.XyzG.y());
233  XyzP = helixO.at(thePath);
234  Tpc2Global.MasterToLocal(XyzP.xyz(),XyzPL.xyz());
235  StThreeVectorD g3 = fgeoOut.momentum();
236  dirP = g3.unit();
237  Tpc2Global.MasterToLocalVect(dirP.xyz(),dirPL.xyz());
238  dU = StThreeVectorD(999.,999.,999.);
239  if (Raft2Tpc) {
240  Raft2Tpc->MasterToLocal(XyzPL.xyz(),XyzPU.xyz());
241  Raft2Tpc->MasterToLocalVect(dirPL.xyz(),dirPU.xyz());
242  dU = XyzPU - Laser.XyzU;
243  }
244  if (Bundle2Tpc) {
245  Bundle2Tpc->MasterToLocal(XyzPL.xyz(),XyzPB.xyz());
246  Bundle2Tpc->MasterToLocalVect(dirPL.xyz(),dirPB.xyz());
247  }
248  if (Mirror2Tpc) {
249  Mirror2Tpc->MasterToLocal(XyzPL.xyz(),XyzPM.xyz());
250  Mirror2Tpc->MasterToLocalVect(dirPL.xyz(),dirPM.xyz());
251  }
252  Flag += Matched();
253 }
254 //______________________________________________________________________________
255 Vertex::Vertex(StPrimaryVertex *vertex) : mType(kUndefinedVtxId), WestEast(0),
256  Xyz(), numberOfDaughter(0){
257  if (vertex) {
258  mType = vertex->type();
259  Xyz = vertex->position();
260  numberOfDaughter = vertex->numberOfDaughters();
261  if (numberOfDaughter > 24) {
262  for (UInt_t i = 0; i < numberOfDaughter; i++) {
263  StTrack *track = vertex->daughter(i);
264  if (! track ) continue;
265  StPtrVecHit hvec = track->detectorInfo()->hits();
266  for (unsigned int j=0; j<hvec.size(); j++) {// hit loop
267  if (hvec[j]->detector() != kTpcId) continue;
268  StTpcHit *tpcHit = static_cast<StTpcHit *> (hvec[j]);
269  Int_t sector = tpcHit->sector();
270  Int_t WE = sector <= 12 ? 1 : 2;
271  if (! WestEast) WestEast = WE;
272  else {
273  if (WE != WestEast) {
274  WestEast = 0;
275  break;
276  }
277  }
278  }
279  }
280  }
281  }
282 }
283 //________________________________________________________________________________
284 Int_t Track::Matched() {
285  Int_t iok = 0;
286  dPhi = -999;
287  dTheta = -999;
288  iok = 1 << 10;
289  if (Flag != 2) return iok;
290  iok = 10;
291  if (Laser.Sector < 1 || Laser.Sector > 24 ||
292  Laser.Mirror < 1 || Laser.Mirror > 7) return iok;
293  iok = 0;
294  Int_t status = 0;
295  if (TMath::Abs(dU.x()) > 0.05) status |= 1 << 3;
296  if (TMath::Abs(dU.y()) > 0.05) status |= 1 << 4;
297  if (thePath < 5 || thePath > 25) status |= 1 << 5;
298  if (mNumberOfFitPointsTpc < 25) status |= 1 << 6;
299  dTheta = Laser.ThetaG-fgeoOut.dipAngle()-TMath::Pi()/2;
300 #if 1
301  if (TMath::Abs(dTheta) > 0.030) status |= 1 << 7;;
302 #endif
303  dPhi = Laser.PhiG - fgeoOut.psi();
304  if (dPhi >= TMath::Pi()) dPhi -= 2*TMath::Pi();
305  if (dPhi <= -TMath::Pi()) dPhi += 2*TMath::Pi();
306  // if (SectorMirror[m][i].sigma > 0 && TMath::Abs(dPhi) > 5*SectorMirror[m][i].sigma) status |= 1 << 8;
307 #if 1
308  if ( TMath::Abs(dPhi) > 0.020) status |= 1 << 9;
309 #endif
310  static const Double_t pTInv0 = 4.78815e-03;
311  static const Double_t DpTInv0 = 9.75313e-03;
312  if (TMath::Abs(fpTInv - pTInv0) > 3.0*DpTInv0) status |= 1 << 10;
313  return iok + 10*status;
314 }
315 //________________________________________________________________________________
316 Hit::Hit(StTpcHit *tpcHit, Int_t trKey) : sector(0),row(0),charge(0),flag(0),usedInFit(0), trackKey(trKey) {
317  if (tpcHit) {
318  // hit = tpcHit;
319  adc = tpcHit->adc();
320  sector = tpcHit->sector();
321  row = tpcHit->padrow();;
322  charge = tpcHit->charge();
323  flag = tpcHit->flag();
324  usedInFit = tpcHit->usedInFit();
325  xyz = tpcHit->position(); // from StTpcHitMover
326  static StTpcCoordinateTransform transform(gStTpcDb);
327  StGlobalCoordinate glob(xyz.x(),xyz.y(),xyz.z());
328  static StTpcLocalSectorCoordinate local;
329  transform(glob,local,tpcHit->sector(),tpcHit->padrow());
330  xyzL = StThreeVectorF(local.position().x(),local.position().y(),local.position().z());
331  static StTpcLocalCoordinate localTpc;
332  transform(glob,localTpc,tpcHit->sector(),tpcHit->padrow());
333  xyzTpcL = StThreeVectorF(localTpc.position().x(),localTpc.position().y(),localTpc.position().z());
334  Double_t xyzs[3];
335  gStTpcDb->SupS2Tpc(tpcHit->sector()).MasterToLocal(localTpc.position().xyz(),xyzs);
336  xyzS = StThreeVectorF(xyzs[0],xyzs[1],xyzs[2]);
337  pad = tpcHit->pad() - St_tpcPadConfigC::instance()->padsPerRow(sector,row)/2 - 1;
338  tbk = tpcHit->timeBucket();
339  }
340 }
pair< double, double > pathLength(double r) const
path length at given r (cylindrical r)
Definition: StHelix.cc:351
C++ STL includes.
Definition: AgUStep.h:47