BOSS 7.0.5
BESIII Offline Software System
Loading...
Searching...
No Matches
Public Member Functions | List of all members
TagSetAlg Class Reference

#include <TagSetAlg.h>

+ Inheritance diagram for TagSetAlg:

Public Member Functions

 TagSetAlg (const std::string &name, ISvcLocator *pSvcLocator)
 
StatusCode initialize ()
 
StatusCode execute ()
 
StatusCode finalize ()
 
unsigned int IntToTag1 (int val1, int val2, int val3, int val4)
 
unsigned int IntToTag2 (int val1, int val2, int val3, int val4, int val5, int val6)
 
 TagSetAlg (const std::string &name, ISvcLocator *pSvcLocator)
 
StatusCode initialize ()
 
StatusCode execute ()
 
StatusCode finalize ()
 
unsigned int IntToTag1 (int val1, int val2, int val3, int val4)
 
unsigned int IntToTag2 (int val1, int val2, int val3, int val4, int val5, int val6)
 

Detailed Description

Definition at line 12 of file Event/TagSetAlg/TagSetAlg-00-00-11/TagSetAlg/TagSetAlg.h.

Constructor & Destructor Documentation

◆ TagSetAlg() [1/2]

TagSetAlg::TagSetAlg ( const std::string &  name,
ISvcLocator *  pSvcLocator 
)

Definition at line 47 of file TagSetAlg.cxx.

47 :
48 Algorithm(name, pSvcLocator) {
49
50}

◆ TagSetAlg() [2/2]

TagSetAlg::TagSetAlg ( const std::string &  name,
ISvcLocator *  pSvcLocator 
)

Member Function Documentation

◆ execute() [1/2]

StatusCode TagSetAlg::execute ( )

Definition at line 65 of file TagSetAlg.cxx.

65 {
66
67 int nGoodCharged = 0;
68 int nGoodChargedp = 0;
69 int nGoodChargedm = 0;
70 int nCharged = 0;
71 int nNeutrk = 0;
72 //int nGoodNeutrk = 0;
73 int nTottrk = 0;
74 int totCharged = 0;
75 int npionp = 0;
76 int npionm = 0;
77 int nprotonp = 0;
78 int nprotonm = 0;
79 int nkaonp = 0;
80 int nkaonm = 0;
81 int nlambda = 0;
82 int nalambda= 0;
83 int nks = 0;
84 int ngamma = 0;
85 int neta = 0;
86 int npi0 = 0;
87 int nmuonp = 0;
88 int nmuonm = 0;
89 int nelectronp = 0;
90 int nelectronm = 0;
91
92
93 Hep3Vector xorigin(0,0,0);
94 if(m_vtxsvc->isVertexValid()){
95 double* dbv = m_vtxsvc->PrimaryVertex();
96 xorigin.setX(dbv[0]);
97 xorigin.setY(dbv[1]);
98 xorigin.setZ(dbv[2]);
99 }
100
101 SmartDataPtr<EvtRecEvent> evtRecEvent(eventSvc(), EventModel::EvtRec::EvtRecEvent);
102 SmartDataPtr<EvtRecTrackCol> evtRecTrkCol(eventSvc(), EventModel::EvtRec::EvtRecTrackCol);
103
104 nCharged = evtRecEvent->totalCharged();
105 nNeutrk = evtRecEvent->totalNeutral();
106 nTottrk = evtRecEvent->totalTracks();
107
108 //std::cout <<"iEvt: "<<iEvt <<" nCharged: "<<nCharged<<std::endl;
109 //cout<<"nNeutrk= "<<nNeutrk<<", nTottrk = "<<nTottrk<<endl;
110
111 // Good charged track selection
112 vector<int> iGood;
113 iGood.clear();
114 for ( Int_t iCharge = 0; iCharge < evtRecEvent->totalCharged(); ++iCharge ) {
115 EvtRecTrackIterator itTrk=evtRecTrkCol->begin() + iCharge;
116 if(!(*itTrk)->isMdcTrackValid()) continue;
117 RecMdcTrack *mdcTrk = (*itTrk)->mdcTrack();
118 double theta = mdcTrk->theta();
119 HepVector a = mdcTrk->helix();
120 HepSymMatrix Ea = mdcTrk->err();
121 HepPoint3D point0(0.,0.,0.); // the initial point for MDC reconstruction
122 HepPoint3D IP(xorigin[0],xorigin[1],xorigin[2]);
123 VFHelix helixip(point0,a,Ea);
124 helixip.pivot(IP);
125 HepVector vecipa = helixip.a();
126 double Rvxy0=fabs(vecipa[0]); //the nearest distance to IP in xy plane
127 double Rvz0=vecipa[3]; //the nearest distance to IP in z direction
128 double Rvphi0=vecipa[1];
129
130 //std::cout<<"Rvxy0: "<<Rvxy0<<" Rvz0: "<<Rvz0<<std::endl;
131 if(fabs(Rvxy0) >= 1.0) continue;
132 if(fabs(Rvz0) >= 10.0) continue;
133 if(fabs(cos(theta))>=0.93) continue;
134 if(mdcTrk->charge() > 0) nGoodChargedp++;
135 if(mdcTrk->charge() < 0) nGoodChargedm++;
136
137 ++nGoodCharged;
138 totCharged +=mdcTrk->charge();
139 iGood.push_back(iCharge);
140 }
141
142 //the DST data will be reformed according to the ' nGoodCharged '
143 m_tagFilterSvc->setTagData0(nGoodCharged);
144 //std::cout <<"nGoodCharged: "<<nGoodCharged<<", nGoodChargedp= "<<nGoodChargedp<<std::endl;
145 //cout<<"nGoodChargedm = "<<nGoodChargedm<<", totCharged = "<<totCharged<<endl;
146
147
148
149 for(int i = evtRecEvent->totalCharged(); i< evtRecEvent->totalTracks(); i++)
150 {
151 EvtRecTrackIterator iTrk=evtRecTrkCol->begin() + i;
152
153 //get photons
154 if(!(*iTrk)->isEmcShowerValid()) continue;
155 RecEmcShower* emcShower = (*iTrk)->emcShower();
156 //if((t_runNo == 10650 && t_evtNo == 2086135)) cout<<emcShower->energy()<<" "<<emcShower->time()<<" "<<abs(cos(emcShower->theta()))<<" "<<emcShower->energy()<<endl;
157 if (((emcShower->energy() > 0.025) &&
158 (emcShower->time() > 0) && (emcShower->time() < 14) &&
159 (((abs(cos(emcShower->theta())) < 0.80) && (emcShower->energy() > 0.025)) ||
160 ((abs(cos(emcShower->theta())) < 0.92) &&
161 (abs(cos(emcShower->theta())) > 0.86) && (emcShower->energy() > 0.050))))
162 ){
163 //|| m_neutralStudies) //xiongxa
164 ngamma++;
165 }
166 }
167
168
169 //cout<<"Ngamma = "<<ngamma<<endl;
170 // get tracks
171 for(int i=0;i<iGood.size();i++){
172 EvtRecTrackIterator iTrk = evtRecTrkCol->begin() + iGood[i];
173 RecMdcKalTrack* mdcKalTrack = (*iTrk)->mdcKalTrack();
174
175
176 // set up pid (but only make very rough pid cuts)
177
178 ParticleID* pid = ParticleID::instance();
179 pid->init();
180 pid->setMethod(pid->methodProbability());
181 pid->setRecTrack(*iTrk);
182 pid->usePidSys(pid->useDedx() |
183 pid->useTof1() |
184 pid->useTof2());
185 pid->identify(pid->onlyPion() |
186 pid->onlyKaon() |
187 pid->onlyProton() |
188 pid->onlyElectron() |
189 pid->onlyMuon());
190 pid->calculate();
191 if(!(pid->IsPidInfoValid())) continue;
192
193 // get pions
194
195 // if ((pid->IsPidInfoValid() && pid->probPion() > 1.0e-5)) //xiongxa
196 if ((pid->probPion() > pid->probKaon()) && (pid->probPion() > pid->probProton())){
197 if (mdcKalTrack->charge() > 0) npionp++;
198 if (mdcKalTrack->charge() < 0) npionm++;
199 }
200
201 // if ((pid->IsPidInfoValid() && pid->probKaon() > 1.0e-5)) //xiongxa
202 if ((pid->probKaon() > pid->probPion()) && (pid->probKaon() > pid->probProton())){
203 if (mdcKalTrack->charge() > 0) nkaonp++;
204 if (mdcKalTrack->charge() < 0) nkaonm++;
205 }
206
207 // if ((pid->IsPidInfoValid() && pid->probProton() > 1.0e-5)) //xiongxa
208 if ((pid->probProton() > pid->probPion()) && (pid->probProton() > pid->probKaon())){
209 if (mdcKalTrack->charge() > 0) nprotonp++;
210 if (mdcKalTrack->charge() < 0) nprotonm++;
211 }
212
213 //if ((pid->IsPidInfoValid() && pid->probElectron() > 1.0e-5) || m_trackStudies)
214 if (pid->probElectron() > 1.0e-5){ //xiongxa
215 if (mdcKalTrack->charge() > 0) nelectronp++;
216 if (mdcKalTrack->charge() < 0) nelectronm++;
217 }
218
219 //if ((pid->IsPidInfoValid() && pid->probMuon() > 1.0e-5) || m_trackStudies)
220 if (pid->probMuon() > 1.0e-5){ //xiongxa
221 if (mdcKalTrack->charge() > 0) nmuonp++;
222 if (mdcKalTrack->charge() < 0) nmuonm++;
223 }
224 }
225 // cout<<"Npionp = "<<npionp<<", npionm = "<<npionm<<endl;
226 // cout<<"Nkaonp = "<<nkaonp<<", nkaonm = "<<nkaonm<<endl;
227 // cout<<"Nprotonp = "<<nprotonp<<", nprotonm= "<<nprotonm<<endl;
228 // cout<<"nelectronp = "<<nelectronp<<", nelectronm = "<<nelectronm<<endl;
229 // cout<<"nmuonp = "<<nmuonp<<", nmuonm = "<<nmuonm<<endl;
230
231
232 // get etas
233 // eta --> gamma gamma
234
235 SmartDataPtr<EvtRecEtaToGGCol> evtRecEtaToGGCol(eventSvc(), EventModel::EvtRec::EvtRecEtaToGGCol);
236 for(EvtRecEtaToGGCol::iterator iEta = evtRecEtaToGGCol->begin();
237 iEta != evtRecEtaToGGCol->end(); iEta++){
238 if ((((*iEta)->chisq() < 2500) &&
239 ((*iEta)->unconMass() > 0.40) &&
240 //((*iEta)->unconMass() < 0.70)) || m_etaStudies){
241 ((*iEta)->unconMass() < 0.70))){ //xiongxa
242 EvtRecTrack* lo = const_cast<EvtRecTrack*>((*iEta)->loEnGamma());
243 RecEmcShower* loShower = lo->emcShower();
244 if (((loShower->energy() > 0.025) &&
245 (loShower->time() > 0) && (loShower->time() < 14) &&
246 (((abs(cos(loShower->theta())) < 0.80) && (loShower->energy() > 0.025)) ||
247 ((abs(cos(loShower->theta())) < 0.92) &&
248 (abs(cos(loShower->theta())) > 0.86) && (loShower->energy() > 0.050))))
249 //|| m_neutralStudies){
250 ){ //xiongxa
251 EvtRecTrack* hi = const_cast<EvtRecTrack*>((*iEta)->hiEnGamma());
252 RecEmcShower* hiShower = hi->emcShower();
253 if (((hiShower->energy() > 0.025) &&
254 (hiShower->time() > 0) && (hiShower->time() < 14) &&
255 (((abs(cos(hiShower->theta())) < 0.80) && (hiShower->energy() > 0.025)) ||
256 ((abs(cos(hiShower->theta())) < 0.92) &&
257 (abs(cos(hiShower->theta())) > 0.86) && (hiShower->energy() > 0.050))))
258 //|| m_neutralStudies){
259 ){ //xiongxa
260 neta++;
261 }
262 }
263 }
264 }
265
266 // get pi0s
267 // pi0 --> gamma gamma
268
269 SmartDataPtr<EvtRecPi0Col> evtRecPi0Col(eventSvc(), EventModel::EvtRec::EvtRecPi0Col);
270 for(EvtRecPi0Col::iterator iPi0 = evtRecPi0Col->begin();
271 iPi0 != evtRecPi0Col->end(); iPi0++){
272 if ((((*iPi0)->chisq() < 2500) &&
273 ((*iPi0)->unconMass() > 0.107) &&
274 //((*iPi0)->unconMass() < 0.163)) || m_pi0Studies){
275 ((*iPi0)->unconMass() < 0.163))){ //xiongxa
276 EvtRecTrack* lo = const_cast<EvtRecTrack*>((*iPi0)->loEnGamma());
277 RecEmcShower* loShower = lo->emcShower();
278 if (((loShower->energy() > 0.025) &&
279 (loShower->time() > 0) && (loShower->time() < 14) &&
280 (((abs(cos(loShower->theta())) < 0.80) && (loShower->energy() > 0.025)) ||
281 ((abs(cos(loShower->theta())) < 0.92) &&
282 (abs(cos(loShower->theta())) > 0.86) && (loShower->energy() > 0.050))))
283 //|| m_neutralStudies){
284 ){ //xiongxa
285 EvtRecTrack* hi = const_cast<EvtRecTrack*>((*iPi0)->hiEnGamma());
286 RecEmcShower* hiShower = hi->emcShower();
287 if (((hiShower->energy() > 0.025) &&
288 (hiShower->time() > 0) && (hiShower->time() < 14) &&
289 (((abs(cos(hiShower->theta())) < 0.80) && (hiShower->energy() > 0.025)) ||
290 ((abs(cos(hiShower->theta())) < 0.92) &&
291 (abs(cos(hiShower->theta())) > 0.86) && (hiShower->energy() > 0.050))))
292 //|| m_neutralStudies){
293 ){ //xiongxa
294 npi0++;
295 }
296 }
297 }
298 }
299 // get kshorts
300 // [[ vertexId is the pdgID ]]
301 // Ks --> pi+ pi- and Lambda --> p+ pi- and ALambda --> p- pi+
302
303 SmartDataPtr<EvtRecVeeVertexCol> evtRecVeeVertexCol(eventSvc(), EventModel::EvtRec::EvtRecVeeVertexCol);
304 for(EvtRecVeeVertexCol::iterator iKs = evtRecVeeVertexCol->begin();
305 iKs != evtRecVeeVertexCol->end(); iKs++){
306 if ((*iKs)->vertexId() == 310){
307 if ((((*iKs)->mass() > 0.471) &&
308 ((*iKs)->mass() < 0.524) &&
309 //((*iKs)->chi2() < 100)) || m_veeStudies){
310 ((*iKs)->chi2() < 100))){ //xiongxa
311 nks++;
312 }
313 }
314 }
315
316 for(EvtRecVeeVertexCol::iterator iL = evtRecVeeVertexCol->begin();
317 iL != evtRecVeeVertexCol->end(); iL++){
318 if ((*iL)->vertexId() == 3122){
319 if ((((*iL)->mass() > 1.100) &&
320 ((*iL)->mass() < 1.130) &&
321 //((*iL)->chi2() < 100)) || m_veeStudies){
322 ((*iL)->chi2() < 100))){ //xiongxa
323 nlambda++;
324 }
325 }
326 if ((*iL)->vertexId() == -3122){
327 if ((((*iL)->mass() > 1.100) &&
328 ((*iL)->mass() < 1.130) &&
329 //((*iL)->chi2() < 100)) || m_veeStudies){
330 ((*iL)->chi2() < 100))){ //xiongxa
331 nalambda++;
332 }
333 }
334 }
335
336
337 //please set the tag data here
338 unsigned int tagdata1 = IntToTag1(nNeutrk, nTottrk, ngamma, npi0);
339 unsigned int tagdata2 = IntToTag2(npionp, npionm, nkaonp, nkaonm, nprotonp, nprotonm);
340 unsigned int tagdata3 = IntToTag2(nlambda, nalambda, nelectronp, nelectronm, nmuonp, nmuonm);
341 unsigned int tagdata4 = IntToTag2(nks, neta, nCharged, nGoodChargedp, nGoodChargedm, totCharged);
342
343
344 m_tagFilterSvc->setTagData1(tagdata1);
345 m_tagFilterSvc->setTagData2(tagdata2);
346 m_tagFilterSvc->setTagData3(tagdata3);
347 m_tagFilterSvc->setTagData4(tagdata4);
348
349 iEvt++;
350}
cos
double cos(const BesAngle a)
Definition: InstallArea/include/MdcGeom/MdcGeom/BesAngle.h:213
DstMdcTrack::err
const HepSymMatrix err() const
Definition: DstMdcTrack.cxx:155
DstMdcTrack::helix
const HepVector helix() const
......
Definition: DstMdcTrack.cxx:151
ITagFilterSvc::setTagData1
virtual void setTagData1(unsigned int t)=0
ITagFilterSvc::setTagData0
virtual void setTagData0(unsigned int t)=0
ITagFilterSvc::setTagData2
virtual void setTagData2(unsigned int t)=0
ITagFilterSvc::setTagData4
virtual void setTagData4(unsigned int t)=0
ITagFilterSvc::setTagData3
virtual void setTagData3(unsigned int t)=0
IVertexDbSvc::isVertexValid
virtual bool isVertexValid()=0
IVertexDbSvc::PrimaryVertex
virtual double * PrimaryVertex()=0
ParticleID::instance
static ParticleID * instance()
Definition: ParticleID.cxx:22
ParticleID::IsPidInfoValid
bool IsPidInfoValid() const
ParticleID::calculate
void calculate()
Definition: ParticleID.cxx:97
ParticleID::init
void init()
Definition: ParticleID.cxx:27
TagSetAlg::IntToTag2
unsigned int IntToTag2(int val1, int val2, int val3, int val4, int val5, int val6)
Definition: TagSetAlg.cxx:384
TagSetAlg::IntToTag1
unsigned int IntToTag1(int val1, int val2, int val3, int val4)
Definition: TagSetAlg.cxx:375

◆ execute() [2/2]

StatusCode TagSetAlg::execute ( )

◆ finalize() [1/2]

StatusCode TagSetAlg::finalize ( )

Definition at line 355 of file TagSetAlg.cxx.

355 {
356
357 MsgStream log(msgSvc(), name());
358 log << MSG::INFO << "in finalize()" << iEvt << "Events processed" << endmsg;
359 return StatusCode::SUCCESS;
360}

◆ finalize() [2/2]

StatusCode TagSetAlg::finalize ( )

◆ initialize() [1/2]

StatusCode TagSetAlg::initialize ( )

Definition at line 53 of file TagSetAlg.cxx.

53 {
54 MsgStream log(msgSvc(), name());
55
56 log << MSG::INFO << "in initialize()" << endmsg;
57 iEvt=0;
58
59 Gaudi::svcLocator()->service("VertexDbSvc", m_vtxsvc);
60 Gaudi::svcLocator()->service("TagFilterSvc", m_tagFilterSvc);
61
62 return StatusCode::SUCCESS;
63}

◆ initialize() [2/2]

StatusCode TagSetAlg::initialize ( )

◆ IntToTag1() [1/2]

unsigned int TagSetAlg::IntToTag1 ( int  val1,
int  val2,
int  val3,
int  val4 
)

Definition at line 375 of file TagSetAlg.cxx.

376{
377 unsigned int res=0;
378 res = (val4 +(val3<<8) +(val2<<16) +(val1<<24));
379
380 //cout<<"int res = "<<res<<endl;
381 return res;
382}

Referenced by execute().

◆ IntToTag1() [2/2]

unsigned int TagSetAlg::IntToTag1 ( int  val1,
int  val2,
int  val3,
int  val4 
)

◆ IntToTag2() [1/2]

unsigned int TagSetAlg::IntToTag2 ( int  val1,
int  val2,
int  val3,
int  val4,
int  val5,
int  val6 
)

Definition at line 384 of file TagSetAlg.cxx.

385{
386 unsigned int res=0;
387 res = (val6 +(val5<<5) +(val4<<10) +(val3<<15) +(val2<<20) +(val1<<26));
388
389 //cout<<"int res = "<<res<<endl;
390 return res;
391}

Referenced by execute().

◆ IntToTag2() [2/2]

unsigned int TagSetAlg::IntToTag2 ( int  val1,
int  val2,
int  val3,
int  val4,
int  val5,
int  val6 
)
The documentation for this class was generated from the following files: