BOSS 6.6.4.p01
BESIII Offline Software System
Loading...
Searching...
No Matches
MrpcRec.cxx
Go to the documentation of this file.
1// Package: MrpcRec
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 "MrpcRec/MrpcRec.h"
29#include "MrpcRec/MrpcTrack.h"
31#include "MrpcRec/MrpcCount.h"
32#include "MrpcRec/MrpcRecTDS.h"
33#include <iostream>
34
37
38using namespace std;
39using namespace Event;
40
41/////////////////////////////////////////////////////////////////////////////
42
45// ITofQCorrSvc* tofQCorrSvc;
46// ITofQElecSvc* tofQElecSvc;
48
49MrpcRec::MrpcRec(const std::string& name, ISvcLocator* pSvcLocator) :
50 Algorithm(name, pSvcLocator)
51{
52 declareProperty( "AcceleratorStatus", m_acceleratorStatus );
53 declareProperty( "MagneticField", m_magneticField );
54 declareProperty( "ForCalibration", m_forCalibration );
55 declareProperty( "Data", m_data );
56 declareProperty( "CalibData", m_calibData );
57 // declareProperty( "CalibDataBarrel", m_calibDataBarrel );
58 declareProperty( "FirstIteration", m_firstIteration );
59 declareProperty( "CheckTrigger", m_checkTrigger );
60 declareProperty( "SaveRootFile4Rec", m_saveRootFile );
61 declareProperty( "PrintOutInfo", m_printOutInfo );
62 declareProperty( "CheckDigi", m_checkDigi );
63 declareProperty( "CheckDigiRaw", m_checkDigiRaw );
64 declareProperty( "CheckDigiExt", m_checkDigiExt );
65 declareProperty( "CheckMcTruth", m_checkMcTruth );
66 declareProperty( "neighborhood", m_neighborhood=6 );//For the time being: 6 is the standardvalue -> It yield the highest reconstruction efficiency!
67}
68
69// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
70
72
73 MsgStream log(msgSvc(), name());
74 log << MSG::INFO << "MrpcRec in initialize()" << endreq;
75
76 //Get TOF Geomertry Service
77 StatusCode scg = service("TofGeomSvc", tofGeomSvc);
78 if (scg== StatusCode::SUCCESS) {
79 log << MSG::INFO << "MrpcRec Get Geometry Service Sucessfully!" << endreq;
80 }else {
81 log << MSG::ERROR << "MrpcRec Get Geometry Service Failed !" << endreq;
82 return StatusCode::FAILURE;
83 }
84
85 //Get TOF Calibtration Service
86 StatusCode scc = service("TofCaliSvc", tofCaliSvc);
87 if (scc == StatusCode::SUCCESS) {
88 log << MSG::INFO << "MrpcRec Get Calibration Service Sucessfully!" << endreq;
89 } else {
90 log << MSG::ERROR << "MrpcRec Get Calibration Service Failed !" << endreq;
91 return StatusCode::FAILURE;
92 }
93
94 //Get TOF Raw Data Provider Service
95 StatusCode scd = service("RawDataProviderSvc", tofDigiSvc);
96 if (scd == StatusCode::SUCCESS) {
97 log << MSG::INFO << "MrpcRec Get Tof Raw Data Service Sucessfully!" << endreq;
98 } else {
99 log << MSG::ERROR << "MrpcRec Get Tof Raw Data Service Failed !" << endreq;
100 return StatusCode::FAILURE;
101 }
102
103 if( m_checkDigi ) {
104 NTuplePtr nt11(ntupleSvc(),"FILE203/digi");
105 NTuplePtr nt12(ntupleSvc(),"FILE203/barrel");
106 NTuplePtr nt13(ntupleSvc(),"FILE203/endcap");
107 NTuplePtr nt14(ntupleSvc(),"FILE203/ext");
108 NTuplePtr nt15(ntupleSvc(),"FILE203/tof");
109 NTuplePtr nt16(ntupleSvc(),"FILE203/bb");
110
111 if( nt11 || nt12 || nt13 || nt14 || nt15 || nt16 ) {
112 m_tuple_digi = nt11;
113 m_tuple_barrel = nt12;
114 m_tuple_endcap = nt13;
115 m_tuple_ext = nt14;
116 m_tuple_tof = nt15;
117 m_tuple_bb = nt16;
118 }
119 else {
120 m_tuple_digi = ntupleSvc()->book("FILE203/digi",CLID_ColumnWiseTuple,"TofRec");
121 m_tuple_barrel = ntupleSvc()->book("FILE203/barrel",CLID_ColumnWiseTuple,"TofRec");
122 m_tuple_endcap = ntupleSvc()->book("FILE203/endcap",CLID_ColumnWiseTuple,"TofRec");
123 m_tuple_ext = ntupleSvc()->book("FILE203/ext",CLID_ColumnWiseTuple,"TofRec");
124 m_tuple_tof = ntupleSvc()->book("FILE203/tof",CLID_ColumnWiseTuple,"TofRec");
125 m_tuple_bb = ntupleSvc()->book("FILE203/bb",CLID_ColumnWiseTuple,"TofRec");
126
127 if( m_tuple_digi || m_tuple_barrel || m_tuple_endcap || m_tuple_ext || m_tuple_tof || m_tuple_bb ) {
128 m_checkdigi_tuple = new MrpcCheckDigi( m_tuple_digi, m_tuple_barrel, m_tuple_endcap, m_tuple_ext, m_tuple_tof, m_tuple_bb );
129 }
130 else {
131 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_digi) <<endmsg;
132 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_barrel) <<endmsg;
133 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_endcap) <<endmsg;
134 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_ext) <<endmsg;
135 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_tof) <<endmsg;
136 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_bb) <<endmsg;
137 return StatusCode::FAILURE;
138 }
139 }
140 }
141
142 if( m_saveRootFile ) {
143 NTuplePtr nt21(ntupleSvc(),"FILE201/trk");
144 NTuplePtr nt22(ntupleSvc(),"FILE201/cbtrk");
145 NTuplePtr nt23(ntupleSvc(),"FILE201/cetrk");
146
147 if( nt21 || nt22 || nt23 ) {
148 m_tuple_trk = nt21;
149 m_tuple_cbtrk = nt22;
150 m_tuple_cetrk = nt23;
151 }
152 else {
153 m_tuple_trk = ntupleSvc()->book("FILE201/trk",CLID_ColumnWiseTuple,"TofRec");
154 m_tuple_cbtrk = ntupleSvc()->book("FILE201/cbtrk",CLID_ColumnWiseTuple,"TofRec");
155 m_tuple_cetrk = ntupleSvc()->book("FILE201/cetrk",CLID_ColumnWiseTuple,"TofRec");
156 if( m_tuple_trk || m_tuple_cbtrk || m_tuple_cetrk ) {
157 m_checkdata_tuple = new MrpcCheckData( m_tuple_trk, m_tuple_cbtrk, m_tuple_cetrk );
158 }
159 else {
160 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_trk) <<endmsg;
161 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_cbtrk) <<endmsg;
162 log << MSG::ERROR <<"Cannot book N-tuple:" << long(m_tuple_cetrk) <<endmsg;
163 return StatusCode::FAILURE;
164 }
165 }
166 }
167
168 if( m_printOutInfo ) { m_printOut = new MrpcCount; }
169
170 return StatusCode::SUCCESS;
171}
172
173StatusCode MrpcRec::beginRun(){
174 MsgStream log(msgSvc(), name());
175 log << MSG::INFO <<"MrpcRec begin_run!"<< endreq;
176 return StatusCode::SUCCESS;
177}
178
179// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
180
181StatusCode MrpcRec::execute() {
182
183 MsgStream log(msgSvc(), name());
184 log << MSG::INFO << "MrpcRec in execute()!" << endreq;
185
186 SmartDataPtr<Event::EventHeader> eventHeader(eventSvc(),"/Event/EventHeader");
187 if( !eventHeader ) {
188 log << MSG::FATAL << "MrpcRec could not find Event Header!" << endreq;
189 return StatusCode::FAILURE;
190 }
191 int run = eventHeader->runNumber();
192 int event = eventHeader->eventNumber();
193 if( ( event % 1000 == 0 ) && m_printOutInfo ) {
194 std::cout << "run:" << run << " event: " << event << std::endl;
195 }
196 log << MSG::INFO << "run= " << run << " event= " << event << endreq;
197
198 MrpcRecTDS tds;
200
201 // Magnetic Field: Colliding data and Cosmic Ray
202 if( m_acceleratorStatus == "Colliding" ) {
203 // Retrieve Event Start Time
204 SmartDataPtr<RecEsTimeCol> estimeCol(eventSvc(),"/Event/Recon/RecEsTimeCol");
205 if( !estimeCol || ( estimeCol->size() == 0 ) ) {
206 log << MSG::WARNING << "MrpcRec Could not find RecEsTimeCol! Run = " << run << " Event = " << event << endreq;
207 return StatusCode::SUCCESS;
208 }
209 RecEsTimeCol::iterator iter_ESTime=estimeCol->begin();
210 double t0 = (*iter_ESTime)->getTest();
211 int t0Stat = (*iter_ESTime)->getStat();
212 log << MSG::INFO << "t0= " << t0 << " t0Stat= " << t0Stat << endreq;
213
214 // Kalman Filter Track
215 SmartDataPtr<RecMdcKalTrackCol> mdcKalTrackCol(eventSvc(),"/Event/Recon/RecMdcKalTrackCol");
216 if( !mdcKalTrackCol ) {
217 log << MSG::WARNING << "No MdcKalTrackCol in TDS! Run = " << run << " Event = " << event << endreq;
218 return StatusCode::SUCCESS;
219 }
220
221 // Tof Get Extrapolated Track
222 SmartDataPtr<RecExtTrackCol> extTrackCol(eventSvc(),"/Event/Recon/RecExtTrackCol");
223 if( !extTrackCol ) {
224 log << MSG::WARNING << "No ExtTrackCol in TDS! Run = " << run << " Event = " << event << endreq;
225 return StatusCode::SUCCESS;
226 }
227 else {
228 if( m_printOutInfo ) { m_printOut->setExtTrackNum( extTrackCol->size() ); }
229 if( m_checkDigi && m_checkDigiExt ) { m_checkdigi_tuple->FillCol( *eventHeader, mdcKalTrackCol, extTrackCol ); }
230 }
231
232 //sunss add 08/9/17
233 if( m_forCalibration ) { // Bhabha Events Selection
234 if( m_printOutInfo ) { m_printOut->addNumber(0); }
235 // if( t0Stat != 111 && t0Stat != 121 ) return StatusCode::SUCCESS;
236 if( t0Stat%10 != 1 ) return StatusCode::SUCCESS;
237 if( m_printOutInfo ) { m_printOut->addNumber(1); }
238
239 if( extTrackCol->size() != 2 ) return StatusCode::SUCCESS;
240 if( m_printOutInfo ) { m_printOut->addNumber(2); }
241
242 SmartDataPtr<RecMdcTrackCol> mdcTrackCol(eventSvc(),"/Event/Recon/RecMdcTrackCol");
243 if( !mdcTrackCol ) {
244 log << MSG::FATAL << "Could NOT find RecMdcTrackCol in TDS! Run = " << run << " Event = " << event << endreq;
245 return StatusCode::SUCCESS;
246 }
247 if( mdcTrackCol->size() != 2 ) return StatusCode::SUCCESS;
248 if( m_printOutInfo ) { m_printOut->addNumber(3); }
249
250 SmartDataPtr<RecEmcShowerCol> emcShowerCol(eventSvc(),"/Event/Recon/RecEmcShowerCol");
251 if( !emcShowerCol ) {
252 log << MSG::FATAL << "Could NOT find EmcRecShowerCol in TDS! Run = " << run << " Event = " << event << endreq;
253 return StatusCode::SUCCESS;
254 }
255 if( m_printOutInfo ) { m_printOut->addNumber(4); }
256
257 if( emcShowerCol->size() < 2 ) return StatusCode::SUCCESS;
258 if( m_printOutInfo ) { m_printOut->addNumber(5); }
259
260 RecMdcTrackCol::iterator iter_mdc1 = mdcTrackCol->begin();
261 RecMdcTrackCol::iterator iter_mdc2 = mdcTrackCol->begin() + 1;
262
263 RecMdcKalTrackCol::iterator iter_kal1 = mdcKalTrackCol->begin();
264 RecMdcKalTrackCol::iterator iter_kal2 = mdcKalTrackCol->begin() + 1;
265
266 RecExtTrackCol::iterator iter_ext1 = extTrackCol->begin();
267 RecExtTrackCol::iterator iter_ext2 = extTrackCol->begin() + 1;
268 Hep3Vector extPos1 = (*iter_ext1)->emcPosition();
269 Hep3Vector extPos2 = (*iter_ext2)->emcPosition();
270
271 RecEmcShowerCol::iterator iter_emc1 = emcShowerCol->begin();
272 RecEmcShowerCol::iterator iter_emc2 = emcShowerCol->begin() + 1;
273 Hep3Vector emcPos1((*iter_emc1)->x(),(*iter_emc1)->y(),(*iter_emc1)->z());
274 Hep3Vector emcPos2((*iter_emc2)->x(),(*iter_emc2)->y(),(*iter_emc2)->z());
275
276 Hep3Vector pep = (*iter_mdc1)->p3();
277 Hep3Vector pem = (*iter_mdc2)->p3();
278 double delta_angle = 180.0 - pep.angle( pem.unit() )*180.0/pi;
279 double delta_phi = abs( (*iter_mdc1)->phi() - (*iter_mdc2)->phi() )*180.0/pi;
280
281 Hep3Vector distant1 = extPos1 - emcPos1;
282 Hep3Vector distant2 = extPos2 - emcPos1;
283 if( distant1.r() > distant2.r() ) {
284 RecEmcShowerCol::iterator iter_tmp = iter_emc1;
285 iter_emc1 = iter_emc2;
286 iter_emc2 = iter_tmp;
287 Hep3Vector emc_tmp = emcPos1;
288 emcPos1 = emcPos2;
289 emcPos2 = emc_tmp;
290 }
291 distant1 = extPos1 - emcPos1;
292 distant2 = extPos2 - emcPos2;
293
294 double e1 = (*iter_emc1)->energy();
295 double e2 = (*iter_emc2)->energy();
296 double etot = 0.0;
297 RecEmcShowerCol::iterator iter_emc = emcShowerCol->begin();
298 for( ; iter_emc != emcShowerCol->end(); iter_emc++ ) {
299 etot += (*iter_emc)->energy();
300 }
301
302 if( m_checkDigi ) { m_checkdigi_tuple->FillCol( *eventHeader, extTrackCol, mdcTrackCol, emcShowerCol, mdcKalTrackCol ); }
303
304 if( ( (*iter_mdc1)->charge() + (*iter_mdc2)->charge() )!= 0 ) return StatusCode::SUCCESS;
305 if( m_printOutInfo ) { m_printOut->addNumber(6); }
306 // if( (*iter_mdc1)->ndof()<23.0 || (*iter_mdc2)->ndof()<23.0 ) return StatusCode::SUCCESS;
307 // if( ((*iter_mdc1)->chi2()/(*iter_mdc1)->ndof())>2.5 || ((*iter_mdc2)->chi2()/(*iter_mdc2)->ndof())>2.5 ) return StatusCode::SUCCESS;
308 if( delta_angle > 10.0 ) return StatusCode::SUCCESS;
309 if( m_printOutInfo ) { m_printOut->addNumber(7); }
310 // if( distant1.r()>6.0 || distant2.r()>6.0 ) return StatusCode::SUCCESS;
311 if( m_calibData == "Bhabha" ) {
312 if( (*iter_kal1)->getStat(1,0)!=0 || (*iter_kal2)->getStat(1,0)!=0 ) return StatusCode::SUCCESS;
313 if( m_printOutInfo ) { m_printOut->addNumber(8); }
314 // Jpsi data taken in 2009
315 if( m_data == "jpsi" ) {
316 if( distant1.r()>6.0 || distant2.r()>6.0 ) return StatusCode::SUCCESS;
317 if( m_printOutInfo ) { m_printOut->addNumber(9); }
318 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;
319 if( m_printOutInfo ) { m_printOut->addNumber(10); }
320 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;
321 if( m_printOutInfo ) { m_printOut->addNumber(11); }
322 if( delta_phi<174.0 || delta_phi>186.0 ) return StatusCode::SUCCESS;
323 if( m_printOutInfo ) { m_printOut->addNumber(12); }
324 if( e1 < 1.1 || e2 < 1.1 ) return StatusCode::SUCCESS;
325 if( m_printOutInfo ) { m_printOut->addNumber(13); }
326 }
327 // Psi prime data taken in 2009
328 else if( m_data == "psip" ) {
329 if( distant1.r()>6.0 || distant2.r()>6.0 ) return StatusCode::SUCCESS;
330 if( m_printOutInfo ) { m_printOut->addNumber(9); }
331 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;
332 if( m_printOutInfo ) { m_printOut->addNumber(10); }
333 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;
334 if( m_printOutInfo ) { m_printOut->addNumber(11); }
335 if( delta_phi<174.0 || delta_phi>183.0 ) return StatusCode::SUCCESS;
336 if( m_printOutInfo ) { m_printOut->addNumber(12); }
337 if( e1 < 1.4 || e2 < 1.4 ) return StatusCode::SUCCESS;
338 if( m_printOutInfo ) { m_printOut->addNumber(13); }
339 }
340 // Psi double prime data taken in first half of 2010
341 else if( m_data == "psipp10" ) {
342 if( distant1.r()>6.0 || distant2.r()>6.0 ) return StatusCode::SUCCESS;
343 if( m_printOutInfo ) { m_printOut->addNumber(9); }
344 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;
345 if( m_printOutInfo ) { m_printOut->addNumber(10); }
346 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;
347 if( m_printOutInfo ) { m_printOut->addNumber(11); }
348 if( delta_phi<174.0 || delta_phi>186.0 ) return StatusCode::SUCCESS;
349 if( m_printOutInfo ) { m_printOut->addNumber(12); }
350 if( e1 < 1.4 || e2 < 1.4 ) return StatusCode::SUCCESS;
351 if( m_printOutInfo ) { m_printOut->addNumber(13); }
352 }
353 // Psi double prime data taken in end of 2010 and 2011
354 else if( m_data == "psipp11" ) {
355 if( distant1.r()>6.0 || distant2.r()>6.0 ) return StatusCode::SUCCESS;
356 if( m_printOutInfo ) { m_printOut->addNumber(9); }
357 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;
358 if( m_printOutInfo ) { m_printOut->addNumber(10); }
359 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;
360 if( m_printOutInfo ) { m_printOut->addNumber(11); }
361 if( delta_phi<174.0 || delta_phi>184.0 ) return StatusCode::SUCCESS;
362 if( m_printOutInfo ) { m_printOut->addNumber(12); }
363 if( e1 < 1.4 || e2 < 1.4 ) return StatusCode::SUCCESS;
364 if( m_printOutInfo ) { m_printOut->addNumber(13); }
365 }
366 // Psi prime data taken in end of 2011 and 2012
367 else if( m_data == "psip12" ) {
368 if( distant1.r()>6.0 || distant2.r()>6.0 ) return StatusCode::SUCCESS;
369 if( m_printOutInfo ) { m_printOut->addNumber(9); }
370 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;
371 if( m_printOutInfo ) { m_printOut->addNumber(10); }
372 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;
373 if( m_printOutInfo ) { m_printOut->addNumber(11); }
374 if( delta_phi<172.0 || delta_phi>188.0 ) return StatusCode::SUCCESS;
375 if( m_printOutInfo ) { m_printOut->addNumber(12); }
376 if( e1 < 1.4 || e2 < 1.4 ) return StatusCode::SUCCESS;
377 if( m_printOutInfo ) { m_printOut->addNumber(13); }
378 }
379 // Jpsi data taken in 2012
380 else if( m_data == "jpsi12" ) {
381 if( distant1.r()>6.0 || distant2.r()>6.0 ) return StatusCode::SUCCESS;
382 if( m_printOutInfo ) { m_printOut->addNumber(9); }
383 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;
384 if( m_printOutInfo ) { m_printOut->addNumber(10); }
385 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;
386 if( m_printOutInfo ) { m_printOut->addNumber(11); }
387 if( delta_phi<172.0 || delta_phi>188.0 ) return StatusCode::SUCCESS;
388 if( m_printOutInfo ) { m_printOut->addNumber(12); }
389 if( e1 < 1.1 || e2 < 1.1 ) return StatusCode::SUCCESS;
390 if( m_printOutInfo ) { m_printOut->addNumber(13); }
391 }
392 if( ( etot - e1 - e2 ) > 0.3 ) return StatusCode::SUCCESS;
393 if( m_printOutInfo ) { m_printOut->addNumber(14); }
394 }
395 else if( m_calibData == "Dimu" ) {
396 if( (*iter_kal1)->getStat(1,1)!=0 || (*iter_kal2)->getStat(1,1)!=0 ) return StatusCode::SUCCESS;
397 if( m_printOutInfo ) { m_printOut->addNumber(8); }
398 if( e1 > 0.5 || e2 > 0.5 )return StatusCode::SUCCESS;
399 if( m_printOutInfo ) { m_printOut->addNumber(9); }
400 }
401
402 }
403 //sunss add 08/9/17
404
405 TofTrackVec* tofTrackVec = new TofTrackVec;
406 RecExtTrackCol::iterator iter_track = extTrackCol->begin();
407 for( ; iter_track < extTrackCol->end(); iter_track++ ) {
408 RecMdcKalTrackCol::iterator iter_kal = mdcKalTrackCol->begin();
409 for( ; iter_kal != mdcKalTrackCol->end(); iter_kal++ ) {
410 if( (*iter_kal)->trackId() == (*iter_track)->trackId() ) break;
411 }
412 int kal[5] = {-1};
413 if( iter_kal != mdcKalTrackCol->end() ) {
414 for( unsigned int i=0; i<5; i++ ) {
415 kal[i] = (*iter_kal)->getStat( 1, i );
416 }
417 }
418 MrpcTrack* tof = new MrpcTrack;
419 tof->setExtTrack( (*iter_track), kal, t0, t0Stat );
420
421 if( tofTrackVec->size()>0 ) {
422 std::vector<MrpcTrack*>::iterator iterExt = tofTrackVec->begin();
423 for( ; iterExt < tofTrackVec->end(); iterExt++ ) {
424 if( (*iterExt)->isNoHit() ) continue;
425 tof->getMultiHit( *iterExt );
426 }
427 }
428
429 tofTrackVec->push_back( tof );
430 }
431
432 if( m_printOutInfo ) {
433 m_printOut->setTrack1Col( tofTrackVec );
434 }
435
436 // Retrieve Tof Digi Data
437 TofDataMap tofDataMap = tofDigiSvc->tofDataMapTof( t0 );
438 if( tofDataMap.empty() ) {
439 log << MSG::WARNING << "No Tof Data Map in RawDataProviderSvc! Run=" << run << " Event=" << event << endreq;
440 }
441
442 if( m_checkDigi && m_checkDigiRaw ) {
443 SmartDataPtr<TofDigiCol> tofDigiCol(eventSvc(),"/Event/Digi/TofDigiCol");
444 if( !tofDigiCol ) {
445 log << MSG::ERROR << "MrpcRec could not find Tof digi! Event = " << event << endreq;
446 }
447 else { m_checkdigi_tuple->FillCol( *eventHeader, tofDigiCol, t0, t0Stat ); }
448
449 m_checkdigi_tuple->FillCol( *eventHeader, tofDataMap, t0, t0Stat );
450 }
451
452 std::vector<int> deadId;
453 if( m_forCalibration ) {
454 for( unsigned int i=0; i<5; i++ ) {
455 int identmp = tofCaliSvc->BrEast(i);
456 if( identmp != 0x2fffffff ) {
457 deadId.push_back( identmp );
458 }
459 identmp = tofCaliSvc->BrWest(i);
460 if( identmp != 0x2fffffff ) {
461 deadId.push_back( identmp );
462 }
463 identmp = tofCaliSvc->Endcap(i);
464 if( identmp != 0x2fffffff ) {
465 deadId.push_back( identmp );
466 }
467 }
468 }
469
470 std::vector<MrpcTrack*>::iterator iter = tofTrackVec->begin();
471 for( ; iter < tofTrackVec->end(); iter++ ) {
472 if( (*iter)->isNoHit() ) continue;
473 (*iter)->setTofData( tofDataMap,m_neighborhood );
474 if( m_printOutInfo ) { m_printOut->setTrack2( (*iter) ); }
475 if( (*iter)->isNoHit() ) continue;
476 (*iter)->match( m_forCalibration, deadId, tofTrackVec );
477 if( m_printOutInfo ) { m_printOut->setTrack3( (*iter) ); }
478 }
479
480 iter = tofTrackVec->begin();
481 for( ; iter < tofTrackVec->end(); iter++ ) {
482
483 (*iter)->setCalibration();
484
485 if( m_checkDigi ) {
486 if( m_checkTrigger ) {
487 // retrieve trigger data for physics analysis
488 SmartDataPtr<TrigData> trigData(eventSvc(), "/Event/Trig/TrigData");
489 if (!trigData) {
490 log << MSG::FATAL << "Could not find Trigger Data for physics analysis" << endreq;
491 m_checkdigi_tuple->Fill_TofTrack( *eventHeader, *iter, t0, t0Stat );
492 }
493 else {
494 m_checkdigi_tuple->Fill_TofTrack( *eventHeader, *iter, t0, t0Stat, trigData );
495 }
496 }
497 else {
498 if( ( run < 0 ) && m_checkMcTruth ) {
499 SmartDataPtr<TofMcHitCol> tofMcCol(eventSvc(),"/Event/MC/TofMcHitCol");
500 SmartDataPtr<McParticleCol> mcParticleCol(eventSvc(),"/Event/MC/McParticleCol");
501 if ( !tofMcCol || !mcParticleCol ) {
502 m_checkdigi_tuple->Fill_TofTrack( *eventHeader, *iter, t0, t0Stat, mdcKalTrackCol );
503 }
504 else {
505 m_checkdigi_tuple->Fill_TofTrack( *eventHeader, *iter, t0, t0Stat, mdcKalTrackCol, tofMcCol, mcParticleCol, m_calibData );
506 }
507 }
508 else {
509 m_checkdigi_tuple->Fill_TofTrack( *eventHeader, *iter, t0, t0Stat, mdcKalTrackCol );
510 }
511 }
512 }
513 }
514
515 tds.RegisterTDS( eventHeader->runNumber(), eventHeader->eventNumber(), tofTrackVec, m_forCalibration, m_calibData );
516
517 clearTofTrackVec( tofTrackVec );
518
519// Check RecTofTrackCol Registered
520 SmartDataPtr<RecTofTrackCol> tofTrackCol(eventSvc(),"/Event/Recon/RecTofTrackCol");
521 if (!tofTrackCol) {
522 log << MSG::FATAL << "MrpcRec could not find RecTofTrackCol!" << endreq;
523 return StatusCode::FAILURE;
524 }
525 else{
526 if( m_saveRootFile ) {
527 m_checkdata_tuple->FillCol( *eventHeader, tofTrackCol, mdcKalTrackCol );
528 }
529 }
530
531 if( m_forCalibration ) {
532 SmartDataPtr<RecBTofCalHitCol> bhitCol(eventSvc(),"/Event/Recon/RecBTofCalHitCol");
533 if (!bhitCol) {
534 log << MSG::WARNING << "MrpcRec could not find RecBTofCalHitCol!" << endreq;
535 }
536 else {
537 if( m_saveRootFile ) {
538 m_checkdata_tuple->FillCol( *eventHeader, bhitCol );
539 }
540 }
541
542 SmartDataPtr<RecETofCalHitCol> ehitCol(eventSvc(),"/Event/Recon/RecETofCalHitCol");
543 if (!ehitCol) {
544 log << MSG::WARNING << "MrpcRec could not find RecETofCalHitCol!" << endreq;
545 }
546 else {
547 if( m_saveRootFile ) {
548 m_checkdata_tuple->FillCol( *eventHeader, ehitCol );
549 }
550 }
551 }
552
553 }
554 else {
555 log << MSG::FATAL << "In MrpcRec: AcceleratorStatus is NOT correct! m_acceleratorStatus = " << m_acceleratorStatus << endreq;
556 return StatusCode::FAILURE;
557 }
558
559 return StatusCode::SUCCESS;
560
561}
562
563// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
564
565StatusCode MrpcRec::finalize() {
566 if( m_printOutInfo ) {
567 if( m_forCalibration ) {
568 m_printOut->finalBhabha( m_calibData );
569 }
570 m_printOut->final();
571 delete m_printOut;
572 }
573 if( m_checkDigi ) delete m_checkdigi_tuple;
574 if( m_saveRootFile ) delete m_checkdata_tuple;
575 return StatusCode::SUCCESS;
576}
577
578
579//----------------------------------------------------------------------------
580
581void MrpcRec::clearTofTrackVec( std::vector<MrpcTrack*>*& tofTrackVec) {
582 if( tofTrackVec ) {
583 std::vector<MrpcTrack*>::iterator iter = tofTrackVec->begin();
584 for( ; iter < tofTrackVec->end(); iter++ ) {
585 delete (*iter);
586 }
587 tofTrackVec->clear();
588 delete tofTrackVec;
589 }
590 return;
591}
Double_t etot
Double_t e1
Double_t e2
double abs(const EvtComplex &c)
Definition: EvtComplex.hh:212
EvtStreamInputIterator< typename Generator::result_type > iter(Generator gen, int N=0)
ITofGeomSvc * tofGeomSvc
Definition: MrpcRec.cxx:43
ITofCaliSvc * tofCaliSvc
Definition: MrpcRec.cxx:44
IRawDataProviderSvc * tofDigiSvc
Definition: MrpcRec.cxx:47
std::vector< MrpcTrack * > TofTrackVec
Definition: MrpcTrack.h:218
INTupleSvc * ntupleSvc()
IMessageSvc * msgSvc()
std::multimap< unsigned int, TofData * > TofDataMap
Definition: TofData.h:244
ITofCaliSvc * tofCaliSvc
virtual TofDataMap & tofDataMapTof(double estime=0)=0
virtual const int BrWest(unsigned int No)=0
virtual const int BrEast(unsigned int No)=0
virtual const int Endcap(unsigned int No)=0
void FillCol(Event::EventHeader &, RecTofTrackCol &, RecMdcKalTrackCol &)
void Fill_TofTrack(Event::EventHeader &, MrpcTrack *&, double, int)
void FillCol(Event::EventHeader &, TofDigiCol &, double, int)
void setExtTrackNum(unsigned int ntrk)
Definition: MrpcCount.cxx:122
void setTrack3(MrpcTrack *&tof)
Definition: MrpcCount.cxx:169
void setTrack2(MrpcTrack *&tof)
Definition: MrpcCount.cxx:152
void final()
Definition: MrpcCount.cxx:229
void setTrack1Col(std::vector< MrpcTrack * > *&tofTrackVec)
Definition: MrpcCount.cxx:143
void addNumber(unsigned int i)
Definition: MrpcCount.cxx:263
void finalBhabha(std::string calibData)
Definition: MrpcCount.cxx:271
StatusCode RegisterNullRecTofTrackCol()
Definition: MrpcRecTDS.cxx:22
StatusCode RegisterTDS(int runNumber, int eventNumber, std::vector< MrpcTrack * > *&tofTrackVec, bool m_forCalibration, std::string m_calibData)
Definition: MrpcRecTDS.cxx:58
StatusCode beginRun()
Definition: MrpcRec.cxx:173
MrpcRec(const std::string &name, ISvcLocator *pSvcLocator)
Definition: MrpcRec.cxx:49
StatusCode execute()
Definition: MrpcRec.cxx:181
StatusCode initialize()
Definition: MrpcRec.cxx:71
StatusCode finalize()
Definition: MrpcRec.cxx:565
void clearTofTrackVec(std::vector< MrpcTrack * > *&tofTrackVec)
Definition: MrpcRec.cxx:581
void getMultiHit(MrpcTrack *&)
Definition: MrpcTrack.cxx:344
void setExtTrack(RecExtTrack *extTrack, int kal[5], double t0, int t0Stat)
Definition: MrpcTrack.cxx:155
Definition: Event.h:21
const float pi
Definition: vector3.h:133