BOSS 7.0.2
BESIII Offline Software System
Loading...
Searching...
No Matches
TrkExtAlg.cxx
Go to the documentation of this file.
1//
2// File: TrkExtAlg.cxx
3// Author: Wang Liangliang
4// Date: 2005.4.4
5//
6#include "GaudiKernel/MsgStream.h"
7#include "GaudiKernel/AlgFactory.h"
8#include "GaudiKernel/ISvcLocator.h"
9#include "GaudiKernel/SmartDataPtr.h"
10#include "GaudiKernel/IDataProviderSvc.h"
11#include "GaudiKernel/IDataManagerSvc.h"
12#include "GaudiKernel/PropertyMgr.h"
13#include "EventModel/EventHeader.h"
14
15#include "CLHEP/Matrix/SymMatrix.h"
16#include "CLHEP/Vector/ThreeVector.h"
17
18#include "ExtEvent/RecExtTrack.h"
19#include "TrkExtAlg/TrkExtAlg.h"
20#include "TrkExtAlg/ExtSteppingAction.h"
21/////////////////////////////////change
22//#include "TrkExtAlg/Helix.h"
23#include "TrackUtil/Helix.h"
24#include "DetVerSvc/IDetVerSvc.h"
25//////////////////////////////////
26#include "TrkExtAlg/ExtMdcTrack.h"
27
28using namespace Event;
29
30string parName[5]={"e","mu","pi","kaon","proton"};
31
32//Constructor
33TrkExtAlg::TrkExtAlg(const string& name, ISvcLocator* pSvcLocator):Algorithm(name, pSvcLocator)
34{
35 myParticleName = "pi";
36 msgFlag = false;
37 myResultFlag = false;
38 myInputTrk= "Kal";
39
40 declareProperty("ParticleName",myParticleName);
41 declareProperty("GeantGeomOptimization",myGeomOptimization=true);
42 declareProperty("MessageFlag",msgFlag);
43 declareProperty("ResultMessageFlag",myResultFlag);
44 declareProperty("BFieldOn",myBFieldOn=true);
45 declareProperty("InputTrk",myInputTrk);
46 declareProperty("UseMucKalFilter",myUseMucKal=true);
47 declareProperty("MucWindow",myMucWindow=6);
48}
49
50//Destructor
51TrkExtAlg::~TrkExtAlg()
52{
53 if(myExtTrack) delete myExtTrack;
54}
55
56//////////////////////
57StatusCode TrkExtAlg::initialize()
58{
59 MsgStream log(msgSvc(), name());
60 log << MSG::INFO << "initialize()" << endreq;
61
62 IDetVerSvc* detVerSvc;
63 StatusCode sc_det = service("DetVerSvc", detVerSvc);
64 if( sc_det.isFailure() ) {
65 log << MSG::FATAL << "can't retrieve DetVerSvc instance" << endreq;
66 return sc_det;
67 }
68
69 m_detVer = detVerSvc->phase();
70 log << MSG::INFO << "** ~~~~~~ ** : retrieved DetectorStage = " << m_detVer << endreq;
71
72 myExtTrack = new Ext_track(msgFlag,myBFieldOn,myGeomOptimization,m_detVer,myUseMucKal,myMucWindow);
73 myExtTrack->Initialization(msgFlag,myBFieldOn,myGeomOptimization,myUseMucKal,myMucWindow);
74 // myFile = new ofstream("ExtData.txt");
75
76 /*
77 //--------- For Ext Test ----------------
78 NTuplePtr nt(ntupleSvc(),"FILE501/ext");
79 if ( nt ) myNtuple = nt;
80 else {
81 myNtuple=ntupleSvc()->book("FILE501/ext",CLID_ColumnWiseTuple,"TrkExt");
82 if(myNtuple) {
83 myNtuple->addItem("charge",myCharge);
84 myNtuple->addItem("emcHitFlag",myEmcHitFlag);
85 myNtuple->addItem("emcHitTheta",myEmcHitTheta);
86 myNtuple->addItem("emcHitPhi",myEmcHitPhi);
87 myNtuple->addItem("emcVolNum",myEmcVolNum);
88 myNtuple->addItem("emcExtTheta",myEmcExtTheta);
89 myNtuple->addItem("emcExtPhi",myEmcExtPhi);
90 myNtuple->addItem("dTheta",myDTheta);
91 myNtuple->addItem("dPhi",myDPhi);
92 }
93 else { // did not manage to book the N tuple....
94 log << MSG::ERROR <<"Cannot book N-tuple:" << long(myNtuple) << endmsg;
95 return StatusCode::FAILURE;
96 }
97 }
98 //-------- end Ext Test -----------------
99 */
100 return StatusCode::SUCCESS;
101}
102
103
104//////////////////////
105StatusCode TrkExtAlg::execute()
106{
107 //cout<<"a new event "<<endl;
108 MsgStream log(msgSvc(), name());
109 log << MSG::INFO << "execute()" << endreq;
110 int eventNumber, runNumber;
111
112 SmartDataPtr<Event::EventHeader> eventHeader(eventSvc(),"/Event/EventHeader");
113 if (!eventHeader)
114 {
115 log << MSG::FATAL << "Could not find Event Header" << endreq;
116 return( StatusCode::FAILURE);
117 }
118 runNumber = eventHeader->runNumber();
119 eventNumber = eventHeader->eventNumber();
120 // if(eventNumber!=5&&eventNumber!=8&&eventNumber!=142) return( StatusCode::SUCCESS);//1900736 //839313
121 if(msgFlag)
122 {
123 cout<<"TrackExt: ******************* Start a event *******************"<<endl;
124 cout<<"run= "<<runNumber<<"; event= "<<eventNumber<<endl;
125 }
126 // cout<<"run= "<<runNumber<<"; event= "<<eventNumber<<endl;
127 /////
128 SmartDataPtr<MucDigiCol> mucDigiCol(eventSvc(),"/Event/Digi/MucDigiCol");
129 if(!mucDigiCol) {
130 log << MSG::FATAL << "Could not find MUC digi" << endreq;
131 return( StatusCode::SUCCESS);
132 }
133 //cout<<"digi col is "<<mucDigiCol->size()<<endl;
134 ////
135 ExtMdcTrack aExtMdcTrack;
136 aExtMdcTrack.SetMsgFlag(msgFlag);
137
138 bool setTrk=false;
139
140 int parID;
141 if(myParticleName=="e") parID=0;
142 else if(myParticleName=="mu") parID=1;
143 else if(myParticleName=="pi") parID=2;
144 else if(myParticleName=="kaon") parID=3;
145 else if(myParticleName=="proton"||myParticleName=="anti_proton") parID=4;
146
147 if(myInputTrk=="Mdc")
148 {
149 SmartDataPtr<RecMdcTrackCol> aMdcTrackCol(eventSvc(),"/Event/Recon/RecMdcTrackCol");
150 if(!aMdcTrackCol)
151 {
152 log << MSG::WARNING << "Can't find RecMdcTrackCol in TDS!" << endreq;
153 return( StatusCode::SUCCESS);
154 }
155 setTrk=aExtMdcTrack.SetMdcRecTrkCol(aMdcTrackCol);
156
157 }
158 else if(myInputTrk=="Kal")
159 {
160 SmartDataPtr<RecMdcKalTrackCol> aMdcKalTrackCol(eventSvc(),"/Event/Recon/RecMdcKalTrackCol");
161 if(!aMdcKalTrackCol)
162 {
163 log << MSG::WARNING << "Can't find RecMdcKalTrackCol in TDS!" << endreq;
164 return( StatusCode::SUCCESS);
165 }
166 setTrk=aExtMdcTrack.SetMdcKalTrkCol(aMdcKalTrackCol);
167 }
168 else
169 {
170 log << MSG::WARNING << "Wrong type of inputTrk:" << myInputTrk << endreq;
171 return( StatusCode::SUCCESS);
172 }
173
174 RecExtTrackCol *aExtTrackCol = new RecExtTrackCol;
175 if(setTrk)
176 {
177 while(aExtMdcTrack.GetOneGoodTrk())
178 {
179 //cout<<"start a track........................"<<endl;
180 RecExtTrack *aExtTrack = new RecExtTrack;
181
182 for(int i=0; i<5; i++)// extrapolate one track with all 5 hypotheses
183 {
184 if(aExtMdcTrack.ReadTrk(i))
185 {
186 aExtTrack->SetParType(i);
187
188 int trackID = aExtMdcTrack.GetTrackID();
189 aExtTrack->SetTrackId(trackID);
190 Hep3Vector position = aExtMdcTrack.GetPosition();
191 Hep3Vector momentum = aExtMdcTrack.GetMomentum();
192 HepSymMatrix error = aExtMdcTrack.GetErrorMatrix();
193 double pathInMDC = aExtMdcTrack.GetTrackLength();
194 double tofInMdc = aExtMdcTrack.GetTrkTof();
195
196 if(msgFlag)
197 {
198 cout<<"Start From:"<<position.x()<<' '<<position.y()<<' '
199 <<position.z()<<endl;
200 cout<<"Start Momentum:"<<momentum.x()<<' '<<momentum.y()<<' '<<momentum.z()<<endl;
201 cout<<"Start Error matrix:"<<error<<endl;
202 cout<<"Path before start:"<< pathInMDC << endl;
203 }
204
205 G4String aParticleName(parName[i]);
206 double charge = aExtMdcTrack.GetParticleCharge();
207 if(!aParticleName.contains("proton"))
208 {
209 if(charge>0) aParticleName += "+";
210 else aParticleName += "-";
211 }
212 else
213 {
214 if(charge>0) aParticleName = "proton";
215 else aParticleName = "anti_proton";
216 }
217
218 if(msgFlag)
219 {
220 cout<<"Charge: "<<charge<<endl;
221 cout<<"Particle: "<<aParticleName<<endl;
222 }
223
224 ExtSteppingAction *extSteppingAction;
225 extSteppingAction = myExtTrack->GetStepAction();
226 extSteppingAction->Set_which_tof_version(m_detVer);//Set the ToFversionnumber. Required for the MRPC
227
228 extSteppingAction->Reset();
229 extSteppingAction->SetMucDigiColPointer(mucDigiCol);
230 extSteppingAction->SetExtTrackPointer(aExtTrack);
231 bool m_trackstatus = false;
232 int trk_startpart = 0;
233 while(!m_trackstatus)
234 {
235
236 trk_startpart++;//just for protection
237 if(trk_startpart>20)
238 {cout<<"-------has modified more than 20 times---------"<<endl;break;}
239 if(myExtTrack->Set(position,momentum,error,aParticleName,pathInMDC,tofInMdc))
240 {
241 myExtTrack->TrackExtrapotation();
242 extSteppingAction->InfmodMuc(position,momentum,error);
243 m_trackstatus = extSteppingAction->TrackStop();
244 }
245 else
246 m_trackstatus = true;
247 }
248 }
249
250 }
251
252 aExtTrack->SetParType(parID);
253
254 if(msgFlag) cout<<"will add aExtTrack!"<<endl;
255 if(aExtTrackCol)
256 {
257 if(aExtTrack) aExtTrackCol->add(aExtTrack);
258 else if(msgFlag) cout<<"No aExtTrack!"<<endl;
259 }
260 else
261 {
262 if(msgFlag) cout<<"No aExtTrackCol!"<<endl;
263 }
264 if(msgFlag) cout<<"add a aExtTrack!"<<endl;
265
266 /*
267 //For Test
268 if(myFile)
269 {
270 (*myFile)<<endPoint.x()<<' '<<endPoint.y()<<' '
271 <<endPoint.z()<<' '<<endMomentum.x()
272 <<' '<<endMomentum.y()<<' '<<endMomentum.z()
273 <<' '<<endErrorMatrix(1,1)<<' '<<endErrorMatrix(2,2)
274 <<' '<<endErrorMatrix(3,3)<<' '<<endErrorMatrix(4,4)
275 <<' '<<endErrorMatrix(5,5)<<' '<<endErrorMatrix(6,6)
276 <<endl;
277 }
278 else {
279 log << MSG::ERROR <<"can't open file" << endreq;
280
281 }
282 */
283 }//while
284 }//if
285
286 //Register ExtTrackCol to TDS.
287 /* ReconEvent *aReconEvent = new ReconEvent();
288 StatusCode sc = eventSvc()->registerObject("/Event/Recon",aReconEvent);
289 if(sc!=StatusCode::SUCCESS) {
290 log << MSG::FATAL << "Could not register ReconEvent" <<endreq;
291 return( StatusCode::FAILURE);
292 }
293 */
294
295 //cout<<"will write in to service "<<endl;
296 SmartIF<IDataManagerSvc> dataManSvc(eventSvc());
297 //IDataManagerSvc *dataManSvc = dynamic_cast<IDataManagerSvc*> (eventSvc());
298
299 DataObject *extTrackCol;
300 eventSvc()->findObject("/Event/Recon/RecExtTrackCol",extTrackCol);
301 if(extTrackCol != NULL) {
302 dataManSvc->clearSubTree("/Event/Recon/RecExtTrackCol");
303 eventSvc()->unregisterObject("/Event/Recon/RecExtTrackCol");
304 }
305
306 // cout<<"write register object "<<endl;
307
308 StatusCode sc = eventSvc()->registerObject("/Event/Recon/RecExtTrackCol", aExtTrackCol);
309 if(sc!=StatusCode::SUCCESS) {
310 log << MSG::FATAL << "Could not register RecExtTrackCol in TDS!" << endreq;
311 return( StatusCode::FAILURE);
312 }
313
314 //Check ExtTrackCol in TDS.
315 //cout<<"Check ExtTrackCol in TDS."<<endl;
316 SmartDataPtr<RecExtTrackCol> aExtTrkCol(eventSvc(),"/Event/Recon/RecExtTrackCol");
317 if (!aExtTrkCol) {
318 log << MSG::FATAL << "Can't find RecExtTrackCol in TDS!" << endreq;
319 return( StatusCode::FAILURE);
320 }
321
322 RecExtTrackCol::iterator iterOfExtTrk;
323 int j=1;
324 //cout<<"aExtTrkCol->size is "<<aExtTrkCol->size()<<" --------------"<<endl;
325 for(iterOfExtTrk = aExtTrkCol->begin();iterOfExtTrk!=aExtTrkCol->end();iterOfExtTrk++)
326 {
327 if(myResultFlag)
328 {
329 for(int i=0; i<5; i++)
330 {
331 //TOF information.
332 cout<<"##########track"<<j<<": "<<"("<<i<<")"<<endl;
333 cout<<"******TOF1:******"<<endl;
334 cout<<"VolumeName: "<<(*iterOfExtTrk)->tof1VolumeName(i)<<"\t"
335 <<"VolumeNumber: "<<(*iterOfExtTrk)->tof1VolumeNumber(i)<<"\t"<<endl
336 <<"Position: "<<(*iterOfExtTrk)->tof1Position(i)<<"\t"
337 <<"Momentum: "<<(*iterOfExtTrk)->tof1Momentum(i)<<"\t"<<endl
338 <<"Error matrix: "<<(*iterOfExtTrk)->tof1ErrorMatrix(i)
339 <<"Error z: "<<(*iterOfExtTrk)->tof1PosSigmaAlongZ(i)<<"\t"
340 <<"Error Tz: "<<(*iterOfExtTrk)->tof1PosSigmaAlongT(i)<<"\t"
341 <<"Error x: "<<(*iterOfExtTrk)->tof1PosSigmaAlongX(i)<<"\t"
342 <<"Error y: "<<(*iterOfExtTrk)->tof1PosSigmaAlongY(i)<<endl
343 <<"Tof: "<<(*iterOfExtTrk)->tof1(i)<<"\t"
344 <<"PathOF: "<<(*iterOfExtTrk)->tof1Path(i)
345 <<endl;
346 cout<<"******TOF2:******"<<endl;
347 cout<<"VolumeName: "<<(*iterOfExtTrk)->tof2VolumeName(i)<<"\t"
348 <<"VolumeNumber: "<<(*iterOfExtTrk)->tof2VolumeNumber(i)<<"\t"<<endl
349 <<"Position: "<<(*iterOfExtTrk)->tof2Position(i)<<"\t"
350 <<"Momentum: "<<(*iterOfExtTrk)->tof2Momentum(i)<<"\t"<<endl
351 <<"Error matrix: "<<(*iterOfExtTrk)->tof2ErrorMatrix(i)
352 <<"Error z: "<<(*iterOfExtTrk)->tof2PosSigmaAlongZ(i)<<"\t"
353 <<"Error Tz: "<<(*iterOfExtTrk)->tof2PosSigmaAlongT(i)<<"\t"
354 <<"Error x: "<<(*iterOfExtTrk)->tof2PosSigmaAlongX(i)<<"\t"
355 <<"Error y: "<<(*iterOfExtTrk)->tof2PosSigmaAlongY(i)<<endl
356 <<"Tof: "<<(*iterOfExtTrk)->tof2(i)<<"\t"
357 <<"PathOF: "<<(*iterOfExtTrk)->tof2Path(i)
358 <<endl;
359
360 //EMC information.
361 cout<<"******EMC:******"<<endl
362 <<"VolumeName: "<<(*iterOfExtTrk)->emcVolumeName(i)<<"\t"
363 <<"VolumeNumber: "<<(*iterOfExtTrk)->emcVolumeNumber(i)<<"\t"<<endl
364 <<"Position: "<<(*iterOfExtTrk)->emcPosition(i)<<"\t"
365 <<"Momentum: "<<(*iterOfExtTrk)->emcMomentum(i)<<"\t"<<endl
366 <<"Error matrix: "<<(*iterOfExtTrk)->emcErrorMatrix(i)
367 <<"Error theta: "<<(*iterOfExtTrk)->emcPosSigmaAlongTheta(i)<<"\t"
368 <<"Error phi: "<<(*iterOfExtTrk)->emcPosSigmaAlongPhi(i)<<"\t"
369 <<"EMC path: "<<(*iterOfExtTrk)->emcPath(i)
370 <<endl;
371
372 //MUC information
373 cout<<"******MUC:******"<<endl
374 <<"VolumeName: "<<(*iterOfExtTrk)->mucVolumeName(i)<<"\t"
375 <<"VolumeNumber: "<<(*iterOfExtTrk)->mucVolumeNumber(i)<<endl
376 <<"Position: "<<(*iterOfExtTrk)->mucPosition(i)<<"\t"
377 <<"Momentum: "<<(*iterOfExtTrk)->mucMomentum(i)<<"\t"<<endl
378 <<"Error matrix: "<<(*iterOfExtTrk)->mucErrorMatrix(i)
379 <<"Error z: "<<(*iterOfExtTrk)->mucPosSigmaAlongZ(i)<<"\t"
380 <<"Error Tz: "<<(*iterOfExtTrk)->mucPosSigmaAlongT(i)<<"\t"
381 <<"Error x: "<<(*iterOfExtTrk)->mucPosSigmaAlongX(i)<<"\t"
382 <<"Error y: "<<(*iterOfExtTrk)->mucPosSigmaAlongY(i)
383 <<endl;
384
385 cout<<"*******MUC KALMANFILTER***********"<<endl;
386 cout<<"Chisq is "<<(*iterOfExtTrk)->MucKalchi2(i)<<endl;
387 cout<<"Nfit is "<<(*iterOfExtTrk)->MucKaldof(i)<<endl;
388 cout<<"chiL "<<(*iterOfExtTrk)->MucKalchi2()<<endl;
389 // cout<<"Pull is "<<(*iterOfExtTrk)->mucKalPull()<<endl;
390
391 // cout<<"Residual is "<<(*iterOfExtTrk)->mucKalResidual()<<endl;
392 //Muc Ext hits information
393 /* ExtMucHitVec aExtMucHitVec = (*iterOfExtTrk)->GetExtMucHitVec();
394 int numOfMucHits = aExtMucHitVec.size();
395 cout<<"******MUC hits:"<<numOfMucHits<<"******"<<endl;
396 for(int j=0;j<numOfMucHits;j++)
397 {
398 cout<<"###Muc Hit "<<j<<":###"<<endl
399 <<"VolumeName: "<<aExtMucHitVec[j].GetVolumeName()<<"\t"
400 <<"VolumeNumber: "<<aExtMucHitVec[j].GetVolumeNumber()<<"\t"<<endl
401 <<"Position: "<<aExtMucHitVec[j].GetPosition()<<"\t"
402 <<"Momentum: "<<aExtMucHitVec[j].GetMomentum()<<"\t"<<endl
403 <<"Error z: "<<aExtMucHitVec[j].GetPosSigmaAlongZ()<<"\t"
404 <<"Error Tz: "<<aExtMucHitVec[j].GetPosSigmaAlongT()<<"\t"
405 <<"Error x: "<<aExtMucHitVec[j].GetPosSigmaAlongX()<<"\t"
406 <<"Error y: "<<aExtMucHitVec[j].GetPosSigmaAlongY()<<"\t"
407 <<endl;
408 }
409 */
410 }
411 }
412 j++;
413
414 } // loop ExtTrkCol
415
416 if(msgFlag) cout<<"****************** End a event! ****************"<<endl<<endl;
417
418 /*
419 //--------- For Ext Test ----------------
420 // Retrieve mc truth
421 SmartDataPtr<McParticleCol> mcParticleCol(eventSvc(),"/Event/MC/McParticleCol");
422 if (!mcParticleCol) {
423 log << MSG::FATAL << "Could not find McParticle" << endreq;
424 return( StatusCode::FAILURE);
425 }
426
427 int numOfTrack = mcParticleCol->size();
428 if(msgFlag) cout<< "numOfMcTrack: " << numOfTrack << endl;
429 if(numOfTrack!=2) return StatusCode::SUCCESS;
430
431 // Retrieve Emc Mc Hits
432 SmartDataPtr<EmcMcHitCol> emcMcHitCol(eventSvc(),"/Event/MC/EmcMcHitCol");
433 if (!emcMcHitCol) {
434 log << MSG::FATAL << "Could not find EMC truth" << endreq;
435 return( StatusCode::FAILURE);
436 }
437
438 McParticleCol::iterator iter_mc = mcParticleCol->begin();
439 EmcMcHitCol::iterator iterEmcBegin = emcMcHitCol->begin();
440 int numOfEmcHits = emcMcHitCol->size();
441
442 // ExtTrackCol *aExtTrackCol = new ExtTrackCol;
443
444 for (int i = 1;iter_mc != mcParticleCol->end(); iter_mc++, i++) {
445 bool flag = (*iter_mc)->primaryParticle();
446 if(!flag) continue;
447 int particleId = (*iter_mc)->particleProperty();
448 double charge = particleId/abs(particleId);
449 unsigned int sFlag = (*iter_mc)->statusFlags();
450 int trackIdx = (*iter_mc)->getTrackIndex();
451 HepPoint3D iPos = (*iter_mc)->initialPosition();
452 HepPoint3D fPos = (*iter_mc)->finalPosition();
453 HepLorentzVector iFMomentum = (*iter_mc)->initialFourMomentum();
454 HepLorentzVector fFMomentum = (*iter_mc)->finalFourMomentum();
455
456
457 bool emcHitFlag = false;
458 double thetaOfEmcHit = 0;
459 double phiOfEmcHit = 0;
460 Hep3Vector posOfEmcHit(0,0,0);
461 //Get Emc Truth
462 for(int j=0;j<numOfEmcHits;j++)
463 {
464 if(trackIdx==(*(iterEmcBegin+j))->getTrackIndex())
465 {
466 emcHitFlag = true;
467 double x = (*(iterEmcBegin+j))->getPositionX();
468 double y = (*(iterEmcBegin+j))->getPositionY();
469 double z = (*(iterEmcBegin+j))->getPositionZ();
470 Hep3Vector vec(x,y,z);
471 posOfEmcHit = vec;
472 thetaOfEmcHit = posOfEmcHit.theta();
473 phiOfEmcHit = posOfEmcHit.phi();
474 break;
475 }
476 }
477
478 ExtSteppingAction *extSteppingAction;
479 extSteppingAction = myExtTrack->GetStepAction();
480
481 G4String aParticleName(myParticleName);
482 if(!aParticleName.contains("proton"))
483 {
484 if(charge>0) aParticleName += "+";
485 else aParticleName += "-";
486 }
487 else
488 {
489 if(charge>0) aParticleName = "proton";
490 else aParticleName = "anti_proton";
491}
492
493ExtTrack *aExtTrack = new ExtTrack;
494aExtTrack->SetTrackID(trackIdx);
495
496Hep3Vector position(iPos);
497Hep3Vector momentum(iFMomentum.x(),iFMomentum.y(),iFMomentum.z());
498HepSymMatrix error(6,0);
499
500if(myExtTrack->Set(position,momentum,error,aParticleName,0))
501{
502 extSteppingAction->SetExtTrackPointer(aExtTrack);
503 myExtTrack->TrackExtrapotation();
504}
505
506Hep3Vector extEmcPos = aExtTrack->GetEmcPosition();
507double volumeNum = aExtTrack->GetEmcVolumeNumber();
508double thetaExt = extEmcPos.theta();
509double phiExt = extEmcPos.phi();
510
511if(myResultFlag)
512{
513 cout<< "*******Mc Track " << i <<" :******"<<endl;
514 cout<< "Parimary particle :" << flag <<" particleId = " << (*iter_mc)->particleProperty() << endl;
515 cout<< "Track Index: " << trackIdx << endl;
516 cout<< "initialPosition: "<< iPos.x() << "," << iPos.y() << "," << iPos.z() << endl;
517 cout<< "initialFourMomentum: " <<iFMomentum.x()<<","<<iFMomentum.y()<<","<<iFMomentum.z() << endl;
518 cout<<"Emc Truth:"<<emcHitFlag<<"!"<<posOfEmcHit<<endl;
519
520 cout<<"Ext Emc: "<<aExtTrack->GetEmcVolumeName()<<aExtTrack->GetEmcPosition()<<endl;
521}
522
523myCharge = charge;
524myEmcHitFlag = (emcHitFlag)? 1.:-1.;
525myEmcHitTheta = thetaOfEmcHit;
526myEmcHitPhi = phiOfEmcHit;
527myEmcVolNum = volumeNum;
528myEmcExtTheta = thetaExt;
529myEmcExtPhi = phiExt;
530myDTheta = myEmcHitTheta-myEmcExtTheta;
531myDPhi = myEmcHitPhi-myEmcExtPhi;
532while(myDTheta<-1*M_PI) myDTheta+=2.0*M_PI;
533while(myDPhi<-1*M_PI) myDPhi+=2.0*M_PI;
534while(myDPhi>M_PI) myDPhi-=2.0*M_PI;
535while(myDTheta>M_PI) myDTheta-=2.0*M_PI;
536myNtuple->write();
537
538if(aExtTrack) delete aExtTrack;
539}
540
541//--------- end Ext Test ----------------
542*/
543
544return StatusCode::SUCCESS;
545}
546
547/////////////////////
548StatusCode TrkExtAlg::finalize()
549{
550 MsgStream log(msgSvc(), name());
551 log << MSG::INFO << "finalize()" << endreq;
552
553 // delete myExtTrack;
554 // myFile->close();
555
556 return StatusCode::SUCCESS;
557}
momentum
**********INTEGER nmxhep !maximum number of particles DOUBLE PRECISION vhep INTEGER jdahep COMMON hepevt $ !serial number $ !number of particles $ !status code $ !particle ident KF $ !parent particles $ !childreen particles $ !four momentum
Definition: BesBdkRc/BesBdkRc-00-00-02/src/fortran/HepEvt.h:23
RecExtTrackCol
ObjectVector< RecExtTrack > RecExtTrackCol
Definition: Event/ExtEvent/ExtEvent-00-00-13/ExtEvent/RecExtTrack.h:209
detVerSvc
IDetVerSvc * detVerSvc
Definition: TofRawDataProvider.cxx:31
parName
string parName[5]
Definition: TrkExtAlg.cxx:30
ExtMdcTrack::GetOneGoodTrk
bool GetOneGoodTrk()
Definition: ExtMdcTrack.cxx:51
ExtMdcTrack::GetParticleCharge
double GetParticleCharge() const
Definition: ExtMdcTrack.cxx:356
ExtMdcTrack::GetTrackLength
double GetTrackLength() const
Definition: ExtMdcTrack.cxx:309
ExtMdcTrack::SetMdcKalTrkCol
bool SetMdcKalTrkCol(RecMdcKalTrackCol *aPointer)
Definition: ExtMdcTrack.cxx:42
ExtMdcTrack::SetMdcRecTrkCol
bool SetMdcRecTrkCol(RecMdcTrackCol *aPointer)
Definition: ExtMdcTrack.cxx:33
ExtMdcTrack::GetPosition
const Hep3Vector GetPosition() const
Definition: ExtMdcTrack.cxx:227
ExtMdcTrack::ReadTrk
bool ReadTrk(int pid)
Definition: ExtMdcTrack.cxx:77
ExtMdcTrack::GetErrorMatrix
const HepSymMatrix GetErrorMatrix() const
Definition: ExtMdcTrack.cxx:251
ExtMdcTrack::GetMomentum
const Hep3Vector GetMomentum() const
Definition: ExtMdcTrack.cxx:242
ExtSteppingAction::InfmodMuc
void InfmodMuc(Hep3Vector &pos, Hep3Vector &mom, HepSymMatrix &err)
Definition: ExtSteppingAction.cxx:1018
Ext_track
Definition: InstallArea/include/TrkExtAlg/TrkExtAlg/Ext_track.h:40
Ext_track::TrackExtrapotation
void TrackExtrapotation()
Definition: Ext_track.cxx:343
Ext_track::Initialization
void Initialization(const bool aMsgFlag, const bool Bfield, const bool GeomOptimization, const bool aUseMucKal, const int aMucWindow)
Definition: Ext_track.cxx:90
Ext_track::GetStepAction
ExtSteppingAction * GetStepAction()
Definition: InstallArea/include/TrkExtAlg/TrkExtAlg/Ext_track.h:90
Ext_track::Set
bool Set(const Hep3Vector &xv3, const Hep3Vector &pv3, const HepSymMatrix &err, const std::string &particleName, const double pathInMDC, const double tofInMdc)
Definition: Ext_track.cxx:235
IDetVerSvc::phase
virtual int phase()=0
TrkExtAlg::finalize
StatusCode finalize()
Definition: TrkExtAlg.cxx:548
TrkExtAlg::execute
StatusCode execute()
Definition: TrkExtAlg.cxx:105
TrkExtAlg::TrkExtAlg
TrkExtAlg(const std::string &name, ISvcLocator *pSvcLocator)
TrkExtAlg::initialize
StatusCode initialize()
Definition: TrkExtAlg.cxx:57