35{
36
37 int nGood=iGood.size();
38 int nGam=iGam.size();
39
40 iGoodtag.clear();
41 iGamtag.clear();
42
43 double mass_bcgg, delE_tag_temp;
44 int m_chargetag, m_chargek1,m_chargek2,m_chargepi1,m_chargepi2;
45 int ik1_temp, ik2_temp, ipi1_temp, ipi2_temp;
46 HepLorentzVector pddd;
47 HepLorentzVector pddd_temp;
48
49 IDataProviderSvc* eventSvc = NULL;
50 Gaudi::svcLocator()->service("EventDataSvc", eventSvc);
52 SmartDataPtr<Event::EventHeader> eventHeader(eventSvc,"/Event/EventHeader");
53
54 int runNo=eventHeader->runNumber();
55 int rec=eventHeader->eventNumber();
56
57 double xecm=2*Ebeam;
58
59 kkpipimd=false;
60 double tagmode=0;
61
62 if((evtRecEvent->totalCharged() < 4)){ return; }
63
65
67 Gaudi::svcLocator()->service("SimplePIDSvc", simple_pid);
68
69 double deltaE_tem = 0.20;
70 int ncount1 = 0;
71
72 Hep3Vector xorigin(0,0,0);
74 Gaudi::svcLocator()->service("VertexDbSvc", vtxsvc);
76 {
79 xorigin.setX(dbv[0]);
80 xorigin.setY(dbv[1]);
81 xorigin.setZ(dbv[2]);
82 }
83
84 double xv=xorigin.x();
85 double yv=xorigin.y();
86 double zv=xorigin.z();
87
89 HepPoint3D IP(xorigin[0],xorigin[1],xorigin[2]);
90
91 for(int i = 0; i < evtRecEvent->totalCharged(); i++) {
93
94 int ik1 = (*itTrk)->trackId();
95
96 if(!(*itTrk)->isMdcKalTrackValid()) continue;
99
100 m_chargek1 = mdcKalTrk1->
charge();
101 if(m_chargek1 != 1) continue;
102
103
106
107 VFHelix helixip3_1(point0,a1,Ea1);
108 helixip3_1.pivot(IP);
109 HepVector vecipa1 = helixip3_1.a();
110
111 double dr1 = fabs(vecipa1[0]);
112 double dz1 = fabs(vecipa1[3]);
113 double costheta1 =
cos(mdcKalTrk1->
theta());
114
115 if ( dr1 >= 1.0) continue;
116 if ( dz1 >= 10.0) continue;
117 if ( fabs(costheta1) >= 0.93) continue;
118
119 if(PID_flag == 5) {
122 }
123
124
126
127
128
129
130 for(int j = 0; j< evtRecEvent->totalCharged();j++) {
132
133 int ik2 = (*itTrk)->trackId();
134 if(ik1==ik2) continue;
135
136 if(!(*itTrk)->isMdcKalTrackValid()) continue;
139
140 m_chargek2=mdcKalTrk2->
charge();
141 if((m_chargek1 + m_chargek2) != 0) continue;
142
143
146 VFHelix helixip3_2(point0,a2,Ea2);
147 helixip3_2.pivot(IP);
148 HepVector vecipa2 = helixip3_2.a();
149
150 double dr2 = fabs(vecipa2[0]);
151 double dz2 = fabs(vecipa2[3]);
152 double costheta2 =
cos(mdcKalTrk2->
theta());
153 if ( dr2 >= 1.0) continue;
154 if ( dz2 >= 10.0) continue;
155 if ( fabs(costheta2) >= 0.93) continue;
156
157 if(PID_flag == 5) {
160 }
161
162
164
165
166
167
168 for(int k = 0; k< evtRecEvent->totalCharged(); k++) {
170
171 int ipi1= (*itTrk)->trackId();
172 if(ipi1==ik1 || ipi1==ik2) continue;
173
174 if(!(*itTrk)->isMdcKalTrackValid()) continue;
177
178 m_chargepi1=mdcKalTrk3->
charge();
179 if(m_chargepi1 != 1) continue;
180
181
183 HepSymMatrix Ea3 = mdcKalTrk3->
getZError();
184 VFHelix helixip3_3(point0,a3,Ea3);
185 helixip3_3.pivot(IP);
186 HepVector vecipa3 = helixip3_3.a();
187
188 double dr3 = fabs(vecipa3[0]);
189 double dz3 = fabs(vecipa3[3]);
190 double costheta3 =
cos(mdcKalTrk3->
theta());
191 if ( dr3 >= 1.0) continue;
192 if ( dz3 >= 10.0) continue;
193 if ( fabs(costheta3) >= 0.93) continue;
194
195 if(PID_flag == 5) {
198 }
199
201
202
203
204
205 for(int l = 0; l< evtRecEvent->totalCharged();l++) {
207
208 int ipi2= (*itTrk)->trackId();
209 if(ipi2==ik1 || ipi2==ik2 || ipi2==ipi1 ) continue;
210
211 if(!(*itTrk)->isMdcKalTrackValid()) continue;
214
215 m_chargepi2 = mdcKalTrk4->
charge();
216 if((m_chargepi1 + m_chargepi2) != 0) continue;
217
218
220 HepSymMatrix Ea4 = mdcKalTrk4->
getZError();
221 VFHelix helixip3_4(point0,a4,Ea4);
222 helixip3_4.pivot(IP);
223 HepVector vecipa4 = helixip3_4.a();
224
225 double dr4 = fabs(vecipa4[0]);
226 double dz4 = fabs(vecipa4[3]);
227 double costheta4 =
cos(mdcKalTrk4->
theta());
228 if ( dr4 >= 1.0) continue;
229 if ( dz4 >= 10.0) continue;
230 if ( fabs(costheta4) >= 0.93) continue;
231
232 if(PID_flag == 5) {
235 }
236
237
239
240
241 HepPoint3D vx(xorigin.x(), xorigin.y(), xorigin.z());
242 HepSymMatrix Evx(3, 0);
243 double bx = 1E+6; Evx[0][0] = bx*bx;
244 double by = 1E+6; Evx[1][1] = by*by;
245 double bz = 1E+6; Evx[2][2] = bz*bz;
247
248
256 if(!vtxfit->
Fit(0))
continue;
258
263
264 HepVector kap_val = HepVector(7,0);
265 HepVector kam_val = HepVector(7,0);
266 HepVector pip_val = HepVector(7,0);
267 HepVector pim_val = HepVector(7,0);
272
273 HepLorentzVector P_KAP(kap_val[0],kap_val[1],kap_val[2],kap_val[3]);
274 HepLorentzVector P_KAM(kam_val[0],kam_val[1],kam_val[2],kam_val[3]);
275 HepLorentzVector P_PIP(pip_val[0],pip_val[1],pip_val[2],pip_val[3]);
276 HepLorentzVector P_PIM(pim_val[0],pim_val[1],pim_val[2],pim_val[3]);
277
278 P_KAP.boost(-0.011,0,0);
279 P_KAM.boost(-0.011,0,0);
280 P_PIP.boost(-0.011,0,0);
281 P_PIM.boost(-0.011,0,0);
282 pddd = P_KAP + P_KAM + P_PIP + P_PIM;
283
284 double pkkpipi=pddd.rho();
285
286 double temp1 = (
ecms/2)*(
ecms/2)-pkkpipi*pkkpipi ;
287 if(temp1<0) temp1 =0;
288 double mass_bc_tem = sqrt(temp1);
289 if(mass_bc_tem < 1.82 || mass_bc_tem > 1.89) continue;
290
291 double delE_tag_tag =
ecms/2-pddd.e();
292
293 if(fabs(delE_tag_tag)<deltaE_tem) {
294 deltaE_tem = fabs(delE_tag_tag);
295 delE_tag_temp = delE_tag_tag;
296 mass_bcgg = mass_bc_tem;
297
298 pddd_temp = pddd;
299
300 ik1_temp=ik1;
301 ik2_temp=ik2;
302 ipi1_temp=ipi1;
303 ipi2_temp=ipi2;
304 ncount1 = 1;
305
306 }
307 }
308 }
309 }
310 }
311
312 if(ncount1 == 1){
313 tagmode=20;
314 if(m_chargetag<0) tagmode=-20;
315 tagmd=tagmode;
316 mass_bc = mass_bcgg;
317 delE_tag = delE_tag_temp;
318 cqtm = 0.0;
319
320 iGoodtag.push_back(ik1_temp);
321 iGoodtag.push_back(ik2_temp);
322 iGoodtag.push_back(ipi1_temp);
323 iGoodtag.push_back(ipi2_temp);
324
325 iGamtag.push_back(9999);
326 iGamtag.push_back(9999);
327 iGamtag.push_back(9999);
328 iGamtag.push_back(9999);
329
330 ptag = pddd_temp;
331
332 kkpipimd=true;
333
334 }
335}
EvtRecTrackCol::iterator EvtRecTrackIterator
double cos(const BesAngle a)
const double theta() const
static void setPidType(PidType pidType)
virtual double probKaon()=0
virtual void preparePID(EvtRecTrack *track)=0
virtual double probPion()=0
virtual bool isVertexValid()=0
virtual double * SigmaPrimaryVertex()=0
virtual double * PrimaryVertex()=0
const HepVector & getZHelix() const
const HepSymMatrix & getZError() const
HepSymMatrix & getZErrorK()
void AddTrack(const int number, const double mass, const RecMdcTrack *trk)
WTrackParameter wtrk(int n) const
void AddVertex(int number, VertexParameter vpar, std::vector< int > lis)
static VertexFit * instance()
void setEvx(const HepSymMatrix &eVx)
void setVx(const HepPoint3D &vx)
_EXTERN_ std::string EvtRecEvent