3#include "GaudiKernel/MsgStream.h"
4#include "GaudiKernel/SmartDataPtr.h"
5#include "GaudiKernel/PropertyMgr.h"
6#include "GaudiKernel/Bootstrap.h"
8#include "GaudiKernel/INTupleSvc.h"
9#include "GaudiKernel/NTuple.h"
10#include "GaudiKernel/ITHistSvc.h"
12#include "CLHEP/Vector/ThreeVector.h"
13#include "CLHEP/Vector/LorentzVector.h"
15#include "EventModel/EventModel.h"
16#include "EventModel/Event.h"
18#include "EvtRecEvent/EvtRecEvent.h"
19#include "EvtRecEvent/EvtRecTrack.h"
20#include "DstEvent/TofHitStatus.h"
21#include "EventModel/EventHeader.h"
24#include "ParticleID/ParticleID.h"
26#include "DQAEvent/DQAEvent.h"
27#include "DQA_EMC/DQA_EMC.h"
29#ifndef ENABLE_BACKWARDS_COMPATIBILITY
32using CLHEP::HepLorentzVector;
38 Algorithm(name, pSvcLocator) {
45 MsgStream log(
msgSvc(), name());
47 log << MSG::INFO <<
"in initialize()" << endmsg;
55 if ( nt ) m_tuple = nt;
57 m_tuple =
ntupleSvc()->book(
"DQAFILE/EMC", CLID_ColumnWiseTuple,
"EMC ntuple");
59 status = m_tuple->addItem(
"ixtal", m_ixtal);
60 status = m_tuple->addItem(
"npart", m_npart);
61 status = m_tuple->addItem(
"ntheta",m_ntheta);
62 status = m_tuple->addItem(
"nphi", m_nphi);
63 status = m_tuple->addItem(
"theta", m_theta);
64 status = m_tuple->addItem(
"phi", m_phi);
65 status = m_tuple->addItem(
"emcX", m_emcX);
66 status = m_tuple->addItem(
"emcY", m_emcY);
67 status = m_tuple->addItem(
"eSeed", m_eSeed);
68 status = m_tuple->addItem(
"e5x5", m_e5x5);
69 status = m_tuple->addItem(
"energy",m_energy);
70 status = m_tuple->addItem(
"time", m_time);
73 log << MSG::ERROR <<
"Can not book N-tuple:" << long(m_tuple) << endreq;
78 if(service(
"THistSvc", m_thistsvc).isFailure()) {
79 log << MSG::ERROR <<
"Couldn't get THistSvc" << endreq;
80 return StatusCode::FAILURE;
86 std::string HistName0, HistName;
87 HistName0=
"/DQAHist/EMC/";
89 sprintf( name,
"EMC_Bhabha_ShowerEneBarrelVsEvent");
90 sprintf(
title,
"EMC Bhabha ShowerEneBarrel vs event");
92 m_HistEnergyB =
new TH1F(name,
title,200,0.,2.0);
93 m_HistEnergyB->GetXaxis()->SetTitle(
"shower energy(GeV)");
94 m_HistEnergyB->GetXaxis()->CenterTitle();
95 m_HistEnergyB->GetYaxis()->SetTitle(
"Number of event");
96 m_HistEnergyB->GetYaxis()->CenterTitle();
98 HistName=HistName0+name;
99 if(m_thistsvc->regHist(HistName, m_HistEnergyB).isFailure()){
100 log << MSG::ERROR <<
"Couldn't register " <<name<< endreq;
103 sprintf( name,
"EMC_Bhabha_ShowerEneEastVsEvent");
104 sprintf(
title,
"EMC Bhabha ShowerEneEast vs event");
105 m_HistEnergyEast =
new TH1F(name,
title,200,0.,2.0);
106 m_HistEnergyEast->GetXaxis()->SetTitle(
"shower energy(GeV)");
107 m_HistEnergyEast->GetXaxis()->CenterTitle();
108 m_HistEnergyEast->GetYaxis()->SetTitle(
"Number of event");
109 m_HistEnergyEast->GetYaxis()->CenterTitle();
111 HistName=HistName0+name;
112 if(m_thistsvc->regHist(HistName, m_HistEnergyEast).isFailure()){
113 log << MSG::ERROR <<
"Couldn't register "<<name << endreq;
117 sprintf( name,
"EMC_Bhabha_ShowerEneWestVsEvent");
118 sprintf(
title,
"EMC Bhabha ShowerEneWest vs event");
119 m_HistEnergyWest =
new TH1F(name,
title,200,0.,2.0);
120 m_HistEnergyWest->GetXaxis()->SetTitle(
"shower energy(GeV)");
121 m_HistEnergyWest->GetXaxis()->CenterTitle();
122 m_HistEnergyWest->GetYaxis()->SetTitle(
"Number of event");
123 m_HistEnergyWest->GetYaxis()->CenterTitle();
126 HistName=HistName0+name;
127 if( m_thistsvc->regHist(HistName, m_HistEnergyWest).isFailure()){
128 log << MSG::ERROR <<
"Couldn't register "<<name << endreq;
132 sprintf( name,
"EMC_Bhabha_ShowerThetaVsvent");
133 sprintf(
title,
"EMC Bhabha ShowerTheta vs event");
134 m_HistTheta =
new TH1F(name,
title,56, 0, 56);
135 m_HistTheta->GetXaxis()->SetTitle(
"shower ID(theta)");
136 m_HistTheta->GetXaxis()->CenterTitle();
137 m_HistTheta->GetYaxis()->SetTitle(
"Number of event");
138 m_HistTheta->GetYaxis()->CenterTitle();
140 HistName=HistName0+name;
141 if( m_thistsvc->regHist(HistName, m_HistTheta).isFailure()){
142 log << MSG::ERROR <<
"Couldn't register" <<name<< endreq;
145 sprintf(name,
"EMC_Bhabha_ShowerCosTheta");
146 sprintf(
title,
"Emc Bhabha Costheta");
147 m_HistCosTheta =
new TH1F(name,
title,200, -1.0, 1.0);
148 m_HistCosTheta->GetXaxis()->SetTitle(
"shower cos(theta)");
149 m_HistCosTheta->GetXaxis()->CenterTitle();
150 m_HistCosTheta->GetYaxis()->SetTitle(
"Number of event");
151 m_HistCosTheta->GetYaxis()->CenterTitle();
153 HistName=HistName0+name;
154 if(m_thistsvc->regHist(HistName, m_HistCosTheta).isFailure()){
155 log << MSG::ERROR <<
"Couldn't register "<<name << endreq;
159 sprintf( name,
"EMC_Bhabha_ShowerPhiBarrelVsEvent");
160 sprintf(
title,
"EMC Bhabha ShowerPhiBarrel vs event");
161 m_HistPhiB =
new TH1F(name,
title,120, 0, 120);
162 m_HistPhiB->GetXaxis()->SetTitle(
"shower ID(phi)");
163 m_HistPhiB->GetXaxis()->CenterTitle();
164 m_HistPhiB->GetYaxis()->SetTitle(
"Number of event");
165 m_HistPhiB->GetYaxis()->CenterTitle();
168 HistName=HistName0+name;
169 if( m_thistsvc->regHist(HistName, m_HistPhiB).isFailure()){
170 log << MSG::ERROR <<
"Couldn't register "<<name << endreq;
173 sprintf( name,
"EMC_Bhabha_ShowerPhiEastVsEvent");
174 sprintf(
title,
"EMC Bhabha ShowerPhiEast vs event");
175 m_HistPhiEast =
new TH1F(name,
title,256, -3.14, 3.14);
176 m_HistPhiEast->GetXaxis()->SetTitle(
"shower phi(radian)");
177 m_HistPhiEast->GetXaxis()->CenterTitle();
178 m_HistPhiEast->GetYaxis()->SetTitle(
"Number of event");
179 m_HistPhiEast->GetYaxis()->CenterTitle();
181 HistName=HistName0+name;
182 if(m_thistsvc->regHist(HistName, m_HistPhiEast).isFailure()){
183 log << MSG::ERROR <<
"Couldn't register"<<name << endreq;
186 sprintf( name,
"EMC_Bhabha_ShowerPhiWestVsEvent");
187 sprintf(
title,
"EMC Bhabha ShowerPhiWest vs event");
188 m_HistPhiWest =
new TH1F(name,
title,256, -3.14, 3.14);
189 m_HistPhiWest->GetXaxis()->SetTitle(
"shower phi(radian)");
190 m_HistPhiWest->GetXaxis()->CenterTitle();
191 m_HistPhiWest->GetYaxis()->SetTitle(
"Number of event");
192 m_HistPhiWest->GetYaxis()->CenterTitle();
194 HistName=HistName0+name;
195 if(m_thistsvc->regHist(HistName, m_HistPhiWest).isFailure()){
196 log << MSG::ERROR <<
"Couldn't register" <<name<<endreq;
199 sprintf( name,
"EMC_Bhabha_ShowerThetaPhi");
200 sprintf(
title,
"EMC Bhabha ShowerThetaPhi");
201 m_ThetaPhi =
new TH2F(name,
"Theta versus Phi",
202 2000, -3.15, 3.15, 2000, 0.1, 3.0);
203 m_ThetaPhi->GetXaxis()->SetTitle(
"shower phi(radian)");
204 m_ThetaPhi->GetXaxis()->CenterTitle();
205 m_ThetaPhi->GetYaxis()->SetTitle(
"shower theta(radian)");
206 m_ThetaPhi->GetYaxis()->CenterTitle();
208 HistName=HistName0+name;
209 if(m_thistsvc->regHist(HistName, m_ThetaPhi).isFailure()){
210 log << MSG::ERROR <<
"Couldn't register" <<name<< endreq;
215 sprintf( name,
"EMC_Bhabha_Time-T0");
216 sprintf(
title,
"EMC Bhabha Time-T0 distribution");
217 m_HistTime =
new TH1F(name,
title,100, -40, 60);
218 m_HistTime->GetXaxis()->SetTitle(
"EmcTime-T0 (50ns)");
219 m_HistTime->GetXaxis()->CenterTitle();
220 m_HistTime->GetYaxis()->SetTitle(
"Number of event");
221 m_HistTime->GetYaxis()->CenterTitle();
223 HistName=HistName0+name;
224 if(m_thistsvc->regHist(HistName, m_HistTime).isFailure()){
225 log << MSG::ERROR <<
"Couldn't register" <<name<<endreq;
228 sprintf( name,
"EMC_Bhabha_ixtal");
229 sprintf(
title,
"EMC Bhabha ixtal distribution");
230 m_HistHitMap =
new TH1F(name,
title,6240, 0, 6240);
231 m_HistHitMap->GetXaxis()->SetTitle(
"ixtalNumber");
232 m_HistHitMap->GetXaxis()->CenterTitle();
233 m_HistHitMap->GetYaxis()->SetTitle(
"Number of event");
234 m_HistHitMap->GetYaxis()->CenterTitle();
236 HistName=HistName0+name;
237 if(m_thistsvc->regHist(HistName, m_HistHitMap).isFailure()){
238 log << MSG::ERROR <<
"Couldn't register" <<name<<endreq;
241 sprintf( name,
"EMC_Bhabha_eSeedvsIxtal ");
242 sprintf(
title,
"EMC Bhabha eSeed vs ixtal");
243 m_eSeedIxtal =
new TH2F(name,
"eSeed:Ixtal",
244 6240, 0, 6240, 2000, 0, 1.8);
245 m_eSeedIxtal->GetXaxis()->SetTitle(
"Ixtal)");
246 m_eSeedIxtal->GetXaxis()->CenterTitle();
247 m_eSeedIxtal->GetYaxis()->SetTitle(
"eSeed(GeV)");
248 m_eSeedIxtal->GetYaxis()->CenterTitle();
250 HistName=HistName0+name;
251 if(m_thistsvc->regHist(HistName, m_eSeedIxtal).isFailure()){
252 log << MSG::ERROR <<
"Couldn't register" <<name<< endreq;
255 sprintf( name,
"EMC_Bhabha_emcX:emcYeast ");
256 sprintf(
title,
"EMC Bhabha emcX vs emcY of east endcap");
257 m_XYeast =
new TH2F(name,
"emcX:emcY",
258 2000, -100, 100, 2000, -100, 100);
259 m_XYeast->GetXaxis()->SetTitle(
"emcX");
260 m_XYeast->GetXaxis()->CenterTitle();
261 m_XYeast->GetYaxis()->SetTitle(
"emcY");
262 m_XYeast->GetYaxis()->CenterTitle();
264 HistName=HistName0+name;
265 if(m_thistsvc->regHist(HistName, m_XYeast).isFailure()){
266 log << MSG::ERROR <<
"Couldn't register" <<name<< endreq;
269 sprintf( name,
"EMC_Bhabha_emcX:emcYwest ");
270 sprintf(
title,
"EMC Bhabha emcX vs emcY of west endcap");
271 m_XYwest =
new TH2F(name,
"emcX:emcY",
272 2000, -100, 100, 2000, -100, 100);
273 m_XYwest->GetXaxis()->SetTitle(
"emcX");
274 m_XYwest->GetXaxis()->CenterTitle();
275 m_XYwest->GetYaxis()->SetTitle(
"emcY");
276 m_XYwest->GetYaxis()->CenterTitle();
278 HistName=HistName0+name;
279 if(m_thistsvc->regHist(HistName, m_XYwest).isFailure()){
280 log << MSG::ERROR <<
"Couldn't register" <<name<< endreq;
285 scCalib = Gaudi::svcLocator() -> service(
"EmcCalibConstSvc", m_emcCalibConstSvc);
286 if( scCalib != StatusCode::SUCCESS){
287 log << MSG::ERROR <<
"can not use EmcCalibConstSvc" << endreq;
291 log << MSG::INFO <<
"successfully return from initialize()" <<endmsg;
292 return StatusCode::SUCCESS;
300 MsgStream log(
msgSvc(), name());
301 log << MSG::INFO <<
"in execute()" << endreq;
303 SmartDataPtr<Event::EventHeader> eventHeader(eventSvc(),
"/Event/EventHeader");
304 m_runNo=eventHeader->runNumber();
305 m_event=eventHeader->eventNumber();
306 log << MSG::DEBUG <<
"run, evtnum = "
313 SmartDataPtr<DQAEvent::DQAEvent> dqaevt(eventSvc(),
"/Event/DQATag");
315 log << MSG::INFO <<
"success get DQAEvent" << endreq;
317 log << MSG::ERROR <<
"Error accessing DQAEvent" << endreq;
318 return StatusCode::FAILURE;
321 log << MSG::DEBUG <<
"event tag = " << dqaevt->EventTag() << endreq;
324 if ( dqaevt->Bhabha() ) {
325 log << MSG::INFO <<
"Bhabha event" << endreq;
329 double eneShower,theta,phi,costheta;
334 std::string HistName;
336 for(
int i = 0; i < evtRecEvent->totalCharged(); i++){
339 log << MSG::DEBUG << i <<
" " << (*itTrk)->partId() <<
" "
340 << (*itTrk)->quality() << endreq;
343 if ( (*itTrk)->partId() != 1 )
continue;
355 if ( mdcTrk->
charge() > 0 ) {
356 log << MSG::DEBUG <<
"is electron" << endreq;
358 log << MSG::DEBUG <<
"is positron" << endreq;
361 if((*itTrk)->isEmcShowerValid()) {
364 eneShower=theShower->
energy();
365 theta = theShower->
theta();
366 phi= theShower->
phi();
382 m_ixtal=m_emcCalibConstSvc->
getIndex(npart,ntheta,nphi);
385 m_emcX = theShower->
x();
386 m_emcY = theShower->
y();
387 m_eSeed = theShower->
eSeed();
388 m_e5x5 = theShower->
e5x5();
389 m_energy = theShower->
energy();
390 m_time = theShower->
time();
396 sprintf(
Name,
"/DQAHist/EMC/EMC_Bhabha_Time-T0");
398 if (m_thistsvc->getHist(HistName, hmom).isSuccess()) {
401 log << MSG::ERROR <<
"Couldn't retrieve" <<HistName<< endreq;
404 sprintf(
Name,
"/DQAHist/EMC/EMC_Bhabha_ixtal");
406 if (m_thistsvc->getHist(HistName, hmom).isSuccess()) {
409 log << MSG::ERROR <<
"Couldn't retrieve" <<HistName<< endreq;
413 m_eSeedIxtal->Fill(m_ixtal,m_eSeed);
420 sprintf(
Name,
"/DQAHist/EMC/EMC_Bhabha_ShowerEneBarrelVsEvent");
422 if (m_thistsvc->getHist(HistName, hmom).isSuccess()) {
423 hmom->Fill(eneShower);
425 log << MSG::ERROR <<
"Couldn't retrieve" <<HistName<< endreq;
428 sprintf(
Name,
"/DQAHist/EMC/EMC_Bhabha_ShowerPhiBarrelVsEvent");
430 if (m_thistsvc->getHist(HistName, hmom ).isSuccess()) {
433 log << MSG::ERROR <<
"Couldn't retrieve" <<HistName<< endreq;
440 sprintf(
Name,
"/DQAHist/EMC/EMC_Bhabha_ShowerEneEastVsEvent");
442 if (m_thistsvc->getHist(HistName, hmom ).isSuccess()) {
443 hmom->Fill( eneShower);
445 log << MSG::ERROR <<
"Couldn't retrieve" <<HistName<< endreq;
448 sprintf(
Name,
"/DQAHist/EMC/EMC_Bhabha_ShowerPhiEastVsEvent");
450 if (m_thistsvc->getHist(HistName, hmom ).isSuccess()) {
453 log << MSG::ERROR <<
"Couldn't retrieve" <<HistName<< endreq;
456 m_XYeast->Fill(m_emcX,m_emcY);
461 sprintf(
Name,
"/DQAHist/EMC/EMC_Bhabha_ShowerEneWestVsEvent");
463 if (m_thistsvc->getHist(HistName, hmom ).isSuccess()) {
464 hmom->Fill( eneShower);
466 log << MSG::ERROR <<
"Couldn't retrieve" <<HistName<< endreq;
469 sprintf(
Name,
"/DQAHist/EMC/EMC_Bhabha_ShowerPhiWestVsEvent");
471 if (m_thistsvc->getHist(HistName, hmom ).isSuccess()) {
474 log << MSG::ERROR <<
"Couldn't retrieve" <<HistName<< endreq;
477 m_XYwest->Fill(m_emcX,m_emcY);
492 sprintf(
Name,
"/DQAHist/EMC/EMC_Bhabha_ShowerThetaVsvent");
494 if (m_thistsvc->getHist(HistName, hmom ).isSuccess()) {
497 log << MSG::ERROR <<
"Couldn't retrieve" <<HistName<< endreq;
500 sprintf(
Name,
"/DQAHist/EMC/EMC_Bhabha_ShowerCosTheta");
502 if (m_thistsvc->getHist(HistName, hmom ).isSuccess()) {
503 hmom->Fill(costheta);
505 log << MSG::ERROR <<
"Couldn't retrieve" <<HistName<< endreq;
516 m_ThetaPhi->Fill(phi, theta);
524 return StatusCode::SUCCESS;
531 MsgStream log(
msgSvc(), name());
532 log << MSG::INFO <<
"in finalize()" << endmsg;
533 return StatusCode::SUCCESS;
HepGeom::Point3D< double > HepPoint3D
EvtRecTrackCol::iterator EvtRecTrackIterator
double cos(const BesAngle a)
DQA_EMC(const std::string &name, ISvcLocator *pSvcLocator)
static unsigned int barrel_ec(const Identifier &id)
Values of different levels (failure returns 0)
static unsigned int theta_module(const Identifier &id)
static unsigned int phi_module(const Identifier &id)
virtual int getIndex(unsigned int PartId, unsigned int ThetaIndex, unsigned int PhiIndex) const =0
RecEmcID getShowerId() const
_EXTERN_ std::string EvtRecEvent
_EXTERN_ std::string EvtRecTrackCol