MpdKalmanTrack.cxx 6.7 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253
  1. /// \class MpdKalmanTrack
  2. ///
  3. /// Kalman filter track object for the MPD detector
  4. /// \author Alexander Zinchenko (LHEP, JINR, Dubna)
  5. #include "MpdKalmanTrack.h"
  6. #include "MpdKalmanFilter.h"
  7. #include <TMatrixD.h>
  8. #include <TMatrixDSym.h>
  9. #include <TObjArray.h>
  10. #include <TVector3.h>
  11. #include <TMath.h>
  12. #include <Riostream.h>
  13. //TpcLheKalmanFilter* TpcLheKalmanTrack::fgKF = 0x0;
  14. //__________________________________________________________________________
  15. MpdKalmanTrack::MpdKalmanTrack()
  16. : TObject(),
  17. fID(0),
  18. fNhits(0),
  19. fTrackDir(kInward),
  20. fTrackType(kBarrel),
  21. fLastLay(0),
  22. fNofWrong(0),
  23. fNode(""),
  24. fNodeNew(""),
  25. fPartID(0),
  26. fPos(0.),
  27. fPosNew(0.),
  28. fPosAtHit(0.),
  29. fChi2(0.),
  30. fChi2Vertex(0.),
  31. fLength(0.),
  32. fLengAtHit(0.),
  33. fParam(new TMatrixD(5,1)),
  34. fParamNew(new TMatrixD(5,1)),
  35. fParamAtHit(new TMatrixD(5,1)),
  36. fCovar(new TMatrixDSym(5)),
  37. fWeight(0x0),
  38. fWeightAtHit(new TMatrixDSym(5)),
  39. fHits(0x0),
  40. fFlag(kOk)
  41. {
  42. /// Default constructor
  43. }
  44. //__________________________________________________________________________
  45. MpdKalmanTrack::MpdKalmanTrack(Double_t pos, TVector3 &vertex)
  46. : TObject(),
  47. fID(0),
  48. fNhits(0),
  49. fTrackDir(kInward),
  50. fTrackType(kBarrel),
  51. fLastLay(0),
  52. fNofWrong(0),
  53. fNode(""),
  54. fNodeNew(""),
  55. fPartID(0),
  56. fPos(pos),
  57. fPosNew(pos),
  58. fPosAtHit(0.),
  59. fChi2(0.),
  60. fChi2Vertex(0.),
  61. fLength(0.),
  62. fLengAtHit(0.),
  63. fParam(new TMatrixD(5,1)),
  64. fParamNew(new TMatrixD(5,1)),
  65. fParamAtHit(new TMatrixD(5,1)),
  66. fCovar(new TMatrixDSym(5)),
  67. fWeight(new TMatrixDSym(5)),
  68. fWeightAtHit(new TMatrixDSym(5)),
  69. fVertex(vertex),
  70. fHits(new TObjArray(70)),
  71. fFlag(kOk)
  72. {
  73. /// Constructor from a track position and vertex
  74. //fHits->SetOwner(kFALSE);
  75. }
  76. //__________________________________________________________________________
  77. MpdKalmanTrack::MpdKalmanTrack (const MpdKalmanTrack& track)
  78. : TObject(track),
  79. fID(track.fID),
  80. fNhits(track.fNhits),
  81. fTrackDir(track.fTrackDir),
  82. fTrackType(track.fTrackType),
  83. fLastLay(track.fLastLay),
  84. fNofWrong(track.fNofWrong),
  85. fNode(track.fNode),
  86. fNodeNew(track.fNodeNew),
  87. fPartID(track.fPartID),
  88. fPos(track.fPos),
  89. fPosNew(track.fPosNew),
  90. fPosAtHit(track.fPosAtHit),
  91. fChi2(track.fChi2),
  92. fChi2Vertex(track.fChi2Vertex),
  93. fLength(track.fLength),
  94. fLengAtHit(track.fLengAtHit),
  95. fParam(new TMatrixD(*(track.fParam))),
  96. //fParamNew(new TMatrixD(*(track.fParamNew))),
  97. fParamAtHit(new TMatrixD(*(track.fParamAtHit))),
  98. fCovar(new TMatrixDSym(*(track.fCovar))),
  99. //fWeight(new TMatrixDSym(*(track.fWeight))),
  100. fWeightAtHit(new TMatrixDSym(*(track.fWeightAtHit))),
  101. fVertex(track.fVertex),
  102. fHits(new TObjArray(70)),
  103. fStepMap(track.fStepMap),
  104. fFlag(track.fFlag)
  105. {
  106. ///copy constructor
  107. if (track.fParamNew) fParamNew = new TMatrixD(*(track.fParamNew));
  108. else fParamNew = new TMatrixD(5,1);
  109. if (track.fWeight) fWeight = new TMatrixDSym(*(track.fWeight));
  110. else fWeight = new TMatrixDSym(5);
  111. if (track.fHits == 0x0) return;
  112. Int_t nHits = track.fHits->GetEntriesFast();
  113. for (Int_t i = 0; i < nHits; ++i) {
  114. fHits->Add(track.fHits->UncheckedAt(i));
  115. }
  116. }
  117. //__________________________________________________________________________
  118. MpdKalmanTrack & MpdKalmanTrack::operator=(const MpdKalmanTrack& track)
  119. {
  120. /// Asignment operator
  121. // check assignement to self
  122. if (this == &track) return *this;
  123. // base class assignement
  124. TObject::operator=(track);
  125. Clear();
  126. fID = track.fID;
  127. fNhits = track.fNhits;
  128. fTrackDir = track.fTrackDir;
  129. fTrackType = track.fTrackType;
  130. fLastLay = track.fLastLay;
  131. fNofWrong = track.fNofWrong;
  132. fNode = track.fNode;
  133. fNodeNew = track.fNodeNew;
  134. fPartID = track.fPartID;
  135. fPos = track.fPos;
  136. fPosNew = track.fPosNew;
  137. fPosAtHit = track.fPosAtHit;
  138. fChi2 = track.fChi2;
  139. fChi2Vertex = track.fChi2Vertex;
  140. fLength = track.fLength;
  141. fLengAtHit = track.fLengAtHit;
  142. fParam = new TMatrixD(*(track.fParam));
  143. fParamAtHit = new TMatrixD(*(track.fParamAtHit));
  144. fCovar = new TMatrixDSym(*(track.fCovar));
  145. fWeightAtHit = new TMatrixDSym(*(track.fWeightAtHit));
  146. fVertex = track.fVertex;
  147. fHits = new TObjArray(70);
  148. fStepMap = track.fStepMap;
  149. fFlag = track.fFlag;
  150. if (track.fParamNew) fParamNew = new TMatrixD(*(track.fParamNew));
  151. else fParamNew = new TMatrixD(5,1);
  152. if (track.fWeight) fWeight = new TMatrixDSym(*(track.fWeight));
  153. else fWeight = new TMatrixDSym(5);
  154. if (track.fHits == 0x0) return *this;
  155. Int_t nHits = track.fHits->GetEntriesFast();
  156. for (Int_t i = 0; i < nHits; ++i) {
  157. fHits->Add(track.fHits->UncheckedAt(i));
  158. }
  159. return *this;
  160. }
  161. //__________________________________________________________________________
  162. MpdKalmanTrack::~MpdKalmanTrack()
  163. {
  164. /// Destructor
  165. delete fParam;
  166. delete fParamNew;
  167. delete fParamAtHit;
  168. fParam = fParamNew = fParamAtHit = NULL;
  169. delete fCovar; fCovar = NULL;
  170. delete fWeight;
  171. delete fWeightAtHit;
  172. fWeight = fWeightAtHit = NULL;
  173. delete fHits; fHits = NULL;
  174. }
  175. //__________________________________________________________________________
  176. void MpdKalmanTrack::Clear()
  177. {
  178. /// Reset track
  179. delete fParam;
  180. delete fParamNew;
  181. delete fParamAtHit;
  182. fParam = fParamNew = fParamAtHit = NULL;
  183. delete fCovar; fCovar = NULL;
  184. delete fWeight;
  185. delete fWeightAtHit;
  186. fWeight = fWeightAtHit = NULL;
  187. delete fHits; fHits = NULL;
  188. fStepMap.clear();
  189. }
  190. //__________________________________________________________________________
  191. TMatrixDSym* MpdKalmanTrack::Weight2Cov()
  192. {
  193. /// Get covariance matrix from weight matrix
  194. *fCovar = *fWeight;
  195. //fCovar->Invert();
  196. Int_t iok = 0;
  197. MpdKalmanFilter::Instance()->MnvertLocal(fCovar->GetMatrixArray(), 5, 5, 5, iok);
  198. return fCovar;
  199. }
  200. //__________________________________________________________________________
  201. void MpdKalmanTrack::ReSetWeight()
  202. {
  203. /// Set weight matrix from covariance
  204. *fWeight = *fCovar;
  205. Int_t iok = 0;
  206. MpdKalmanFilter::Instance()->MnvertLocal(fWeight->GetMatrixArray(), 5, 5, 5, iok);
  207. }
  208. //__________________________________________________________________________
  209. Double_t MpdKalmanTrack::DCA() const
  210. {
  211. /// Return signed track DCA
  212. if (fPosNew < 1.e-6) return fPosNew;
  213. Float_t phiTr = GetParam(2);
  214. TVector3 vecTr(TMath::Cos(phiTr), TMath::Sin(phiTr), 0.0);
  215. Float_t phi = GetParam(0) / fPosNew;
  216. TVector3 vec(TMath::Cos(phi), TMath::Sin(phi), 0.0);
  217. return fPosNew * TMath::Sign(1.0,vec.Cross(vecTr).Z());
  218. }
  219. //__________________________________________________________________________
  220. void MpdKalmanTrack::PrintNode(Bool_t nNew) const
  221. {
  222. /// Print node name (for debugging)
  223. if (nNew) printf("%s\n",fNodeNew.Data());
  224. else printf("%s\n",fNode.Data());
  225. }
  226. ///__________________________________________________________________________
  227. ClassImp(MpdKalmanTrack)