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

#include <TwoGamma.h>

+ Inheritance diagram for TwoGamma:

Public Member Functions

 TwoGamma (const std::string &name, ISvcLocator *pSvcLocator)
 
StatusCode initialize ()
 
StatusCode execute ()
 
StatusCode finalize ()
 

Detailed Description

Definition at line 7 of file TwoGamma.h.

Constructor & Destructor Documentation

◆ TwoGamma()

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

Definition at line 47 of file TwoGamma.cxx.

47 :
48 Algorithm(name, pSvcLocator)
49{
50 // m_tuple1 = 0;
51 for(int i = 0; i < 10; i++)
52 {
53 m_pass[i] = 0;
54 }
55 m_event = 0;
56
57 Ndata1 = 0;
58 Ndata2 = 0;
59 m_runNo = 0;
60
61 //Declare the properties
62 declareProperty("CmsEnergy", m_ecms = 3.686);
63
64 declareProperty("max1", m_max1 = 1.2);
65 declareProperty("max2", m_max2 = 0.8);
66 declareProperty("costheta", m_costheta = 0.8);
67 declareProperty("dphi1", m_dphi1 = 2.5);
68 declareProperty("dphi2", m_dphi2 = 5);
69 declareProperty("eff", m_eff = 0.07);
70 declareProperty("sec", m_sec = 184.1);
71}

Member Function Documentation

◆ execute()

StatusCode TwoGamma::execute ( )

Definition at line 136 of file TwoGamma.cxx.

137{
138 StatusCode sc=StatusCode::SUCCESS;
139 m_event++;
140
141 MsgStream log(msgSvc(),name());
142 log<<MSG::INFO<<"in execute()"<<endreq;
143 SmartDataPtr<Event::EventHeader> eventHeader(eventSvc(),"/Event/EventHeader");
144 SmartDataPtr<EvtRecEvent> evtRecEvent(eventSvc(),EventModel::EvtRec::EvtRecEvent);
145 SmartDataPtr<EvtRecTrackCol> evtRecTrkCol(eventSvc(),EventModel::EvtRec::EvtRecTrackCol);
146
147 m_run = eventHeader->runNumber();
148 m_rec = eventHeader->eventNumber();
149 int runNo = m_run;
150 int event = m_rec;
151 int time = eventHeader->time();
152 if(m_event > 1 && runNo != m_runNo) return sc;
153 m_runNo = runNo;
154 m_time = time;
155
156 if(m_rec%1000==0)cout<<"Run "<<m_run<<" Event "<<m_rec<<endl;
157
158 HepLorentzVector p4psip(0.011*m_ecms, 0.0, 0.005, m_ecms);
159 double psipBeta = (p4psip.vect()).mag()/(p4psip.e());
160
161 m_pass[0]+=1;
162
163 int NCharge = evtRecEvent->totalCharged();
164 if(NCharge != 0)return sc;
165
166 m_pass[1]+=1;
167
168 double Emax1=0.0;
169 double Emax2=0.0;
170 int Imax[2];
171
172 if(((evtRecEvent->totalTracks() - evtRecEvent->totalCharged() )<2) || ((evtRecEvent->totalTracks() - evtRecEvent->totalCharged() )>15)) return sc;
173 m_pass[2]+=1;
174
175 HepLorentzVector p4Gam_1;
176 HepLorentzVector p4Gam_2;
177
178 for(int i=evtRecEvent->totalCharged(); i<evtRecEvent->totalTracks(); i++)
179 {
180 EvtRecTrackIterator itTrk=evtRecTrkCol->begin() + i;
181 if(!(*itTrk)->isEmcShowerValid()) continue;
182 RecEmcShower *emcTrk = (*itTrk)->emcShower();
183
184 HepLorentzVector p4Gam;
185 double Phi = emcTrk->phi();
186 double Theta = emcTrk->theta();
187 double Ener=emcTrk->energy();
188
189 p4Gam.setPx(Ener * sin(Theta) * cos(Phi));
190 p4Gam.setPy(Ener * sin(Theta) * sin(Phi));
191 p4Gam.setPz(Ener * cos(Theta));
192 p4Gam.setE(Ener);
193 p4Gam.boost(-0.011, 0 ,-0.005 * 1.0 / 3.686);
194
195 if(Ener>Emax2)
196 {
197 Emax2=Ener;
198 Imax[1]=i;
199 p4Gam_2 = p4Gam;
200 }
201 if(Ener>Emax1)
202 {
203 Emax2=Emax1;
204 p4Gam_2 = p4Gam_1;
205 Imax[1]=Imax[0];
206 Emax1=Ener;
207 p4Gam_1 = p4Gam;
208 Imax[0]=i;
209 }
210 }
211 m_e1_lab = Emax1;
212 m_e2_lab = Emax2;
213 m_e_lab = Emax1 + Emax2;
214 log << MSG::INFO << "Emax1 = " << Emax1 <<"Emax2= "<<Emax2<< endreq;
215 if(Emax1 < m_max1 || Emax2 < m_max2) return sc;
216 m_pass[3]+=1;
217
218 //P4 in Lab
219 double emcphi[2],emctht[2];
220 for(int i = 0;i < 2; i++)
221 {
222 if(i>=evtRecTrkCol->size()) break;
223 EvtRecTrackIterator itTrk=evtRecTrkCol->begin() + Imax[i];
224 if(!(*itTrk)->isEmcShowerValid()) continue;
225 RecEmcShower *emcTrk = (*itTrk)->emcShower();
226 emcphi[i]=emcTrk->phi();
227 emctht[i]=emcTrk->theta();
228 }
229 double dltphi_lab = (fabs(emcphi[0] - emcphi[1]) - pai) * 180.0 / pai;
230 double dlttheta_lab = (fabs(emctht[0] + emctht[1]) - pai) * 180.0 / pai;
231 m_costheta1_lab = cos(emctht[0]);
232 m_costheta2_lab = cos(emctht[1]);
233 m_phi1_lab = emcphi[0] * 180.0 / pai;
234 m_phi2_lab = emcphi[1] * 180.0 / pai;
235 m_dlttheta_lab = dlttheta_lab;
236 m_dltphi_lab = dltphi_lab;
237
238 if(fabs(m_costheta1_lab) > m_costheta || fabs(m_costheta2_lab) > m_costheta) return sc;
239 m_pass[4]+=1;
240 //P4 in Lab
241
242 //P4 in Cms
243 double px1 = p4Gam_1.px();
244 double py1 = p4Gam_1.py();
245 double pz1 = p4Gam_1.pz();
246 double pxy1 = sqrt(px1 * px1 + py1 * py1);
247 double e1 = p4Gam_1.e();
248
249 double px2 = p4Gam_2.px();
250 double py2 = p4Gam_2.py();
251 double pz2 = p4Gam_2.pz();
252 double pxy2 = sqrt(px2 * px2 + py2 * py2);
253 double e2 = p4Gam_2.e();
254
255
256 double phi1 = 0;
257 double phi2 = 0;
258 if(atan(py1 * 1.0 / px1) > 0){
259 if(px1 > 0)phi1 = atan(py1 * 1.0 / px1);
260 else phi1 = -(pai - atan(py1 * 1.0 / px1));
261 }
262 else{
263 if(px1 > 0)phi1 = atan(py1 * 1.0 / px1);
264 else phi1 = pai + atan(py1 * 1.0 / px1);
265 }
266
267 if(atan(py2 * 1.0 / px2) > 0){
268 if(px2 > 0)phi2 = atan(py2 * 1.0 / px2);
269 else phi2 = -(pai - atan(py2 * 1.0 / px2));
270 }
271 else{
272 if(px2 > 0)phi2 = atan(py2 * 1.0 / px2);
273 else phi2 = pai + atan(py2 * 1.0 / px2);
274 }
275
276 double dltphi = (fabs(phi1 - phi2) - pai) * 180.0 / pai;
277
278 double theta1 = 0;
279 double theta2 = 0;
280
281 if(pz1 > 0) theta1 = atan(pxy1 * 1.0 / pz1);
282 else theta1 = (pai + atan(pxy1 * 1.0 / pz1));
283
284 if(pz2 > 0) theta2 = atan(pxy2 * 1.0 / pz2);
285 else theta2 = (pai + atan(pxy2 * 1.0 / pz2));
286
287 double dlttheta = (theta1 + theta2 - pai) * 180.0 / pai;
288
289 double xBoost = p4Gam_1.px() + p4Gam_2.px();
290 double yBoost = p4Gam_1.py() + p4Gam_2.py();
291 double zBoost = p4Gam_1.pz() + p4Gam_2.pz();
292 m_xBoost = xBoost;
293 m_yBoost = yBoost;
294 m_zBoost = zBoost;
295
296 m_costheta1 = cos(theta1);
297 m_costheta2 = cos(theta2);
298 m_phi1 = phi1 * 180.0 / pai;
299 m_phi2 = phi2 * 180.0 / pai;
300 m_dlttheta = dlttheta;
301 m_dltphi = dltphi;
302 m_e1 = e1;
303 m_e2 = e2;
304 m_e = e1 + e2;
305 //P4 in Cms
306
307 if(fabs(m_dltphi) < m_dphi1) Ndata1++;
308 if(fabs(m_dltphi) > m_dphi1 && fabs(m_dltphi) < m_dphi2) Ndata2++;
309
310
311 m_tuple1->write();
312// runNo = 0;
313
314 return StatusCode::SUCCESS;
315}
pai
const double pai
Definition BbEmc.cxx:48
sin
double sin(const BesAngle a)
Definition BesAngle.h:210
cos
double cos(const BesAngle a)
Definition BesAngle.h:213
Phi
double Phi(RecMdcKalTrack *trk)
Definition CalibEventSelect.cxx:99
runNo
int runNo
Definition DQA_TO_DB.cxx:12
phi2
Double_t phi2
Definition DataBase/tau_mode.c:7
dltphi
Double_t dltphi
Definition DataBase/tau_mode.c:7
phi1
Double_t phi1
Definition DataBase/tau_mode.c:7
time
Double_t time
Definition DataBase/tau_mode.c:7
e1
Double_t e1
Definition DataBase/tau_mode.c:7
e2
Double_t e2
Definition DataBase/tau_mode.c:7
EvtRecTrackIterator
EvtRecTrackCol::iterator EvtRecTrackIterator
Definition EvtRecTrack.h:111
msgSvc
IMessageSvc * msgSvc()
Definition ServiceAccessor.h:13
theta2
@ theta2
Definition TrkKalDeriv.h:24
theta1
@ theta1
Definition TrkKalDeriv.h:24
DstEmcShower::theta
double theta() const
Definition DstEmcShower.h:38
DstEmcShower::phi
double phi() const
Definition DstEmcShower.h:39
DstEmcShower::energy
double energy() const
Definition DstEmcShower.h:45
EventModel::EvtRec::EvtRecEvent
_EXTERN_ std::string EvtRecEvent
Definition EventModel.h:116
EventModel::EvtRec::EvtRecTrackCol
_EXTERN_ std::string EvtRecTrackCol
Definition EventModel.h:117

◆ finalize()

StatusCode TwoGamma::finalize ( )

Definition at line 318 of file TwoGamma.cxx.

319{
320 char head[100];
321 char foot[100] = ".txt";
322 sprintf(head,"OffLineLum_%04d",m_runNo);
323 strcat(head,foot);
324 ofstream fout(head,ios::out);
325 fout<<"==============================================="<<endl;
326 fout<<"DESIGNED For OffLine Luminosity BY 2Gam Events"<<endl<<endl;
327 fout<<"MANY THANKS to Prof. C.Z Yuan & Z.Y Wang"<<endl<<endl;
328 fout<<" 2009/05"<<endl;
329 // fout<<" Charmonium Group"<<endl;
330 fout<<"==============================================="<<endl<<endl<<endl;
331
332 double lum = (Ndata1 - Ndata2) * 1.0 / (m_eff * m_sec);
333 fout<<" Table List "<<endl<<endl;
334 fout<<"-----------------------------------------------"<<endl;
335 fout<<"Firstly Energy Gamma "<<m_max1<<" Gev"<<endl;
336 fout<<"Secondly Energy Gamma "<<m_max2<<" Gev"<<endl;
337 fout<<"Solid Angle Acceptance "<<m_costheta<<endl;
338 fout<<"QED XSection "<<m_sec<<" NB"<<endl;
339 fout<<"Monte Carto Efficiency "<<m_eff * 100<<"%"<<endl<<endl;
340 fout<<"runNo Luminosity "<<endl<<endl;
341 fout<<m_runNo<<" "<<lum<<" nb^(-1)"<<endl;
342 fout<<"-----------------------------------------------"<<endl;
343 fout.close();
344
345 MsgStream log(msgSvc(), name());
346 cout<< "in finalize()" <<endl;
347 cout<< "all events " <<m_pass[0]<<endl;
348 cout<< "NCharged==0 " <<m_pass[1]<<endl;
349 cout<< "Number of Gam " <<m_pass[2]<<endl;
350 cout<< "N events " <<m_pass[3]<<endl;
351 cout<< "N events 0.8 " <<m_pass[4]<<endl;
352 return StatusCode::SUCCESS;
353}
sprintf
sprintf(cut,"kal_costheta0_em>-0.93&&kal_costheta0_em<0.93&&kal_pxy0_em>=0.05+%d*0.1&&kal_pxy0_em<0.15+%d*0.1&&NGch>=2", j, j)
lum
Double_t lum[10]
Definition DataBase/tau_mode.c:36

◆ initialize()

StatusCode TwoGamma::initialize ( )

Definition at line 74 of file TwoGamma.cxx.

75{
76 MsgStream log(msgSvc(), name());
77
78 log << MSG::INFO << "in initialize()" << endmsg;
79
80 StatusCode status;
81
82 NTuplePtr nt1(ntupleSvc(), "FILE1/DiGam");
83 if ( nt1 ) m_tuple1 = nt1;
84 else
85 {
86 m_tuple1 = ntupleSvc()->book ("FILE1/DiGam", CLID_ColumnWiseTuple, "DiGam N-Tuple signal");
87 if ( m_tuple1 )
88 {
89 status = m_tuple1->addItem ("run", m_run );
90 status = m_tuple1->addItem ("rec", m_rec );
91 status = m_tuple1->addItem ("time", m_time );
92
93 status = m_tuple1->addItem ("ngood", m_ngood);
94 status = m_tuple1->addItem ("nchrg", m_nchrg);
95
96 status = m_tuple1->addItem ("Cms_e1", m_e1);
97 status = m_tuple1->addItem ("Cms_e2", m_e2);
98 status = m_tuple1->addItem ("Cms_e", m_e);
99 status = m_tuple1->addItem ("Cms_costheta1", m_costheta1);
100 status = m_tuple1->addItem ("Cms_costheta2", m_costheta2);
101 status = m_tuple1->addItem ("Cms_dltphi", m_dltphi);
102 status = m_tuple1->addItem ("Cms_dltphi_1", m_dltphi_1);
103 status = m_tuple1->addItem ("Cms_dlttheta", m_dlttheta);
104 status = m_tuple1->addItem ("Cms_phi1", m_phi1);
105 status = m_tuple1->addItem ("Cms_phi2", m_phi2);
106
107 status = m_tuple1->addItem ("Lab_e1", m_e1_lab);
108 status = m_tuple1->addItem ("Lab_e2", m_e2_lab);
109 status = m_tuple1->addItem ("Lab_e", m_e_lab);
110 status = m_tuple1->addItem ("Lab_costheta1", m_costheta1_lab);
111 status = m_tuple1->addItem ("Lab_costheta2", m_costheta2_lab);
112 status = m_tuple1->addItem ("Lab_dltphi", m_dltphi_lab);
113 status = m_tuple1->addItem ("Lab_dlttheta", m_dlttheta_lab);
114 status = m_tuple1->addItem ("Lab_phi1", m_phi1_lab);
115 status = m_tuple1->addItem ("Lab_phi2", m_phi2_lab);
116
117 status = m_tuple1->addItem ("xBoost", m_xBoost);
118 status = m_tuple1->addItem ("yBoost", m_yBoost);
119 status = m_tuple1->addItem ("zBoost", m_zBoost);
120 }
121 else
122 {
123 log << MSG::ERROR << "Cannot book N-tuple1:"<<long(m_tuple1)<<endmsg;
124 return StatusCode::FAILURE;
125 }
126 }
127 Ndata1 = 0;
128 Ndata2 = 0;
129 m_runNo = 0;
130
131 log << MSG::INFO << "successfully return from initialize()" <<endmsg;
132 return StatusCode::SUCCESS;
133}
ntupleSvc
INTupleSvc * ntupleSvc()
Definition ServiceAccessor.h:55
The documentation for this class was generated from the following files: