FairRoot/PandaRoot
PndFts2.cxx
Go to the documentation of this file.
1 // -------------------------------------------------------------------------
2 // ----- PndFts source file -----
3 // -------------------------------------------------------------------------
4 
5 #include "PndFts2.h"
6 
7 #include "PndGeoFtsPar.h"
8 #include "PndGeoFts.h"
9 #include "PndDetectorList.h"
10 #include "PndStack.h"
11 #include "PndFtsMapCreator.h"
12 
13 #include "FairRun.h"
14 #include "FairGeoInterface.h"
15 #include "FairGeoLoader.h"
16 #include "FairGeoNode.h"
17 #include "FairGeoRootBuilder.h"
18 #include "FairRootManager.h"
19 #include "FairRuntimeDb.h"
20 #include "FairVolume.h"
21 
22 #include "TClonesArray.h"
23 #include "TLorentzVector.h"
24 #include "TParticle.h"
25 #include "TVirtualMC.h"
26 #include "TGeoMatrix.h"
27 #include "TObjArray.h"
28 #include "TGeoTube.h"
29 #include "TGeoMedium.h"
30 #include "TGeoVolume.h"
31 
32 #include <iostream>
33 
34 using std::cout;
35 using std::endl;
36 using std::string;
37 
38 // TODO: read this from geant initialization
39 #define innerStrawDiameter 1.
40 //#define redefineLambdaChargedDecay 0
41 
42 // ----- Default constructor -------------------------------------------
44  : fTrackID(0), fVolumeID(0), fPos(0,0,0,0), fPosIn(0,0,0,0), fPosOut(0,0,0,0), fPosInLocal(0,0,0,0), fPosOutLocal(0,0,0,0),
45  fMomIn(0,0,0,0), fMomOut(0,0,0,0), fTime(0), fLength(0), fELoss(0), fMass(0), fIsInitialized(kFALSE), fPosIndex(0),
46  fFtsCollection(0), fpostot(0,0,0,0), fpostotin(0,0,0,0), fpostotout(0,0,0,0), fPassNodes(new TObjArray()), valid(kFALSE), fGeoType(0)
47 {
48  fFtsCollection = new TClonesArray("PndFtsPoint");
49  fVerboseLevel = 0;
50 }
51 // -------------------------------------------------------------------------
52 
53 
54 
55 // ----- Standard constructor ------------------------------------------
56 PndFts2::PndFts2(const char* name, Bool_t active)
57  : FairDetector(name, active), fTrackID(0), fVolumeID(0), fPos(0,0,0,0), fPosIn(0,0,0,0), fPosOut(0,0,0,0), fPosInLocal(0,0,0,0), fPosOutLocal(0,0,0,0),
58  fMomIn(0,0,0,0), fMomOut(0,0,0,0), fTime(0), fLength(0), fELoss(0), fMass(0), fIsInitialized(kFALSE), fPosIndex(0),
59  fFtsCollection(0), fpostot(0,0,0,0), fpostotin(0,0,0,0), fpostotout(0,0,0,0), fPassNodes(new TObjArray()), valid(kFALSE), fGeoType(0)
60 {
61  fFtsCollection = new TClonesArray("PndFtsPoint");
62  fVerboseLevel = 0;
63 }
64 // -------------------------------------------------------------------------
65 
66 
67 
68 
69 // ----- Destructor ----------------------------------------------------
71 {
72  if (fFtsCollection)
73  {
74  fFtsCollection->Delete();
75  delete fFtsCollection;
76  }
77 }
78 // -------------------------------------------------------------------------
79 
80 
81 // ----- Private method GetSquaredDistanceFromWire -----------------------
83 {
84  TLorentzVector
85  entryPosition;
86 
87  float
88  positionInMother[3],
89  positionInStraw[3];
90 
91  gMC->TrackPosition(entryPosition);
92  positionInMother[0] = entryPosition.X();
93  positionInMother[1] = entryPosition.Y();
94  positionInMother[2] = entryPosition.Z();
95  gMC->Gmtod(positionInMother, positionInStraw, 1);
96 
97  return positionInStraw[0] * positionInStraw[0] +
98  positionInStraw[1] * positionInStraw[1];
99 }
100 // -------------------------------------------------------------------------
101 
102 bool PndFts2::Split(string &aDest, string &aSrc, char aDelim)
103 {
104  if(aSrc.empty())
105  return false;
106 
107  string::size_type
108  pos = aSrc.find(aDelim);
109 
110  aDest = aSrc.substr(0, pos);
111 
112  if(pos != string::npos)
113  aSrc = aSrc.substr(pos + 1);
114  else
115  aSrc = "";
116 
117  return true;
118 }
119 
120 string PndFts2::GetStringPart(string &aSrc, Int_t part, char aDelim)
121 {
122  string
123  retval = "",
124  sub;
125 
126  int
127  counter = 0;
128 
129  while(Split(sub, aSrc, aDelim))
130  {
131  if (counter == part)
132  {
133  retval = sub;
134  break;
135  }
136  counter++;
137  }
138 
139  return retval;
140 }
141 
142 // ----- Public method ProcessHits --------------------------------------
143 Bool_t PndFts2::ProcessHits(FairVolume* vol)
144 {
145 
146 
147  //std::cout<<"ProcessHit PndFts####################################################"<<std::endl;
148 
149  //TParticle* particle = gMC->GetStack()->GetCurrentTrack(); //[R.K. 01/2017] unused variable?
150  //TGeoMedium *medium = (TGeoMedium*) vol->getGeoNode()->getRootVolume()->GetMedium(); //[R.K. 01/2017] unused variable?
151  //Double_t epsil = medium->GetParam(6); //[R.K. 01/2017] unused variable?
152 
153  TString vol_name(gMC->CurrentVolName());
154  TGeoHMatrix M;
155  gMC->GetTransformation(gMC->CurrentVolPath(),M);
156  TString name(gMC->CurrentVolName());
157 
158 
159  if (gMC->TrackCharge() != 0.)
160  {
161 
162  if ( gMC->IsTrackEntering() )
163  {
164  valid = kTRUE;
165 
166  // Set parameters at entrance of volume. Reset ELoss.
167  fELoss = 0.;
168  fTime = gMC->TrackTime() * 1.0e09;
169  fLength = gMC->TrackLength();
170  gMC->TrackPosition(fPos);
171  gMC->TrackMomentum(fMomIn);
172  gMC->TrackPosition(fpostotin);// da cancellare
173  Double_t globalPos[3] = {0., 0., 0.}; // stt1 modified
174  Double_t localPos[3] = {0., 0., 0.}; // stt1 modified
175 
176  globalPos[0] = fPos.X();
177  globalPos[1] = fPos.Y();
178  globalPos[2] = fPos.Z();
179 
180  gMC->Gmtod(globalPos, localPos, 1);
181  fPosInLocal.SetXYZM(localPos[0], localPos[1], localPos[2], 0.0);
182 
183  }
184 
185  // Sum energy loss for all steps in the active volume
186  fELoss += gMC->Edep();
187 
188  //std::cout<<"bool valid = "<<valid<<std::endl;
189 
190  // Create PndFtsPoint at exit of active volume -- but not into the wire
191  if ( gMC->IsTrackExiting() && valid==kTRUE )
192  {
193  //std::cout<<"sono nel ciclo"<<std::endl;
194  valid=kFALSE;
195  fTrackID = gMC->GetStack()->GetCurrentTrackNumber();
196  fVolumeID = kFTS;//vol->getMCid();
197 
198  if(fTrackID !=0){
199  if (fVerboseLevel>2) std::cout<<"test Vol----------"<<vol->getMCid()<<std::endl;
200  if (fVerboseLevel>2) std::cout<<"fTrackID-----"<<fTrackID<<std::endl;
201  }
202 
203  fMass = gMC->TrackMass(); // mass (GeV)
204  gMC->TrackPosition(fPosOut);
205  gMC->TrackMomentum(fMomOut);
206  gMC->TrackPosition(fpostotout);// da cancellare
207  Double_t globalPos[3] = {0., 0., 0.}; // stt1 modified
208  Double_t localPos[3] = {0., 0., 0.}; // stt1 modified
209 
210  gMC->Gdtom(localPos, globalPos, 1);
211 
212  fPos.SetXYZM(globalPos[0], globalPos[1], globalPos[2], 0.0);
213 
214  globalPos[0] = fPosOut.X();
215  globalPos[1] = fPosOut.Y();
216  globalPos[2] = fPosOut.Z();
217 
218  gMC->Gmtod(globalPos, localPos, 1);
219  fPosOutLocal.SetXYZM(localPos[0], localPos[1], localPos[2], 0.0);
220  /*
221  // string basename("stt1tube");
222  string basename;
223  TString volumename;
224  string
225  hashmark("#"),
226  volName,
227  fullName,
228  number,
229  specialname,
230  volPath(gMC->CurrentVolPath()),
231  volPath2(gMC->CurrentVolPath());
232 
233  volumename = volPath;
234  volName = GetStringPart(volPath, 2, '/');
235  number = GetStringPart(volName, 1, '_');
236 
237  Int_t start =volPath2.find("fts",0);
238  specialname = volPath2.substr(start,volPath2.find("_",start)-start);
239 
240  if(volumename.Contains("fts01")) basename = "fts01tube";
241  if(volumename.Contains("fts02")) basename = "fts02tube";
242  if(volumename.Contains("fts31")) basename = "fts31tube";
243  if(volumename.Contains("fts32")) basename = "fts32tube";
244  if(volumename.Contains("fts33")) basename = "fts33tube";
245  if(volumename.Contains("fts34")) basename = "fts34tube";
246  if(volumename.Contains("fts35")) basename = "fts35tube";
247  if(volumename.Contains("fts36")) basename = "fts36tube";
248  if(volumename.Contains("fts37")) basename = "fts37tube";
249  if(volumename.Contains("fts38")) basename = "fts38tube";
250  if(volumename.Contains("fts05")) basename = "fts05tube";
251  if(volumename.Contains("fts06")) basename = "fts06tube";
252 
253  fullName = basename + hashmark + number;
254 
255  FairGeoNode
256  *volnode = dynamic_cast<FairGeoNode*> (fPassNodes->FindObject(fullName.c_str()));
257  if (!volnode)
258  {
259  cout << "-I- PndFts: No volume " << fullName.c_str() << " found in geometry container." << endl;
260  return kFALSE;
261  }
262 
263 
264  if(fTrackID!=0){
265  if (fVerboseLevel>2) std::cout<<"befor Mapper------"<<std::endl;
266  }
267 
268  fpostot.SetXYZM((fpostotin.X() + fpostotout.X())/2., (fpostotin.Y() + fpostotout.Y())/2.,(fpostotin.Z() + fpostotout.Z())/2.,0.0);
269 
270 
271  //CHECK map creator-------------------------------------------------
272  PndFtsMapCreator *mapper = new PndFtsMapCreator(fGeoType);
273  Int_t tubeID=0;
274  //testTubeID=tube number as in the geometry file....
275  //we need to calculate the real number in order to have a different
276  //number for each tube (also for up and down short tubes)
277  Int_t testTubeID = mapper->GetTubeIDFromPath(gMC->CurrentVolPath());
278  Int_t chamberID=mapper->GetChamberIDFromPath(gMC->CurrentVolPath());
279  Int_t layerID=mapper->GetLayerID(chamberID,testTubeID,gMC->CurrentVolPath());
280  //std::cout<<"LAYERID==========="<<layerID<<std::endl;
281  //tubeID=at each tube corresponds only one id number.....
282  tubeID = mapper->GetTubeIDTot(chamberID,layerID,testTubeID,gMC->CurrentVolPath());
283  //std::cout<<"PndFts.cxx. chamber, tube ID = "<<chamberID<<" "<<tubeID<<" "<<testTubeID<<std::endl;
285 */
286 fpostot.SetXYZM((fpostotin.X() + fpostotout.X())/2., (fpostotin.Y() + fpostotout.Y())/2.,(fpostotin.Z() + fpostotout.Z())/2.,0.0);
287  Int_t tubeID = 0, chamberID = 0, layerID = 0;
288  AddHit(fTrackID, fVolumeID, tubeID, chamberID,layerID,
289  TVector3(fpostot.X(), fpostot.Y(), fpostot.Z()),
290  TVector3(fPosInLocal.X(), fPosInLocal.Y(), fPosInLocal.Z()),
291  TVector3(fPosOutLocal.X(), fPosOutLocal.Y(), fPosOutLocal.Z()),
292  TVector3(fMomIn.Px(), fMomIn.Py(), fMomIn.Pz()),
293  TVector3(fMomOut.Px(), fMomOut.Py(), fMomOut.Pz()),
294  fTime, fLength, fELoss, fMass);
295 
296  if (fTrackID !=0){
297  if (fVerboseLevel>2) std::cout<<"AddHit PndFts.cxx= "<<fTrackID<<std::endl;
298  }
299 
300 
301  //Increment number of stt points for TParticle
302  PndStack* stack = (PndStack*) gMC->GetStack();
303  stack->AddPoint(kFTS);
304  ResetParameters();
305 
306 
307 
308  }
309 
310  }
311 
312  return kTRUE;
313 }
314 // -------------------------------------------------------------------------
315 
316 
317 
318 // ----- Public method EndOfEvent --------------------------------------
320 {
321  if (fVerboseLevel)
322  Print();
323  fFtsCollection->Delete();
324  fPosIndex = 0;
325 }
326 // -------------------------------------------------------------------------
327 
328 
329 
330 // ----- Public method Register ----------------------------------------
332 {
333  FairRootManager::Instance()->Register("FTSPoint", "Fts", fFtsCollection, kTRUE);
334 }
335 // -------------------------------------------------------------------------
336 
337 
338 
339 // ----- Public method GetCollection -----------------------------------
340 TClonesArray* PndFts2::GetCollection(Int_t iColl) const
341 {
342  if (iColl == 0)
343  return fFtsCollection;
344  else
345  return NULL;
346 }
347 // -------------------------------------------------------------------------
348 
349 
350 
351 // ----- Public method Print -------------------------------------------
352 void PndFts2::Print() const
353 {
354  Int_t
355  nHits = fFtsCollection->GetEntriesFast();
356 
357  cout << "-I- PndFts: " << nHits << " points registered in this event." << endl;
358 
359  if (fVerboseLevel>1)
360  for (Int_t i=0; i<nHits; i++)
361  (*fFtsCollection)[i]->Print();
362 }
363 // -------------------------------------------------------------------------
364 
365 
366 
367 // ----- Public method Reset -------------------------------------------
369 {
370  fFtsCollection->Delete();
371  ResetParameters();
372 }
373 // -------------------------------------------------------------------------
374 
375 
376 
377 // ----- Public method CopyClones --------------------------------------
378 void PndFts2::CopyClones(TClonesArray* cl1, TClonesArray* cl2, Int_t offset)
379 {
380  Int_t
381  nEntries = cl1->GetEntriesFast();
382 
383  cout << "-I- PndFts: " << nEntries << " entries to add." << endl;
384 
385  TClonesArray& clref = *cl2;
386 
388  *oldpoint = NULL;
389  for (Int_t i=0; i<nEntries; i++)
390  {
391  oldpoint = (PndFtsPoint*) cl1->At(i);
392 
393  Int_t
394  index = oldpoint->GetTrackID() + offset;
395 
396  oldpoint->SetTrackID(index);
397  new (clref[fPosIndex]) PndFtsPoint(*oldpoint);
398  fPosIndex++;
399  }
400  cout << "-I- PndFts: " << cl2->GetEntriesFast() << " merged entries."
401  << endl;
402 }
403 // -------------------------------------------------------------------------
404 
406  cout << " -I- Initializing PndFts()" <<endl;
408 
409 }
410 
411 
412 
413 // ----- Public method ConstructGeometry -------------------------------
415 {
416  FairGeoLoader* geoLoad = FairGeoLoader::Instance();
417  FairGeoInterface* geoFace = geoLoad->getGeoInterface();
418  PndGeoFts* Geo = new PndGeoFts();
419  Geo->setGeomFile(GetGeometryFileName());
420  geoFace->addGeoModule(Geo);
421 
422  Bool_t rc = geoFace->readSet(Geo);
423  if (rc) Geo->create(geoLoad->getGeoBuilder());
424 
425  TList* volList = Geo->getListOfVolumes();
426 
427  // store geo parameter
428  FairRun *fRun = FairRun::Instance();
429  FairRuntimeDb *rtdb= FairRun::Instance()->GetRuntimeDb();
430  PndGeoFtsPar* par=(PndGeoFtsPar*)(rtdb->getContainer("PndGeoFtsPar"));
431  TObjArray *fSensNodes = par->GetGeoSensitiveNodes();
433 
434  TListIter iter(volList);
435  FairGeoNode* node = NULL;
436  FairGeoVolume *aVol=NULL;
437 
438 
439  while( (node = (FairGeoNode*)iter.Next()) ) {
440  aVol = dynamic_cast<FairGeoVolume*> ( node );
441  if ( node->isSensitive() ) {
442  fSensNodes->AddLast( aVol );
443  }else{
444  fPassNodes->AddLast( aVol );
445  }
446  }
447 
448  //std::cout<<"BBBB->"<<fSensNodes->GetEntries()<<std::endl;
449  //std::cout<<"cccc->"<<fPassNodes->GetEntries()<<std::endl;
450 
451  fGeoType = 1; // CHECK
453  par->SetTubeInRad(0.5); // cm
454  par->SetTubeOutRad(0.001); // cm
455 
456 
457  par->setChanged();
458  par->setInputVersion(fRun->GetRunId(),1);
459  ProcessNodes ( volList );
460 
461 }
462 // -------------------------------------------------------------------------
463 
464 
465 
466 
467 // ----- Private method AddHit -----------------------------------------
468 PndFtsPoint* PndFts2::AddHit(Int_t trackID, Int_t detID, Int_t tubeID, Int_t chamberID, Int_t layerID,
469  TVector3 pos, TVector3 posInLocal, TVector3 posOutLocal,
470  TVector3 momIn, TVector3 momOut,
471  Double_t time, Double_t length, Double_t eLoss, Double_t mass) // da cancellare postot
472 {
473  TClonesArray&
474  clref = *fFtsCollection;
475 
476  Int_t
477  size = clref.GetEntriesFast();
478 
479  PndFtsPoint *pointnew = new(clref[size]) PndFtsPoint(trackID, detID, tubeID, chamberID, layerID, pos, posInLocal, posOutLocal,momIn, momOut, time, length, eLoss, mass);
480 
481  pointnew->SetTubeID(tubeID);
482  pointnew->SetChamberID(chamberID);
483  pointnew->SetLayerID(layerID);
484 
485  return pointnew;
486 }
487 // -------------------------------------------------------------------------
488 
489 
490 
492 
TVector3 pos
virtual void CopyClones(TClonesArray *cl1, TClonesArray *cl2, Int_t offset)
Definition: PndFts2.cxx:378
TLorentzVector fMomOut
momentum
Definition: PndFts2.h:119
Double_t fTime
momentum
Definition: PndFts2.h:120
Bool_t valid
Definition: PndFts2.h:134
PndFtsPoint * AddHit(Int_t trackID, Int_t detID, Int_t tubeID, Int_t chamberID, Int_t layerID, TVector3 pos, TVector3 posInLocal, TVector3 posOutLocal, TVector3 momIn, TVector3 momOut, Double_t time, Double_t length, Double_t eLoss, Double_t mass)
Definition: PndFts2.cxx:468
FairGeoLoader * geoLoad
Int_t i
Definition: run_full.C:25
TLorentzVector fpostot
Hit collection.
Definition: PndFts2.h:128
TObjArray * fPassNodes
Definition: PndFts2.h:132
Double_t fMass
energy loss
Definition: PndFts2.h:123
Double_t par[3]
void SetLayerID(Int_t layerid)
Definition: PndFtsPoint.h:88
TVector3 offset(2, 0, 0)
void SetTubeInRad(Double_t inrad)
Definition: PndGeoFtsPar.h:35
TLorentzVector fpostotout
Definition: PndFts2.h:130
float GetSquaredDistanceFromWire()
Definition: PndFts2.cxx:82
bool Split(std::string &aDest, std::string &aSrc, char aDelim)
Definition: PndFts2.cxx:102
void AddPoint(DetectorId iDet)
Definition: PndStack.cxx:408
Double_t fLength
time
Definition: PndFts2.h:121
PndFts2()
Definition: PndFts2.cxx:43
FairRunAna * fRun
Definition: hit_dirc.C:58
void SetChamberID(Int_t chamberid)
Definition: PndFtsPoint.h:87
Double_t fELoss
length
Definition: PndFts2.h:122
Int_t fPosIndex
Definition: PndFts2.h:126
TClonesArray * fFtsCollection
Definition: PndFts2.h:127
virtual void Reset()
Definition: PndFts2.cxx:368
int counter
Definition: ZeeAnalysis.C:59
int nHits
Definition: RiemannTest.C:16
Double_t
void SetTubeOutRad(Double_t outrad)
Definition: PndGeoFtsPar.h:36
virtual Bool_t ProcessHits(FairVolume *vol=0)
Definition: PndFts2.cxx:143
void Initialize()
Definition: PndFts2.cxx:405
FairRuntimeDb * rtdb
Definition: hit_dirc.C:66
TLorentzVector fpostotin
Definition: PndFts2.h:129
void ResetParameters()
Definition: PndFts2.h:173
virtual void Print() const
Definition: PndFts2.cxx:352
TString name
Int_t fTrackID
Definition: PndFts2.h:111
TLorentzVector fPosOutLocal
entry position in straw frame
Definition: PndFts2.h:117
TObjArray * GetGeoPassiveNodes()
Definition: PndGeoFtsPar.h:31
TLorentzVector fPosInLocal
exit position in global frame
Definition: PndFts2.h:116
TLorentzVector fPos
volume id
Definition: PndFts2.h:113
TLorentzVector fPosOut
entry position in global frame
Definition: PndFts2.h:115
void SetGeometryType(Int_t geoType)
Definition: PndGeoFtsPar.h:34
TLorentzVector fMomIn
exit position in straw frame
Definition: PndFts2.h:118
virtual void Register()
Definition: PndFts2.cxx:331
virtual ~PndFts2()
Definition: PndFts2.cxx:70
Int_t fGeoType
Definition: PndFts2.h:135
ClassImp(PndAnaContFact)
virtual TClonesArray * GetCollection(Int_t iColl) const
Definition: PndFts2.cxx:340
Int_t fVolumeID
track index
Definition: PndFts2.h:112
FairGeoInterface * geoFace
TObjArray * GetGeoSensitiveNodes()
Definition: PndGeoFtsPar.h:30
void SetTubeID(Int_t tubeid)
Definition: PndFtsPoint.h:85
Mvd Initialize()
virtual void ConstructGeometry()
Definition: PndFts2.cxx:414
std::string GetStringPart(std::string &aSrc, Int_t part, char aDelim)
Definition: PndFts2.cxx:120
virtual void EndOfEvent()
Definition: PndFts2.cxx:319