star-travex
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
TStiKalmanTrackNode.cxx
Go to the documentation of this file.
2 
3 #include <boost/regex.hpp>
4 #include <algorithm>
5 
6 #include "StarClassLibrary/StThreeVector.hh"
7 #include "Sti/StiPlacement.h"
9 
10 
12 
13 
15  fStiTrackNode(nullptr),
16  fTrack(nullptr), fValid(false), fIsInsideVolume(-1),
17  fPosition(),
18  fError(),
19  fPositionLocal(),
20  fProjPositionLocal(),
21  fProjError(),
22  fTrackP(), fEnergyLosses(-1), fNodeRadius(0), fNodeCenterRefAngle(0), fNodeMaterialDensity(0),
23  fNodeTrackLength(0),
24  fNodeRelRadLength(0), fVolumeName(""), fAcceptedStiHit(nullptr), fClosestStiHit(nullptr),
25  fCandidateStiHits()
26 {
27 }
28 
29 
30 TStiKalmanTrackNode::TStiKalmanTrackNode(const StiKalmanTrackNode &stiKTN, TStiKalmanTrack* const parent) : TObject(),
31  fStiTrackNode(&stiKTN),
32  fTrack(parent), fValid(stiKTN.isValid()), fIsInsideVolume(stiKTN.inside(1+2+4)),
33  fPosition(stiKTN.x_g(), stiKTN.y_g(), stiKTN.z_g()),
34  fError(stiKTN.fitErrs()._cXX, stiKTN.fitErrs()._cYY, stiKTN.fitErrs()._cZZ),
35  fPositionLocal(stiKTN.x(), stiKTN.y(), stiKTN.z()),
36  fProjPositionLocal(),
37  fProjError(),
38  fTrackP(stiKTN.getGlobalMomentum().x(), stiKTN.getGlobalMomentum().y(), stiKTN.getGlobalMomentum().z()),
39  fEnergyLosses(-1), fNodeRadius(0), fNodeCenterRefAngle(0), fNodeMaterialDensity(0),
40  fNodeTrackLength(stiKTN.getTrackLength()),
41  fNodeRelRadLength(0), fVolumeName(""), fAcceptedStiHit(nullptr), fClosestStiHit(nullptr),
42  fCandidateStiHits()
43 {
44  fEnergyLosses = stiKTN.getEnergyLosses() * 1e6; // Get losses in volume material and convert GeV to keV
45  fNodeRelRadLength = stiKTN.getRelRadLength();
46 
47  const StiDetector* stiKTNDet = stiKTN.getDetector();
48 
49  if (stiKTNDet)
50  {
51  fVolumeName = stiKTNDet->getName();
52  StiPlacement* stiPlacement = stiKTNDet->getPlacement();
53  assert(stiPlacement);
54  fNodeRadius = stiPlacement->getLayerRadius();
55  fNodeCenterRefAngle = stiPlacement->getCenterRefAngle();
56 
57  StiMaterial* stiMaterial = stiKTNDet->getMaterial();
58  assert(stiMaterial);
59 
60  fNodeMaterialDensity = stiMaterial->getDensity();
61  }
62 
63  if (stiKTN.getHit())
64  {
65  auto resultPair = fTrack->AddToParentEvent( TStiHit(*stiKTN.getHit()) );
66  // Save the pointer to the hit in the parent event
67  fAcceptedStiHit = &(*resultPair.first);
68  }
69 
70  // Note (as of now) this info is available for nodes with r < kRMinTpc [= 55 cm]
71  // see node->saveInfo() in StiKalmanTrackFinder.cxx
72  const StiNodeInf* prefitKTNParams = stiKTN.getInfo();
73 
74  if (prefitKTNParams) {
75  fProjPositionLocal.SetXYZ(prefitKTNParams->mPP.x(), prefitKTNParams->mPP.y(), prefitKTNParams->mPP.z());
76  fProjError.SetXYZ(sqrt(prefitKTNParams->mPE._cXX), sqrt(prefitKTNParams->mPE._cYY), sqrt(prefitKTNParams->mPE._cZZ));
77  }
78 }
79 
80 
81 TVector3 TStiKalmanTrackNode::CalcPullToHit(const TStiHit& hit) const
82 {
83  double pullX = (hit.GetPositionLocal().X() - fProjPositionLocal.X()) / fProjError.X();
84  double pullY = (hit.GetPositionLocal().Y() - fProjPositionLocal.Y()) / fProjError.Y();
85  double pullZ = (hit.GetPositionLocal().Z() - fProjPositionLocal.Z()) / fProjError.Z();
86 
87  return TVector3(pullX, pullY, pullZ);
88 }
89 
90 
96 {
97  if (!fClosestStiHit) return TVector3(DBL_MAX, DBL_MAX, DBL_MAX);
98 
100 }
101 
102 
106 double TStiKalmanTrackNode::CalcChi2(const TStiHit& hit) const
107 {
108  if (!fStiTrackNode || !hit.GetStiHit())
109  return -1;
110 
111  return fStiTrackNode->evaluateChi2Info(hit.GetStiHit());
112 }
113 
114 
115 bool TStiKalmanTrackNode::MatchedVolName(const std::string & pattern) const
116 {
117  if (fVolumeName.empty()) return true;
118 
119  boost::regex r(pattern);
120  bool matched = boost::regex_match(fVolumeName, r);
121 
122  //if (matched)
123  // Info("MatchedVolName", "Volume [%s] matched pattern [%s]", fVolumeName.c_str(), pattern.c_str());
124 
125  return matched;
126 }
127 
128 
129 bool TStiKalmanTrackNode::MatchedVolName(const std::set<std::string> & patterns) const
130 {
131  if (fVolumeName.empty() || patterns.empty())
132  return true;
133 
134  std::set<std::string>::const_iterator iPattern = patterns.begin();
135 
136  for( ; iPattern != patterns.end(); ++iPattern )
137  {
138  //Info("MatchedVolName", "Looking for pattern [%s] in fVolumeName [%s]", (*iPattern).c_str(), fVolumeName.c_str());
139  if ( MatchedVolName(*iPattern) )
140  return true;
141  }
142 
143  return false;
144 }
145 
146 
147 void TStiKalmanTrackNode::Print(Option_t *opt) const
148 {
149  fPosition.Print();
150  fTrackP.Print();
151  printf("fEnergyLosses: %f\n", fEnergyLosses);
152  printf("fNodeRadius: %f\n", fNodeRadius);
153  printf("fVolumeName: %s\n", fVolumeName.c_str());
154 }
155 
156 
158 {
159  return lhs.fNodeRadius < rhs.fNodeRadius;
160 }
161 
162 
163 std::set<const TStiHit*> TStiKalmanTrackNode::GetCandidateHits() const
164 {
165  std::set<const TStiHit*> candidateHits;
166 
167  std::transform(fCandidateStiHits.begin(), fCandidateStiHits.end(),
168  std::inserter(candidateHits, candidateHits.begin()), TStiHitProxy::GetBareStiHit);
169 
170  return candidateHits;
171 }
172 
173 
178 void TStiKalmanTrackNode::FindClosestHit(const std::set<TStiHit>& stiHits) const
179 {
180  TVector3 distVec;
181  double min_dist = DBL_MAX;
182 
183  for (const auto& hit : stiHits)
184  {
185  distVec = GetPosition() - hit.GetPosition();
186 
187  double dist = distVec.Mag();
188  if (dist < min_dist) {
189  min_dist = dist;
190  fClosestStiHit = &hit;
191  }
192  }
193 }
194 
195 
209 void TStiKalmanTrackNode::FindCandidateHits(const std::set<TStiHit>& stiHits) const
210 {
211  TVector3 distVec; // Can be a static variable to avoid ROOT's memory allocation and overhead...
212  double min_dist = DBL_MAX;
213 
214  for (const auto& hit : stiHits)
215  {
216  // We require the hit to be associated with the same detector volume as
217  // the track node itself
218  if (fVolumeName != hit.GetVolumeName())
219  continue;
220 
221  distVec = GetPositionLocal() - hit.GetPositionLocal();
222 
223  if (fabs(distVec.Y()) < 5*fProjError.Y() &&
224  fabs(distVec.Z()) < 5*fProjError.Z() )
225  {
226  fCandidateStiHits.insert( TStiHitProxy(hit, *this) );
227 
228  double dist = distVec.Mag();
229  if (dist < min_dist) {
230  min_dist = dist;
231  fClosestStiHit = &hit;
232  }
233  }
234  }
235 }
236 
237 
238 void TStiKalmanTrackNode::FindCandidateHitsByChi2(const std::set<TStiHit>& stiHits) const
239 {
240  for (const auto& hit : stiHits)
241  {
242  if (fVolumeName != hit.GetVolumeName()) continue;
243 
244  if (CalcChi2(hit) < 20)
245  {
246  fCandidateStiHits.insert(TStiHitProxy(hit, *this));
247  }
248  }
249 }
virtual void Print(Option_t *opt="") const
TVector3 fProjError
The projection error to the node before the re-fit.
static const TStiHit * GetBareStiHit(const TStiHitProxy &hitProxy)
Definition: TStiHitProxy.h:22
float fNodeRadius
The nominal radius of the Sti volume associated with this node.
const TVector3 & GetPositionLocal() const
Definition: TStiHit.h:20
TVector3 fPosition
Global coordinates of the final (post re-fit) track node position.
void FindCandidateHitsByChi2(const std::set< TStiHit > &stiHits) const
float fNodeCenterRefAngle
Angle to the center of the Sti volume associated with this node.
std::string fVolumeName
Name of Sti volume.
const TVector3 & GetPosition() const
const StiHit * GetStiHit() const
Definition: TStiHit.h:17
TVector3 CalcPullToHit(const TStiHit &hit) const
TVector3 CalcPullClosestHit() const
Calculates and returns uncorrelated three components which can be used in the pull distribution...
const StiKalmanTrackNode * fStiTrackNode
Transient pointer to original StiKalmanTrackNode to access non-persistent info.
const TStiHit * fAcceptedStiHit
Pointer to the hit associated with this node by the reconstruction algorithm, if any.
void FindClosestHit(const std::set< TStiHit > &stiHits) const
Traverses the provided collection of hits and sets the internal pointer to the hit closest to this tr...
bool operator<(const TStiKalmanTrackNode &lhs, const TStiKalmanTrackNode &rhs)
bool MatchedVolName(const std::string &pattern) const
std::set< TStiHitProxy > fCandidateStiHits
Collection of hits in some proximity of mean track projection.
const TStiHit * fClosestStiHit
Pointer to the hit closest to this node if any.
void FindCandidateHits(const std::set< TStiHit > &stiHits) const
In the provided collection of hits this class method finds all hits within a 5-sigma vicinity around ...
TVector3 fProjPositionLocal
Local coordinates of the projected (pre re-fit) track node position.
ClassImp(TStiKalmanTrackNode) TStiKalmanTrackNode
float fNodeRelRadLength
Relative radiation length.
const TVector3 & GetPositionLocal() const
float fEnergyLosses
Energy lost in the volume.
float fNodeMaterialDensity
Density of the material of this node/volume.
std::pair< std::set< TStiHit >::iterator, bool > AddToParentEvent(const TStiHit &stiHit)
double CalcChi2(const TStiHit &hit) const
Calculate chi2 for a user provided hit using the original Sti methods.
TStiKalmanTrack * fTrack
transient member
TVector3 fTrackP
Track momentum vector in global CS.
std::set< const TStiHit * > GetCandidateHits() const