BOSS 7.1.1
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)
 

Detailed Description

Definition at line 12 of file TagSetAlg.h.

Constructor & Destructor Documentation

◆ TagSetAlg()

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

Definition at line 48 of file TagSetAlg.cxx.

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

Member Function Documentation

◆ execute()

StatusCode TagSetAlg::execute ( )

Definition at line 67 of file TagSetAlg.cxx.

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

◆ finalize()

StatusCode TagSetAlg::finalize ( )

Definition at line 356 of file TagSetAlg.cxx.

356 {
357
358 MsgStream log(msgSvc(), name());
359 log << MSG::INFO << "in finalize()" << iEvt << "Events processed" << endmsg;
360 return StatusCode::SUCCESS;
361}
msgSvc
IMessageSvc * msgSvc()
Definition ServiceAccessor.h:13

◆ initialize()

StatusCode TagSetAlg::initialize ( )

Definition at line 54 of file TagSetAlg.cxx.

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

◆ IntToTag1()

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

Definition at line 376 of file TagSetAlg.cxx.

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

Referenced by execute().

◆ IntToTag2()

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

Definition at line 385 of file TagSetAlg.cxx.

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

Referenced by execute().

The documentation for this class was generated from the following files: