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