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