FairRoot/PandaRoot
Public Types | Public Member Functions | Public Attributes | Protected Member Functions | Protected Attributes | List of all members
PndFtsHoughTrackerTask Class Reference

Interface between PandaRoot (data input and output) and PndFtsHoughTrackFinder (implementation of PR algorithm). More...

#include <PndFtsHoughTrackerTask.h>

Inheritance diagram for PndFtsHoughTrackerTask:
PndPersistencyTask PndFtsHoughTrackerTaskQA

Public Types

enum  DebugLevelMultiplicators {
  kHoughSpaces =2, kAllFoundPeaksTogether =3, kEachFoundPeakSeparately =5, kMcTruthPeaksExclusively =7,
  kMcTruthPeaksProjected =11, kHitCurvesExclusively =13, kHitCurvesProjected =17
}
 

Public Member Functions

 PndFtsHoughTrackerTask (Int_t verbose=0, Bool_t persistence=kTRUE)
 Constructor with flags. Can also be used as standard constructor. More...
 
 ~PndFtsHoughTrackerTask ()
 Destructor. More...
 
virtual void SetParContainers ()
 Loads the parameter container from the runtime database. More...
 
virtual InitStatus Init ()
 Initialization of task at the beginning of a run. More...
 
virtual InitStatus ReInit ()
 ReInitiliazation of task when the runID changes. More...
 
virtual void Exec (Option_t *opt)
 Executed for each event. More...
 
virtual void FinishEvent ()
 When is this executed? After each event? More...
 
virtual void Finish ()
 Writes output to root file, I guess. Called at the end of the run. More...
 
void SetVerbose (Int_t verbose)
 
void SetPersistence (Bool_t val)
 
void SetSaveDebugInfo (Int_t saveDebugInfo)
 
Int_t GetVerbose () const
 
Int_t GetSaveDebugInfo () const
 Returns the verbosity level. More...
 
UInt_t GetEventNr () const
 Returns the save debug flag. More...
 
Int_t GetNFtsHits () const
 Returns the event number. More...
 
const PndFtsHitGetFtsHit (UInt_t hitId) const
 Returns pointer to the hit with index hitId in the FTS hit array. More...
 
Int_t getMcTruthIdForHitId (UInt_t hitId) const
 
const PndFtsTubeGetFtsTube (const PndFtsHit *const myHit) const
 Returns pointer to the FTS tube corresponding to input FTS hit. More...
 
const TVector3 GetFtsHitPosErrors (const PndFtsHit *const ftsHit) const
 Returns the position error (based on FTS straw geometry) for the hit with index hitId in the FTS hit array. More...
 
const TMatrixT< Double_tGetFtsHitCovMatrix (const PndFtsHit *const ftsHit) const
 Returns the position covariance matrix (based on FTS straw geometry) for the hit with index hitId in the FTS hit array. More...
 
Int_t getFtsBranchId () const
 Returns detector Id of FTS. Try not to use it. More...
 
TClonesArray * getFtsHitArrayPtr () const
 Returns pointer to the hit array in which FTS hits are saved as PndFtsHit. Try not to use it. More...
 
FairField * getMagneticFieldPtr () const
 
void SetPersistency (Bool_t val=kTRUE)
 
Bool_t GetPersistency ()
 

Public Attributes

FairLogger * fLogger
 Returns pointer to the B field. More...
 

Protected Member Functions

void throwError (const TString s) const
 For error reporting. More...
 
void CheckForDuplicateFtsHits ()
 
 PndFtsHoughTrackerTask (const PndFtsHoughTrackerTask &)
 
PndFtsHoughTrackerTask operator= (const PndFtsHoughTrackerTask &)
 
 ClassDef (PndFtsHoughTrackerTask, 1)
 

Protected Attributes

Int_t fFtsBranchId
 Detector Id of FTS. More...
 
TClonesArray * fFtsHitArray
 Input array of PndFtsHit. More...
 
TClonesArray * fFtsMcPoints
 Input array of McPoints. More...
 
PndGeoFtsParfFtsParameters
 Needed for FTS map creator. More...
 
TClonesArray * fFtsTubeArray
 Input array of PndFtsTube (map of FTS tubes). More...
 
FairField * fField
 For B field access. More...
 
TString fTracksArrayName
 Branch name where to store the Track candidates. More...
 
TClonesArray * fTrackCands
 Array of found track candidates in PndTrackCand (for output) More...
 
TClonesArray * fTracks
 Array of found tracks in PndTrack (for output) More...
 
Int_t fSaveDebugInfo
 Debug information will be created if >0. More...
 
UInt_t fEventNr
 Event number for debugging purposes. More...
 

Detailed Description

Interface between PandaRoot (data input and output) and PndFtsHoughTrackFinder (implementation of PR algorithm).

Author
Martin J. Galuska <martin [dot] j [dot] galuska (at) physik [dot] uni (minus) giessen [dot] de>

This task provides functionality / data for the FTS PR classes by passing a pointer to itself as an argument in the object's constructor.

This class was originally modeled after mvd/MvdTracking/PndMvdRiemannTrackFinderTask (among many others)

TODO Make this task work for time-based simulation as well [have a look at PndMvdRiemannTrackFinderTask::FillHitArray()]

Created: 18.06.2013

Definition at line 49 of file PndFtsHoughTrackerTask.h.

Member Enumeration Documentation

Constructor & Destructor Documentation

PndFtsHoughTrackerTask::PndFtsHoughTrackerTask ( Int_t  verbose = 0,
Bool_t  persistence = kTRUE 
)

Constructor with flags. Can also be used as standard constructor.

Parameters
[in]verboseVerbosity level: 0 least, higher -> more output.
[in]persistencekFALSE does not write track candidates from PR to output root file.
[in]saveDebugInfokTRUE will write internal representation of track(let) candidates to output root file.

Definition at line 78 of file PndFtsHoughTrackerTask.cxx.

References fVerbose, and PndPersistencyTask::SetPersistency().

79 : PndPersistencyTask("PndFtsHoughTrackerTask", verbose),
80  fLogger(FairLogger::GetLogger()),
81  fFtsBranchId(0),
82  fFtsHitArray(0),
83  fFtsMcPoints(0),
84  fFtsParameters(0),
85  fFtsTubeArray(0),
86  fField(0),
87  fTracksArrayName("FTSTrkHough"),
88  fTrackCands(0),
89  fTracks(0),
90  fSaveDebugInfo(kFALSE),
91  fEventNr(0)
92  // fOutFile(0),
93  // fHoughTrackCands(0),
94 {
95  if(3<fVerbose) std::cout << "PndFtsHoughTrackerTask is the tracker ptr " << this << '\n';
96  SetPersistency(persistence);
97 }
int fVerbose
Definition: poormantracks.C:24
UInt_t fEventNr
Event number for debugging purposes.
FairField * fField
For B field access.
#define verbose
void SetPersistency(Bool_t val=kTRUE)
TClonesArray * fFtsHitArray
Input array of PndFtsHit.
TClonesArray * fTracks
Array of found tracks in PndTrack (for output)
Int_t fFtsBranchId
Detector Id of FTS.
TClonesArray * fFtsMcPoints
Input array of McPoints.
TString fTracksArrayName
Branch name where to store the Track candidates.
FairLogger * fLogger
Returns pointer to the B field.
PndGeoFtsPar * fFtsParameters
Needed for FTS map creator.
Int_t fSaveDebugInfo
Debug information will be created if &gt;0.
TClonesArray * fFtsTubeArray
Input array of PndFtsTube (map of FTS tubes).
TClonesArray * fTrackCands
Array of found track candidates in PndTrackCand (for output)
PndFtsHoughTrackerTask::~PndFtsHoughTrackerTask ( )

Destructor.

Definition at line 100 of file PndFtsHoughTrackerTask.cxx.

References fLogger, fTrackCands, fTracks, and fVerbose.

101 {
102  if(fVerbose>3) fLogger->Info(MESSAGE_ORIGIN,"Destructor of PndFtsHoughTrackerTask");
103  // TODO Is that correct?!?
104 // delete fHoughTrackCands;
105  delete fTrackCands;
106  delete fTracks;
107  // fOutFile->Close();
108 }
int fVerbose
Definition: poormantracks.C:24
TClonesArray * fTracks
Array of found tracks in PndTrack (for output)
FairLogger * fLogger
Returns pointer to the B field.
TClonesArray * fTrackCands
Array of found track candidates in PndTrackCand (for output)
PndFtsHoughTrackerTask::PndFtsHoughTrackerTask ( const PndFtsHoughTrackerTask )
protected

Member Function Documentation

void PndFtsHoughTrackerTask::CheckForDuplicateFtsHits ( )
inlineprotected

Definition at line 227 of file PndFtsHoughTrackerTask.h.

References GetEventNr(), GetFtsHit(), GetNFtsHits(), and PndFtsHit::GetTubeID().

Referenced by Exec().

228  {
229  for (int iHit1 = 0; iHit1 < GetNFtsHits(); ++iHit1)
230  {
231  const PndFtsHit *const myHit1 = GetFtsHit(iHit1);
232 
233  const Int_t tubeIdHit1 = myHit1->GetTubeID();
234 
235  for (int iHit2 = iHit1+1; iHit2 < GetNFtsHits(); ++iHit2)
236  {
237  const PndFtsHit *const myHit2 = GetFtsHit(iHit2);
238 
239  const Int_t tubeIdHit2 = myHit2->GetTubeID();
240 
241  if ( tubeIdHit1 == tubeIdHit2 ) std::cout << "Event " << GetEventNr() << "tubeID1="<<tubeIdHit1<<" tubeID2="<<tubeIdHit2 <<": HitIdx " << iHit1 << " and HitIdx " << iHit2 << " are duplicate!\n";
242  }
243 
244  } // for loop over all hits
245  };
Int_t GetTubeID() const
Definition: PndFtsHit.h:70
Int_t GetNFtsHits() const
Returns the event number.
const PndFtsHit * GetFtsHit(UInt_t hitId) const
Returns pointer to the hit with index hitId in the FTS hit array.
UInt_t GetEventNr() const
Returns the save debug flag.
PndFtsHoughTrackerTask::ClassDef ( PndFtsHoughTrackerTask  ,
 
)
protected
void PndFtsHoughTrackerTask::Exec ( Option_t *  opt)
virtual

Executed for each event.

Parameters
optNot used.

Reimplemented in PndFtsHoughTrackerTaskQA.

Definition at line 294 of file PndFtsHoughTrackerTask.cxx.

References PndTrackCand::CalcTimeStamp(), CheckForDuplicateFtsHits(), fEventNr, PndFtsHoughTrackFinderQA::FindTracks(), fTrackCands, fTracks, fVerbose, PndFtsHoughTrackFinder::GetHoughTrack(), PndFtsHoughTrackFinder::GetPndTrack(), PndFtsHoughTrackFinder::GetPndTrackCand(), PndFtsHoughTrackFinder::NTracks(), PndTrack::Print(), PndFtsHoughTrackCand::Print(), PndTrackCand::Print(), PndTrack::SetTrackCandRef(), and trackFinder.

295 {
296  if(1<fVerbose) Info("Exec","Exec of PndFtsHoughTrackerTask on event %i", fEventNr);
297 
298  // Reset output array
299  if ( ! fTrackCands ) Fatal("Exec", "No track cand array");
300 
301  fTrackCands->Delete();
302  fTracks->Delete();
303 // fHoughTrackCands->Delete();
304  //fHoughSpaces->Delete();
305 
307 
308  if(3<fVerbose) std::cout << "PndFtsHoughTrackFinder::Exec tracker ptr " << this << '\n';
310  // trackFinder.SetMinPeakHeightZxLineParabola(4);
311  // trackFinder.SetMinPeakHeightZxParabola(6);
312  // trackFinder.SetMinPeakHeightZxParabolaLine(4);
313  // trackFinder.SetMinPeakHeightZyLine(4);
314  trackFinder.FindTracks();
315 
316 
317  // have a look at PndMvdRiemannTrackFinderTask line 161 for reference
318  // store the found tracks as PndTrack and PndTrackCand
319  for (Int_t iFoundTrack = 0; iFoundTrack < trackFinder.NTracks(); ++iFoundTrack){
320 
321 // // for debug output get PndFtsHoughTrackCand
322 // // FIXME this caused trouble, but I don't need it
323 // if (0<fSaveDebugInfo) {
324 // PndFtsHoughTrackCand* myHoughCand = new ((*fHoughTrackCands)[iFoundTrack])PndFtsHoughTrackCand(trackFinder.GetHoughTrack(iFoundTrack));
325 // }
326 
327  // get output as PndTrackCand and store into TCA
328  PndTrackCand* myCand = new ((*fTrackCands)[iFoundTrack])PndTrackCand(trackFinder.GetPndTrackCand(iFoundTrack));
329  if (1<fVerbose)
330  {
331  std::cout << "Track " << iFoundTrack << '\n';
332  std::cout << "Links: ";
333  ((FairMultiLinkedData*) myCand)->Print();
334  std::cout << '\n';
335  }
336 
337  myCand->CalcTimeStamp(); // TODO Why is this needed?
338  if (1<fVerbose) trackFinder.GetHoughTrack(iFoundTrack).Print();
339 
340 
341  PndTrack* myTrack = new ((*fTracks)[iFoundTrack])PndTrack(trackFinder.GetPndTrack(iFoundTrack)); // TODO Some parameters are missing
342  if (myCand->GetTimeStamp() == 0){
343  myTrack->SetTimeStamp(0.0001 * (iFoundTrack+1));
344  myCand->SetTimeStamp(0.0001 * (iFoundTrack+1));
345 
346  } else {
347  myTrack->SetTimeStamp(myCand->GetTimeStamp());
348  myTrack->SetTimeStampError(myCand->GetTimeStampError());
349  }
350  myTrack->SetLink(FairLink("FTSHoughTrackCand", iFoundTrack));
351  myTrack->SetTrackCandRef(myCand);
352 
353 
354  if (1<fVerbose) {
355  std::cout << iFoundTrack << ": ";
356  myTrack->Print();
357  }
358  }
359  fTrackCands->Sort();
360  fTracks->Sort();
361 
362 
363 
364  if(3<fVerbose) Info("Exec","End eventloop.");
365  ++fEventNr;
366 }
int fVerbose
Definition: poormantracks.C:24
UInt_t fEventNr
Event number for debugging purposes.
void Print() const
void SetTrackCandRef(PndTrackCand *candPointer)
Definition: PndTrack.h:44
TClonesArray * fTracks
Array of found tracks in PndTrack (for output)
void CalcTimeStamp()
PndLheTrackFinder * trackFinder
void Print()
Definition: PndTrack.cxx:39
Implementation of the QA for the Hough transform based FTS PR. Use this for parameter optimization...
TClonesArray * fTrackCands
Array of found track candidates in PndTrackCand (for output)
void PndFtsHoughTrackerTask::Finish ( )
virtual

Writes output to root file, I guess. Called at the end of the run.

Reimplemented in PndFtsHoughTrackerTaskQA.

Definition at line 379 of file PndFtsHoughTrackerTask.cxx.

References fLogger, fTrackCands, and fVerbose.

380 {
381  if(3<fVerbose) Info("Finish","Finish of PndFtsHoughTrackerTask");
382  // Get a handle from the IO manager
383  FairRootManager* ioman = FairRootManager::Instance();
384  if ( ! ioman ) {
385  fLogger->Fatal(MESSAGE_ORIGIN,"RootManager not instantiated, return!");
386  }
387  ioman->Write();
388  if(3<fVerbose) Info("Finish","Found %i tracks.",fTrackCands->GetEntriesFast());
389 }
int fVerbose
Definition: poormantracks.C:24
FairLogger * fLogger
Returns pointer to the B field.
TClonesArray * fTrackCands
Array of found track candidates in PndTrackCand (for output)
void PndFtsHoughTrackerTask::FinishEvent ( )
virtual

When is this executed? After each event?

Reimplemented in PndFtsHoughTrackerTaskQA.

Definition at line 370 of file PndFtsHoughTrackerTask.cxx.

References fTrackCands, and fTracks.

371 {
372  fTrackCands->Delete();
373  fTracks->Delete();
374 // fHoughTrackCands->Delete();
375 }
TClonesArray * fTracks
Array of found tracks in PndTrack (for output)
TClonesArray * fTrackCands
Array of found track candidates in PndTrackCand (for output)
UInt_t PndFtsHoughTrackerTask::GetEventNr ( ) const
inline
Int_t PndFtsHoughTrackerTask::getFtsBranchId ( ) const
inline

Returns detector Id of FTS. Try not to use it.

See Also
GetFtsHit
GetNHits
Returns
Detector Id of FTS.

Definition at line 152 of file PndFtsHoughTrackerTask.h.

References fFtsBranchId.

152 { return fFtsBranchId; };
Int_t fFtsBranchId
Detector Id of FTS.
const PndFtsHit* PndFtsHoughTrackerTask::GetFtsHit ( UInt_t  hitId) const
inline

Returns pointer to the hit with index hitId in the FTS hit array.

Parameters
[in]hitIdIndex (in FTS hit array) of the hit which should be returned.
Returns
Pointer to hit with index hitId in FTS hit array.

Definition at line 110 of file PndFtsHoughTrackerTask.h.

References fFtsHitArray, GetNFtsHits(), and throwError().

Referenced by PndFtsHoughSpace::AddHitToHS(), CheckForDuplicateFtsHits(), PndFtsHoughSpace::filterInputHits(), PndFtsHoughSpace::getHitFromHS(), getMcTruthIdForHitId(), and PndFtsHoughTrackCand::getTrackParPForHit().

110  {
111  if ( hitId >= (UInt_t)GetNFtsHits() ) throwError("GetFtsHit: hitId is too large.");
112  // TClonesArray *ftsHitArray= (TClonesArray *)FairRootManager::Instance()->GetObject("FTSHit");
113  const PndFtsHit *myHit = (PndFtsHit*) fFtsHitArray->At(hitId);
114  if (0 == myHit) throwError("GetFtsHit was not able to get the hit.");
115  return myHit;
116  };
TClonesArray * fFtsHitArray
Input array of PndFtsHit.
void throwError(const TString s) const
For error reporting.
Int_t GetNFtsHits() const
Returns the event number.
TClonesArray* PndFtsHoughTrackerTask::getFtsHitArrayPtr ( ) const
inline

Returns pointer to the hit array in which FTS hits are saved as PndFtsHit. Try not to use it.

See Also
GetFtsHit
GetNHits
Returns
Pointer to the hit array in which hits are saved as PndFtsHit.

Definition at line 159 of file PndFtsHoughTrackerTask.h.

References fFtsHitArray.

159 { return fFtsHitArray; };
TClonesArray * fFtsHitArray
Input array of PndFtsHit.
const TMatrixT< Double_t > PndFtsHoughTrackerTask::GetFtsHitCovMatrix ( const PndFtsHit *const  ftsHit) const

Returns the position covariance matrix (based on FTS straw geometry) for the hit with index hitId in the FTS hit array.

Parameters
[in]ftsHitpointer to hit for which the error should be returned.
Returns
Error in cm.

Definition at line 266 of file PndFtsHoughTrackerTask.cxx.

References Double_t, GetFtsTube(), PndFtsTube::GetHalfLength(), PndFtsTube::GetRadIn(), and PndFtsTube::GetRotationMatrix().

267 {
268  const PndFtsTube *const tube = GetFtsTube(ftsHit);
269 
270  const Double_t sizeSigmaCoeff = 1.5; // TODO Check value
271  const Double_t rhoError = tube->GetRadIn()/sizeSigmaCoeff;
272  const Double_t zError = tube->GetHalfLength()/sizeSigmaCoeff; // TODO might need additional factor
273  TMatrixT<Double_t> rotationMatrix = tube->GetRotationMatrix();
274 
275  TMatrixT<Double_t> unrotatedCovMatrix(3,3);
276  // initialize with 0
277  for (Int_t firstIdx=0; firstIdx < 3; ++firstIdx){
278  for (Int_t secondIdx=0; secondIdx < 3; ++secondIdx){
279  unrotatedCovMatrix[firstIdx][secondIdx] = 0;
280  }
281  }
282  unrotatedCovMatrix[0][0] = pow(rhoError, 2);
283  unrotatedCovMatrix[1][1] = pow(rhoError, 2);
284  unrotatedCovMatrix[2][2] = pow(zError, 2);
285 
286  TMatrixT<Double_t> rotatedCovMatrix = rotationMatrix*unrotatedCovMatrix;
287  rotatedCovMatrix *= rotationMatrix.Transpose(rotationMatrix);
288 
289  return rotatedCovMatrix;
290 }
Double_t
TMatrixT< Double_t > GetRotationMatrix() const
Definition: PndFtsTube.cxx:71
Double_t GetHalfLength() const
Definition: PndFtsTube.cxx:80
Double_t GetRadIn() const
Definition: PndFtsTube.cxx:74
const PndFtsTube * GetFtsTube(const PndFtsHit *const myHit) const
Returns pointer to the FTS tube corresponding to input FTS hit.
const TVector3 PndFtsHoughTrackerTask::GetFtsHitPosErrors ( const PndFtsHit *const  ftsHit) const

Returns the position error (based on FTS straw geometry) for the hit with index hitId in the FTS hit array.

Parameters
[in]pointerto hit for which the error should be returned.
Returns
Error in cm.

Definition at line 253 of file PndFtsHoughTrackerTask.cxx.

References Double_t, GetFtsTube(), PndFtsTube::GetHalfLength(), and PndFtsTube::GetRadIn().

Referenced by PndFtsHoughTrackCand::getTrackParPForHit().

254 {
255  const PndFtsTube *const tube = GetFtsTube(ftsHit);
256 
257  const Double_t sizeSigmaCoeff = 1.5; // TODO Check value
258  const Double_t rhoError = tube->GetRadIn()/sizeSigmaCoeff;
259  const Double_t zError = tube->GetHalfLength()/sizeSigmaCoeff; // TODO might need additional factor
260 
261  TVector3 hitPosErrors(rhoError,rhoError,zError);
262 
263  return hitPosErrors;
264 }
Double_t
Double_t GetHalfLength() const
Definition: PndFtsTube.cxx:80
Double_t GetRadIn() const
Definition: PndFtsTube.cxx:74
const PndFtsTube * GetFtsTube(const PndFtsHit *const myHit) const
Returns pointer to the FTS tube corresponding to input FTS hit.
const PndFtsTube* PndFtsHoughTrackerTask::GetFtsTube ( const PndFtsHit *const  myHit) const
inline

Returns pointer to the FTS tube corresponding to input FTS hit.

Parameters
[in]myHit,:FTS hit for which the tube should be returned.
Returns
Pointer to tube corresponding to myHit.

Definition at line 131 of file PndFtsHoughTrackerTask.h.

References fFtsTubeArray, and PndFtsHit::GetTubeID().

Referenced by PndFtsHoughSpace::CalculateHitPosFromIntersectionsWithZxTrackModel(), GetFtsHitCovMatrix(), and GetFtsHitPosErrors().

131  {
132  Int_t tubeID = myHit->GetTubeID();
133  const PndFtsTube *tube = (PndFtsTube*) fFtsTubeArray->At(tubeID);
134  return tube;
135  }
Int_t GetTubeID() const
Definition: PndFtsHit.h:70
TClonesArray * fFtsTubeArray
Input array of PndFtsTube (map of FTS tubes).
FairField* PndFtsHoughTrackerTask::getMagneticFieldPtr ( ) const
inline

Definition at line 166 of file PndFtsHoughTrackerTask.h.

References fField.

166 { return fField; };
FairField * fField
For B field access.
Int_t PndFtsHoughTrackerTask::getMcTruthIdForHitId ( UInt_t  hitId) const
inline

Definition at line 117 of file PndFtsHoughTrackerTask.h.

References fFtsMcPoints, GetFtsHit(), and throwError().

Referenced by PndFtsHoughSpace::WriteHistoOfAllPathsForEachMcTruthTrack().

117  {
118  const PndFtsHit *const ftsHit = GetFtsHit(hitId);
119  Int_t mcPointId=ftsHit->GetRefIndex();
120  if(0 > mcPointId) throwError("getMcTruthIdForHitId: negative mcPointId.");
121  FairMCPoint* myPoint = (FairMCPoint*)(fFtsMcPoints->At(mcPointId));
122  if(0==myPoint) throwError("getMcTruthIdForHitId: Could not get point belonging to hit.");
123  Int_t mcTrackId = myPoint->GetTrackID();
124  if(mcTrackId<0) throwError("getMcTruthIdForHitId: negative mcTrackId.");
125  return mcTrackId;
126  }
void throwError(const TString s) const
For error reporting.
TClonesArray * fFtsMcPoints
Input array of McPoints.
const PndFtsHit * GetFtsHit(UInt_t hitId) const
Returns pointer to the hit with index hitId in the FTS hit array.
Int_t PndFtsHoughTrackerTask::GetNFtsHits ( ) const
inline

Returns the event number.

Returns the number of FTS hits.

Returns
Number of FTS hits.

Definition at line 105 of file PndFtsHoughTrackerTask.h.

References fFtsHitArray.

Referenced by CheckForDuplicateFtsHits(), PndFtsHoughSpace::filterInputHits(), PndFtsHoughTrackFinder::FindTracks(), and GetFtsHit().

105 { return fFtsHitArray->GetEntriesFast(); };
TClonesArray * fFtsHitArray
Input array of PndFtsHit.
Bool_t PndPersistencyTask::GetPersistency ( )
inlineinherited

Definition at line 32 of file PndPersistencyTask.h.

References PndPersistencyTask::fPersistency.

Referenced by PndLmdPixelHitProducerFast::GetPersistance(), PndMdtDigitization::Init(), PndMdtHitProducerIdeal::Init(), PndMdtClusterTask::Init(), PndFtsHitProducerRealFast::Init(), PndDiscTaskReconstruction::Init(), PndRichHitProducer::Init(), PndSttHitProducerRealFast::Init(), PndSttHelixHitProducer::Init(), PndDiscTaskPID::Init(), PndIdealTrackFinder::Init(), PndSttMvdGemTracking::Init(), PndMdtTrkProducer::Init(), PndFtsHitProducerRealFull::Init(), PndLmdPixelClusterTask::Init(), PndSttHitProducerRealFull::Init(), PndLmdStripClusterTask::Init(), PndEmcApdHitProducer::Init(), PndMissingPzCleanerTask::Init(), PndEmcMakeRecoHit::Init(), PndEmcMakeClusterOnline::Init(), PndTrackSmearTask::Init(), PndEmcFWEndcapTimebasedWaveforms::Init(), PndSttHitProducerIdeal::Init(), PndEmcFWEndcapDigi::Init(), PndFtsHitProducerIdeal::Init(), PndEmcMakeCluster::Init(), PndMdtPointsToWaveform::Init(), PndDiscTaskDigitization::Init(), PndEmcMakeDigi::Init(), PndSdsTimeWalkCorrTask::Init(), PndLmdPixelHitProducerFast::Init(), PndDrcHitFinder::Init(), PndRichHitFinder::Init(), PndEmcMakeCorr::Init(), PndFtofHitProducerIdeal::Init(), PndEmcHitsToWaveform::Init(), PndSciTDigiTask::Init(), PndDrcHitProducerIdeal::Init(), PndSdsHitProducerIdeal::Init(), PndSciTHitProducerIdeal::Init(), PndEmcHitProducer::Init(), PndRecoMultiKalmanTask2::Init(), PndDrcHitProducerReal::Init(), PndDskFLGHitProducerIdeal::Init(), PndEmcTmpWaveformToDigi::Init(), PndDrcDigiTask::Init(), PndEmcWaveformToDigi::Init(), PndSttMatchTracks::Init(), PndEmcWaveformToCalibratedDigi::Init(), PndTrkTracking2::Init(), PndSttFindTracks::Init(), PndEmcMultiWaveformToCalibratedDigi::Init(), PndDrcTimeDigiTask::Init(), PndRecoKalmanTask2::Init(), PndEmcExpClusterSplitter::Init(), PndSdsNoiseProducer::Init(), Init(), PndEmcPhiBumpSplitter::Init(), PndSdsHybridHitProducer::Init(), PndSdsIdealRecoTask::Init(), PndRecoMultiKalmanTask::Init(), PndSdsIdealClusterTask::Init(), PndRecoKalmanTask::Init(), PndSdsStripHitProducerDif::Init(), PndSdsStripHitProducer::Init(), PndGemDigitize::Init(), PndGemFindHits::Init(), PndSdsPixelClusterTask::Init(), PndSdsStripClusterTask::Init(), PndMvdGemTrackFinderOnHits::Init(), PndBarrelTrackFinder::Init(), PndEmcFullDigiTask::PndEmcFullDigiTask(), PndEmcMakeBump::PndEmcMakeBump(), PndUnassignedHitsTask::RegisterBranches(), PndMvdClusterTask::SetPersistance(), PndMvdDigiTask::SetPersistance(), PndEmcMakeBump::SetStorageOfData(), and PndEmcFullDigiTask::StoreDigi().

32 { return fPersistency; }
Int_t PndFtsHoughTrackerTask::GetSaveDebugInfo ( ) const
inline
Int_t PndFtsHoughTrackerTask::GetVerbose ( ) const
inline
InitStatus PndFtsHoughTrackerTask::Init ( )
virtual

Initialization of task at the beginning of a run.

Reimplemented in PndFtsHoughTrackerTaskQA.

Definition at line 124 of file PndFtsHoughTrackerTask.cxx.

References fField, fFtsBranchId, fFtsHitArray, fFtsMcPoints, fFtsParameters, fFtsTubeArray, PndFtsMapCreator::FillTubeArray(), fLogger, fTrackCands, fTracks, fTracksArrayName, fVerbose, and PndPersistencyTask::GetPersistency().

Referenced by PndFtsHoughTrackerTaskQA::Init().

125 {
126  if(fVerbose>3) Info("Init","Initilization of PndFtsHoughTrackerTask");
127 
128  // Get a handle from the IO manager
129  FairRootManager* ioman = FairRootManager::Instance();
130  if ( ! ioman ) {
131  fLogger->Fatal(MESSAGE_ORIGIN,"RootManager not instantiated, return!");
132  return kFATAL;
133  }
134 
135  // Get a pointer to the previous already existing data level
136  /*
137  <InputDataLevel> = (TClonesArray*) ioman->GetObject("InputDataLevelName");
138  if ( ! <InputLevel> ) {
139  fLogger->Error(MESSAGE_ORIGIN,"No InputDataLevelName array!\n PndFtsHoughTrackerTask will be inactive");
140  return kERROR;
141  }
142  */
143 
144  // Create the TClonesArray for the output data and register
145  // it in the IO manager
146  /*
147  <OutputDataLevel> = new TClonesArray("OutputDataLevelName", 100);
148  ioman->Register("OutputDataLevelName","OutputDataLevelName",<OutputDataLevel>,kTRUE);
149  */
150 
151  // Do whatever else is needed at the initilization stage
152  // Create histograms to be filled
153  // initialize variables
154 
155 
156 
157 
158  // FTS Hits
159  fFtsHitArray= (TClonesArray *)ioman->GetObject("FTSHit");
160  if ( ! fFtsHitArray ) {
161  fLogger->Info(MESSAGE_ORIGIN,"No FTSHit array!");
162  return kERROR;
163  }
164  fFtsMcPoints = dynamic_cast<TClonesArray *> (ioman->GetObject("FTSPoint"));
165  if ( ! fFtsMcPoints ) {
166  fLogger->Info(MESSAGE_ORIGIN,"No McPoints array!");
167  return kERROR;
168  }
169 
170  // FTS Branch
171  fFtsBranchId = ioman->GetBranchId("FTSHit");
172 
173  // FTS Tube Array
174  PndFtsMapCreator mapperFts(fFtsParameters);
175  fFtsTubeArray = mapperFts.FillTubeArray();
176 
177 
178 
179  // B field
180  if(fVerbose>3) Info("Init","Try to get B field.");
181  fField = FairRunAna::Instance()->GetField();
182  if ( ! fField ) {
183  if(fVerbose>3) fLogger->Info(MESSAGE_ORIGIN,"No fField array!");
184  return kERROR;
185  }
186 
187  // Debugging
188  // if (fSaveDebugInfo){
189  // InitOutFileForDebugging();
190  // }
191  // FIXME this caused trouble, but I don't need it
192 // fHoughTrackCands = new TClonesArray("PndFtsHoughTrackCand");
193 // ioman->Register("FTSTrkDebugCand", "HoughTrackCand", fHoughTrackCands, fSaveDebugInfo);
194 
195  //fHoughSpaces = new TClonesArray("PndFtsHoughSpace");
196  //ioman->Register("FTSTrkDebugHS", "HoughSpaces", fHoughSpaces, fSaveDebugInfo);
197 
198 
199  // Output
200  fTrackCands = new TClonesArray("PndTrackCand");
201  fTracks = new TClonesArray("PndTrack");
202  ioman->Register(fTracksArrayName,"FTSTrk", fTracks, GetPersistency()); // for PndTrack
203  ioman->Register(fTracksArrayName+"Cand","FTSTrk", fTrackCands, GetPersistency()); // for PndTrackCand // TODO Is that correct, should it not be FTSTrkCand or something?
204 
205  if(3<fVerbose) Info("Register","Done.");
206 
207  return kSUCCESS;
208 
209 }
int fVerbose
Definition: poormantracks.C:24
FairField * fField
For B field access.
TClonesArray * fFtsHitArray
Input array of PndFtsHit.
TClonesArray * fTracks
Array of found tracks in PndTrack (for output)
Int_t fFtsBranchId
Detector Id of FTS.
TClonesArray * fFtsMcPoints
Input array of McPoints.
TString fTracksArrayName
Branch name where to store the Track candidates.
FairLogger * fLogger
Returns pointer to the B field.
PndGeoFtsPar * fFtsParameters
Needed for FTS map creator.
TClonesArray * fFtsTubeArray
Input array of PndFtsTube (map of FTS tubes).
TClonesArray * fTrackCands
Array of found track candidates in PndTrackCand (for output)
PndFtsHoughTrackerTask PndFtsHoughTrackerTask::operator= ( const PndFtsHoughTrackerTask )
protected
InitStatus PndFtsHoughTrackerTask::ReInit ( )
virtual

ReInitiliazation of task when the runID changes.

Reimplemented in PndFtsHoughTrackerTaskQA.

Definition at line 245 of file PndFtsHoughTrackerTask.cxx.

References fLogger, and fVerbose.

Referenced by PndFtsHoughTrackerTaskQA::ReInit().

246 {
247  InitStatus stat=kSUCCESS;
248  if(3<fVerbose) fLogger->Info(MESSAGE_ORIGIN,"Re- Initilization of PndFtsHoughTrackerTask");
249  return stat;
250 }
int fVerbose
Definition: poormantracks.C:24
FairLogger * fLogger
Returns pointer to the B field.
void PndFtsHoughTrackerTask::SetParContainers ( )
virtual

Loads the parameter container from the runtime database.

Definition at line 114 of file PndFtsHoughTrackerTask.cxx.

References fFtsParameters, fLogger, fVerbose, and rtdb.

115 {
116  if(fVerbose>3) fLogger->Info(MESSAGE_ORIGIN,"SetParContainers of PndFtsHoughTrackerTask");
117 
118  // FTS parameters
119  FairRuntimeDb *rtdb= FairRun::Instance()->GetRuntimeDb();
120  fFtsParameters=(PndGeoFtsPar*)(rtdb->getContainer("PndGeoFtsPar"));
121 }
int fVerbose
Definition: poormantracks.C:24
FairRuntimeDb * rtdb
Definition: hit_dirc.C:66
FairLogger * fLogger
Returns pointer to the B field.
PndGeoFtsPar * fFtsParameters
Needed for FTS map creator.
void PndFtsHoughTrackerTask::SetPersistence ( Bool_t  val)
inline

Definition at line 86 of file PndFtsHoughTrackerTask.h.

References PndPersistencyTask::SetPersistency().

86 { SetPersistency(val);};
Double_t val[nBoxes][nFEBox]
Definition: createCalib.C:11
void SetPersistency(Bool_t val=kTRUE)
void PndPersistencyTask::SetPersistency ( Bool_t  val = kTRUE)
inlineinherited

Definition at line 31 of file PndPersistencyTask.h.

References PndPersistencyTask::fPersistency, and val.

Referenced by barrelTrackFinder(), digi_complete(), digi_complete_newSTT(), digiOnly_complete(), PndBarrelTrackFinder::PndBarrelTrackFinder(), PndCATracking::PndCATracking(), PndDrcHitFinder::PndDrcHitFinder(), PndEmc2DLocMaxFinder::PndEmc2DLocMaxFinder(), PndEmcExpClusterSplitter::PndEmcExpClusterSplitter(), PndEmcFullDigiTask::PndEmcFullDigiTask(), PndEmcFWEndcapDigi::PndEmcFWEndcapDigi(), PndEmcFWEndcapTimebasedWaveforms::PndEmcFWEndcapTimebasedWaveforms(), PndEmcHitProducer::PndEmcHitProducer(), PndEmcHitsToWaveform::PndEmcHitsToWaveform(), PndEmcMakeBump::PndEmcMakeBump(), PndEmcMakeCluster::PndEmcMakeCluster(), PndEmcMakeClusterOnline::PndEmcMakeClusterOnline(), PndEmcMakeDigi::PndEmcMakeDigi(), PndEmcMakeRecoHit::PndEmcMakeRecoHit(), PndEmcMultiWaveformToCalibratedDigi::PndEmcMultiWaveformToCalibratedDigi(), PndEmcPhiBumpSplitter::PndEmcPhiBumpSplitter(), PndEmcTmpWaveformToDigi::PndEmcTmpWaveformToDigi(), PndEmcWaveformToCalibratedDigi::PndEmcWaveformToCalibratedDigi(), PndEmcWaveformToDigi::PndEmcWaveformToDigi(), PndFtofHitProducerIdeal::PndFtofHitProducerIdeal(), PndFtsCATracking::PndFtsCATracking(), PndFtsHitProducerIdeal::PndFtsHitProducerIdeal(), PndFtsHitProducerRealFast::PndFtsHitProducerRealFast(), PndFtsHitProducerRealFull::PndFtsHitProducerRealFull(), PndFtsHoughTrackerTask(), PndGemDigitize::PndGemDigitize(), PndGemFindHits::PndGemFindHits(), PndIdealTrackFinder::PndIdealTrackFinder(), PndLmdPixelClusterTask::PndLmdPixelClusterTask(), PndLmdPixelHitProducerFast::PndLmdPixelHitProducerFast(), PndMdtClusterTask::PndMdtClusterTask(), PndMdtDigitization::PndMdtDigitization(), PndMdtHitProducerIdeal::PndMdtHitProducerIdeal(), PndMdtPointsToWaveform::PndMdtPointsToWaveform(), PndMdtTrkProducer::PndMdtTrkProducer(), PndMissingPzCleanerTask::PndMissingPzCleanerTask(), PndMvdGemTrackFinderOnHits::PndMvdGemTrackFinderOnHits(), PndMvdHitProducerIdeal::PndMvdHitProducerIdeal(), PndMvdPixelClusterTask::PndMvdPixelClusterTask(), PndMvdTimeWalkCorrTask::PndMvdTimeWalkCorrTask(), PndMvdToPix4ClusterTask::PndMvdToPix4ClusterTask(), PndRecoKalmanTask::PndRecoKalmanTask(), PndRecoKalmanTask2::PndRecoKalmanTask2(), PndRecoMultiKalmanTask::PndRecoMultiKalmanTask(), PndRecoMultiKalmanTask2::PndRecoMultiKalmanTask2(), PndRichHitFinder::PndRichHitFinder(), PndRichHitProducer::PndRichHitProducer(), PndSciTDigiTask::PndSciTDigiTask(), PndSciTHitProducerIdeal::PndSciTHitProducerIdeal(), PndSdsHitProducerIdeal::PndSdsHitProducerIdeal(), PndSdsHybridHitProducer::PndSdsHybridHitProducer(), PndSdsIdealClusterTask::PndSdsIdealClusterTask(), PndSdsIdealRecoTask::PndSdsIdealRecoTask(), PndSdsNoiseProducer::PndSdsNoiseProducer(), PndSdsPixelClusterTask::PndSdsPixelClusterTask(), PndSdsStripClusterTask::PndSdsStripClusterTask(), PndSdsStripHitProducer::PndSdsStripHitProducer(), PndSdsTimeWalkCorrTask::PndSdsTimeWalkCorrTask(), PndSttFindTracks::PndSttFindTracks(), PndSttHelixHitProducer::PndSttHelixHitProducer(), PndSttHitProducerIdeal::PndSttHitProducerIdeal(), PndSttHitProducerRealFast::PndSttHitProducerRealFast(), PndSttHitProducerRealFull::PndSttHitProducerRealFull(), PndSttMatchTracks::PndSttMatchTracks(), PndSttMvdGemTracking::PndSttMvdGemTracking(), PndTrackSmearTask::PndTrackSmearTask(), PndTrkTracking2::PndTrkTracking2(), reco(), reco_complete(), reco_complete_gf2(), reco_complete_newSTT(), reco_complete_sec(), recoideal_complete(), PndMvdClusterTask::SetPersistance(), PndMvdDigiTask::SetPersistance(), PndLmdPixelHitProducerFast::SetPersistance(), PndSdsHitProducerIdeal::SetPersistance(), PndSttMvdGemTracking::SetPersistenc(), PndMdtClusterTask::SetPersistence(), PndSttHelixHitProducer::SetPersistence(), PndMissingPzCleanerTask::SetPersistence(), PndFtsHitProducerRealFast::SetPersistence(), PndFtsHitProducerRealFull::SetPersistence(), PndSttHitProducerIdeal::SetPersistence(), PndSttHitProducerRealFull::SetPersistence(), PndSttHitProducerRealFast::SetPersistence(), PndFtsHitProducerIdeal::SetPersistence(), PndTrackSmearTask::SetPersistence(), PndSciTHitProducerIdeal::SetPersistence(), PndIdealTrackFinder::SetPersistence(), PndSttMatchTracks::SetPersistence(), PndSttFindTracks::SetPersistence(), SetPersistence(), PndTrkTracking2::SetPersistence(), PndEmcMakeRecoHit::SetStorageOfData(), PndEmcMakeClusterOnline::SetStorageOfData(), PndEmcFWEndcapDigi::SetStorageOfData(), PndEmcFWEndcapTimebasedWaveforms::SetStorageOfData(), PndEmcMakeDigi::SetStorageOfData(), PndMdtPointsToWaveform::SetStorageOfData(), PndEmc2DLocMaxFinder::SetStorageOfData(), PndEmcMakeCluster::SetStorageOfData(), PndEmcHitsToWaveform::SetStorageOfData(), PndEmcMakeBump::SetStorageOfData(), PndEmcTmpWaveformToDigi::SetStorageOfData(), PndEmcWaveformToDigi::SetStorageOfData(), PndEmcWaveformToCalibratedDigi::SetStorageOfData(), PndEmcMultiWaveformToCalibratedDigi::SetStorageOfData(), PndEmcExpClusterSplitter::SetStorageOfData(), PndEmcPhiBumpSplitter::SetStorageOfData(), standard_tracking(), and PndEmcFullDigiTask::StoreDigi().

31 { fPersistency = val; }
Double_t val[nBoxes][nFEBox]
Definition: createCalib.C:11
void PndFtsHoughTrackerTask::SetSaveDebugInfo ( Int_t  saveDebugInfo)
inline

Definition at line 87 of file PndFtsHoughTrackerTask.h.

References fSaveDebugInfo.

Referenced by reco().

87 { fSaveDebugInfo = saveDebugInfo;};
Int_t fSaveDebugInfo
Debug information will be created if &gt;0.
void PndFtsHoughTrackerTask::SetVerbose ( Int_t  verbose)
inline

Definition at line 85 of file PndFtsHoughTrackerTask.h.

References fVerbose, and verbose.

Referenced by reco(), and reco_fts().

85 { fVerbose = verbose;};
int fVerbose
Definition: poormantracks.C:24
#define verbose
void PndFtsHoughTrackerTask::throwError ( const TString  s) const
inlineprotected

For error reporting.

Definition at line 225 of file PndFtsHoughTrackerTask.h.

Referenced by GetFtsHit(), and getMcTruthIdForHitId().

225 { throw std::runtime_error(s.Data()); };
TLorentzVector s
Definition: Pnd2DStar.C:50

Member Data Documentation

UInt_t PndFtsHoughTrackerTask::fEventNr
protected

Event number for debugging purposes.

Definition at line 223 of file PndFtsHoughTrackerTask.h.

Referenced by PndFtsHoughTrackerTaskQA::Exec(), Exec(), and GetEventNr().

FairField* PndFtsHoughTrackerTask::fField
protected

For B field access.

Definition at line 202 of file PndFtsHoughTrackerTask.h.

Referenced by getMagneticFieldPtr(), and Init().

Int_t PndFtsHoughTrackerTask::fFtsBranchId
protected

Detector Id of FTS.

Definition at line 185 of file PndFtsHoughTrackerTask.h.

Referenced by getFtsBranchId(), and Init().

TClonesArray* PndFtsHoughTrackerTask::fFtsHitArray
protected

Input array of PndFtsHit.

Definition at line 186 of file PndFtsHoughTrackerTask.h.

Referenced by GetFtsHit(), getFtsHitArrayPtr(), GetNFtsHits(), and Init().

TClonesArray* PndFtsHoughTrackerTask::fFtsMcPoints
protected

Input array of McPoints.

Definition at line 187 of file PndFtsHoughTrackerTask.h.

Referenced by getMcTruthIdForHitId(), and Init().

PndGeoFtsPar* PndFtsHoughTrackerTask::fFtsParameters
protected

Needed for FTS map creator.

I don't really know what that does...

Definition at line 194 of file PndFtsHoughTrackerTask.h.

Referenced by Init(), and SetParContainers().

TClonesArray* PndFtsHoughTrackerTask::fFtsTubeArray
protected

Input array of PndFtsTube (map of FTS tubes).

Is filled by map creator

Definition at line 199 of file PndFtsHoughTrackerTask.h.

Referenced by GetFtsTube(), and Init().

FairLogger* PndFtsHoughTrackerTask::fLogger
Int_t PndFtsHoughTrackerTask::fSaveDebugInfo
protected

Debug information will be created if >0.

Definition at line 222 of file PndFtsHoughTrackerTask.h.

Referenced by GetSaveDebugInfo(), and SetSaveDebugInfo().

TClonesArray* PndFtsHoughTrackerTask::fTrackCands
protected

Array of found track candidates in PndTrackCand (for output)

Definition at line 209 of file PndFtsHoughTrackerTask.h.

Referenced by Exec(), PndFtsHoughTrackerTaskQA::Finish(), Finish(), FinishEvent(), Init(), and ~PndFtsHoughTrackerTask().

TClonesArray* PndFtsHoughTrackerTask::fTracks
protected

Array of found tracks in PndTrack (for output)

Definition at line 210 of file PndFtsHoughTrackerTask.h.

Referenced by Exec(), FinishEvent(), Init(), and ~PndFtsHoughTrackerTask().

TString PndFtsHoughTrackerTask::fTracksArrayName
protected

Branch name where to store the Track candidates.

Definition at line 208 of file PndFtsHoughTrackerTask.h.

Referenced by Init().


The documentation for this class was generated from the following files: