BOSS 7.1.0
BESIII Offline Software System
Loading...
Searching...
No Matches
TofRec.cxx
Go to the documentation of this file.
1// Package: TofRec
2// BESIII Tof Reconstruction Algorithm.
3// Created by Sun Shengsen (EPC IHEP)
4//
5#include "GaudiKernel/MsgStream.h"
6#include "GaudiKernel/AlgFactory.h"
7#include "GaudiKernel/ISvcLocator.h"
8#include "GaudiKernel/SmartDataPtr.h"
9#include "GaudiKernel/SmartDataLocator.h"
10#include "GaudiKernel/IDataProviderSvc.h"
11#include "GaudiKernel/PropertyMgr.h"
15#include "McTruth/McParticle.h"
16#include "McTruth/TofMcHit.h"
17#include "McTruth/McParticle.h"
28#include "TofRec/TofRec.h"
29#include "TofRec/TofTrack.h"
30#include "TofRec/TofCheckDigi.h"
31#include "TofRec/TofCount.h"
32#include "TofRec/TofRecTDS.h"
33#include <iostream>
34
37
39
40using namespace std;
41using namespace Event;
42
43/////////////////////////////////////////////////////////////////////////////
44
47// ITofQCorrSvc* tofQCorrSvc;
48// ITofQElecSvc* tofQElecSvc;
50
51DECLARE_COMPONENT(TofRec)
52TofRec::TofRec(const std::string& name, ISvcLocator* pSvcLocator) :
53 Algorithm(name, pSvcLocator)
54{
55 declareProperty( "AcceleratorStatus", m_acceleratorStatus );
56 declareProperty( "MagneticField", m_magneticField );
57 declareProperty( "ForCalibration", m_forCalibration );
58 declareProperty( "Data", m_data );
59 declareProperty( "CalibData", m_calibData );
60 // declareProperty( "CalibDataBarrel", m_calibDataBarrel );
61 declareProperty( "FirstIteration", m_firstIteration );
62 declareProperty( "CheckTrigger", m_checkTrigger );
63 declareProperty( "SaveRootFile4Rec", m_saveRootFile );
64 declareProperty( "PrintOutInfo", m_printOutInfo );
65 declareProperty( "CheckDigi", m_checkDigi );
66 declareProperty( "CheckDigiRaw", m_checkDigiRaw );
67 declareProperty( "CheckDigiExt", m_checkDigiExt );
68 declareProperty( "CheckMcTruth", m_checkMcTruth );
69}
70
71// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
72
73StatusCode TofRec::initialize(){
74
75 MsgStream log(msgSvc(), name());
76 log << MSG::INFO << "TofRec in initialize()" << endreq;
77
78 //Get TOF Geomertry Service
79 StatusCode scg = service("TofGeomSvc", tofGeomSvc);
80 if (scg== StatusCode::SUCCESS) {
81 log << MSG::INFO << "TofRec Get Geometry Service Sucessfully!" << endreq;
82 }else {
83 log << MSG::ERROR << "TofRec Get Geometry Service Failed !" << endreq;
84 return StatusCode::FAILURE;
85 }
86
87 //Get TOF Calibtration Service
88 StatusCode scc = service("TofCaliSvc", tofCaliSvc);
89 if (scc == StatusCode::SUCCESS) {
90 log << MSG::INFO << "TofRec Get Calibration Service Sucessfully!" << endreq;
91 } else {
92 log << MSG::ERROR << "TofRec Get Calibration Service Failed !" << endreq;
93 return StatusCode::FAILURE;
94 }
95
96 //Get TOF Raw Data Provider Service
97 StatusCode scd = service("RawDataProviderSvc", tofDigiSvc);
98 if (scd == StatusCode::SUCCESS) {
99 log << MSG::INFO << "TofRec Get Tof Raw Data Service Sucessfully!" << endreq;
100 } else {
101 log << MSG::ERROR << "TofRec Get Tof Raw Data Service Failed !" << endreq;
102 return StatusCode::FAILURE;
103 }
104
105 if( m_checkDigi ) {
106 NTuplePtr nt11(ntupleSvc(),"FILE203/digi");
107 NTuplePtr nt12(ntupleSvc(),"FILE203/barrel");
108 NTuplePtr nt13(ntupleSvc(),"FILE203/endcap");
109 NTuplePtr nt14(ntupleSvc(),"FILE203/mrpc");
110 NTuplePtr nt15(ntupleSvc(),"FILE203/ext");
111 NTuplePtr nt16(ntupleSvc(),"FILE203/tof");
112 NTuplePtr nt17(ntupleSvc(),"FILE203/bb");
113
114 if( nt11 || nt12 || nt13 || nt14 || nt15 || nt16 || nt17) {
115 m_tuple_digi = nt11;
116 m_tuple_barrel = nt12;
117 m_tuple_endcap = nt13;
118 m_tuple_endcap = nt14;
119 m_tuple_ext = nt15;
120 m_tuple_tof = nt16;
121 m_tuple_bb = nt17;
122 }
123 else {
124 m_tuple_digi = ntupleSvc()->book("FILE203/digi",CLID_ColumnWiseTuple,"TofRec");
125 m_tuple_barrel = ntupleSvc()->book("FILE203/barrel",CLID_ColumnWiseTuple,"TofRec");
126 m_tuple_endcap = ntupleSvc()->book("FILE203/endcap",CLID_ColumnWiseTuple,"TofRec");
127 m_tuple_mrpc = ntupleSvc()->book("FILE203/mrpc",CLID_ColumnWiseTuple,"TofRec");
128 m_tuple_ext = ntupleSvc()->book("FILE203/ext",CLID_ColumnWiseTuple,"TofRec");
129 m_tuple_tof = ntupleSvc()->book("FILE203/tof",CLID_ColumnWiseTuple,"TofRec");
130 m_tuple_bb = ntupleSvc()->book("FILE203/bb",CLID_ColumnWiseTuple,"TofRec");
131
132 if( m_tuple_digi || m_tuple_barrel || m_tuple_endcap || m_tuple_mrpc || m_tuple_ext || m_tuple_tof || m_tuple_bb ) {
133 m_checkdigi_tuple = new TofCheckDigi( m_tuple_digi, m_tuple_barrel, m_tuple_endcap, m_tuple_mrpc, m_tuple_ext, m_tuple_tof, m_tuple_bb );
134 }
135 else {
136 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_digi) <<endmsg;
137 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_barrel) <<endmsg;
138 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_endcap) <<endmsg;
139 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_mrpc) <<endmsg;
140 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_ext) <<endmsg;
141 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_tof) <<endmsg;
142 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_bb) <<endmsg;
143 return StatusCode::FAILURE;
144 }
145 }
146 }
147
148 if( m_saveRootFile ) {
149 NTuplePtr nt21(ntupleSvc(),"FILE201/trk");
150 NTuplePtr nt22(ntupleSvc(),"FILE201/cbtrk");
151 NTuplePtr nt23(ntupleSvc(),"FILE201/cetrk");
152 NTuplePtr nt24(ntupleSvc(),"FILE201/cetftrk");
153
154 if( nt21 || nt22 || nt23 || nt24 ) {
155 m_tuple_trk = nt21;
156 m_tuple_cbtrk = nt22;
157 m_tuple_cetrk = nt23;
158 m_tuple_cetftrk = nt24;
159 }
160 else {
161 m_tuple_trk = ntupleSvc()->book("FILE201/trk",CLID_ColumnWiseTuple,"TofRec");
162 m_tuple_cbtrk = ntupleSvc()->book("FILE201/cbtrk",CLID_ColumnWiseTuple,"TofRec");
163 m_tuple_cetrk = ntupleSvc()->book("FILE201/cetrk",CLID_ColumnWiseTuple,"TofRec");
164 m_tuple_cetftrk = ntupleSvc()->book("FILE201/cetftrk",CLID_ColumnWiseTuple,"TofRec");
165 if( m_tuple_trk || m_tuple_cbtrk || m_tuple_cetrk || m_tuple_cetftrk) {
166 m_checkdata_tuple = new TofCheckData( m_tuple_trk, m_tuple_cbtrk, m_tuple_cetrk, m_tuple_cetftrk );
167 }
168 else {
169 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_trk) <<endmsg;
170 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_cbtrk) <<endmsg;
171 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_cetrk) <<endmsg;
172 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_cetftrk) <<endmsg;
173 return StatusCode::FAILURE;
174 }
175 }
176 }
177
178 if( m_printOutInfo ) { m_printOut = new TofCount; }
179
180 return StatusCode::SUCCESS;
181}
182
183StatusCode TofRec::beginRun(){
184 MsgStream log(msgSvc(), name());
185 log << MSG::INFO <<"TofRec begin_run!"<< endreq;
186 return StatusCode::SUCCESS;
187}
188
189// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
190
191StatusCode TofRec::execute() {
192
193 MsgStream log(msgSvc(), name());
194 log << MSG::INFO << "TofRec in execute()!" << endreq;
195
196 SmartDataPtr<Event::EventHeader> eventHeader(eventSvc(),"/Event/EventHeader");
197 if( !eventHeader ) {
198 log << MSG::FATAL << "TofRec could not find Event Header!" << endreq;
199 return StatusCode::FAILURE;
200 }
201 int run = eventHeader->runNumber();
202 int event = eventHeader->eventNumber();
203 if( ( event % 1000 == 0 ) && m_printOutInfo ) {
204 std::cout << "run:" << run << " event: " << event << std::endl;
205 }
206 log << MSG::INFO << "run= " << run << " event= " << event << endreq;
207 if( tofCaliSvc->chooseConstants(run, event) == StatusCode:: FAILURE ) { return StatusCode:: FAILURE; }
208
209 TofRecTDS tds;
211
212 // Magnetic Field: Colliding data and Cosmic Ray
213 if( m_acceleratorStatus == "Colliding" ) {
214 // Retrieve Event Start Time
215 SmartDataPtr<RecEsTimeCol> estimeCol(eventSvc(),"/Event/Recon/RecEsTimeCol");
216 if( !estimeCol || ( estimeCol->size() == 0 ) ) {
217 log << MSG::WARNING << "TofRec Could not find RecEsTimeCol! Run = " << run << " Event = " << event << endreq;
218 return StatusCode::SUCCESS;
219 }
220 RecEsTimeCol::iterator iter_ESTime=estimeCol->begin();
221 double t0 = (*iter_ESTime)->getTest();
222 int t0Stat = (*iter_ESTime)->getStat();
223 log << MSG::INFO << "t0= " << t0 << " t0Stat= " << t0Stat << endreq;
224
225 // Read Offset of barrel and endcap for Event Start Time
226 IDataProviderSvc* m_pCalibDataSvc;
227 StatusCode sc = service("CalibDataSvc", m_pCalibDataSvc, true);
228 std::string fullpath = "/Calib/EsTimeCal";
229 SmartDataPtr<CalibData::EsTimeCalibData> Est( m_pCalibDataSvc, fullpath );
230
231 if( ( event % 1000 == 0 ) && m_printOutInfo ) {
232 if( !Est ) {
233 std::cout << " t0 = " << t0 << " t0stat = " << t0Stat << std::endl;
234 }
235 else {
236 std::cout << " t0 = " << t0 << " t0stat = " << t0Stat << " bunch = " << Est->getBunchTime() << " offset_b = " << Est->getToffsetb() << " offset_e = " << Est->getToffsete() << std::endl;
237 }
238 }
239
240
241 // Kalman Filter Track
242 SmartDataPtr<RecMdcKalTrackCol> mdcKalTrackCol(eventSvc(),"/Event/Recon/RecMdcKalTrackCol");
243 if( !mdcKalTrackCol ) {
244 log << MSG::WARNING << "No MdcKalTrackCol in TDS! Run = " << run << " Event = " << event << endreq;
245 return StatusCode::SUCCESS;
246 }
247 SmartDataPtr<RecMdcTrackCol> mdcTrackCol(eventSvc(),"/Event/Recon/RecMdcTrackCol");
248 if( !mdcTrackCol ) {
249 log << MSG::FATAL << "Could NOT find RecMdcTrackCol in TDS! Run = " << run << " Event = " << event << endreq;
250 return StatusCode::SUCCESS;
251 }
252
253 // Tof Get Extrapolated Track
254 SmartDataPtr<RecExtTrackCol> extTrackCol(eventSvc(),"/Event/Recon/RecExtTrackCol");
255 if( !extTrackCol ) {
256 log << MSG::WARNING << "No ExtTrackCol in TDS! Run = " << run << " Event = " << event << endreq;
257 return StatusCode::SUCCESS;
258 }
259 else {
260 if( m_printOutInfo ) { m_printOut->setExtTrackNum( extTrackCol->size() ); }
261 if( m_checkDigi && m_checkDigiExt ) { m_checkdigi_tuple->FillCol( *eventHeader, mdcTrackCol, mdcKalTrackCol, extTrackCol ); }
262 }
263
264 //sunss add 08/9/17
265 if( m_forCalibration ) { // Bhabha Events Selection
266 if( m_printOutInfo ) { m_printOut->addNumber(0); }
267 // if( t0Stat != 111 && t0Stat != 121 ) return StatusCode::SUCCESS;
268 if( t0Stat%10 != 1 ) return StatusCode::SUCCESS;
269 if( m_printOutInfo ) { m_printOut->addNumber(1); }
270
271 if( extTrackCol->size() != 2 ) return StatusCode::SUCCESS;
272 if( m_printOutInfo ) { m_printOut->addNumber(2); }
273
274 if( mdcTrackCol->size() != 2 ) return StatusCode::SUCCESS;
275 if( m_printOutInfo ) { m_printOut->addNumber(3); }
276
277 SmartDataPtr<RecEmcShowerCol> emcShowerCol(eventSvc(),"/Event/Recon/RecEmcShowerCol");
278 if( !emcShowerCol ) {
279 log << MSG::FATAL << "Could NOT find EmcRecShowerCol in TDS! Run = " << run << " Event = " << event << endreq;
280 return StatusCode::SUCCESS;
281 }
282 if( m_printOutInfo ) { m_printOut->addNumber(4); }
283
284 if( emcShowerCol->size() < 2 ) return StatusCode::SUCCESS;
285 if( m_printOutInfo ) { m_printOut->addNumber(5); }
286
287 RecMdcTrackCol::iterator iter_mdc1 = mdcTrackCol->begin();
288 RecMdcTrackCol::iterator iter_mdc2 = mdcTrackCol->begin() + 1;
289
290 RecMdcKalTrackCol::iterator iter_kal1 = mdcKalTrackCol->begin();
291 RecMdcKalTrackCol::iterator iter_kal2 = mdcKalTrackCol->begin() + 1;
292
293 RecExtTrackCol::iterator iter_ext1 = extTrackCol->begin();
294 RecExtTrackCol::iterator iter_ext2 = extTrackCol->begin() + 1;
295 Hep3Vector extPos1 = (*iter_ext1)->emcPosition();
296 Hep3Vector extPos2 = (*iter_ext2)->emcPosition();
297
298 RecEmcShowerCol::iterator iter_emc1 = emcShowerCol->begin();
299 RecEmcShowerCol::iterator iter_emc2 = emcShowerCol->begin() + 1;
300 Hep3Vector emcPos1((*iter_emc1)->x(),(*iter_emc1)->y(),(*iter_emc1)->z());
301 Hep3Vector emcPos2((*iter_emc2)->x(),(*iter_emc2)->y(),(*iter_emc2)->z());
302
303 Hep3Vector pep = (*iter_mdc1)->p3();
304 Hep3Vector pem = (*iter_mdc2)->p3();
305 double delta_angle = 180.0 - pep.angle( pem.unit() )*180.0/pi;
306 double delta_phi = abs( (*iter_mdc1)->phi() - (*iter_mdc2)->phi() )*180.0/pi;
307
308 Hep3Vector distant1 = extPos1 - emcPos1;
309 Hep3Vector distant2 = extPos2 - emcPos1;
310 if( distant1.r() > distant2.r() ) {
311 RecEmcShowerCol::iterator iter_tmp = iter_emc1;
312 iter_emc1 = iter_emc2;
313 iter_emc2 = iter_tmp;
314 Hep3Vector emc_tmp = emcPos1;
315 emcPos1 = emcPos2;
316 emcPos2 = emc_tmp;
317 }
318 distant1 = extPos1 - emcPos1;
319 distant2 = extPos2 - emcPos2;
320
321 double p1 = (*iter_mdc1)->p();
322 double p2 = (*iter_mdc2)->p();
323 double e1 = (*iter_emc1)->energy();
324 double e2 = (*iter_emc2)->energy();
325 double etot = 0.0;
326 RecEmcShowerCol::iterator iter_emc = emcShowerCol->begin();
327 for( ; iter_emc != emcShowerCol->end(); iter_emc++ ) {
328 etot += (*iter_emc)->energy();
329 }
330
331 if( m_checkDigi ) { m_checkdigi_tuple->FillCol( *eventHeader, extTrackCol, mdcTrackCol, emcShowerCol, mdcKalTrackCol ); }
332
333 if( ( (*iter_mdc1)->charge() + (*iter_mdc2)->charge() )!= 0 ) return StatusCode::SUCCESS;
334 if( m_printOutInfo ) { m_printOut->addNumber(6); }
335
336 if( delta_angle > 10.0 ) return StatusCode::SUCCESS;
337 if( m_printOutInfo ) { m_printOut->addNumber(7); }
338
339 if( (*iter_kal1)->getStat(1,0)!=0 || (*iter_kal2)->getStat(1,0)!=0 ) return StatusCode::SUCCESS;
340 if( m_printOutInfo ) { m_printOut->addNumber(8); }
341
342 if( distant1.r()>6.0 || distant2.r()>6.0 ) return StatusCode::SUCCESS;
343 if( m_printOutInfo ) { m_printOut->addNumber(9); }
344
345 // Jpsi data taken in 2009
346 if( m_data == "jpsi09" ) {
347 if( (*iter_mdc1)->x()<-0.2 || (*iter_mdc1)->x()>0.4 || (*iter_mdc1)->y()<-0.5 || (*iter_mdc1)->y()>0.1 || abs((*iter_mdc1)->z())>4.0 ) return StatusCode::SUCCESS;
348 if( m_printOutInfo ) { m_printOut->addNumber(10); }
349 if( (*iter_mdc2)->x()<-0.2 || (*iter_mdc2)->x()>0.4 || (*iter_mdc2)->y()<-0.5 || (*iter_mdc2)->y()>0.1 || abs((*iter_mdc2)->z())>4.0 ) return StatusCode::SUCCESS;
350 if( m_printOutInfo ) { m_printOut->addNumber(11); }
351 if( delta_phi<174.0 || delta_phi>186.0 ) return StatusCode::SUCCESS;
352 if( m_printOutInfo ) { m_printOut->addNumber(12); }
353 if( m_calibData == "Bhabha" ) {
354 if( e1 < 1.1 || e2 < 1.1 ) return StatusCode::SUCCESS;
355 if( m_printOutInfo ) { m_printOut->addNumber(13); }
356 }
357 }
358 // Psi prime data taken in 2009
359 else if( m_data == "psip09" ) {
360 if( (*iter_mdc1)->x()<-0.2 || (*iter_mdc1)->x()>0.4 || (*iter_mdc1)->y()<-0.5 || (*iter_mdc1)->y()>0.1 || abs((*iter_mdc1)->z())>4.0 ) return StatusCode::SUCCESS;
361 if( m_printOutInfo ) { m_printOut->addNumber(10); }
362 if( (*iter_mdc2)->x()<-0.2 || (*iter_mdc2)->x()>0.4 || (*iter_mdc2)->y()<-0.5 || (*iter_mdc2)->y()>0.1 || abs((*iter_mdc2)->z())>4.0 ) return StatusCode::SUCCESS;
363 if( m_printOutInfo ) { m_printOut->addNumber(11); }
364 if( delta_phi<174.0 || delta_phi>183.0 ) return StatusCode::SUCCESS;
365 if( m_printOutInfo ) { m_printOut->addNumber(12); }
366 if( m_calibData == "Bhabha" ) {
367 if( e1 < 1.4 || e2 < 1.4 ) return StatusCode::SUCCESS;
368 if( m_printOutInfo ) { m_printOut->addNumber(13); }
369 }
370 }
371 // Psi double prime data taken in first half of 2010
372 else if( m_data == "psipp10" ) {
373 if( (*iter_mdc1)->x()<-0.2 || (*iter_mdc1)->x()>1.2 || (*iter_mdc1)->y()<-0.9 || (*iter_mdc1)->y()>0.5 || abs((*iter_mdc1)->z())>6.0 ) return StatusCode::SUCCESS;
374 if( m_printOutInfo ) { m_printOut->addNumber(10); }
375 if( (*iter_mdc2)->x()<-0.2 || (*iter_mdc2)->x()>1.2 || (*iter_mdc2)->y()<-0.9 || (*iter_mdc2)->y()>0.5 || abs((*iter_mdc2)->z())>6.0 ) return StatusCode::SUCCESS;
376 if( m_printOutInfo ) { m_printOut->addNumber(11); }
377 if( delta_phi<174.0 || delta_phi>186.0 ) return StatusCode::SUCCESS;
378 if( m_printOutInfo ) { m_printOut->addNumber(12); }
379 if( m_calibData == "Bhabha" ) {
380 if( e1 < 1.4 || e2 < 1.4 ) return StatusCode::SUCCESS;
381 if( m_printOutInfo ) { m_printOut->addNumber(13); }
382 }
383 }
384 // Psi double prime data taken in end of 2010 and 2011
385 else if( m_data == "psipp11" ) {
386 if( (*iter_mdc1)->x()<-0.15 || (*iter_mdc1)->x()>0.3 || (*iter_mdc1)->y()<-0.3 || (*iter_mdc1)->y()>0.15 || abs((*iter_mdc1)->z())>6.0 ) return StatusCode::SUCCESS;
387 if( m_printOutInfo ) { m_printOut->addNumber(10); }
388 if( (*iter_mdc2)->x()<-0.15 || (*iter_mdc2)->x()>0.3 || (*iter_mdc2)->y()<-0.3 || (*iter_mdc2)->y()>0.15 || abs((*iter_mdc2)->z())>6.0 ) return StatusCode::SUCCESS;
389 if( m_printOutInfo ) { m_printOut->addNumber(11); }
390 if( delta_phi<174.0 || delta_phi>184.0 ) return StatusCode::SUCCESS;
391 if( m_printOutInfo ) { m_printOut->addNumber(12); }
392 if( m_calibData == "Bhabha" ) {
393 if( e1 < 1.4 || e2 < 1.4 ) return StatusCode::SUCCESS;
394 if( m_printOutInfo ) { m_printOut->addNumber(13); }
395 }
396 }
397 // Psi prime data taken in end of 2011 and 2012
398 else if( m_data == "psip12" ) {
399 if( (*iter_mdc1)->x()<-0.25 || (*iter_mdc1)->x()>0.3 || (*iter_mdc1)->y()<-0.3 || (*iter_mdc1)->y()>0.15 || abs((*iter_mdc1)->z())>6.0 ) return StatusCode::SUCCESS;
400 if( m_printOutInfo ) { m_printOut->addNumber(10); }
401 if( (*iter_mdc2)->x()<-0.25 || (*iter_mdc2)->x()>0.3 || (*iter_mdc2)->y()<-0.3 || (*iter_mdc2)->y()>0.15 || abs((*iter_mdc2)->z())>6.0 ) return StatusCode::SUCCESS;
402 if( m_printOutInfo ) { m_printOut->addNumber(11); }
403 if( delta_phi<172.0 || delta_phi>188.0 ) return StatusCode::SUCCESS;
404 if( m_printOutInfo ) { m_printOut->addNumber(12); }
405 if( m_calibData == "Bhabha" ) {
406 if( e1 < 1.4 || e2 < 1.4 ) return StatusCode::SUCCESS;
407 if( m_printOutInfo ) { m_printOut->addNumber(13); }
408 }
409 }
410 // Jpsi data taken in 2012
411 else if( m_data == "jpsi12" ) {
412 if( (*iter_mdc1)->x()<-0.2 || (*iter_mdc1)->x()>0.4 || (*iter_mdc1)->y()<-0.4 || (*iter_mdc1)->y()>0.2 || abs((*iter_mdc1)->z())>4.0 ) return StatusCode::SUCCESS;
413 if( m_printOutInfo ) { m_printOut->addNumber(10); }
414 if( (*iter_mdc2)->x()<-0.2 || (*iter_mdc2)->x()>0.4 || (*iter_mdc2)->y()<-0.4 || (*iter_mdc2)->y()>0.2 || abs((*iter_mdc2)->z())>4.0 ) return StatusCode::SUCCESS;
415 if( m_printOutInfo ) { m_printOut->addNumber(11); }
416 if( delta_phi<172.0 || delta_phi>188.0 ) return StatusCode::SUCCESS;
417 if( m_printOutInfo ) { m_printOut->addNumber(12); }
418 if( m_calibData == "Bhabha" ) {
419 if( e1 < 1.1 || e2 < 1.1 ) return StatusCode::SUCCESS;
420 if( m_printOutInfo ) { m_printOut->addNumber(13); }
421 }
422 }
423 // Y(4260) and Y(4360) taken in 2013
424 else if( m_data == "psi13" ) {
425 if( (*iter_mdc1)->x()<-0.15 || (*iter_mdc1)->x()>0.35 || (*iter_mdc1)->y()<-0.35 || (*iter_mdc1)->y()>0.15 || abs((*iter_mdc1)->z())>4.0 ) return StatusCode::SUCCESS;
426 if( m_printOutInfo ) { m_printOut->addNumber(10); }
427 if( (*iter_mdc2)->x()<-0.15 || (*iter_mdc2)->x()>0.35 || (*iter_mdc2)->y()<-0.35 || (*iter_mdc2)->y()>0.15 || abs((*iter_mdc2)->z())>4.0 ) return StatusCode::SUCCESS;
428 if( m_printOutInfo ) { m_printOut->addNumber(11); }
429 if( delta_phi<172.0 || delta_phi>188.0 ) return StatusCode::SUCCESS;
430 if( m_printOutInfo ) { m_printOut->addNumber(12); }
431 if( m_calibData == "Bhabha" ) {
432 if( e1 < 1.5 || e2 < 1.5 ) return StatusCode::SUCCESS;
433 if( m_printOutInfo ) { m_printOut->addNumber(13); }
434 }
435 }
436 else if( m_data == "rxyz14" ) {
437 if( (*iter_mdc1)->x()<-0.15 || (*iter_mdc1)->x()>0.35 || (*iter_mdc1)->y()<-0.35 || (*iter_mdc1)->y()>0.15 || abs((*iter_mdc1)->z())>4.0 ) return StatusCode::SUCCESS;
438 if( m_printOutInfo ) { m_printOut->addNumber(10); }
439 if( (*iter_mdc2)->x()<-0.15 || (*iter_mdc2)->x()>0.35 || (*iter_mdc2)->y()<-0.35 || (*iter_mdc2)->y()>0.15 || abs((*iter_mdc2)->z())>4.0 ) return StatusCode::SUCCESS;
440 if( m_printOutInfo ) { m_printOut->addNumber(11); }
441 if( delta_phi<175.0 || delta_phi>185.0 ) return StatusCode::SUCCESS;
442 if( m_printOutInfo ) { m_printOut->addNumber(12); }
443 if( m_calibData == "Bhabha" ) {
444 if( e1/p1 < 0.75 || e2/p2 < 0.75 ) return StatusCode::SUCCESS;
445 if( m_printOutInfo ) { m_printOut->addNumber(13); }
446 }
447 }
448 else if( m_data == "r15" ) {
449 if( (*iter_mdc1)->x()<-0.15 || (*iter_mdc1)->x()>0.35 || (*iter_mdc1)->y()<-0.35 || (*iter_mdc1)->y()>0.15 || abs((*iter_mdc1)->z())>4.0 ) return StatusCode::SUCCESS;
450 if( m_printOutInfo ) { m_printOut->addNumber(10); }
451 if( (*iter_mdc2)->x()<-0.15 || (*iter_mdc2)->x()>0.35 || (*iter_mdc2)->y()<-0.35 || (*iter_mdc2)->y()>0.15 || abs((*iter_mdc2)->z())>4.0 ) return StatusCode::SUCCESS;
452 if( m_printOutInfo ) { m_printOut->addNumber(11); }
453 if( delta_phi<175.0 || delta_phi>185.0 ) return StatusCode::SUCCESS;
454 if( m_printOutInfo ) { m_printOut->addNumber(12); }
455 if( m_calibData == "Bhabha" ) {
456 if( e1/p1 < 0.75 || e2/p2 < 0.75 ) return StatusCode::SUCCESS;
457 if( m_printOutInfo ) { m_printOut->addNumber(13); }
458 }
459 }
460 else if( m_data == "data16" ) {
461 if( (*iter_mdc1)->x()<-0.15 || (*iter_mdc1)->x()>0.35 || (*iter_mdc1)->y()<-0.35 || (*iter_mdc1)->y()>0.2 || abs((*iter_mdc1)->z())>4.0 ) return StatusCode::SUCCESS;
462 if( m_printOutInfo ) { m_printOut->addNumber(10); }
463 if( (*iter_mdc2)->x()<-0.15 || (*iter_mdc2)->x()>0.35 || (*iter_mdc2)->y()<-0.35 || (*iter_mdc2)->y()>0.2 || abs((*iter_mdc2)->z())>4.0 ) return StatusCode::SUCCESS;
464 if( m_printOutInfo ) { m_printOut->addNumber(11); }
465 if( delta_phi<170.0 || delta_phi>190.0 ) return StatusCode::SUCCESS;
466 if( m_printOutInfo ) { m_printOut->addNumber(12); }
467 if( m_calibData == "Bhabha" ) {
468 if( e1/p1 < 0.75 || e2/p2 < 0.75 ) return StatusCode::SUCCESS;
469 if( m_printOutInfo ) { m_printOut->addNumber(13); }
470 }
471 }
472 else if( m_data == "data17" ) {
473 if( (*iter_mdc1)->x()<-0.15 || (*iter_mdc1)->x()>0.35 || (*iter_mdc1)->y()<-0.3 || (*iter_mdc1)->y()>0.2 || (*iter_mdc1)->z()<-3.5 || (*iter_mdc1)->z()>4.5 ) return StatusCode::SUCCESS;
474 if( m_printOutInfo ) { m_printOut->addNumber(10); }
475 if( (*iter_mdc2)->x()<-0.15 || (*iter_mdc2)->x()>0.35 || (*iter_mdc2)->y()<-0.3 || (*iter_mdc2)->y()>0.2 || (*iter_mdc2)->z()<-3.5 || (*iter_mdc2)->z()>4.5 ) return StatusCode::SUCCESS;
476 if( m_printOutInfo ) { m_printOut->addNumber(11); }
477 if( delta_phi<170.0 || delta_phi>190.0 ) return StatusCode::SUCCESS;
478 if( m_printOutInfo ) { m_printOut->addNumber(12); }
479 if( m_calibData == "Bhabha" ) {
480 if( e1/p1 < 0.75 || e2/p2 < 0.75 ) return StatusCode::SUCCESS;
481 if( m_printOutInfo ) { m_printOut->addNumber(13); }
482 }
483 }
484 else if( m_data == "jpsi18" ) {
485 if( (*iter_mdc1)->x()<-0.1 || (*iter_mdc1)->x()>0.4 || (*iter_mdc1)->y()<-0.3 || (*iter_mdc1)->y()>0.2 || (*iter_mdc1)->z()<-3.5 || (*iter_mdc1)->z()>4.0 ) return StatusCode::SUCCESS;
486 if( m_printOutInfo ) { m_printOut->addNumber(10); }
487 if( (*iter_mdc2)->x()<-0.1 || (*iter_mdc2)->x()>0.4 || (*iter_mdc2)->y()<-0.3 || (*iter_mdc2)->y()>0.2 || (*iter_mdc2)->z()<-3.5 || (*iter_mdc2)->z()>4.0 ) return StatusCode::SUCCESS;
488 if( m_printOutInfo ) { m_printOut->addNumber(11); }
489 if( delta_phi<170.0 || delta_phi>190.0 ) return StatusCode::SUCCESS;
490 if( m_printOutInfo ) { m_printOut->addNumber(12); }
491 if( m_calibData == "Bhabha" ) {
492 if( e1/p1 < 0.75 || e2/p2 < 0.75 ) return StatusCode::SUCCESS;
493 if( m_printOutInfo ) { m_printOut->addNumber(13); }
494 }
495 }
496 else if( m_data == "jpsi19" ) {
497 if( (*iter_mdc1)->x()<-0.1 || (*iter_mdc1)->x()>0.4 || (*iter_mdc1)->y()<-0.3 || (*iter_mdc1)->y()>0.2 || (*iter_mdc1)->z()<-3.5 || (*iter_mdc1)->z()>4.0 ) return StatusCode::SUCCESS;
498 if( m_printOutInfo ) { m_printOut->addNumber(10); }
499 if( (*iter_mdc2)->x()<-0.1 || (*iter_mdc2)->x()>0.4 || (*iter_mdc2)->y()<-0.3 || (*iter_mdc2)->y()>0.2 || (*iter_mdc2)->z()<-3.5 || (*iter_mdc2)->z()>4.0 ) return StatusCode::SUCCESS;
500 if( m_printOutInfo ) { m_printOut->addNumber(11); }
501 if( delta_phi<170.0 || delta_phi>190.0 ) return StatusCode::SUCCESS;
502 if( m_printOutInfo ) { m_printOut->addNumber(12); }
503 if( m_calibData == "Bhabha" ) {
504 if( e1/p1 < 0.75 || e2/p2 < 0.75 ) return StatusCode::SUCCESS;
505 if( m_printOutInfo ) { m_printOut->addNumber(13); }
506 }
507 }
508
509 if( m_calibData == "Bhabha" ) {
510 if( ( etot - e1 - e2 ) > 0.3 ) return StatusCode::SUCCESS;
511 if( m_printOutInfo ) { m_printOut->addNumber(14); }
512 }
513 else if( m_calibData == "Dimu" ) {
514 if( e1 > 0.5 || e2 > 0.5 )return StatusCode::SUCCESS;
515 if( m_printOutInfo ) { m_printOut->addNumber(13); }
516 }
517 }
518 //sunss add 08/9/17
519
520 TofTrackVec* tofTrackVec = new TofTrackVec;
521 RecExtTrackCol::iterator iter_track = extTrackCol->begin();
522 for( ; iter_track < extTrackCol->end(); iter_track++ ) {
523 RecMdcTrackCol::iterator iter_mdc = mdcTrackCol->begin();
524 for( ; iter_mdc != mdcTrackCol->end(); iter_mdc++ ) {
525 if( (*iter_mdc)->trackId() == (*iter_track)->trackId() ) break;
526 }
527 double costheta = cos( (*iter_mdc)->theta() );
528 RecMdcKalTrackCol::iterator iter_kal = mdcKalTrackCol->begin();
529 for( ; iter_kal != mdcKalTrackCol->end(); iter_kal++ ) {
530 if( (*iter_kal)->trackId() == (*iter_track)->trackId() ) break;
531 }
532 double p[5] = {-1.0};
533 int kal[5] = {-1};
534 if( iter_kal != mdcKalTrackCol->end() ) {
535 for( unsigned int i=0; i<5; i++ ) {
536 if( i==0 ) { (*iter_kal)->setPidType(RecMdcKalTrack::electron); }
537 else if( i==1 ) { (*iter_kal)->setPidType(RecMdcKalTrack::muon); }
538 else if( i==2 ) { (*iter_kal)->setPidType(RecMdcKalTrack::pion); }
539 else if( i==3 ) { (*iter_kal)->setPidType(RecMdcKalTrack::kaon); }
540 else if( i==4 ) { (*iter_kal)->setPidType(RecMdcKalTrack::proton); }
541 p[i] = (*iter_kal)->p3().mag();
542 kal[i] = (*iter_kal)->getStat(0,i);
543 }
544 }
545 TofTrack* tof = new TofTrack( run, event );
546 tof->setExtTrack( (*iter_track), costheta, p, kal, t0, t0Stat );
547
548 if( tofTrackVec->size()>0 ) {
549 std::vector<TofTrack*>::iterator iterExt = tofTrackVec->begin();
550 for( ; iterExt < tofTrackVec->end(); iterExt++ ) {
551 if( (*iterExt)->isNoHit() ) continue;
552 tof->getMultiHit( *iterExt );
553 }
554 }
555
556 tofTrackVec->push_back( tof );
557 }
558
559 if( m_printOutInfo ) {
560 m_printOut->setTrack1Col( tofTrackVec );
561 }
562
563 // Retrieve Tof Digi Data
564 TofDataMap tofDataMap = tofDigiSvc->tofDataMapTof( t0 );
565 if( tofDataMap.empty() ) {
566 log << MSG::WARNING << "No Tof Data Map in RawDataProviderSvc! Run=" << run << " Event=" << event << endreq;
567 }
568
569 if( m_checkDigi && m_checkDigiRaw ) {
570 SmartDataPtr<TofDigiCol> tofDigiCol(eventSvc(),"/Event/Digi/TofDigiCol");
571 if( !tofDigiCol ) {
572 log << MSG::ERROR << "TofRec could not find Tof digi! Event = " << event << endreq;
573 }
574 else { m_checkdigi_tuple->FillCol( *eventHeader, tofDigiCol, t0, t0Stat ); }
575
576 m_checkdigi_tuple->FillCol( *eventHeader, tofDataMap, t0, t0Stat );
577 }
578
579 std::vector<int> deadId;
580 if( m_forCalibration ) {
581 for( unsigned int i=0; i<5; i++ ) {
582 int identmp = tofCaliSvc->BrEast(i);
583 if( identmp != 0x2fffffff ) {
584 deadId.push_back( identmp );
585 }
586 identmp = tofCaliSvc->BrWest(i);
587 if( identmp != 0x2fffffff ) {
588 deadId.push_back( identmp );
589 }
590 identmp = tofCaliSvc->Endcap(i);
591 if( identmp != 0x2fffffff ) {
592 deadId.push_back( identmp );
593 }
594 }
595 }
596
597 std::vector<TofTrack*>::iterator iter = tofTrackVec->begin();
598 for( ; iter < tofTrackVec->end(); iter++ ) {
599 if( (*iter)->isNoHit() ) continue;
600 (*iter)->setTofData( tofDataMap );
601 if( m_printOutInfo ) { m_printOut->setTrack2( (*iter) ); }
602 if( (*iter)->isNoHit() ) continue;
603 (*iter)->match( m_forCalibration, deadId, tofTrackVec );
604 if( m_printOutInfo ) { m_printOut->setTrack3( (*iter) ); }
605 }
606
607 iter = tofTrackVec->begin();
608 for( ; iter < tofTrackVec->end(); iter++ ) {
609
610 (*iter)->setCalibration();
611
612 if( m_checkDigi ) {
613 if( m_checkTrigger ) {
614 // retrieve trigger data for physics analysis
615 SmartDataPtr<TrigData> trigData(eventSvc(), "/Event/Trig/TrigData");
616 if (!trigData) {
617 log << MSG::FATAL << "Could not find Trigger Data for physics analysis" << endreq;
618 m_checkdigi_tuple->Fill_TofTrack( *eventHeader, *iter, t0, t0Stat );
619 }
620 else {
621 m_checkdigi_tuple->Fill_TofTrack( *eventHeader, *iter, t0, t0Stat, trigData );
622 }
623 }
624 else {
625 if( ( run < 0 ) && m_checkMcTruth ) {
626 SmartDataPtr<TofMcHitCol> tofMcCol(eventSvc(),"/Event/MC/TofMcHitCol");
627 SmartDataPtr<McParticleCol> mcParticleCol(eventSvc(),"/Event/MC/McParticleCol");
628 if ( !tofMcCol || !mcParticleCol ) {
629 m_checkdigi_tuple->Fill_TofTrack( *eventHeader, *iter, t0, t0Stat, mdcKalTrackCol );
630 }
631 else {
632 m_checkdigi_tuple->Fill_TofTrack( *eventHeader, *iter, t0, t0Stat, mdcKalTrackCol, tofMcCol, mcParticleCol, m_calibData );
633 }
634 }
635 else {
636 m_checkdigi_tuple->Fill_TofTrack( *eventHeader, *iter, t0, t0Stat, mdcKalTrackCol );
637 }
638 }
639 }
640 }
641
642 tds.RegisterTDS( eventHeader->runNumber(), eventHeader->eventNumber(), tofTrackVec, m_forCalibration, m_calibData );
643
644 clearTofTrackVec( tofTrackVec );
645
646// Check RecTofTrackCol Registered
647 SmartDataPtr<RecTofTrackCol> tofTrackCol(eventSvc(),"/Event/Recon/RecTofTrackCol");
648 if (!tofTrackCol) {
649 log << MSG::FATAL << "TofRec could not find RecTofTrackCol!" << endreq;
650 return StatusCode::FAILURE;
651 }
652 else{
653 if( m_saveRootFile ) {
654 m_checkdata_tuple->FillCol( *eventHeader, tofTrackCol, mdcKalTrackCol );
655 }
656 }
657
658 if( m_forCalibration ) {
659 SmartDataPtr<RecBTofCalHitCol> bhitCol(eventSvc(),"/Event/Recon/RecBTofCalHitCol");
660 if (!bhitCol) {
661 log << MSG::WARNING << "TofRec could not find RecBTofCalHitCol!" << endreq;
662 }
663 else {
664 if( m_saveRootFile ) {
665 m_checkdata_tuple->FillCol( *eventHeader, bhitCol );
666 }
667 }
668
669 SmartDataPtr<RecETofCalHitCol> ehitCol(eventSvc(),"/Event/Recon/RecETofCalHitCol");
670 if (!ehitCol) {
671 log << MSG::WARNING << "TofRec could not find RecETofCalHitCol!" << endreq;
672 }
673 else {
674 if( m_saveRootFile ) {
675 m_checkdata_tuple->FillCol( *eventHeader, ehitCol );
676 }
677 }
678 }
679
680 }
681 else {
682 log << MSG::FATAL << "In TofRec: AcceleratorStatus is NOT correct! m_acceleratorStatus = " << m_acceleratorStatus << endreq;
683 return StatusCode::FAILURE;
684 }
685
686 return StatusCode::SUCCESS;
687
688}
689
690// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
691
692StatusCode TofRec::finalize() {
693 if( m_printOutInfo ) {
694 if( m_forCalibration ) {
695 m_printOut->finalBhabha( m_calibData );
696 }
697 m_printOut->final();
698 delete m_printOut;
699 }
700 if( m_checkDigi ) delete m_checkdigi_tuple;
701 if( m_saveRootFile ) delete m_checkdata_tuple;
702 return StatusCode::SUCCESS;
703}
704
705
706//----------------------------------------------------------------------------
707
708void TofRec::clearTofTrackVec( std::vector<TofTrack*>*& tofTrackVec) {
709 if( tofTrackVec ) {
710 std::vector<TofTrack*>::iterator iter = tofTrackVec->begin();
711 for( ; iter < tofTrackVec->end(); iter++ ) {
712 delete (*iter);
713 }
714 tofTrackVec->clear();
715 delete tofTrackVec;
716 }
717 return;
718}
double cos(const BesAngle a)
Definition: BesAngle.h:213
Double_t etot
Double_t e1
Double_t e2
EvtStreamInputIterator< typename Generator::result_type > iter(Generator gen, int N=0)
INTupleSvc * ntupleSvc()
IMessageSvc * msgSvc()
std::multimap< unsigned int, TofData * > TofDataMap
Definition: TofData.h:244
IRawDataProviderSvc * tofDigiSvc
ITofCaliSvc * tofCaliSvc
ITofGeomSvc * tofGeomSvc
Definition: TofRec.cxx:45
ITofCaliSvc * tofCaliSvc
Definition: TofRec.cxx:46
IRawDataProviderSvc * tofDigiSvc
Definition: TofRec.cxx:49
std::vector< TofTrack * > TofTrackVec
Definition: TofTrack.h:220
virtual TofDataMap & tofDataMapTof(double estime=0)=0
virtual const int BrWest(unsigned int No)=0
virtual const int BrEast(unsigned int No)=0
virtual StatusCode chooseConstants(int run, int event)=0
virtual const int Endcap(unsigned int No)=0
void FillCol(Event::EventHeader &, RecTofTrackCol &, RecMdcKalTrackCol &)
void Fill_TofTrack(Event::EventHeader &, TofTrack *&, double, int)
void FillCol(Event::EventHeader &, TofDigiCol &, double, int)
void finalBhabha(std::string calibData)
Definition: TofCount.cxx:215
void setExtTrackNum(unsigned int ntrk)
Definition: TofCount.cxx:92
void final()
Definition: TofCount.cxx:182
void setTrack1Col(std::vector< TofTrack * > *&tofTrackVec)
Definition: TofCount.cxx:111
void addNumber(unsigned int i)
Definition: TofCount.cxx:207
void setTrack3(TofTrack *&tof)
Definition: TofCount.cxx:135
void setTrack2(TofTrack *&tof)
Definition: TofCount.cxx:120
StatusCode RegisterTDS(int runNumber, int eventNumber, std::vector< TofTrack * > *&tofTrackVec, bool m_forCalibration, std::string m_calibData)
Definition: TofRecTDS.cxx:58
StatusCode RegisterNullRecTofTrackCol()
Definition: TofRecTDS.cxx:22
Definition: TofRec.h:17
StatusCode initialize()
Definition: TofRec.cxx:73
StatusCode finalize()
Definition: TofRec.cxx:692
StatusCode execute()
Definition: TofRec.cxx:191
void clearTofTrackVec(std::vector< TofTrack * > *&tofTrackVec)
Definition: TofRec.cxx:708
StatusCode beginRun()
Definition: TofRec.cxx:183
void getMultiHit(TofTrack *&)
Definition: TofTrack.cxx:426
void setExtTrack(RecExtTrack *extTrack, double costheta, double p[5], int kal[5], double t0, int t0Stat)
Definition: TofTrack.cxx:157
Definition: Event.h:21
float costheta
double * p2
Definition: qcdloop1.h:76
const float pi
Definition: vector3.h:133