CGEM BOSS 6.6.5.g
BESIII Offline Software System
Loading...
Searching...
No Matches
KalmanKinematicFit.h
Go to the documentation of this file.
1#ifndef VertexFit_KalmanKinematicFit_H
2#define VertexFit_KalmanKinematicFit_H
3
4#include <vector>
9#include "TGraph2D.h"
10
12
13 public:
14 // constructor & destructor
15
18 //
19 // Resonance Constraints
20 //
21 void AddResonance(int number, double mres, std::vector<int> tlis);
22 void AddResonance(int number, double mres, int n1);
23 void AddResonance(int number, double mres, int n1, int n2);
24 void AddResonance(int number, double mres, int n1, int n2, int n3);
25 void AddResonance(int number, double mres, int n1, int n2, int n3, int n4);
26 void AddResonance(int number, double mres, int n1, int n2, int n3, int n4,
27 int n5);
28 void AddResonance(int number, double mres, int n1, int n2, int n3, int n4,
29 int n5, int n6);
30 void AddResonance(int number, double mres, int n1, int n2, int n3, int n4,
31 int n5, int n6, int n7);
32 void AddResonance(int number, double mres, int n1, int n2, int n3, int n4,
33 int n5, int n6, int n7, int n8);
34 void AddResonance(int number, double mres, int n1, int n2, int n3, int n4,
35 int n5, int n6, int n7, int n8, int n9);
36 void AddResonance(int number, double mres, int n1, int n2, int n3, int n4,
37 int n5, int n6, int n7, int n8, int n9,
38 int n10);
39 void AddResonance(int number, double mres, int n1, int n2, int n3, int n4,
40 int n5, int n6, int n7, int n8, int n9,
41 int n10, int n11);
42 void AddResonance(int number, double mres, int n1, int n2, int n3, int n4,
43 int n5, int n6, int n7, int n8, int n9,
44 int n10, int n11, int n12);
45 //
46 // Total Energy Constraints
47 //
48 void AddTotalEnergy(int number, double etot, std::vector<int> lis);
49 void AddTotalEnergy(int number, double etot, int n1);
50 void AddTotalEnergy(int number, double etot, int n1, int n2);
51 void AddTotalEnergy(int number, double etot, int n1, int n2, int n3);
52 void AddTotalEnergy(int number, double etot, int n1, int n2, int n3, int n4);
53 void AddTotalEnergy(int number, double etot, int n1, int n2, int n3, int n4,
54 int n5);
55 void AddTotalEnergy(int number, double etot, int n1, int n2, int n3, int n4,
56 int n5, int n6);
57 void AddTotalEnergy(int number, double etot, int n1, int n2, int n3, int n4,
58 int n5, int n6, int n7);
59 void AddTotalEnergy(int number, double etot, int n1, int n2, int n3, int n4,
60 int n5, int n6, int n7, int n8);
61 void AddTotalEnergy(int number, double etot, int n1, int n2, int n3, int n4,
62 int n5, int n6, int n7, int n8, int n9);
63 void AddTotalEnergy(int number, double etot, int n1, int n2, int n3, int n4,
64 int n5, int n6, int n7, int n8, int n9, int n10);
65 void AddTotalEnergy(int number, double etot, int n1, int n2, int n3, int n4,
66 int n5, int n6, int n7, int n8, int n9, int n10, int n11);
67 void AddTotalEnergy(int number, double etot, int n1, int n2, int n3, int n4,
68 int n5, int n6, int n7, int n8, int n9,
69 int n10, int n11, int n12);
70
71 //
72 // Total Momentum Constraints
73 //
74 void AddTotalMomentum(int number, double ptot, std::vector<int> lis);
75 void AddTotalMomentum(int number, double ptot, int n1);
76 void AddTotalMomentum(int number, double ptot, int n1, int n2);
77 void AddTotalMomentum(int number, double ptot, int n1, int n2, int n3);
78 void AddTotalMomentum(int number, double ptot, int n1, int n2, int n3, int n4);
79 void AddTotalMomentum(int number, double ptot, int n1, int n2, int n3, int n4,
80 int n5);
81 void AddTotalMomentum(int number, double ptot, int n1, int n2, int n3, int n4,
82 int n5, int n6);
83 void AddTotalMomentum(int number, double ptot, int n1, int n2, int n3, int n4,
84 int n5, int n6, int n7);
85 void AddTotalMomentum(int number, double ptot, int n1, int n2, int n3, int n4,
86 int n5, int n6, int n7, int n8);
87 void AddTotalMomentum(int number, double ptot, int n1, int n2, int n3, int n4,
88 int n5, int n6, int n7, int n8, int n9);
89 void AddTotalMomentum(int number, double ptot, int n1, int n2, int n3, int n4,
90 int n5, int n6, int n7, int n8, int n9,
91 int n10);
92 void AddTotalMomentum(int number, double ptot, int n1, int n2, int n3, int n4,
93 int n5, int n6, int n7, int n8, int n9,
94 int n10, int n11);
95 void AddTotalMomentum(int number, double ptot, int n1, int n2, int n3, int n4,
96 int n5, int n6, int n7, int n8, int n9,
97 int n10, int n11, int n12);
98 //
99 // Three Momentum Constraints
100 //
101 void AddThreeMomentum(int number, Hep3Vector p3);
102 //
103 // Four Momentum Constraints
104 //
105 void AddFourMomentum(int number, HepLorentzVector p4);
106 void AddFourMomentum(int number, double etot);
107
108 //
109 // Equal Mass Constraints
110 //
111
112 void AddEqualMass(int number, std::vector<int> tlis1, std::vector<int> tlis2);
113
114 //
115 // Position Constraints
116 //
117 // void AddPosition(int number, HepPoint3D Xbeam, HepSymMatrix Vbeam,std::vector<int> tlis_V);
118
119 //
120 // Build Virtual Particle
121 //
122 void BuildVirtualParticle(int number);
123
124 // initialization called before fit
125
126 void init();
127 //
128 // set iteration number and chisq cut
129 //
130 void setFlag(const bool flag = 1) {m_flag = flag;}
131 void setIterNumber(const int niter = 5) {m_niter = niter;}
132 void setChisqCut(const double chicut = 200, const double chiter=0.05) {m_chicut = chicut;m_chiter=chiter;}
133 //
134 // set error matrix parameters
135 //
136 void setEspread (const double espread = 0.0009) {m_espread = espread;}
137 void setCollideangle (const double collideangle = 11e-3) {m_collideangle = collideangle;}
138 void setDynamicerror (const bool dynamicerror = 1) {m_dynamicerror = dynamicerror;}
139 void setTgraph ( TGraph2D* graph2d) {m_graph2d = graph2d;}
140 //
141 // fit function
142 //
143 bool Fit();
144 bool Fit(int n);
145 //
146 // Fit Results
147 //
148 // chisq of fit
149 //
150 double chisq() const {return m_chi;}
151 double chisq(int n) const {return m_chisq[n];}
152 // updated WTrack parameter after kinematic fit
153 // HepLorentzVector pfit(int n) const {return wTrackInfit(n).p();}
154 HepLorentzVector pfit(int n) const {return p4Infit(n);}
155 // added by yanl 05.10.07
156 // HepLorentzVector pfit1(int n) const {return wTrackOrigin(n).p();}
157 HepLorentzVector pfit1(int n) {return p4Origin(n);}
158 HepVector xfit() {return m_q.sub(1,3);}
159
160 WTrackParameter origin(int n) const {return wTrackOrigin(n);}
161 WTrackParameter infit(int n) const {return wTrackInfit(n);}
162
163 HepVector pull(int n) ;
164// WTrackParameter wVirtualTrack(int n) const {return m_virtual_wtrk[n];}
165 // error Matrix parameters
166 double espread() const {return m_espread;}
167 double collideangle() const {return m_collideangle;}
168 bool dynamicerror() const {return m_dynamicerror;}
169 // cpu time
170 HepVector cpu() const {return m_cpu;}
171 HepSymMatrix getCOrigin(int i) const;
172 HepSymMatrix getCInfit(int i) const;
173
174 WTrackParameter wVirtualTrack(int n) const {return m_virtual_wtrk[n];}
175
176private:
177 //bulid virtual particle
178 std::vector<WTrackParameter> m_virtual_wtrk;
179
180 private:
181 void updateConstraints(KinematicConstraints kc);
182 void fit();
183 // void covMatrix();
184 void upCovmtx();
185 void upTrkpar();
186 void clearDDT();
187 void fit(int n);
188 void covMatrix(int n);
189 void gda();
190 private:
191 std::vector<KinematicConstraints> m_kc;
192 std::vector<double> m_chisq;
193 double m_chi;
194 HepSymMatrix m_Vm;
195 HepMatrix m_A;
196 void setA(int ic, int itk, const HepMatrix &p) {m_A.sub(ic+1, itk, p);}
197 HepMatrix m_AT;
198 void setAT(int itk, int ic, const HepMatrix &p) { m_AT.sub(itk, ic+1, p);}
199 HepVector m_G;
200 HepSymMatrix m_W;
201 HepMatrix m_KP;
202
203 // virtual particle
204 HepMatrix m_B;
205 void setB(int ic, int itk, const HepMatrix &p) {m_B.sub(ic+1, itk, p);}
206 HepMatrix m_BT;
207 void setBT(int itk, int ic, const HepMatrix &p) { m_BT.sub(itk, ic+1, p);}
208
209
210 HepMatrix m_KQ;
211 int m_nc;
212 int m_nktrk;
213
214
215 HepVector m_p0;
216 HepVector m_p;
217 HepSymMatrix m_C0;
218 HepSymMatrix m_C;
219 HepVector pOrigin(int i) const ;
220 HepLorentzVector p4Origin(int i) const { HepVector p(4, 0); p = pOrigin(i); return HepLorentzVector(p[0], p[1], p[2], p[3]);}
221 HepVector pInfit(int i) const ;
222 HepLorentzVector p4Infit(int i) const { HepVector p(4, 0); p = pInfit(i); return HepLorentzVector(p[0], p[1], p[2], p[3]); }
223
224
225 void setPOrigin(int i, const HepVector &p) { m_p0.sub(i, p);}
226 void setPInfit(int i, const HepVector &p) {m_p.sub(i, p);}
227 void setCOrigin(int i, const HepSymMatrix &D) {m_C0.sub(i, D);}
228 void setCInfit(int i, const HepSymMatrix &D) {m_C.sub(i,D);}
229
230 HepVector m_q0;
231 HepVector m_q;
232 HepSymMatrix m_D0;
233 HepSymMatrix m_D;
234 HepSymMatrix m_D0inv;
235 HepSymMatrix m_Dinv;
236
237 void setQOrigin(int i, const HepVector &q) { m_q0.sub(i, q);}
238 void setQInfit(int i, const HepVector &q) {m_q.sub(i, q);}
239 void setDOrigin(int i, const HepSymMatrix &D) {m_D0.sub(i, D);}
240 void setDInfit(int i, const HepSymMatrix &D) {m_D.sub(i,D);}
241 void setDOriginInv(int i, const HepSymMatrix &Dinv) {m_D0inv.sub(i,Dinv);}
242
243
244
245 private:
247 static KalmanKinematicFit * m_pointer;
248 private:
249 int m_niter;
250 bool m_flag;
251 double m_chicut;
252 double m_chiter;
253 private:
254 double m_espread;
255 double m_collideangle;
256 private:
257 HepVector m_cpu;
258 private:
259 bool m_dynamicerror;
260 private:
261 static const int NTRKPAR;
262 static const int NKFPAR;
263 static const int Resonance;
264 static const int TotalEnergy;
265 static const int TotalMomentum;
266 static const int ThreeMomentum;
267 static const int FourMomentum;
268 static const int EqualMass;
269 static const int Position;
270 private:
271 TGraph2D* m_graph2d;
272
273
274};
275#endif
276
const Int_t n
Double_t etot
****INTEGER imax DOUBLE PRECISION m_pi *DOUBLE PRECISION m_amfin DOUBLE PRECISION m_Chfin DOUBLE PRECISION m_Xenph DOUBLE PRECISION m_sinw2 DOUBLE PRECISION m_GFermi DOUBLE PRECISION m_MfinMin DOUBLE PRECISION m_ta2 INTEGER m_out INTEGER m_KeyFSR INTEGER m_KeyQCD *COMMON c_Semalib $ !copy of input $ !CMS energy $ !beam mass $ !final mass $ !beam charge $ !final charge $ !smallest final mass $ !Z mass $ !Z width $ !EW mixing angle $ !Gmu Fermi $ alphaQED at q
Definition: KKsem.h:33
NTuple::Item< double > m_p
Definition: MdcHistItem.h:75
NTuple::Item< double > m_q
Definition: MdcHistItem.h:84
********INTEGER modcns REAL m_C
Definition: PseuMar.h:13
int n2
Definition: SD0Tag.cxx:55
int n1
Definition: SD0Tag.cxx:54
double espread() const
double chisq(int n) const
void AddTotalEnergy(int number, double etot, std::vector< int > lis)
void setIterNumber(const int niter=5)
void AddFourMomentum(int number, HepLorentzVector p4)
void setDynamicerror(const bool dynamicerror=1)
void setChisqCut(const double chicut=200, const double chiter=0.05)
WTrackParameter origin(int n) const
void AddThreeMomentum(int number, Hep3Vector p3)
void setCollideangle(const double collideangle=11e-3)
void setTgraph(TGraph2D *graph2d)
HepLorentzVector pfit(int n) const
WTrackParameter infit(int n) const
void BuildVirtualParticle(int number)
HepVector cpu() const
HepLorentzVector pfit1(int n)
WTrackParameter wVirtualTrack(int n) const
HepSymMatrix getCInfit(int i) const
void AddResonance(int number, double mres, std::vector< int > tlis)
void setEspread(const double espread=0.0009)
void AddTotalMomentum(int number, double ptot, std::vector< int > lis)
HepSymMatrix getCOrigin(int i) const
static KalmanKinematicFit * instance()
double collideangle() const
void setFlag(const bool flag=1)
void AddEqualMass(int number, std::vector< int > tlis1, std::vector< int > tlis2)
std::vector< WTrackParameter > wTrackInfit() const
Definition: TrackPool.h:73
std::vector< WTrackParameter > wTrackOrigin() const
Definition: TrackPool.h:72