4#include "GaudiKernel/MsgStream.h"
5#include "GaudiKernel/AlgFactory.h"
6#include "GaudiKernel/ISvcLocator.h"
7#include "GaudiKernel/SmartDataPtr.h"
8#include "GaudiKernel/IDataProviderSvc.h"
9#include "GaudiKernel/PropertyMgr.h"
11#include "GaudiKernel/Bootstrap.h"
12#include "GaudiKernel/ISvcLocator.h"
23#include "GaudiKernel/INTupleSvc.h"
24#include "GaudiKernel/NTuple.h"
25#include "GaudiKernel/Bootstrap.h"
26#include "GaudiKernel/IHistogramSvc.h"
27#include "CLHEP/Vector/ThreeVector.h"
28#include "CLHEP/Vector/LorentzVector.h"
29#include "CLHEP/Vector/TwoVector.h"
31using CLHEP::Hep3Vector;
32using CLHEP::Hep2Vector;
33using CLHEP::HepLorentzVector;
34#include "CLHEP/Geometry/Point3D.h"
44#ifndef ENABLE_BACKWARDS_COMPATIBILITY
59using CLHEP::HepLorentzVector;
62const double mpi = 0.13957;
63const double mk = 0.493677;
64const double xmass[5] = {0.000511, 0.105658, 0.139570,0.493677, 0.938272};
65const double velc = 299.792458;
66typedef std::vector<int>
Vint;
67typedef std::vector<HepLorentzVector>
Vp4;
69static int counter[13]={0,0,0,0,0,0,0,0,0,0,0,0,0};
75 Algorithm(name, pSvcLocator) {
78 declareProperty(
"writentuple",m_writentuple =
true );
79 declareProperty(
"useVertexDB", m_useVertexDB =
false );
80 declareProperty(
"ecms",m_ecms = 3.097);
81 declareProperty(
"beamangle",m_beamangle = 0.022);
82 declareProperty(
"Vr0cut", m_vr0cut=1.0);
83 declareProperty(
"Vz0cut", m_vz0cut=10.0);
84 declareProperty(
"Coscut", m_coscut=0.93);
86 declareProperty(
"EnergyThreshold", m_energyThreshold=0.05);
87 declareProperty(
"GammaPhiCut", m_gammaPhiCut=20.0);
88 declareProperty(
"GammaThetaCut", m_gammaThetaCut=20.0);
89 declareProperty(
"GammaTrkCut", m_gammaTrkCut=20.0);
90 declareProperty(
"GammaTLCut", m_gammatlCut=0);
91 declareProperty(
"GammaTHCut", m_gammathCut=14);
94 declareProperty (
"acoll_e_cut", m_acoll_e_cut=6.);
95 declareProperty (
"acopl_e_cut", m_acopl_e_cut=6.);
96 declareProperty (
"poeb_e_cut", m_poeb_e_cut=0.5);
97 declareProperty (
"dtof_e_cut", m_dtof_e_cut=4.);
98 declareProperty (
"eoeb_e_cut", m_eoeb_e_cut=0.4);
99 declareProperty (
"etotal_e_cut", m_etotal_e_cut=0.8);
100 declareProperty (
"tpoeb_e_cut", m_tpoeb_e_cut=0.95);
101 declareProperty (
"tptotal_e_cut", m_tptotal_e_cut=0.16);
102 declareProperty (
"tetotal_e_cut", m_tetotal_e_cut=0.65);
106 declareProperty (
"m_useEMConly", m_useEMConly=
false);
107 declareProperty (
"m_usePID", m_usePID=
false);
108 declareProperty (
"m_useMDC", m_useMDC=
true);
109 declareProperty (
"m_useDEDX", m_useDEDX=
false);
110 declareProperty (
"m_useTOF", m_useTOF=
false);
111 declareProperty (
"m_useEMC", m_useEMC=
true);
112 declareProperty (
"m_useMUC", m_useMUC=
false);
114 declareProperty (
"m_use_acolle", m_use_acolle =
false);
115 declareProperty (
"m_use_acople",m_use_acople =
false);
116 declareProperty (
"m_use_eoeb",m_use_eoeb =
false);
117 declareProperty (
"m_use_deltatof", m_use_deltatof =
false);
118 declareProperty (
"m_use_poeb", m_use_poeb=
false);
119 declareProperty (
"m_use_mucinfo",m_use_mucinfo=
false);
120 declareProperty (
"m_use_ne",m_use_ne =
false);
127 MsgStream log(
msgSvc(), name());
129 log << MSG::INFO <<
"in initialize()" << endmsg;
131 status = service(
"THistSvc", m_thistsvc);
132 if(status.isFailure() ){
133 log << MSG::INFO <<
"Unable to retrieve pointer to THistSvc" << endreq;
138 e_z0 =
new TH1F(
"e_z0",
"e_z0",100,0,50);
139 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/e_z0", e_z0);
140 e_r0 =
new TH1F(
"e_r0",
"e_r0",100,0,10);
141 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/e_r0", e_r0);
143 m_ee_mass =
new TH1F(
"ee_mass",
"ee_mass", 500, m_ecms-0.5, m_ecms+0.5 );
144 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_mass", m_ee_mass);
145 m_ee_acoll =
new TH1F(
"ee_acoll",
"ee_acoll", 60, 0, 6 );
146 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_acoll", m_ee_acoll);
147 m_ee_eop_ep =
new TH1F(
"ee_eop_ep",
"ee_eop_ep", 100,0.4,1.4 );
148 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_eop_ep", m_ee_eop_ep);
149 m_ee_eop_em =
new TH1F(
"ee_eop_em",
"ee_eop_em", 100,0.4,1.4 );
150 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_eop_em", m_ee_eop_em);
151 m_ee_costheta_ep =
new TH1F(
"ee_costheta_ep",
"ee_costheta_ep", 100,-1,1 );
152 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_costheta_ep", m_ee_costheta_ep);
153 m_ee_costheta_em =
new TH1F(
"ee_costheta_em",
"ee_costheta_em", 100,-1,1 );
154 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_costheta_em", m_ee_costheta_em);
156 m_ee_phi_ep =
new TH1F(
"ee_phi_ep",
"ee_phi_ep", 120,-3.2,3.2 );
157 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_phi_ep", m_ee_phi_ep);
158 m_ee_phi_em =
new TH1F(
"ee_phi_em",
"ee_phi_em", 120,-3.2,3.2 );
159 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_phi_em", m_ee_phi_em);
161 m_ee_nneu =
new TH1I(
"ee_nneu",
"ee_nneu", 5,0,5 );
162 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_nneu", m_ee_nneu);
166 m_ee_eemc_ep=
new TH1F(
"ee_eemc_ep",
"ee_eemc_ep",100,m_ecms/2.-0.8,m_ecms/2.+0.3);
167 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_eemc_ep", m_ee_eemc_ep);
168 m_ee_eemc_em=
new TH1F(
"ee_eemc_em",
"ee_eemc_em",100,m_ecms/2.-0.8,m_ecms/2.+0.3);
169 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_eemc_em", m_ee_eemc_em);
170 m_ee_x_ep=
new TH1F(
"ee_x_ep",
"ee_x_ep",100,-1.0,1.0);
171 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_x_ep", m_ee_x_ep);
172 m_ee_y_ep=
new TH1F(
"ee_y_ep",
"ee_y_ep",100,-1.0,1.0);
173 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_y_ep", m_ee_y_ep);
174 m_ee_z_ep=
new TH1F(
"ee_z_ep",
"ee_z_ep",100,-10.0,10.0);
175 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_z_ep", m_ee_z_ep);
176 m_ee_x_em=
new TH1F(
"ee_x_em",
"ee_x_em",100,-1.0,1.0);
177 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_x_em", m_ee_x_em);
178 m_ee_y_em=
new TH1F(
"ee_y_em",
"ee_y_em",100,-1.0,1.0);
179 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_y_em", m_ee_y_em);
180 m_ee_z_em=
new TH1F(
"ee_z_em",
"ee_z_em",100,-10.0,10.0);
181 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_z_em", m_ee_z_em);
183 m_ee_px_ep=
new TH1F(
"ee_px_ep",
"ee_px_ep",200,-2.2,2.2);
184 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_px_ep", m_ee_px_ep);
185 m_ee_py_ep=
new TH1F(
"ee_py_ep",
"ee_py_ep",200,-2.2,2.2);
186 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_py_ep", m_ee_py_ep);
187 m_ee_pz_ep=
new TH1F(
"ee_pz_ep",
"ee_pz_ep",200,-2.,2.5);
188 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_pz_ep", m_ee_pz_ep);
189 m_ee_p_ep=
new TH1F(
"ee_p_ep",
"ee_p_ep",100, m_ecms/2.-0.8 , m_ecms/2.+0.5 );
190 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_p_ep", m_ee_p_ep);
191 m_ee_px_em=
new TH1F(
"ee_px_em",
"ee_px_em",100,-2.2,2.2);
192 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_px_em", m_ee_px_em);
193 m_ee_py_em=
new TH1F(
"ee_py_em",
"ee_py_em",100,-2.2,2.2);
194 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_py_em", m_ee_py_em);
195 m_ee_pz_em=
new TH1F(
"ee_pz_em",
"ee_pz_em",100,-2.5,2.);
196 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_pz_em", m_ee_pz_em);
197 m_ee_p_em=
new TH1F(
"ee_p_em",
"ee_p_em",100,m_ecms/2.-0.8,m_ecms/2.+0.5);
198 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_p_em", m_ee_p_em);
199 m_ee_deltatof=
new TH1F(
"ee_deltatof",
"ee_deltatof",50,0.0,10.0);
200 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_deltatof", m_ee_deltatof);
202 m_ee_pidchidedx_ep=
new TH1F(
"ee_pidchidedx_ep",
"ee_pidchidedx_ep",160,-4,4);
203 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_pidchidedx_ep", m_ee_pidchidedx_ep);
204 m_ee_pidchidedx_em=
new TH1F(
"ee_pidchidedx_em",
"ee_pidchidedx_em",160,-4,4);
205 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_pidchidedx_em", m_ee_pidchidedx_em);
206 m_ee_pidchitof1_ep=
new TH1F(
"ee_pidchitof1_ep",
"ee_pidchitof1_ep",160,-4,4);
207 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_pidchitof1_ep", m_ee_pidchitof1_ep);
208 m_ee_pidchitof1_em=
new TH1F(
"ee_pidchitof1_em",
"ee_pidchitof1_em",160,-4,4);
209 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_pidchitof1_em", m_ee_pidchitof1_em);
210 m_ee_pidchitof2_ep=
new TH1F(
"ee_pidchitof2_ep",
"ee_pidchitof2_ep",160,-4,4);
211 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_pidchitof2_ep", m_ee_pidchitof2_ep);
212 m_ee_pidchitof2_em=
new TH1F(
"ee_pidchitof2_em",
"ee_pidchitof2_em",160,-4,4);
213 status = m_thistsvc->regHist(
"/DQAHist/Bhabha/ee_pidchitof2_em", m_ee_pidchitof2_em);
217 NTuplePtr nt1(
ntupleSvc(),
"DQAFILE/Bhabha");
218 if ( nt1 ) m_tuple1 = nt1;
220 m_tuple1 =
ntupleSvc()->book (
"DQAFILE/Bhabha", CLID_ColumnWiseTuple,
"N-Tuple");
222 status = m_tuple1->addItem (
"run", m_run);
223 status = m_tuple1->addItem (
"rec", m_rec);
224 status = m_tuple1->addItem (
"Nchrg", m_ncharg);
225 status = m_tuple1->addItem (
"Nneu", m_nneu,0,40);
226 status = m_tuple1->addItem (
"NGch", m_ngch, 0, 40);
227 status = m_tuple1->addItem (
"NGam", m_nGam);
230 status = m_tuple1->addItem (
"bhabhatag", m_bhabhatag);
232 status = m_tuple1->addItem (
"acoll", m_acoll);
233 status = m_tuple1->addItem (
"acopl", m_acopl);
234 status = m_tuple1->addItem (
"deltatof", m_deltatof);
235 status = m_tuple1->addItem (
"eop1", m_eop1);
236 status = m_tuple1->addItem (
"eop2", m_eop2);
237 status = m_tuple1->addItem (
"eoeb1", m_eoeb1);
238 status = m_tuple1->addItem (
"eoeb2", m_eoeb2);
239 status = m_tuple1->addItem (
"poeb1", m_poeb1);
240 status = m_tuple1->addItem (
"poeb2", m_poeb2);
241 status = m_tuple1->addItem (
"etoeb1", m_etoeb1);
242 status = m_tuple1->addItem (
"etoeb2", m_etoeb2);
243 status = m_tuple1->addItem (
"mucinfo1", m_mucinfo1);
244 status = m_tuple1->addItem (
"mucinfo2", m_mucinfo2);
248 status = m_tuple1->addIndexedItem (
"delang",m_nneu, m_delang);
249 status = m_tuple1->addIndexedItem (
"delphi",m_nneu, m_delphi);
250 status = m_tuple1->addIndexedItem (
"delthe",m_nneu, m_delthe);
251 status = m_tuple1->addIndexedItem (
"npart",m_nneu, m_npart);
252 status = m_tuple1->addIndexedItem (
"nemchits",m_nneu, m_nemchits);
253 status = m_tuple1->addIndexedItem (
"module",m_nneu, m_module);
254 status = m_tuple1->addIndexedItem (
"x",m_nneu, m_x);
255 status = m_tuple1->addIndexedItem (
"y",m_nneu, m_y);
256 status = m_tuple1->addIndexedItem (
"z",m_nneu, m_z);
257 status = m_tuple1->addIndexedItem (
"px",m_nneu, m_px);
258 status = m_tuple1->addIndexedItem (
"py",m_nneu, m_py);
259 status = m_tuple1->addIndexedItem (
"pz",m_nneu, m_pz);
260 status = m_tuple1->addIndexedItem (
"theta",m_nneu, m_theta);
261 status = m_tuple1->addIndexedItem (
"phi",m_nneu, m_phi);
262 status = m_tuple1->addIndexedItem (
"dx",m_nneu, m_dx);
263 status = m_tuple1->addIndexedItem (
"dy",m_nneu, m_dy);
264 status = m_tuple1->addIndexedItem (
"dz",m_nneu, m_dz);
265 status = m_tuple1->addIndexedItem (
"dtheta",m_nneu, m_dtheta);
266 status = m_tuple1->addIndexedItem (
"dphi",m_nneu, m_dphi);
267 status = m_tuple1->addIndexedItem (
"energy",m_nneu, m_energy);
268 status = m_tuple1->addIndexedItem (
"dE",m_nneu, m_dE);
269 status = m_tuple1->addIndexedItem (
"eSeed",m_nneu, m_eSeed);
270 status = m_tuple1->addIndexedItem (
"nSeed",m_nneu, m_nSeed);
271 status = m_tuple1->addIndexedItem (
"e3x3",m_nneu, m_e3x3);
272 status = m_tuple1->addIndexedItem (
"e5x5",m_nneu, m_e5x5);
273 status = m_tuple1->addIndexedItem (
"secondMoment",m_nneu, m_secondMoment);
274 status = m_tuple1->addIndexedItem (
"latMoment",m_nneu, m_latMoment);
275 status = m_tuple1->addIndexedItem (
"a20Moment",m_nneu, m_a20Moment);
276 status = m_tuple1->addIndexedItem (
"a42Moment",m_nneu, m_a42Moment);
277 status = m_tuple1->addIndexedItem (
"getTime",m_nneu, m_getTime);
278 status = m_tuple1->addIndexedItem (
"getEAll",m_nneu, m_getEAll);
282 status = m_tuple1->addIndexedItem(
"charge", m_ngch, m_charge);
283 status = m_tuple1->addIndexedItem (
"vx", m_ngch, m_vx0);
284 status = m_tuple1->addIndexedItem (
"vy", m_ngch, m_vy0);
285 status = m_tuple1->addIndexedItem (
"vz", m_ngch, m_vz0);
286 status = m_tuple1->addIndexedItem (
"r0", m_ngch, m_vr0);
288 status = m_tuple1->addIndexedItem (
"px", m_ngch, m_px) ;
289 status = m_tuple1->addIndexedItem (
"py", m_ngch, m_py) ;
290 status = m_tuple1->addIndexedItem (
"pz", m_ngch, m_pz) ;
291 status = m_tuple1->addIndexedItem (
"p", m_ngch, m_p) ;
295 status = m_tuple1->addIndexedItem (
"kal_vx", m_ngch, m_kal_vx0);
296 status = m_tuple1->addIndexedItem (
"kal_vy", m_ngch, m_kal_vy0);
297 status = m_tuple1->addIndexedItem (
"kal_vz", m_ngch, m_kal_vz0);
300 status = m_tuple1->addIndexedItem (
"kal_px", m_ngch, m_kal_px) ;
301 status = m_tuple1->addIndexedItem (
"kal_py", m_ngch, m_kal_py) ;
302 status = m_tuple1->addIndexedItem (
"kal_pz", m_ngch, m_kal_pz) ;
303 status = m_tuple1->addIndexedItem (
"kal_p", m_ngch, m_kal_p) ;
306 status = m_tuple1->addIndexedItem (
"probPH" , m_ngch, m_probPH) ;
307 status = m_tuple1->addIndexedItem (
"normPH" , m_ngch, m_normPH) ;
308 status = m_tuple1->addIndexedItem (
"chie" , m_ngch, m_chie) ;
309 status = m_tuple1->addIndexedItem (
"chimu" , m_ngch, m_chimu) ;
310 status = m_tuple1->addIndexedItem (
"chipi" , m_ngch, m_chipi) ;
311 status = m_tuple1->addIndexedItem (
"chik" , m_ngch, m_chik) ;
312 status = m_tuple1->addIndexedItem (
"chip" , m_ngch, m_chip) ;
313 status = m_tuple1->addIndexedItem (
"ghit" , m_ngch, m_ghit) ;
314 status = m_tuple1->addIndexedItem (
"thit" , m_ngch, m_thit) ;
316 status = m_tuple1->addIndexedItem (
"e_emc" , m_ngch, m_e_emc) ;
317 status = m_tuple1->addIndexedItem (
"phi_emc" , m_ngch, m_phi_emc) ;
318 status = m_tuple1->addIndexedItem (
"theta_emc" , m_ngch, m_theta_emc) ;
320 status = m_tuple1->addIndexedItem (
"nhit_muc" , m_ngch, m_nhit_muc) ;
321 status = m_tuple1->addIndexedItem (
"nlay_muc" , m_ngch, m_nlay_muc) ;
322 status = m_tuple1->addIndexedItem (
"t_btof" , m_ngch, m_t_btof );
323 status = m_tuple1->addIndexedItem (
"t_etof" , m_ngch, m_t_etof );
324 status = m_tuple1->addIndexedItem (
"qual_etof" , m_ngch, m_qual_etof );
325 status = m_tuple1->addIndexedItem (
"tof_etof" , m_ngch, m_tof_etof );
326 status = m_tuple1->addIndexedItem (
"te_etof" , m_ngch, m_te_etof );
327 status = m_tuple1->addIndexedItem (
"tmu_etof" , m_ngch, m_tmu_etof );
328 status = m_tuple1->addIndexedItem (
"tpi_etof" , m_ngch, m_tpi_etof );
329 status = m_tuple1->addIndexedItem (
"tk_etof" , m_ngch, m_tk_etof );
330 status = m_tuple1->addIndexedItem (
"tp_etof" , m_ngch, m_tp_etof );
332 status = m_tuple1->addIndexedItem (
"qual_btof1", m_ngch, m_qual_btof1 );
333 status = m_tuple1->addIndexedItem (
"tof_btof1" , m_ngch, m_tof_btof1 );
334 status = m_tuple1->addIndexedItem (
"te_btof1" , m_ngch, m_te_btof1 );
335 status = m_tuple1->addIndexedItem (
"tmu_btof1" , m_ngch, m_tmu_btof1 );
336 status = m_tuple1->addIndexedItem (
"tpi_btof1" , m_ngch, m_tpi_btof1 );
337 status = m_tuple1->addIndexedItem (
"tk_btof1" , m_ngch, m_tk_btof1 );
338 status = m_tuple1->addIndexedItem (
"tp_btof1" , m_ngch, m_tp_btof1 );
340 status = m_tuple1->addIndexedItem (
"qual_btof2", m_ngch, m_qual_btof2 );
341 status = m_tuple1->addIndexedItem (
"tof_btof2" , m_ngch, m_tof_btof2 );
342 status = m_tuple1->addIndexedItem (
"te_btof2" , m_ngch, m_te_btof2 );
343 status = m_tuple1->addIndexedItem (
"tmu_btof2" , m_ngch, m_tmu_btof2 );
344 status = m_tuple1->addIndexedItem (
"tpi_btof2" , m_ngch, m_tpi_btof2 );
345 status = m_tuple1->addIndexedItem (
"tk_btof2" , m_ngch, m_tk_btof2 );
346 status = m_tuple1->addIndexedItem (
"tp_btof2" , m_ngch, m_tp_btof2 );
347 status = m_tuple1->addIndexedItem (
"pidcode" , m_ngch, m_pidcode);
348 status = m_tuple1->addIndexedItem (
"pidprob" , m_ngch, m_pidprob);
349 status = m_tuple1->addIndexedItem (
"pidchiDedx" , m_ngch, m_pidchiDedx);
350 status = m_tuple1->addIndexedItem (
"pidchiTof1" , m_ngch, m_pidchiTof1);
351 status = m_tuple1->addIndexedItem (
"pidchiTof2" , m_ngch, m_pidchiTof2);
354 status = m_tuple1->addItem (
"dedx_GoodHits_ep" , m_dedx_goodhits_ep);
355 status = m_tuple1->addItem (
"dedx_chi_ep" , m_dedx_chiep);
356 status = m_tuple1->addItem (
"dedx_GoodHits_em" , m_dedx_goodhits_em);
357 status = m_tuple1->addItem (
"dedx_chi_em" , m_dedx_chiem);
359 status = m_tuple1->addItem (
"px_cms_ep", m_px_cms_ep);
360 status = m_tuple1->addItem (
"py_cms_ep", m_py_cms_ep);
361 status = m_tuple1->addItem (
"pz_cms_ep", m_pz_cms_ep);
362 status = m_tuple1->addItem (
"e_cms_ep", m_e_cms_ep);
363 status = m_tuple1->addItem (
"p_cms_ep", m_p_cms_ep);
365 status = m_tuple1->addItem (
"cos_ep", m_cos_ep);
366 status = m_tuple1->addItem (
"kal_p_ep", m_kal_p_ep);
367 status = m_tuple1->addItem (
"kal_px_ep", m_kal_px_ep);
368 status = m_tuple1->addItem (
"kal_py_ep", m_kal_py_ep);
369 status = m_tuple1->addItem (
"kal_pz_ep", m_kal_pz_ep);
371 status = m_tuple1->addItem (
"px_cms_em", m_px_cms_em);
372 status = m_tuple1->addItem (
"py_cms_em", m_py_cms_em);
373 status = m_tuple1->addItem (
"pz_cms_em", m_pz_cms_em);
374 status = m_tuple1->addItem (
"e_cms_em", m_e_cms_em);
375 status = m_tuple1->addItem (
"p_cms_em", m_p_cms_em);
377 status = m_tuple1->addItem (
"cos_em", m_cos_em);
378 status = m_tuple1->addItem (
"kal_p_em", m_kal_p_em);
379 status = m_tuple1->addItem (
"kal_px_em", m_kal_px_em);
380 status = m_tuple1->addItem (
"kal_py_em", m_kal_py_em);
381 status = m_tuple1->addItem (
"kal_pz_em", m_kal_pz_em);
383 status = m_tuple1->addItem (
"mass_ee", m_mass_ee);
384 status = m_tuple1->addItem (
"px_ee", m_px_ee);
385 status = m_tuple1->addItem (
"py_ee", m_py_ee);
386 status = m_tuple1->addItem (
"pz_ee", m_pz_ee);
387 status = m_tuple1->addItem (
"e_ee", m_e_ee);
388 status = m_tuple1->addItem (
"p_ee", m_p_ee);
390 status = m_tuple1->addItem (
"nep", m_nep );
391 status = m_tuple1->addItem (
"nem", m_nem );
396 log << MSG::ERROR <<
" Cannot book N-tuple:" << long(m_tuple1) << endmsg;
397 return StatusCode::FAILURE;
405 log << MSG::INFO <<
"successfully return from initialize()" <<endmsg;
406 return StatusCode::SUCCESS;
412 const double beamEnergy = m_ecms/2.;
413 const HepLorentzVector
p_cms(m_ecms*
sin(m_beamangle*0.5),0.0,0.0,m_ecms);
414 const Hep3Vector
u_cms = -
p_cms.boostVector();
416 MsgStream log(
msgSvc(), name());
417 log << MSG::INFO <<
"in execute()" << endreq;
419 setFilterPassed(
false);
421 SmartDataPtr<Event::EventHeader> eventHeader(eventSvc(),
"/Event/EventHeader");
424 log << MSG::FATAL <<
"Could not find Event Header" << endreq;
425 return StatusCode::SUCCESS;
428 m_run = eventHeader->runNumber();
429 m_rec = eventHeader->eventNumber();
436 log << MSG::FATAL <<
"Could not find EvtRecEvent" << endreq;
437 return StatusCode::SUCCESS;
439 log << MSG::INFO <<
"ncharg, nneu, tottks = "
440 << evtRecEvent->totalCharged() <<
" , "
441 << evtRecEvent->totalNeutral() <<
" , "
442 << evtRecEvent->totalTracks() <<endreq;
444 m_ncharg = evtRecEvent->totalCharged();
445 m_nneu = evtRecEvent->totalNeutral();
449 HepSymMatrix Evx(3, 0);
452 Gaudi::svcLocator()->service(
"VertexDbSvc", vtxsvc);
459 Evx[0][0]=vv[0]*vv[0];
460 Evx[0][1]=vv[0]*vv[1];
461 Evx[1][1]=vv[1]*vv[1];
462 Evx[1][2]=vv[1]*vv[2];
463 Evx[2][2]=vv[2]*vv[2];
470 log << MSG::FATAL <<
"Could not find EvtRecTrackCol" << endreq;
471 return StatusCode::SUCCESS;
478 for(
int i = 0; i < evtRecEvent->totalCharged(); i++){
480 if(!(*itTrk)->isMdcTrackValid())
continue;
481 if(!(*itTrk)->isMdcKalTrackValid())
continue;
484 double pch=mdcTrk->
p();
485 double x0=mdcTrk->
x();
486 double y0=mdcTrk->
y();
487 double z0=mdcTrk->
z();
488 double phi0=mdcTrk->
helix(1);
492 double Rxy=(x0-xv)*
cos(phi0)+(y0-yv)*
sin(phi0);
499 if(fabs(z0) >= m_vz0cut)
continue;
501 if(fabs(Rxy) >= m_vr0cut)
continue;
504 nCharge += mdcTrk->
charge();
508 int nGood = iGood.size();
510 log << MSG::DEBUG <<
"ngood, totcharge = " << nGood <<
" , " << nCharge << endreq;
512 if((nGood != 2)||(nCharge!=0)){
513 return StatusCode::SUCCESS;
530 Vint ipip, ipim, iep,iem,imup,imum;
541 for(
int i = 0; i < m_ngch; i++) {
566 HepLorentzVector
ptrk;
567 ptrk.setPx(mdcTrk->
px()) ;
568 ptrk.setPy(mdcTrk->
py()) ;
569 ptrk.setPz(mdcTrk->
pz()) ;
570 double p3 =
ptrk.mag() ;
573 m_pidprob[i]=pid->
prob(0);
574 m_pidchiDedx[i]=pid->
chiDedx(0);
575 m_pidchiTof1[i]=pid->
chiTof1(0);
576 m_pidchiTof2[i]=pid->
chiTof2(0);
577 if(mdcTrk->
charge() > 0) {
578 iep.push_back(iGood[i]);
580 if (mdcTrk->
charge() < 0) {
581 iem.push_back(iGood[i]);
597 for(
int i = evtRecEvent->totalCharged(); i< evtRecEvent->totalTracks(); i++) {
598 if(i>=evtRecTrkCol->size())
break;
600 if(!(*itTrk)->isEmcShowerValid())
continue;
602 Hep3Vector emcpos(emcTrk->
x(), emcTrk->
y(), emcTrk->
z());
607 int module=emcTrk->
module();
608 double x = emcTrk->
x();
609 double y = emcTrk->
y();
610 double z = emcTrk->
z();
611 double dx = emcTrk->
dx();
612 double dy = emcTrk->
dy();
613 double dth = emcTrk->
dtheta();
614 double dph = emcTrk->
dphi();
615 double dz = emcTrk->
dz();
617 double dE = emcTrk->
dE();
618 double eSeed = emcTrk->
eSeed();
619 double e3x3 = emcTrk->
e3x3();
620 double e5x5 = emcTrk->
e5x5();
623 double getTime = emcTrk->
time();
624 double getEAll = emcTrk->
getEAll();
630 m_nemchits[iphoton]=
n;
631 m_npart[iphoton]=npart;
632 m_module[iphoton]=module;
633 m_theta[iphoton]=EmcPos.theta();
634 m_phi[iphoton]=EmcPos.phi();
641 m_dtheta[iphoton]=dth;
645 m_eSeed[iphoton]=eSeed;
646 m_nSeed[iphoton]=nseed;
647 m_e3x3[iphoton]=e3x3;
648 m_e5x5[iphoton]=e5x5;
649 m_secondMoment[iphoton]=secondMoment;
650 m_latMoment[iphoton]=latMoment;
651 m_getTime[iphoton]=getTime;
652 m_getEAll[iphoton]=getEAll;
653 m_a20Moment[iphoton]=a20Moment;
654 m_a42Moment[iphoton]=a42Moment;
661 for(
int j = 0; j < nGood; j++) {
663 if (!(*jtTrk)->isMdcTrackValid())
continue;
665 double jtcharge = jtmdcTrk->
charge();
666 if(!(*jtTrk)->isExtTrackValid())
continue;
671 double angd = extpos.angle(emcpos);
672 double thed = extpos.theta() - emcpos.theta();
673 double phid = extpos.deltaPhi(emcpos);
674 thed = fmod(thed+CLHEP::twopi+CLHEP::twopi+
pi, CLHEP::twopi) - CLHEP::pi;
675 phid = fmod(phid+CLHEP::twopi+CLHEP::twopi+
pi, CLHEP::twopi) - CLHEP::pi;
677 if(fabs(thed) < fabs(dthe)) dthe = thed;
678 if(fabs(phid) < fabs(dphi)) dphi = phid;
679 if(angd < dang) dang = angd;
688 dthe = dthe * 180 / (CLHEP::pi);
689 dphi = dphi * 180 / (CLHEP::pi);
690 dang = dang * 180 / (CLHEP::pi);
691 double eraw = emcTrk->
energy();
693 double phi = emcTrk->
phi();
694 double the = emcTrk->
theta();
696 m_delphi[iphoton]=dphi;
697 m_delthe[iphoton]=dthe;
698 m_delang[iphoton]=dang;
701 if (fc_theta < 0.80 ) {
702 if ( eraw <= m_energyThreshold/2)
continue;
703 }
else if (fc_theta > 0.86 && fc_theta < 0.92 ) {
704 if (eraw <= m_energyThreshold )
continue;
708 if(getTime>m_gammathCut||getTime<m_gammatlCut)
continue;
711 if(dang< m_gammaTrkCut)
continue;
715 if(iphoton>=40)
return StatusCode::SUCCESS;
718 if(iphoton>0) counter[4]++;
719 int nGam = iGam.size();
730 for(
int i = 0; i < m_nGam; i++) {
732 if(!(*itTrk)->isEmcShowerValid())
continue;
734 double eraw = emcTrk->
energy();
735 double phi = emcTrk->
phi();
736 double the = emcTrk->
theta();
737 HepLorentzVector
ptrk;
738 ex_gam+=eraw*
sin(the)*
cos(phi);
739 ey_gam+=eraw*
sin(the)*
sin(phi);
740 ez_gam+=eraw*
cos(the);
741 et_gam+=eraw*
sin(the);
768 for(
int i = 0; i < m_ngch; i++ ){
772 if(!(*itTrk)->isMdcTrackValid())
continue;
773 if(!(*itTrk)->isMdcKalTrackValid())
continue;
778 m_charge[i] = mdcTrk->
charge();
779 m_vx0[i] = mdcTrk->
x();
780 m_vy0[i] = mdcTrk->
y();
781 m_vz0[i] = mdcTrk->
z();
783 m_px[i] = mdcTrk->
px();
784 m_py[i] = mdcTrk->
py();
785 m_pz[i] = mdcTrk->
pz();
786 m_p[i] = mdcTrk->
p();
799 m_kal_vx0[i]= mdcKalTrk->
dr()*
cos(mdcKalTrk->
fi0());
801 m_kal_vy0[i] = mdcKalTrk->
dr()*
sin(mdcKalTrk->
fi0());
803 m_kal_vz0[i]= mdcKalTrk->
dz();
806 m_kal_px[i] = mdcKalTrk->
px();
807 m_kal_py[i] = mdcKalTrk->
py();
808 m_kal_pz[i] = mdcKalTrk->
pz();
810 t_pxy2 = m_kal_px[i]*m_kal_px[i] + m_kal_py[i]*m_kal_py[i];
811 m_kal_p[i] = sqrt(t_pxy2 + m_kal_pz[i]*m_kal_pz[i]);
812 double ptrk = m_kal_p[i];
816 pt_had+=sqrt(t_pxy2);
818 e_had+=sqrt(m_kal_p[i]*m_kal_p[i]+mdcKalTrk->
mass()*mdcKalTrk->
mass());
821 if((*itTrk)->isMdcDedxValid()) {
824 m_probPH[i]= dedxTrk->
probPH();
825 m_normPH[i]= dedxTrk->
normPH();
827 m_chie[i] = dedxTrk->
chiE();
828 m_chimu[i] = dedxTrk->
chiMu();
829 m_chipi[i] = dedxTrk->
chiPi();
830 m_chik[i] = dedxTrk->
chiK();
831 m_chip[i] = dedxTrk->
chiP();
836 if((*itTrk)->isEmcShowerValid()) {
838 m_e_emc[i] = emcTrk->
energy();
839 m_phi_emc[i] = emcTrk->
phi();
840 m_theta_emc[i] = emcTrk->
theta();
846 if((*itTrk)->isMucTrackValid()){
848 m_nhit_muc[i] = mucTrk->
numHits() ;
852 if((*itTrk)->isTofTrackValid()) {
854 SmartRefVector<RecTofTrack> tofTrkCol = (*itTrk)->tofTrack();
855 SmartRefVector<RecTofTrack>::iterator iter_tof = tofTrkCol.begin();
857 for(;iter_tof != tofTrkCol.end(); iter_tof++ ) {
861 status->
setStatus((*iter_tof)->status());
863 if( (status->
is_cluster()) ) m_t_etof[i] = (*iter_tof)->tof();
864 if( !(status->
is_counter()) ){
if(status)
delete status;
continue; }
865 if( status->
layer()!=0 ) {
if(status)
delete status;
continue;}
867 double path=(*iter_tof)->path();
868 double tof = (*iter_tof)->tof();
869 double ph = (*iter_tof)->ph();
870 double rhit = (*iter_tof)->zrhit();
871 double qual = 0.0 + (*iter_tof)->quality();
872 double cntr = 0.0 + (*iter_tof)->tofID();
874 for(
int j = 0; j < 5; j++) {
876 double beta = gb/sqrt(1+gb*gb);
877 texp[j] = path /beta/
velc;
879 m_qual_etof[i] = qual;
880 m_tof_etof[i] = tof ;
883 if( (status->
is_cluster()) ) m_t_btof[i] = (*iter_tof)->tof();
884 if( !(status->
is_counter()) ){
if(status)
delete status;
continue;}
885 if(status->
layer()==1){
886 double path=(*iter_tof)->path();
887 double tof = (*iter_tof)->tof();
888 double ph = (*iter_tof)->ph();
889 double rhit = (*iter_tof)->zrhit();
890 double qual = 0.0 + (*iter_tof)->quality();
891 double cntr = 0.0 + (*iter_tof)->tofID();
893 for(
int j = 0; j < 5; j++) {
895 double beta = gb/sqrt(1+gb*gb);
896 texp[j] = path /beta/
velc;
898 m_qual_btof1[i] = qual;
899 m_tof_btof1[i] = tof ;
902 if(status->
layer()==2){
903 double path=(*iter_tof)->path();
904 double tof = (*iter_tof)->tof();
905 double ph = (*iter_tof)->ph();
906 double rhit = (*iter_tof)->zrhit();
907 double qual = 0.0 + (*iter_tof)->quality();
908 double cntr = 0.0 + (*iter_tof)->tofID();
910 for(
int j = 0; j < 5; j++) {
912 double beta = gb/sqrt(1+gb*gb);
913 texp[j] = path /beta/
velc;
915 m_qual_btof2[i] = qual;
916 m_tof_btof2[i] = tof ;
919 if(status)
delete status;
930 if(m_ngch != 2 || nCharge != 0 )
return StatusCode::SUCCESS;
936 HepLorentzVector p41e,p42e,p4le;
937 Hep3Vector p31e,p32e,p3le;
938 HepLorentzVector p41m,p42m,p4lm;
939 Hep3Vector p31m,p32m,p3lm;
940 HepLorentzVector p41h,p42h,p4lh;
941 Hep3Vector p31h,p32h,p3lh;
958 double ex1, ey1, ez1, epx1,epy1, epz1,epp1,pidchidedx1,pidchitof11,pidchitof21,eoeb1,exoeb1,eyoeb1,ezoeb1,etoeb1,kalpp,cmsp;
959 double ex2, ey2, ez2, epx2,epy2, epz2,epp2,pidchidedx2,pidchitof12,pidchitof22,eoeb2,exoeb2,eyoeb2,ezoeb2,etoeb2,kalpm,cmsm;
962 for(
int i = 0; i < m_ngch; i++ ){
964 itTrk1= evtRecTrkCol->begin() + iGood[i];
965 mdcKalTrk1 = (*itTrk1)->mdcKalTrack();
967 if(!(*itTrk1)->isMdcDedxValid())
continue;
971 m_dedx_chiep=dedxTrk1->
chiE();
981 m_px_cms_ep=p41e.px();
982 m_py_cms_ep=p41e.py();
983 m_pz_cms_ep=p41e.pz();
985 m_p_cms_ep=sqrt(p41e.px()*p41e.px()+p41e.py()*p41e.py()+p41e.pz()*p41e.pz());
988 m_kal_p_ep=m_kal_p[i];
1004 pidchidedx1=m_pidchiDedx[i];
1005 pidchitof11=m_pidchiTof1[i];
1006 pidchitof21=m_pidchiTof2[i];
1008 eoeb1=m_e_emc[i]/beamEnergy;
1010 if(p41e.rho()>0) eope1=m_e_emc[i]/p41e.rho();
1013 exoeb1=m_e_emc[i]*
sin(m_theta_emc[i])*
cos(m_phi_emc[i])/beamEnergy;
1014 eyoeb1=m_e_emc[i]*
sin(m_theta_emc[i])*
sin(m_phi_emc[i])/beamEnergy;
1015 ezoeb1=m_e_emc[i]*
cos(m_theta_emc[i])/beamEnergy;
1016 etoeb1=m_e_emc[i]*
sin(m_theta_emc[i])/beamEnergy;
1018 if(m_nhit_muc[i]>=2&&m_nlay_muc[i]>=2) mucinfo1=1;
1025 itTrk2= evtRecTrkCol->begin() + iGood[i];
1026 mdcKalTrk2 = (*itTrk2)->mdcKalTrack();
1029 if(!(*itTrk2)->isMdcDedxValid())
continue;
1033 m_dedx_chiem=dedxTrk2->
chiE();
1042 m_px_cms_em=p42e.px();
1043 m_py_cms_em=p42e.py();
1044 m_pz_cms_em=p42e.pz();
1045 m_e_cms_em=p42e.e();
1046 m_p_cms_em=sqrt(p42e.px()*p42e.px()+p42e.py()*p42e.py()+p42e.pz()*p42e.pz());
1048 m_kal_p_em=m_kal_p[i];
1064 pidchidedx2=m_pidchiDedx[i];
1065 pidchitof12=m_pidchiTof1[i];
1066 pidchitof22=m_pidchiTof2[i];
1068 eoeb2=m_e_emc[i]/beamEnergy;
1071 if(p42e.rho()>0) eope2=m_e_emc[i]/p42e.rho();
1073 exoeb2=m_e_emc[i]*
sin(m_theta_emc[i])*
cos(m_phi_emc[i])/beamEnergy;
1074 eyoeb2=m_e_emc[i]*
sin(m_theta_emc[i])*
sin(m_phi_emc[i])/beamEnergy;
1075 ezoeb2=m_e_emc[i]*
cos(m_theta_emc[i])/beamEnergy;
1076 etoeb2=m_e_emc[i]*
sin(m_theta_emc[i])/beamEnergy;
1078 if(m_nhit_muc[i]>=2&&m_nlay_muc[i]>=2) mucinfo2=1;
1085 p4le=( e01 > e02 ) ?p41e:p42e;
1086 p4lm=( e01 > e02 ) ?p41m:p42m;
1087 p3le=( e01 > e02 ) ?p31e:p32e;
1088 p3lm=( e01 > e02 ) ?p31m:p32m;
1090 double acolle=180.-p31e.angle(p32e)* 180.0 / CLHEP::pi;
1091 double acople=180.- (p31e.perpPart()).angle(p32e.perpPart ())* 180.0 / CLHEP::pi;
1092 double poeb1e=p41e.rho()/beamEnergy;
1093 double poeb2e=p42e.rho()/beamEnergy;
1095 int ilarge=( e01 > e02 ) ?iip:iim;
1097 double eoebl=m_e_emc[ilarge]/beamEnergy;
1098 if(p4le.rho()>0)eopl=m_e_emc[ilarge]/p4le.rho();
1099 double exoebl= m_e_emc[ilarge]*
sin(m_theta_emc[ilarge])*
cos(m_phi_emc[ilarge])/beamEnergy;
1100 double eyoebl= m_e_emc[ilarge]*
sin(m_theta_emc[ilarge])*
sin(m_phi_emc[ilarge])/beamEnergy;
1101 double ezoebl=m_e_emc[ilarge]*
cos(m_theta_emc[ilarge])/beamEnergy;
1102 double etoebl=m_e_emc[ilarge]*
sin(m_theta_emc[ilarge])/beamEnergy;
1103 int mucinfol=(m_nhit_muc[ilarge]>=2&&m_nlay_muc[ilarge]>=2) ? 1 : 0;
1105 int pidel=( e01 > e02 ) ? m_nep : m_nem;
1107 if(m_t_btof[iip]*m_t_btof[iim]!=0) deltatof=fabs(m_t_btof[iip]-m_t_btof[iim]);
1108 if(m_t_etof[iip]*m_t_etof[iim]!=0) deltatof=fabs(m_t_etof[iip]-m_t_etof[iim]);
1112 if (m_use_acolle&&acolle>m_acoll_e_cut)
return StatusCode::SUCCESS;
1117 if (m_use_acople&&acople>m_acopl_e_cut)
return StatusCode::SUCCESS;
1121 if (m_use_eoeb&&sqrt((eoeb1-1)*(eoeb1-1)+(eoeb2-1)*(eoeb2-1))>m_tetotal_e_cut)
return StatusCode::SUCCESS;
1127 if (m_use_deltatof&&m_useTOF&&(deltatof>m_dtof_e_cut))
return StatusCode::SUCCESS;
1132 if (m_use_poeb&&poeb1e<m_tpoeb_e_cut&&poeb2e<m_tpoeb_e_cut&&(sqrt((poeb1e-1)*(poeb1e-1)+(poeb2e-1)*(poeb2e-1))>m_tptotal_e_cut))
return StatusCode::SUCCESS;
1137 if (m_use_mucinfo&&m_useMUC&&(mucinfo1!=0&&mucinfo2!=0))
return StatusCode::SUCCESS;
1142 if (m_use_ne&&m_usePID&&(m_nep!=1||m_nem!=1))
return StatusCode::SUCCESS;
1155 m_cos_ep=p41e.cosTheta ();
1156 m_cos_em=p42e.cosTheta ();
1157 m_mass_ee=(p41e+p42e).m();
1158 m_px_ee=(p41e+p42e).px();
1159 m_py_ee=(p41e+p42e).py();
1160 m_pz_ee=(p41e+p42e).pz();
1161 m_e_ee=(p41e+p42e).e();
1162 m_p_ee=(p41e+p42e).rho();
1164 m_deltatof=deltatof;
1171 m_mucinfo1=mucinfo1;
1172 m_mucinfo2=mucinfo2;
1185 (*itTrk1)->tagElectron();
1186 (*itTrk2)->tagElectron();
1192 (*itTrk1)->setQuality(0);
1193 (*itTrk2)->setQuality(0);
1198 setFilterPassed(
true);
1200 m_ee_mass->Fill((p41e+p42e).m());
1201 m_ee_acoll->Fill(acolle);
1202 m_ee_eop_ep->Fill(eope1);
1203 m_ee_eop_em->Fill(eope2);
1204 m_ee_costheta_ep->Fill(p41e.cosTheta ());
1205 m_ee_costheta_em->Fill(p42e.cosTheta ());
1207 m_ee_phi_ep->Fill(p41e.phi ());
1208 m_ee_phi_em->Fill(p42e.phi ());
1210 m_ee_nneu->Fill(m_nGam);
1213 m_ee_eemc_ep->Fill(e01);
1214 m_ee_eemc_em->Fill(e02);
1216 m_ee_x_ep->Fill(ex1);
1217 m_ee_y_ep->Fill(ey1);
1218 m_ee_z_ep->Fill(ez1);
1220 m_ee_x_em->Fill(ex2);
1221 m_ee_y_em->Fill(ey2);
1222 m_ee_z_em->Fill(ez2);
1225 m_ee_px_ep->Fill(epx1);
1226 m_ee_py_ep->Fill(epy1);
1227 m_ee_pz_ep->Fill(epz1);
1228 m_ee_p_ep->Fill(epp1);
1230 m_ee_px_em->Fill(epx2);
1231 m_ee_py_em->Fill(epy2);
1232 m_ee_pz_em->Fill(epz2);
1233 m_ee_p_em->Fill(epp2);
1235 m_ee_deltatof->Fill(deltatof);
1237 m_ee_pidchidedx_ep->Fill(pidchidedx1);
1238 m_ee_pidchidedx_em->Fill(pidchidedx2);
1239 m_ee_pidchitof1_ep->Fill(pidchitof11);
1240 m_ee_pidchitof1_em->Fill(pidchitof12);
1241 m_ee_pidchitof2_ep->Fill(pidchitof21);
1242 m_ee_pidchitof2_em->Fill(pidchitof22);
1244 if(m_writentuple) m_tuple1 -> write();
1246 return StatusCode::SUCCESS;
1252 MsgStream log(
msgSvc(), name());
1253 log << MSG::INFO <<
"in finalize()" << endmsg;
1257 for(
int i=0;i<13;i++){
1258 p[i]=(counter[i]*0.1/(counter[0]*0.1))*100;
1261 cout<<
"**************************************************************************************"<<endl<<endl;
1262 cout<<
"Cuts of Bhabha Selection"<<
'\t'<<
'\t'<<
'\t'<<
"Numbers"<<
'\t'<<
'\t'<<
'\t'<<
"Ratio"<<endl<<endl;
1264 cout<<
"Total Number Before All Cuts"<<
'\t'<<
'\t'<<
'\t'<<counter[0]<<
'\t'<<
'\t'<<
'\t'<<p[0]<<
"%"<<endl;
1265 cout<<
"Finish Good Charged Track Selection"<<
'\t'<<
'\t'<<counter[1]<<
'\t'<<
'\t'<<
'\t'<<p[1]<<
"%"<<endl;
1266 cout<<
"Finish Particle ID"<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<counter[2]<<
'\t'<<
'\t'<<
'\t'<<p[2]<<
"%"<<endl;
1267 cout<<
"Finish Good Neutral Track Selection"<<
'\t'<<
'\t'<<counter[3]<<
'\t'<<
'\t'<<
'\t'<<p[3]<<
"%"<<endl;
1268 cout<<
"Good Neutral Track Not Zero"<<
'\t'<<
'\t'<<
'\t'<<counter[4]<<endl;
1269 cout<<
"Finish Check Good Charged Track's Info."<<
'\t'<<
'\t'<<counter[5]<<
'\t'<<
'\t'<<
'\t'<<p[5]<<
"%"<<endl;
1270 if(m_use_acolle) cout<<
"After Acolle Cut"<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<counter[6]<<
'\t'<<
'\t'<<
'\t'<<p[6]<<
"%"<<endl;
1271 else cout<<
"No Acolle Cut"<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<
"NULL"<<
'\t'<<
'\t'<<
'\t'<<
"NULL"<<endl;
1272 if(m_use_acople) cout<<
"After Acople Cut"<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<counter[7]<<
'\t'<<
'\t'<<
'\t'<<p[7]<<
"%"<<endl;
1273 else cout<<
"No Acople Cut"<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<
"NULL"<<
'\t'<<
'\t'<<
'\t'<<
"NULL"<<endl;
1274 if(m_use_eoeb) cout<<
"After Eoeb Cut"<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<counter[8]<<
'\t'<<
'\t'<<
'\t'<<p[8]<<
"%"<<endl;
1275 else cout<<
"No Eoeb Cut"<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<
"NULL"<<
'\t'<<
'\t'<<
'\t'<<
"NULL"<<endl;
1276 if(m_use_deltatof) cout<<
"After Deltatof Cut"<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<counter[9]<<
'\t'<<
'\t'<<
'\t'<<p[9]<<
"%"<<endl;
1277 else cout<<
"No Deltatof Cut"<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<
"NULL"<<
'\t'<<
'\t'<<
'\t'<<
"NULL"<<endl;
1278 if(m_use_poeb) cout<<
"After Poeb Cut"<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<counter[10]<<
'\t'<<
'\t'<<
'\t'<<p[10]<<
"%"<<endl;
1279 else cout<<
"No Poeb Cut"<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<
"NULL"<<
'\t'<<
'\t'<<
'\t'<<
"NULL"<<endl;
1280 if(m_use_mucinfo) cout<<
"After Mucinfo Cut"<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<counter[11]<<
'\t'<<
'\t'<<
'\t'<<p[11]<<
"%"<<endl;
1281 else cout<<
"No Mucinfo Cut"<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<
"NULL"<<
'\t'<<
'\t'<<
'\t'<<
"NULL"<<endl;
1282 if(m_use_ne) cout<<
"After PID_ne Cut"<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<counter[12]<<
'\t'<<
'\t'<<
'\t'<<p[12]<<
"%"<<endl<<endl;
1283 else cout<<
"No PID_ne Cut"<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<
'\t'<<
"NULL"<<
'\t'<<
'\t'<<
'\t'<<
"NULL"<<endl;
1285 cout<<
"**************************************************************************************"<<endl<<endl<<endl;
1288 return StatusCode::SUCCESS;
double sin(const BesAngle a)
double cos(const BesAngle a)
const HepLorentzVector p_cms(0.034067, 0.0, 0.0, 3.097)
HepGeom::Point3D< double > HepPoint3D
std::vector< HepLorentzVector > Vp4
EvtRecTrackCol::iterator EvtRecTrackIterator
************Class m_ypar INTEGER m_KeyWgt INTEGER m_KeyIHVP INTEGER m_KeyGPS INTEGER m_IsBeamPolarized INTEGER m_EvtGenInterface DOUBLE PRECISION m_Emin DOUBLE PRECISION m_sphot DOUBLE PRECISION m_Xenph DOUBLE PRECISION m_q2 DOUBLE PRECISION m_PolBeam2 DOUBLE PRECISION m_xErrPb *COMMON c_KK2f $ !CMS energy average $ !Spin Polarization vector first beam $ !Spin Polarization vector second beam $ !Beam energy spread[GeV] $ !minimum hadronization energy[GeV] $ !input READ never touch them !$ !debug facility $ !maximum weight $ !inverse alfaQED $ !minimum real photon energy
DQASelBhabha(const std::string &name, ISvcLocator *pSvcLocator)
double secondMoment() const
const Hep3Vector emcPosition() const
const int emcVolumeNumber() const
const double dz(void) const
const double dr(void) const
static void setPidType(PidType pidType)
const double fi0(void) const
const double mass() const
const HepVector helix() const
......
static unsigned int barrel_ec(const Identifier &id)
Values of different levels (failure returns 0)
virtual bool isVertexValid()=0
virtual double * SigmaPrimaryVertex()=0
virtual double * PrimaryVertex()=0
int methodProbability() const
void setChiMinCut(const double chi=4)
void setRecTrack(EvtRecTrack *trk)
void setMethod(const int method)
double chiTof2(int n) const
void identify(const int pidcase)
double probElectron() const
void usePidSys(const int pidsys)
static ParticleID * instance()
bool IsPidInfoValid() const
double chiTof1(int n) const
double probProton() const
double chiDedx(int n) const
RecEmcID getShowerId() const
RecEmcEnergy getEAll() const
HepSymMatrix & getZErrorE()
unsigned int layer() const
void setStatus(unsigned int status)
HepLorentzVector p() const
_EXTERN_ std::string EvtRecEvent
_EXTERN_ std::string EvtRecTrackCol