116 {
117 StatusCode sc = StatusCode::SUCCESS;
118
119 MsgStream log(
msgSvc(), name() );
120 log << MSG::INFO << "in execute()" << endmsg;
121
122
123
124
125 setFilterPassed( false );
126
127 m_pass[0] += 1;
128
129 SmartDataPtr<Event::EventHeader> eventHeader( eventSvc(), "/Event/EventHeader" );
132
133 Vint iGood, ikm, ikp, ipip, ipim, iGam;
134 iGood.clear();
135 ikm.clear();
136 ikp.clear();
137 ipip.clear();
138 ipim.clear();
139 iGam.clear();
140
141 Vp4 ppip, ppim, pkm, pkp;
142 ppip.clear();
143 ppim.clear();
144 pkm.clear();
145 pkp.clear();
146
147 int TotCharge = 0;
148 for ( int i = 0; i < evtRecEvent->totalCharged(); i++ )
149 {
151 if ( !( *itTrk )->isMdcTrackValid() ) continue;
152 RecMdcTrack* mdcTrk = ( *itTrk )->mdcTrack();
153 if ( fabs( mdcTrk->
z() ) >= m_vz0cut )
continue;
154 if ( mdcTrk->
r() >= m_vr0cut )
continue;
155 iGood.push_back( i );
156 TotCharge += mdcTrk->
charge();
157 }
158
159
160
161 int nGood = iGood.size();
162
163
164
165
166
167 if ( ( nGood < 2 ) || ( TotCharge != 0 ) ) return sc;
168
169 m_pass[1] += 1;
170
171
172
173
175 for ( int i = 0; i < nGood; i++ )
176 {
178
185
186
189
190 if ( !( *itTrk )->isMdcKalTrackValid() ) continue;
191 RecMdcKalTrack* mdcKalTrk = ( *itTrk )->mdcKalTrack();
193 {
195 HepLorentzVector ptrk;
196 ptrk.setPx( mdcKalTrk->
px() );
197 ptrk.setPy( mdcKalTrk->
py() );
198 ptrk.setPz( mdcKalTrk->
pz() );
199 double p3 = ptrk.mag();
200 ptrk.setE( sqrt( p3 * p3 +
mk *
mk ) );
201 if ( mdcKalTrk->
charge() > 0 )
202 {
203 ikp.push_back( iGood[i] );
204 pkp.push_back( ptrk );
205 }
206 else
207 {
208 ikm.push_back( iGood[i] );
209 pkm.push_back( ptrk );
210 }
211 }
212
214 {
216 HepLorentzVector ptrk;
217 ptrk.setPx( mdcKalTrk->
px() );
218 ptrk.setPy( mdcKalTrk->
py() );
219 ptrk.setPz( mdcKalTrk->
pz() );
220 double p3 = ptrk.mag();
221 ptrk.setE( sqrt( p3 * p3 +
mpi *
mpi ) );
222 if ( mdcKalTrk->
charge() > 0 )
223 {
224 ipip.push_back( iGood[i] );
225 ppip.push_back( ptrk );
226 }
227 else
228 {
229 ipim.push_back( iGood[i] );
230 ppim.push_back( ptrk );
231 }
232 }
233 }
234 int npip = ipip.size();
235 int npim = ipim.size();
236 int nkm = ikm.size();
237 int nkp = ikp.size();
238
239 m_nkm = nkm;
240 m_nkp = nkp;
241 m_npip = npip;
242 m_npim = npim;
243 m_ncharge = nGood;
244 m_pass[3] += 1;
245
246 if ( npip < 1 && npim < 1 ) return sc;
247 if ( nkp < 1 && nkm < 1 ) return sc;
248
249 m_pass[4] += 1;
250
251
252
253
254
255 HepLorentzVector pkstar, pkstar1, ppi1, ppi2, pTot;
257 ikstar.clear();
258
259 double difchi0 = 99999.0;
260 int ixpim = -1;
261 int ixkp = -1;
262
263 for ( int i = 0; i < npim; i++ )
264 {
265 for ( int j = 0; j < nkp; j++ )
266 {
267 pkstar = ppim[i] + pkp[j];
268 double difchi = fabs( pkstar.m() -
mkstar0 );
269 if ( difchi < difchi0 )
270 {
271 difchi0 = difchi;
272 ixpim = i;
273 ixkp = j;
274 }
275 }
276 }
277
278 m_difchikp = difchi0;
279
280 if ( ixpim != -1 ) m_kstarkp = ( ppim[ixpim] + pkp[ixkp] ).m();
281
282 double difchi1 = 99999.0;
283 int ixpip = -1;
284 int ixkm = -1;
285
286 for ( int ii = 0; ii < npip; ii++ )
287 {
288 for ( int jj = 0; jj < nkm; jj++ )
289 {
290 pkstar1 = ppip[ii] + pkm[jj];
291 double difchi2 = fabs( pkstar1.m() -
mkstar0 );
292 if ( difchi2 < difchi1 )
293 {
294 difchi1 = difchi2;
295 ixpip = ii;
296 ixkm = jj;
297 }
298 }
299 }
300
301 m_difchikm = difchi1;
302
303 if ( ixpip != -1 ) m_kstarkm = ( ppip[ixpip] + pkm[ixkm] ).m();
304
305
306 if ( ixpip == -1 && ixpim == -1 ) return sc;
307
308 if ( difchi0 < difchi1 ) { pTot = ppim[ixpim] + pkp[ixkp]; }
309 else { pTot = ppip[ixpip] + pkm[ixkm]; }
310 m_mkstar = pTot.m();
311 m_pkstar = pTot.rho();
312
313
314 TH1* h( 0 );
315 if ( m_thsvc->getHist( "/DQAHist/InclKstar/InclKstar_mass", h ).isSuccess() )
316 { h->Fill( m_mkstar ); }
317 else { log << MSG::ERROR << "Couldn't retrieve InclKstar_mass" << endmsg; }
318
319 m_tuple2->write().ignore();
320
321
322
323 if ( ixpim != -1 )
324 {
325
326
327 ( *( evtRecTrkCol->begin() + ipim[ixpim] ) )->tagPion();
328 ( *( evtRecTrkCol->begin() + ikp[ixkp] ) )->tagKaon();
329 ( *( evtRecTrkCol->begin() + ipim[ixpim] ) )->setQuality( 3 );
330 ( *( evtRecTrkCol->begin() + ikp[ixkp] ) )->setQuality( 3 );
331 }
332 else
333 {
334
335
336 ( *( evtRecTrkCol->begin() + ipip[ixpip] ) )->tagPion();
337 ( *( evtRecTrkCol->begin() + ikm[ixkm] ) )->tagKaon();
338 ( *( evtRecTrkCol->begin() + ipip[ixpip] ) )->setQuality( 3 );
339 ( *( evtRecTrkCol->begin() + ikm[ixkm] ) )->setQuality( 3 );
340 }
341
342
343
344 setFilterPassed( true );
345
346 return StatusCode::SUCCESS;
347}
EvtRecTrackCol::iterator EvtRecTrackIterator
std::vector< HepLorentzVector > Vp4
static void setPidType(PidType pidType)
int methodProbability() const
void setChiMinCut(const double chi=4)
void setRecTrack(EvtRecTrack *trk)
void setMethod(const int method)
void identify(const int pidcase)
void usePidSys(const int pidsys)
static ParticleID * instance()
bool IsPidInfoValid() const
_EXTERN_ std::string EvtRecEvent
_EXTERN_ std::string EvtRecTrackCol