From 52ac3921df0f46f12eb36f5063d50b612fa4fc5a Mon Sep 17 00:00:00 2001
From: Chengdong Fu <fucd@ihep.ac.cn>
Date: Fri, 24 Jul 2020 19:07:24 +0800
Subject: [PATCH] add FTD/SIT support
---
.../src/SiliconTrackingAlg.cpp | 145 ++++++++----------
.../SiliconTracking/src/SiliconTrackingAlg.h | 15 +-
Service/TrackSystemSvc/CMakeLists.txt | 5 +-
Service/TrackSystemSvc/src/HelixTrack.cc | 15 +-
Service/TrackSystemSvc/src/MarlinKalTest.cc | 18 +--
.../TrackSystemSvc/src/MarlinKalTestTrack.cc | 43 ++----
Service/TrackSystemSvc/src/MarlinTrkUtils.cc | 122 ++++++---------
7 files changed, 144 insertions(+), 219 deletions(-)
diff --git a/Reconstruction/SiliconTracking/src/SiliconTrackingAlg.cpp b/Reconstruction/SiliconTracking/src/SiliconTrackingAlg.cpp
index 156d9afd..74ed5f22 100644
--- a/Reconstruction/SiliconTracking/src/SiliconTrackingAlg.cpp
+++ b/Reconstruction/SiliconTracking/src/SiliconTrackingAlg.cpp
@@ -2,6 +2,7 @@
#include "GearSvc/IGearSvc.h"
#include "EventSeeder/IEventSeeder.h"
#include "TrackSystemSvc/ITrackSystemSvc.h"
+#include "DataHelper/Navagation.h"
#include "edm4hep/MCParticle.h"
#include "edm4hep/TrackerHit.h"
//#include "edm4hep/TrackerHitPlane.h"
@@ -83,7 +84,8 @@ SiliconTrackingAlg::SiliconTrackingAlg(const std::string& name, ISvcLocator* svc
declareProperty("FTDPixelHitCollection", _inFTDPixelColHdl, "Handle of the Input FTD TrackerHits collection");
declareProperty("FTDSpacePointCollection", _inFTDSpacePointColHdl, "Handle of the Input FTD SpacePoints collection");
declareProperty("SITHitCollection", _inSITColHdl, "Handle of the Input SIT TrackerHits collection");
-
+ declareProperty("SITRawHitCollection", _inSITRawColHdl, "Handle of the SIT SpacePoints raw hit collection collection");
+ declareProperty("FTDRawHitCollection", _inFTDRawColHdl, "Handle of the FTD SpacePoints raw hit collection collection");
// Output Collections
// ^^^^^^^^^^^^^^^^^^
declareProperty("SiTrackCollection", _outColHdl, "Handle of the SiTrack output collection");
@@ -92,7 +94,6 @@ SiliconTrackingAlg::SiliconTrackingAlg(const std::string& name, ISvcLocator* svc
// ^^^^^^^^^^^^^^^^^^^
_output_track_col_quality = _output_track_col_quality_GOOD;
-
}
@@ -120,8 +121,8 @@ StatusCode SiliconTrackingAlg::initialize() {
_trksystem->setOption( IMarlinTrkSystem::CFG::useQMS, _MSOn ) ;
_trksystem->setOption( IMarlinTrkSystem::CFG::usedEdx, _ElossOn) ;
_trksystem->setOption( IMarlinTrkSystem::CFG::useSmoothing, _SmoothOn) ;
- _trksystem->init() ;
- std::cout << "fucd ==============" << _trksystem << std::endl;
+ _trksystem->init() ;
+ debug() << "TrackSystem builded at address=" << _trksystem << endmsg;
#ifdef MARLINTRK_DIAGNOSTICS_ON
void * dcv = _trksystem->getDiagnositicsPointer();
@@ -155,7 +156,7 @@ StatusCode SiliconTrackingAlg::initialize() {
}
StatusCode SiliconTrackingAlg::execute(){
-
+ Navagation::Instance()->Initialize();
//_current_event = evt;
//_allHits.reserve(1000);
@@ -177,7 +178,6 @@ StatusCode SiliconTrackingAlg::execute(){
//debug() << "Processing Run[" << runNo << "]::Event[" << evtNo << "]" << endmsg;
_trackImplVec.reserve(100);
- _allHits.reserve(1000);
int successVTX = InitialiseVTX();
int successFTD = 0;
@@ -240,7 +240,6 @@ StatusCode SiliconTrackingAlg::execute(){
//edm4hep::LCRelationCollection* relCol = nullptr;
auto trkCol = _outColHdl.createAndPut();
//auto relCol = _outRelColHdl.createAndPut();
- //std::cout << "fucd------------------" << std::endl;
/*
LCCollectionVec * trkCol = new LCCollectionVec(LCIO::TRACK);
// if we want to point back to the hits we need to set the flag
@@ -252,7 +251,7 @@ StatusCode SiliconTrackingAlg::execute(){
*/
//FinalRefit(trkCol, relCol);
FinalRefit(trkCol);
- //std::cout << "fucd------------------" << std::endl;
+ debug() << "FinalRefit finish" << endmsg;
// set the quality of the output collection
switch (_output_track_col_quality) {
@@ -402,7 +401,6 @@ void SiliconTrackingAlg::CleanUp() {
}
}
_trackImplVec.clear();
- _allHits.clear();
}
int SiliconTrackingAlg::InitialiseFTD() {
@@ -451,7 +449,7 @@ int SiliconTrackingAlg::InitialiseFTD() {
// V must be the global z axis
if( fabs(V.dot(Z)) > eps ) {
error() << "SiliconTrackingAlg: VXD Hit measurment vectors V is not in the global X-Y plane. \n\n exit(1) called from file " << __FILE__ << " and line " << __LINE__ << endmsg;
- exit(1);
+ exit(1);
}
if( fabs(U.dot(Z)) > eps ) {
@@ -499,7 +497,7 @@ int SiliconTrackingAlg::InitialiseFTD() {
}
if (layer >= _nlayersFTD) {
- error() << "SiliconTrackingAlg => fatal error in FTD : layer is outside allowed range : " << layer << " number of layers = " << _nlayersFTD << endmsg;
+ fatal() << "SiliconTrackingAlg => fatal error in FTD : layer is outside allowed range : " << layer << " number of layers = " << _nlayersFTD << endmsg;
exit(1);
}
@@ -528,7 +526,19 @@ int SiliconTrackingAlg::InitialiseFTD() {
success = 0;
}
+ const edm4hep::TrackerHitCollection* rawHitCol = nullptr;
if(hitFTDSpacePointCol){
+ try{
+ rawHitCol = _inFTDRawColHdl.get();
+ }
+ catch ( GaudiException &e ) {
+ fatal() << "Collection " << _inFTDRawColHdl.fullKey() << " is unavailable in event " << _nEvt << endmsg;
+ success = 0;
+ }
+ }
+
+ if(hitFTDSpacePointCol&&rawHitCol){
+ Navagation::Instance()->AddTrackerHitCollection(rawHitCol);
//LCCollection * hitCollection = evt->getCollection(_FTDSpacePointCollection.c_str());
//_colNamesTrackerHits[hitCollection] = _FTDSpacePointCollection;
@@ -625,12 +635,10 @@ int SiliconTrackingAlg::InitialiseFTD() {
}
int SiliconTrackingAlg::InitialiseVTX() {
- //std::cout << "fucd================" << std::endl;
_nTotalVTXHits = 0;
_nTotalSITHits = 0;
_sectors.clear();
_sectors.resize(_nLayers*_nDivisionsInPhi*_nDivisionsInTheta);
- //std::cout << "fucd================" << std::endl;
int success = 1;
// Reading out VTX Hits Collection
//^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
@@ -647,16 +655,13 @@ int SiliconTrackingAlg::InitialiseVTX() {
//_colNamesTrackerHits[hitCollection] = _VTXHitCollection;
//_colTrackerHits.push_back(hitCollection);
- //std::cout << "fucd================1" << std::endl;
int nelem = hitVTXCol->size();
- //std::cout << "fucd================2" << std::endl;
debug() << "Number of VTX hits = " << nelem << endmsg;
_nTotalVTXHits = nelem;
for (int ielem=0; ielem<nelem; ++ielem) {
//for(auto hit : *hitVTXCol){
edm4hep::TrackerHit hit = hitVTXCol->at(ielem);
- //_allHits.push_back(hit);
//gear::Vector3D U(1.0,hit->getU()[1],hit->getU()[0],gear::Vector3D::spherical);
//gear::Vector3D V(1.0,hit->getV()[1],hit->getV()[0],gear::Vector3D::spherical);
gear::Vector3D U(1.0,hit.getCovMatrix()[1],hit.getCovMatrix()[0],gear::Vector3D::spherical);
@@ -674,11 +679,8 @@ int SiliconTrackingAlg::InitialiseVTX() {
error() << "SiliconTrackingAlg: VXD Hit measurment vectors U is not in the global X-Y plane. \n\n exit(1) called from file " << __FILE__ << " and line " << __LINE__ << endmsg;
exit(1);
}
- //std::cout << "fucd: " << &hit << " " << &_allHits.back() << std::endl;
- //TrackerHitExtended * hitExt = new TrackerHitExtended( &_allHits.back() );
TrackerHitExtended * hitExt = new TrackerHitExtended(hit);
- std::cout << "Saved TrackerHit pointer in TrackerHitExtended " << ielem << ": " << hitExt->getTrackerHit() << std::endl;
- //std::cout << (&_allHits.back())->getPosition()[0] << " " << hit.getCovMatrix()[2] << " " << hit.getCovMatrix()[5] << std::endl;
+ debug() << "Saved TrackerHit pointer in TrackerHitExtended " << ielem << ": " << hitExt->getTrackerHit() << std::endl;
// SJA:FIXME: just use planar res for now
hitExt->setResolutionRPhi(hit.getCovMatrix()[2]);
@@ -731,12 +733,21 @@ int SiliconTrackingAlg::InitialiseVTX() {
debug() << "Collection " << _inSITColHdl.fullKey() << " is unavailable in event " << _nEvt << endmsg;
success = 0;
}
+
+ const edm4hep::TrackerHitCollection* rawHitCol = nullptr;
+ if(hitSITCol){
+ try{
+ rawHitCol = _inSITRawColHdl.get();
+ }
+ catch ( GaudiException &e ) {
+ warning() << "Collection " << _inSITRawColHdl.fullKey() << " is unavailable in event " << _nEvt << ", if SIT is Space Point, it needed " << endmsg;
+ //success = 0;
+ }
+ }
+
if(hitSITCol){
- //LCCollection *hitCollection = evt->getCollection(_SITHitCollection.c_str());
-
- //_colNamesTrackerHits[hitCollection] = _SITHitCollection;
- //_colTrackerHits.push_back(hitCollection);
-
+ if(rawHitCol) Navagation::Instance()->AddTrackerHitCollection(rawHitCol);
+ //debug() << "SITHitCollection pointer = " << hitSITCol << endmsg;
int nelem = hitSITCol->size();
debug() << "Number of SIT hits = " << nelem << endmsg;
@@ -774,13 +785,11 @@ int SiliconTrackingAlg::InitialiseVTX() {
error() << "SiliconTrackingAlg => fatal error in SIT : layer is outside allowed range : " << layer << endmsg;
exit(1);
}
-
+
// first check that we have not been given 1D hits by mistake, as they won't work here
if ( UTIL::BitSet32( trkhit.getType() )[ UTIL::ILDTrkHitTypeBit::ONE_DIMENSIONAL ] ) {
-
- error() << "SiliconTrackingAlg: SIT Hit cannot be of type UTIL::ILDTrkHitTypeBit::ONE_DIMENSIONAL COMPOSITE SPACEPOINTS must be use instead. \n\n exit(1) called from file " << __FILE__ << " and line " << __LINE__ << endmsg;
- exit(1);
-
+ fatal() << "SiliconTrackingAlg: SIT Hit cannot be of type UTIL::ILDTrkHitTypeBit::ONE_DIMENSIONAL COMPOSITE SPACEPOINTS must be use instead. \n\n exit(1) called from file " << __FILE__ << " and line " << __LINE__ << endmsg;
+ exit(1);
}
// most likely case: COMPOSITE_SPACEPOINT hits formed from stereo strip hits
else if ( UTIL::BitSet32( trkhit.getType() )[ UTIL::ILDTrkHitTypeBit::COMPOSITE_SPACEPOINT ] ) {
@@ -1042,11 +1051,9 @@ void SiliconTrackingAlg::ProcessOneSector(int iPhi, int iTheta) {
for (int iInner=0;iInner<nHitsInner;iInner++) { // loop over hits in the inner sector
TrackerHitExtended * innerHit = hitVecInner[iInner];
HelixClass helix;
- //std::cout <<"fucd++++++++++++++++++++++1" << std::endl;
- // test fit to triplet
+ // test fit to triplet
TrackExtended * trackAR = TestTriplet(outerHit,middleHit,innerHit,helix);
- //std::cout <<"fucd++++++++++++++++++++++2" << std::endl;
- if ( trackAR != NULL ) {
+ if ( trackAR != NULL ) {
int nHits = BuildTrack(outerHit,middleHit,innerHit,helix,nLR[2],
iPhiLowInner,iPhiUpInner,
iThetaLowInner,iThetaUpInner,trackAR);
@@ -1055,7 +1062,6 @@ void SiliconTrackingAlg::ProcessOneSector(int iPhi, int iTheta) {
counter ++ ;
}
- //std::cout <<"fucd++++++++++++++++++++++3" << std::endl;
} // endloop over hits in the inner sector
} // endloop over hits in the middle sector
} // endloop over hits in the outer sector
@@ -1078,13 +1084,11 @@ TrackExtended * SiliconTrackingAlg::TestTriplet(TrackerHitExtended * outerHit,
/*
Methods checks if the triplet of hits satisfies helix hypothesis
*/
- //std::cout << "fucd================1" << std::endl;
// get the tracks already associated with the triplet
TrackExtendedVec& trackOuterVec = outerHit->getTrackExtendedVec();
TrackExtendedVec& trackMiddleVec = middleHit->getTrackExtendedVec();
TrackExtendedVec& trackInnerVec = innerHit->getTrackExtendedVec();
- //std::cout << "fucd================2" << std::endl;
// check if all the hits are already assigned to a track
if ( (!trackOuterVec.empty()) && (!trackMiddleVec.empty()) && (!trackInnerVec.empty())) {
@@ -1124,13 +1128,10 @@ TrackExtended * SiliconTrackingAlg::TestTriplet(TrackerHitExtended * outerHit,
}// for outer
}// for middle
}// if all vectors are not empty
- //std::cout << "fucd================3" << std::endl;
-
// float dZ = FastTripletCheck(innerHit, middleHit, outerHit);
// if (fabs(dZ) > _minDistCutAttach)
// return trackAR;
-
// increase triplet count
++_ntriplets;
@@ -1146,11 +1147,7 @@ TrackExtended * SiliconTrackingAlg::TestTriplet(TrackerHitExtended * outerHit,
float par[5];
float epar[15];
- /*fucd
- for(int i=0;i<_allHits.size();i++){
- std::cout << &_allHits.at(i) << std::endl;
- }
- */
+
// first hit
xh[0] = outerHit->getTrackerHit()->getPosition()[0];
yh[0] = outerHit->getTrackerHit()->getPosition()[1];
@@ -1455,13 +1452,10 @@ int SiliconTrackingAlg::BuildTrack(TrackerHitExtended * outerHit,
for (int layer = innerLayer-1; layer>=0; layer--) { // loop over remaining layers
float distMin = 1.0e+20;
TrackerHitExtended * assignedhit = NULL;
- //std::cout << "fucd---------------------1" << std::endl;
// loop over phi in the Inner region
for (int ipInner=iPhiLow; ipInner<iPhiUp+1;ipInner++) {
-
// loop over theta in the Inner region
for (int itInner=iThetaLow; itInner<iThetaUp+1;itInner++) {
-
int iPhiInner = ipInner;
// catch wrap-around
@@ -1506,12 +1500,9 @@ int SiliconTrackingAlg::BuildTrack(TrackerHitExtended * outerHit,
} // endloop over hits in the Inner sector
} // endloop over theta in the Inner region
} // endloop over phi in the Inner region
- //std::cout << "fucd---------------------2" << std::endl;
// check if closest hit fulfills the min distance cut
if (distMin < _minDistCutAttach) {
-
// if yes try to include it in the fit
-
TrackerHitExtendedVec& hvec = trackAR->getTrackerHitExtendedVec();
int nHits = int(hvec.size());
double * xh = new double[nHits+1];
@@ -1552,7 +1543,7 @@ int SiliconTrackingAlg::BuildTrack(TrackerHitExtended * outerHit,
float chi2RPhi;
float chi2Z;
-// std::cout << "######## number of hits to fit with _fastfitter = " << NPT << endmsg;
+ //debug() << "######## number of hits to fit with _fastfitter = " << NPT << endmsg;
_fastfitter->fastHelixFit(NPT, xh, yh, rh, ph, wrh, zh, wzh,iopt, par, epar, chi2RPhi, chi2Z);
par[3] = par[3]*par[0]/fabs(par[0]);
@@ -1602,14 +1593,10 @@ int SiliconTrackingAlg::BuildTrack(TrackerHitExtended * outerHit,
}
} // endloop over remaining layers
- //std::cout << "fucd---------------------3" << std::endl;
TrackerHitExtendedVec& hvec = trackAR->getTrackerHitExtendedVec();
int nTotalHits = int(hvec.size());
-
-// std::cout << "######## number of hits to return = " << nTotalHits << endmsg;
-
+ debug() << "######## number of hits to return = " << nTotalHits << endmsg;
return nTotalHits;
-
}
@@ -1822,7 +1809,7 @@ void SiliconTrackingAlg::CreateTrack(TrackExtended * trackAR ) {
trackOld->setD0(d0);
trackOld->setZ0(z0);
- // std::cout << "Split track found " << d0 << " " << z0 << endmsg;
+ //debug() << "Split track found " << d0 << " " << z0 << endmsg;
// killeb: In the original SiliconTrackingAlg this was in the NOT simple helix branch.
// The rest of the code uses the simple helix branch, where ndf_D is never set.
@@ -2559,14 +2546,14 @@ void SiliconTrackingAlg::FinalRefit(edm4hep::TrackCollection* trk_col) {
float pxTot = 0.;
float pyTot = 0.;
float pzTot = 0.;
- std::cout << "fucd============" << nTracks << std::endl;
+ debug() << "Total " << nTracks << " candidate tracks will be dealed" << endmsg;
for (int iTrk=0;iTrk<nTracks;++iTrk) {
TrackExtended * trackAR = _trackImplVec[iTrk];
TrackerHitExtendedVec& hitVec = trackAR->getTrackerHitExtendedVec();
int nHits = int(hitVec.size());
- std::cout << "fucd-------------" << iTrk << ": " << nHits << std::endl;
+ //debug() << " Track " << iTrk << " has hits: " << nHits << endmsg;
if( nHits >= _minimalHits) {
// int * lh = new int[nHits];
std::vector<int> lh;
@@ -2586,7 +2573,6 @@ void SiliconTrackingAlg::FinalRefit(edm4hep::TrackCollection* trk_col) {
helix->Initialize_Canonical(phi0, d0, z0, omega,
tanlambda, _bField);
-
// get the point of closest approach to the reference point
// here it is implicitly assumed that the reference point is the origin
float Pos[3];
@@ -2594,7 +2580,6 @@ void SiliconTrackingAlg::FinalRefit(edm4hep::TrackCollection* trk_col) {
Pos[1] = d0*cos(phi0);
Pos[2] = z0;
- std::cout << "fucd------------------1" << std::endl;
// at this point is is possible to have hits from the same layer ...
// so a check is made to ensure that the hit with the smallest distance to the
// current helix hypothosis is used, the other hit has lh set to 0
@@ -2689,6 +2674,7 @@ void SiliconTrackingAlg::FinalRefit(edm4hep::TrackCollection* trk_col) {
// check if the hit has been rejected as being on the same layer and further from the helix lh==0
if (lh[i] == 1) {
edm4hep::TrackerHit * trkHit = hitVec[i]->getTrackerHit();
+ debug() << "TrackerHit " << i << " address = " << trkHit << endmsg;
nFit++;
if(trkHit) {
trkHits.push_back(trkHit);
@@ -2707,12 +2693,10 @@ void SiliconTrackingAlg::FinalRefit(edm4hep::TrackCollection* trk_col) {
debug() << "SiliconTrackingAlg::FinalRefit: Cannot fit less than 3 hits. Number of hits = " << trkHits.size() << endmsg;
continue ;
}
- std::cout << "fucd------------------2" << std::endl;
//TrackImpl* Track = new TrackImpl ;
auto track = trk_col->create();
//fucd
//edm4hep::Track track;// = new edm4hep::Track;
- std::cout << "fucd------------------3" << std::endl;
// setup initial dummy covariance matrix
//std::vector<float> covMatrix;
//covMatrix.resize(15);
@@ -2739,7 +2723,7 @@ void SiliconTrackingAlg::FinalRefit(edm4hep::TrackCollection* trk_col) {
}
sort(r2_values.begin(),r2_values.end());
- //std::cout << "fucd------------------3" << std::endl;
+
trkHits.clear();
trkHits.reserve(r2_values.size());
@@ -2782,18 +2766,16 @@ void SiliconTrackingAlg::FinalRefit(edm4hep::TrackCollection* trk_col) {
}
int status = 0;
- std::cout << "fucd------------------3" << std::endl;
+ debug() << "call createFinalisedLCIOTrack now" << endmsg;
try {
status = MarlinTrk::createFinalisedLCIOTrack(marlinTrk, trkHits, &track, fit_backwards, covMatrix, _bField, _maxChi2PerHit);
} catch (...) {
-
- // delete Track;
+ // delete Track;
// delete marlinTrk;
- error() << "MarlinTrk::createFinalisedLCIOTrack " << endmsg;
- throw ;
-
+ error() << "MarlinTrk::createFinalisedLCIOTrack fail" << endmsg;
+ //throw ;
}
- std::cout << "fucd------------------4" << std::endl;
+ debug() << "createFinalisedLCIOTrack finish" << endmsg;
/*
#ifdef MARLINTRK_DIAGNOSTICS_ON
if ( status != IMarlinTrack::success && _runMarlinTrkDiagnostics ) {
@@ -2812,10 +2794,12 @@ void SiliconTrackingAlg::FinalRefit(edm4hep::TrackCollection* trk_col) {
marlinTrk->getHitsInFit(hits_in_fit);
for ( unsigned ihit = 0; ihit < hits_in_fit.size(); ++ihit) {
- all_hits.push_back(hits_in_fit[ihit].first);
+ debug() << "Hit address=" << hits_in_fit[ihit].first << endmsg;
+ edm4hep::TrackerHit* trk = hits_in_fit[ihit].first;
+ all_hits.push_back(trk);//hits_in_fit[ihit].first);
}
- UTIL::BitField64 cellID_encoder( lcio::ILDCellID0::encoder_string ) ;
+ UTIL::BitField64 cellID_encoder( UTIL::ILDCellID0::encoder_string ) ;
MarlinTrk::addHitNumbersToTrack(&track, all_hits, true, cellID_encoder);
@@ -2844,18 +2828,15 @@ void SiliconTrackingAlg::FinalRefit(edm4hep::TrackCollection* trk_col) {
//if (nhits_in_ftd > 0) Track->setTypeBit( lcio::ILDDetID::FTD ) ;
//if (nhits_in_sit > 0) Track->setTypeBit( lcio::ILDDetID::SIT ) ;
-
-
if( status != IMarlinTrack::success ) {
-
//delete track;
- debug() << "SiliconTrackingAlg::FinalRefit: Track fit failed with error code " << status << " track dropped. Number of hits = "<< trkHits.size() << endmsg;
+ debug() << "FinalRefit: Track fit failed with error code " << status << " track dropped. Number of hits = "<< trkHits.size() << endmsg;
continue ;
}
if( track.getNdf() < 0) {
//delete track;
- debug() << "SiliconTrackingAlg::FinalRefit: Track fit returns " << track.getNdf() << " degress of freedom track dropped. Number of hits = "<< trkHits.size() << endmsg;
+ debug() << "FinalRefit: Track fit returns " << track.getNdf() << " degress of freedom track dropped. Number of hits = "<< trkHits.size() << endmsg;
//delete track;
continue ;
}
@@ -2908,18 +2889,14 @@ void SiliconTrackingAlg::FinalRefit(edm4hep::TrackCollection* trk_col) {
}
}
- debug() << "SiliconTrackingAlg -> run " << _nRun
- << " event " << _nEvt << endmsg;
+ debug() << " -> run " << _nRun << " event " << _nEvt << endmsg;
debug() << "Number of Si tracks = " << nSiSegments << endmsg;
debug() << "Total 4-momentum of Si Tracks : E = " << std::setprecision(7) << eTot
<< " Px = " << pxTot
<< " Py = " << pyTot
<< " Pz = " << pzTot << endmsg;
-
-
}
-
StatusCode SiliconTrackingAlg::setupGearGeom(){
auto _gear = service<IGearSvc>("GearSvc");
if ( !_gear ) {
diff --git a/Reconstruction/SiliconTracking/src/SiliconTrackingAlg.h b/Reconstruction/SiliconTracking/src/SiliconTrackingAlg.h
index e5cfca0c..adb20296 100644
--- a/Reconstruction/SiliconTracking/src/SiliconTrackingAlg.h
+++ b/Reconstruction/SiliconTracking/src/SiliconTrackingAlg.h
@@ -10,7 +10,7 @@
#include "edm4hep/SimTrackerHitCollection.h"
#include "edm4hep/TrackerHitCollection.h"
#include "edm4hep/TrackCollection.h"
-//#include "edm4hep/LCRelationCollection.h"
+#include "edm4hep/MCRecoTrackerAssociationCollection.h"
//#include "lcio.h"
#include <string>
@@ -195,7 +195,7 @@ class SiliconTrackingAlg : public GaudiAlgorithm {
virtual StatusCode finalize() ;
-protected:
+ protected:
int _nRun ;
int _nEvt ;
@@ -214,7 +214,8 @@ protected:
//std::string _MarlinTrkDiagnosticsName;
typedef std::vector<int> IntVec;
- Gaudi::Property<IntVec> _Combinations{this, "LayerCombinations", {8,6,5, 8,6,4, 8,6,3, 8,6,2, 8,5,3, 8,5,2, 8,4,3, 8,4,2, 6,5,3, 6,5,2, 6,4,3, 6,4,2, 6,3,1, 6,3,0, 6,2,1, 6,2,0,
+ Gaudi::Property<IntVec> _Combinations{this, "LayerCombinations",
+ {8,6,5, 8,6,4, 8,6,3, 8,6,2, 8,5,3, 8,5,2, 8,4,3, 8,4,2, 6,5,3, 6,5,2, 6,4,3, 6,4,2, 6,3,1, 6,3,0, 6,2,1, 6,2,0,
5,3,1, 5,3,0, 5,2,1, 5,2,0, 4,3,1, 4,3,0, 4,2,1, 4,2,0}};
Gaudi::Property<IntVec> _CombinationsFTD{this, "LayerCombinationsFTD", {4,3,2, 4,3,1, 4,3,0, 4,2,1, 4,2,0, 4,1,0, 3,2,1, 3,2,0, 3,1,0, 2,1,0,
9,8,7, 9,8,6, 9,8,5, 9,7,6, 9,7,5, 9,6,5, 8,7,6, 8,7,5, 8,6,5, 7,6,5}};
@@ -287,12 +288,12 @@ protected:
// Input collections
DataHandle<edm4hep::EventHeaderCollection> _headerColHdl{"EventHeaderCol", Gaudi::DataHandle::Reader, this};
DataHandle<edm4hep::MCParticleCollection> _inMCColHdl{"MCParticle", Gaudi::DataHandle::Reader, this};
- //DataHandle<edm4hep::TrackerHitPlaneCollection> _inVTXColHdl{"VXDCollection", Gaudi::DataHandle::Reader, this};
- //DataHandle<edm4hep::TrackerHitPlaneCollection> _inFTDPixelColHdl{"FTDPixelCollection", Gaudi::DataHandle::Reader, this};
DataHandle<edm4hep::TrackerHitCollection> _inVTXColHdl{"VXDTrackerHits", Gaudi::DataHandle::Reader, this};
DataHandle<edm4hep::TrackerHitCollection> _inFTDPixelColHdl{"FTDPixelTrackerHits", Gaudi::DataHandle::Reader, this};
DataHandle<edm4hep::TrackerHitCollection> _inFTDSpacePointColHdl{"FTDSpacePoints", Gaudi::DataHandle::Reader, this};
- DataHandle<edm4hep::TrackerHitCollection> _inSITColHdl{"SITTrackerHits", Gaudi::DataHandle::Reader, this};
+ DataHandle<edm4hep::TrackerHitCollection> _inSITColHdl{"SITSpacePoints", Gaudi::DataHandle::Reader, this};
+ DataHandle<edm4hep::TrackerHitCollection> _inSITRawColHdl{"SITTrackerHits", Gaudi::DataHandle::Reader, this};
+ DataHandle<edm4hep::TrackerHitCollection> _inFTDRawColHdl{"FTDStripTrackerHits", Gaudi::DataHandle::Reader, this};
// Output collections
DataHandle<edm4hep::TrackCollection> _outColHdl{"SiTracks", Gaudi::DataHandle::Writer, this};
//DataHandle<edm4hep::LCRelationCollection> _outRelColHdl{"TrackerHitRelations", Gaudi::DataHandle::Reader, this};
@@ -419,8 +420,6 @@ protected:
static const int _output_track_col_quality_GOOD;
static const int _output_track_col_quality_FAIR;
static const int _output_track_col_quality_POOR;
-
- std::vector<edm4hep::TrackerHit> _allHits;
} ;
#endif
diff --git a/Service/TrackSystemSvc/CMakeLists.txt b/Service/TrackSystemSvc/CMakeLists.txt
index f441f576..b44c4a82 100644
--- a/Service/TrackSystemSvc/CMakeLists.txt
+++ b/Service/TrackSystemSvc/CMakeLists.txt
@@ -9,12 +9,11 @@ find_package(KalTest REQUIRED)
find_package(KalDet REQUIRED)
#find_package(streamlog REQUIRED)
-gaudi_depends_on_subdirs(Service/GearSvc Detector/DetInterface)
+gaudi_depends_on_subdirs(Service/GearSvc Detector/DetInterface Utilities/DataHelper)
set(TrackSystemSvc_srcs src/*.cpp)
set(TrackSystemSvcLib_srcs src/*.cc)
-
gaudi_install_headers(TrackSystemSvc)
message( "${KalDet_INCLUDE_DIRS}" )
@@ -23,7 +22,7 @@ message( "${KalDet_LIBRARIES}" )
gaudi_add_library(TrackSystemSvcLib ${TrackSystemSvcLib_srcs}
PUBLIC_HEADERS TrackSystemSvc
INCLUDE_DIRS GaudiKernel ROOT CLHEP gear ${LCIO_INCLUDE_DIRS} ${KalDet_INCLUDE_DIRS} ${KalTest_INCLUDE_DIRS} ${EDM4HEP_INCLUDE_DIRS}
- LINK_LIBRARIES GaudiKernel ROOT CLHEP $ENV{GEAR}/lib/libgear.so $ENV{GEAR}/lib/libgearsurf.so ${LCIO_LIBRARIES} ${KalDet_LIBRARIES} ${KalTest_LIBRARIES} -Wl,--no-as-needed
+ LINK_LIBRARIES DataHelperLib GaudiKernel ROOT CLHEP $ENV{GEAR}/lib/libgear.so $ENV{GEAR}/lib/libgearsurf.so ${LCIO_LIBRARIES} ${KalDet_LIBRARIES} ${KalTest_LIBRARIES} -Wl,--no-as-needed
EDM4HEP::edm4hep EDM4HEP::edm4hepDict
-Wl,--as-needed
)
diff --git a/Service/TrackSystemSvc/src/HelixTrack.cc b/Service/TrackSystemSvc/src/HelixTrack.cc
index 5934082d..8901f503 100644
--- a/Service/TrackSystemSvc/src/HelixTrack.cc
+++ b/Service/TrackSystemSvc/src/HelixTrack.cc
@@ -17,16 +17,13 @@ HelixTrack::HelixTrack( const edm4hep::Vector3d& x1, const edm4hep::Vector3d& x2
TVector3 p2( x2[0], x2[1], x2[2] );
TVector3 p3( x3[0], x3[1], x3[2] );
/*
- streamlog_out(DEBUG2) << "HelixTrack::HelixTrack Create from hits: \n "
- << "P1 x = " << p1.x() << " y = " << p1.y() << " z = " << p1.z() << " r = " << p1.Perp() << "\n "
- << "P2 x = " << p2.x() << " y = " << p2.y() << " z = " << p2.z() << " r = " << p2.Perp() << "\n "
- << "P3 x = " << p3.x() << " y = " << p3.y() << " z = " << p3.z() << " r = " << p3.Perp() << "\n "
- << "Bz = " << Bz << " direction = " << direction
- << std::endl;
+ std::cout << "debug: " << "HelixTrack::HelixTrack Create from hits: \n "
+ << "P1 x = " << p1.x() << " y = " << p1.y() << " z = " << p1.z() << " r = " << p1.Perp() << "\n "
+ << "P2 x = " << p2.x() << " y = " << p2.y() << " z = " << p2.z() << " r = " << p2.Perp() << "\n "
+ << "P3 x = " << p3.x() << " y = " << p3.y() << " z = " << p3.z() << " r = " << p3.Perp() << "\n "
+ << "Bz = " << Bz << " direction = " << direction
+ << std::endl;
*/
-
-
-
THelicalTrack* helicalTrack;
helicalTrack = new THelicalTrack( p1, p2, p3, Bz, direction );
diff --git a/Service/TrackSystemSvc/src/MarlinKalTest.cc b/Service/TrackSystemSvc/src/MarlinKalTest.cc
index ee0cd734..b6e32021 100644
--- a/Service/TrackSystemSvc/src/MarlinKalTest.cc
+++ b/Service/TrackSystemSvc/src/MarlinKalTest.cc
@@ -111,10 +111,8 @@ namespace MarlinTrk{
_det->Install( *supportdet ) ;
}
catch( gear::UnknownParameterException& e){
-
- //streamlog_out( ERROR ) << "MarlinKalTest - Support Material missing in gear file: Cannot proceed as propagations and extrapolations for cannonical track states are impossible: exit(1) called" << std::endl ;
+ std::cout << "Error: " << "MarlinKalTest - Support Material missing in gear file: Cannot proceed as propagations and extrapolations for cannonical track states are impossible: exit(1) called" << std::endl ;
exit(1);
-
}
try{
@@ -124,7 +122,7 @@ namespace MarlinTrk{
_det->Install( *vxddet ) ;
}
catch( gear::UnknownParameterException& e){
- //streamlog_out( MESSAGE ) << " MarlinKalTest - VXD missing in gear file: VXD Material Not Built " << std::endl ;
+ std::cout << "Warning: " << " MarlinKalTest - VXD missing in gear file: VXD Material Not Built " << std::endl ;
}
@@ -148,7 +146,7 @@ namespace MarlinTrk{
_det->Install( *sitdet ) ;
}
catch( gear::UnknownParameterException& e){
- //streamlog_out( MESSAGE ) << " MarlinKalTest - Simple Cylinder Based SIT missing in gear file: Simple Cylinder Based SIT Not Built " << std::endl ;
+ std::cout << "Warning: " << " MarlinKalTest - Simple Cylinder Based SIT missing in gear file: Simple Cylinder Based SIT Not Built " << std::endl ;
}
}
@@ -158,8 +156,8 @@ namespace MarlinTrk{
this->storeActiveMeasurementModuleIDs(setdet);
_det->Install( *setdet ) ;
}
- catch( gear::UnknownParameterException& e){
- //streamlog_out( MESSAGE ) << " MarlinKalTest - SET missing in gear file: SET Not Built " << std::endl ;
+ catch( gear::UnknownParameterException& e){
+ std::cout << "Warning: " << " MarlinKalTest - SET missing in gear file: SET Not Built " << std::endl ;
}
@@ -172,7 +170,7 @@ namespace MarlinTrk{
FTD_found = true ;
}
catch( gear::UnknownParameterException& e){
- //streamlog_out( MESSAGE ) << " MarlinKalTest - Petal Based FTD missing in gear file: Petal Based FTD Not Built " << std::endl ;
+ std::cout << "Warning: " << " MarlinKalTest - Petal Based FTD missing in gear file: Petal Based FTD Not Built " << std::endl ;
}
if( ! FTD_found ){
@@ -183,7 +181,7 @@ namespace MarlinTrk{
_det->Install( *ftddet ) ;
}
catch( gear::UnknownParameterException& e){
- //streamlog_out( MESSAGE ) << " MarlinKalTest - Simple Disc Based FTD missing in gear file: Simple Disc Based FTD Not Built " << std::endl ;
+ std::cout << "Warning: " << " MarlinKalTest - Simple Disc Based FTD missing in gear file: Simple Disc Based FTD Not Built " << std::endl ;
}
}
@@ -194,7 +192,7 @@ namespace MarlinTrk{
_det->Install( *tpcdet ) ;
}
catch( gear::UnknownParameterException& e){
- //streamlog_out( MESSAGE ) << " MarlinKalTest - TPC missing in gear file: TPC Not Built " << std::endl ;
+ std::cout << "Warning: " << " MarlinKalTest - TPC missing in gear file: TPC Not Built " << std::endl ;
}
}
diff --git a/Service/TrackSystemSvc/src/MarlinKalTestTrack.cc b/Service/TrackSystemSvc/src/MarlinKalTestTrack.cc
index 473f5ac6..a98029db 100644
--- a/Service/TrackSystemSvc/src/MarlinKalTestTrack.cc
+++ b/Service/TrackSystemSvc/src/MarlinKalTestTrack.cc
@@ -116,22 +116,21 @@ namespace MarlinTrk {
}
int MarlinKalTestTrack::addHit( edm4hep::TrackerHit* trkhit, const ILDVMeasLayer* ml) {
-
- std::cout << "MarlinKalTestTrack::addHit: trkhit = " << trkhit->id() << " addr: " << trkhit << " ml = " << ml << std::endl ;
-
+ //std::cout << "MarlinKalTestTrack::addHit: trkhit = " << trkhit->id() << " addr: " << trkhit << " ml = " << ml << std::endl ;
if( trkhit && ml ) {
//if(ml){
return this->addHit( trkhit, ml->ConvertLCIOTrkHit(trkhit), ml) ;
}
else {
+ std::cout << "MarlinKalTestTrack::addHit: trkhit = " << trkhit->id() << " addr: " << trkhit << " ml = " << ml << std::endl ;
//streamlog_out( ERROR ) << " MarlinKalTestTrack::addHit - bad inputs " << trkhit << " ml : " << ml << std::endl ;
- return bad_intputs ;
+ return bad_intputs ;
}
return bad_intputs ;
}
int MarlinKalTestTrack::addHit( edm4hep::TrackerHit* trkhit, ILDVTrackHit* kalhit, const ILDVMeasLayer* ml) {
- std::cout << "MarlinKalTestTrack::addHit: trkhit = " << trkhit->id() << " ILDVTrackHit: " << kalhit << " ml = " << ml << std::endl ;
+ //std::cout << "MarlinKalTestTrack::addHit: trkhit = " << trkhit->id() << " ILDVTrackHit: " << kalhit << " ml = " << ml << std::endl ;
if( kalhit && ml ) {
//if(ml){
_kalhits->Add(kalhit ) ; // Add hit and set surface found
@@ -139,36 +138,26 @@ namespace MarlinTrk {
// _kaltest_hits_to_lcio_hits[kalhit] = trkhit ; // add hit to map relating kaltest and lcio hits
}
else {
- delete kalhit;
+ //std::cout << "MarlinKalTestTrack::addHit: trkhit = " << trkhit->id() << " ILDVTrackHit: " << kalhit << " ml = " << ml << std::endl ;
+ if(kalhit) delete kalhit;
return bad_intputs ;
}
- /*
- streamlog_out(DEBUG1) << "MarlinKalTestTrack::addHit: hit added "
- << "number of hits for track = " << _kalhits->GetEntries()
- << std::endl ;
- */
+ //std::cout << "debug: " << "MarlinKalTestTrack::addHit: hit added number of hits for track = " << _kalhits->GetEntries() << std::endl ;
return success ;
-
}
int MarlinKalTestTrack::initialise( bool fitDirection ) {;
-
-
//SJA:FIXME: check here if the track is already initialised, and for now don't allow it to be re-initialised
// if the track is going to be re-initialised then we would need to do it directly on the first site
if ( _initialised ) {
-
throw MarlinTrk::Exception("Track fit already initialised");
-
}
if (_kalhits->GetEntries() < 3) {
-
- //streamlog_out( ERROR) << "<<<<<< MarlinKalTestTrack::initialise: Shortage of Hits! nhits = "
- //<< _kalhits->GetEntries() << " >>>>>>>" << std::endl;
+ std::cout << "Error: <<<<<< MarlinKalTestTrack::initialise: Shortage of Hits! nhits = "
+ << _kalhits->GetEntries() << " >>>>>>>" << std::endl;
return error ;
-
}
_fitDirection = fitDirection ;
@@ -206,7 +195,7 @@ namespace MarlinTrk {
pDummyHit = (new ILDPlanarStripHit(*static_cast<ILDPlanarStripHit*>( startingHit )));
}
else {
- //streamlog_out( ERROR) << "<<<<<<<<< MarlinKalTestTrack::initialise: dynamic_cast failed for hit type >>>>>>>" << std::endl;
+ std::cout << "Error: <<<<<<<<< MarlinKalTestTrack::initialise: dynamic_cast failed for hit type >>>>>>>" << std::endl;
return error ;
}
@@ -972,8 +961,10 @@ namespace MarlinTrk {
int MarlinKalTestTrack::getHitsInFit( std::vector<std::pair<edm4hep::TrackerHit*, double> >& hits ) {
-
+ //std::cout << "debug: _hit_chi2_values address= " << &_hit_chi2_values << " " << &(*(_hit_chi2_values.begin())) << " want to copy to hits address=" << &hits << std::endl;
std::copy( _hit_chi2_values.begin() , _hit_chi2_values.end() , std::back_inserter( hits ) ) ;
+ //hits.resize(_hit_chi2_values.size());
+ //std::copy( _hit_chi2_values.begin() , _hit_chi2_values.end() , hits.begin());
// this needs more thought. What about when the hits are added using addAndFit?
@@ -994,8 +985,7 @@ namespace MarlinTrk {
int MarlinKalTestTrack::getOutliers( std::vector<std::pair<edm4hep::TrackerHit*, double> >& hits ) {
std::copy( _outlier_chi2_values.begin() , _outlier_chi2_values.end() , std::back_inserter( hits ) ) ;
-
-
+
// this needs more thought. What about when the hits are added using addAndFit?
// // need to check the order so that we can return the list ordered in time
// // as they will be added to _hit_chi2_values in the order of fitting
@@ -1006,13 +996,10 @@ namespace MarlinTrk {
// } else {
// std::copy( _outlier_chi2_values.begin() , _outlier_chi2_values.end() , std::back_inserter( hits ) ) ;
// }
-
-
+
return success ;
-
}
-
int MarlinKalTestTrack::getNDF( int& ndf ){
if( _initialised == false ) {
diff --git a/Service/TrackSystemSvc/src/MarlinTrkUtils.cc b/Service/TrackSystemSvc/src/MarlinTrkUtils.cc
index f77df79e..cf8ae82b 100644
--- a/Service/TrackSystemSvc/src/MarlinTrkUtils.cc
+++ b/Service/TrackSystemSvc/src/MarlinTrkUtils.cc
@@ -1,4 +1,3 @@
-
#include "TrackSystemSvc/MarlinTrkUtils.h"
#include <vector>
@@ -7,7 +6,8 @@
#include "TrackSystemSvc/IMarlinTrack.h"
#include "TrackSystemSvc/HelixTrack.h"
-#include "lcio.h"
+#include "DataHelper/Navagation.h"
+//#include "lcio.h"
//#include <IMPL/TrackImpl.h>
//#include <IMPL/TrackStateImpl.h>
//#include <EVENT/TrackerHit.h>
@@ -101,21 +101,19 @@ namespace MarlinTrk {
if ( hit_list.empty() ) return IMarlinTrack::bad_intputs ;
if( track == 0 ){
- throw EVENT::Exception( std::string("MarlinTrk::finaliseLCIOTrack: TrackImpl == NULL ") ) ;
+ throw std::runtime_error( "MarlinTrk::finaliseLCIOTrack: TrackImpl == NULL " ) ;
}
int return_error = 0;
-
-
///////////////////////////////////////////////////////
// produce prefit parameters
///////////////////////////////////////////////////////
edm4hep::TrackState pre_fit ;
- std::cout << "fucd:=====================" << std::endl;
+ //std::cout << "debug:=====================before createPrefit" << std::endl;
return_error = createPrefit(hit_list, &pre_fit, bfield_z, fit_backwards);
- std::cout << "fucd:=====================" << return_error << std::endl;
+ //std::cout << "debug:=====================after createPrefit return code=" << return_error << std::endl;
pre_fit.covMatrix = initial_cov_for_prefit;
///////////////////////////////////////////////////////
@@ -127,12 +125,9 @@ namespace MarlinTrk {
return_error = createFinalisedLCIOTrack( marlinTrk, hit_list, track, fit_backwards, &pre_fit, bfield_z, maxChi2Increment);
} else {
- std::cout << "MarlinTrk::createFinalisedLCIOTrack : Prefit failed error = " << return_error << std::endl;
+ //std::cout << "MarlinTrk::createFinalisedLCIOTrack : Prefit failed error = " << return_error << std::endl;
}
-
-
return return_error;
-
}
int createFinalisedLCIOTrack( IMarlinTrack* marlinTrk, std::vector<edm4hep::TrackerHit*>& hit_list, edm4hep::Track* track, bool fit_backwards, edm4hep::TrackState* pre_fit, float bfield_z, double maxChi2Increment){
@@ -144,11 +139,11 @@ namespace MarlinTrk {
if ( hit_list.empty() ) return IMarlinTrack::bad_intputs ;
if( track == 0 ){
- throw EVENT::Exception( std::string("MarlinTrk::finaliseLCIOTrack: TrackImpl == NULL ") ) ;
+ throw std::runtime_error("MarlinTrk::finaliseLCIOTrack: TrackImpl == NULL ");
}
if( pre_fit == 0 ){
- throw EVENT::Exception( std::string("MarlinTrk::finaliseLCIOTrack: TrackStateImpl == NULL ") ) ;
+ throw std::runtime_error("MarlinTrk::finaliseLCIOTrack: TrackStateImpl == NULL ");
}
@@ -156,10 +151,9 @@ namespace MarlinTrk {
if( fit_status != IMarlinTrack::success ){
- std::cout << "MarlinTrk::createFinalisedLCIOTrack fit failed: fit_status = " << fit_status << std::endl;
+ //std::cout << "MarlinTrk::createFinalisedLCIOTrack fit failed: fit_status = " << fit_status << std::endl;
return fit_status;
-
}
int error = finaliseLCIOTrack(marlinTrk, track, hit_list);
@@ -179,11 +173,11 @@ namespace MarlinTrk {
if ( hit_list.empty() ) return IMarlinTrack::bad_intputs;
if( marlinTrk == 0 ){
- throw EVENT::Exception( std::string("MarlinTrk::createFit: IMarlinTrack == NULL ") ) ;
+ throw std::runtime_error("MarlinTrk::createFit: IMarlinTrack == NULL ");
}
if( pre_fit == 0 ){
- throw EVENT::Exception( std::string("MarlinTrk::createFit: TrackStateImpl == NULL ") ) ;
+ throw std::runtime_error("MarlinTrk::createFit: TrackStateImpl == NULL ");
}
int return_error = 0;
@@ -212,7 +206,7 @@ namespace MarlinTrk {
std::vector<edm4hep::TrackerHit*>::iterator it = hit_list.begin();
// start by trying to add the hits to the track we want to finally use.
- std::cout << "MarlinTrk::createFit Start Fit: AddHits: number of hits to fit " << hit_list.size() << std::endl;
+ //std::cout << "MarlinTrk::createFit Start Fit: AddHits: number of hits to fit " << hit_list.size() << std::endl;
std::vector<edm4hep::TrackerHit*> added_hits;
unsigned int ndof_added = 0;
@@ -221,31 +215,23 @@ namespace MarlinTrk {
edm4hep::TrackerHit* trkHit = *it;
bool isSuccessful = false;
- std::cout << "debug TrackerHit pointer " << trkHit << std::endl;
- if( UTIL::BitSet32( trkHit->getType() )[ UTIL::ILDTrkHitTypeBit::COMPOSITE_SPACEPOINT ] ){ //it is a composite spacepoint
-
+ //std::cout << "debug: TrackerHit pointer " << trkHit << std::endl;
+ if( UTIL::BitSet32( trkHit->getType() )[ UTIL::ILDTrkHitTypeBit::COMPOSITE_SPACEPOINT ] ){ //it is a composite spacepoint
//Split it up and add both hits to the MarlinTrk
//const EVENT::LCObjectVec rawObjects = trkHit->getRawHits();
- std::cout << "space point is not still valid! pelease wait updating..." <<std::endl;
- exit(1);
- /*
+ //std::cout << "space point is not still valid! pelease wait updating..." <<std::endl;
+ //exit(1);
int nRawHit = trkHit->rawHits_size();
-
for( unsigned k=0; k< nRawHit; k++ ){
-
- edm4hep::TrackerHit* rawHit = dynamic_cast< edm4hep::TrackerHit* >(trkHit->getRawHits(k));
-
- if( marlinTrk->addHit( rawHit ) == IMarlinTrack::success ){
-
- isSuccessful = true; //if at least one hit from the spacepoint gets added
+ edm4hep::TrackerHit* rawHit = Navagation::Instance()->GetTrackerHit(trkHit->getRawHits(k));
+ if( marlinTrk->addHit( rawHit ) == IMarlinTrack::success ){
+ isSuccessful = true; //if at least one hit from the spacepoint gets added
++ndof_added;
// streamlog_out(DEBUG4) << "MarlinTrk::createFit ndof_added = " << ndof_added << std::endl;
}
}
- */
}
else { // normal non composite hit
-
if (marlinTrk->addHit( trkHit ) == IMarlinTrack::success ) {
isSuccessful = true;
ndof_added += 2;
@@ -257,7 +243,7 @@ namespace MarlinTrk {
added_hits.push_back(trkHit);
}
else{
- std::cout << "Hit " << it - hit_list.begin() << " Dropped " << std::endl;
+ //std::cout << "MarlinTrkUtils::createFit Hit " << it - hit_list.begin() << " Dropped " << std::endl;
}
}
@@ -286,7 +272,7 @@ namespace MarlinTrk {
// try fit and return error
///////////////////////////////////////////////////////
int status = marlinTrk->fit(maxChi2Increment);
- std::cout << "fucd===================" << status << std::endl;
+ //std::cout << "debug:===================createFit " << status << std::endl;
return status;
@@ -440,47 +426,35 @@ namespace MarlinTrk {
edm4hep::TrackerHit* trkHit = hit_list[ihit];
if( UTIL::BitSet32( trkHit->getType() )[ UTIL::ILDTrkHitTypeBit::COMPOSITE_SPACEPOINT ] ){ //it is a composite spacepoint
- std::cout << "Error: space point is not still valid! pelease wait updating..." <<std::endl;
- exit(1);
- /*
- // get strip hits
+ //std::cout << "Error: space point is not still valid! pelease wait updating..." <<std::endl;
+ //exit(1);
+ // get strip hits
int nRawHit = trkHit->rawHits_size();
-
for( unsigned k=0; k< nRawHit; k++ ){
-
- edm4hep::TrackerHit* rawHit = dynamic_cast< edm4hep::TrackerHit* >(trkHit->getRawHits(k));
-
- bool is_outlier = false;
-
- // here we loop over outliers as this will be faster than looping over the used hits
+ edm4hep::TrackerHit* rawHit = Navagation::Instance()->GetTrackerHit(trkHit->getRawHits(k));
+ bool is_outlier = false;
+ // here we loop over outliers as this will be faster than looping over the used hits
for ( unsigned ohit = 0; ohit < outliers.size(); ++ohit) {
-
- if ( rawHit == outliers[ohit].first ) {
+ if ( rawHit->id() == outliers[ohit].first->id() ) {
is_outlier = true;
break; // break out of loop over outliers
}
}
-
if (is_outlier == false) {
used_hits.push_back(hit_list[ihit]);
- track->addTrackerHit(*used_hits.back());
+ track->addToTrackerHits(*used_hits.back());
break; // break out of loop over rawObjects
}
}
- */
} else {
-
bool is_outlier = false;
-
// here we loop over outliers as this will be faster than looping over the used hits
for ( unsigned ohit = 0; ohit < outliers.size(); ++ohit) {
-
if ( trkHit == outliers[ohit].first ) {
is_outlier = true;
break; // break out of loop over outliers
}
}
-
if (is_outlier == false) {
used_hits.push_back(hit_list[ihit]);
track->addToTrackerHits(*used_hits.back());
@@ -693,36 +667,34 @@ namespace MarlinTrk {
// check inputs
///////////////////////////////////////////////////////
if( track == 0 ){
- throw EVENT::Exception( std::string("MarlinTrk::addHitsToTrack: TrackImpl == NULL ") ) ;
+ throw std::runtime_error( std::string("MarlinTrk::addHitsToTrack: TrackImpl == NULL ") ) ;
}
-
std::map<int, int> hitNumbers;
- hitNumbers[lcio::ILDDetID::VXD] = 0;
- hitNumbers[lcio::ILDDetID::SIT] = 0;
- hitNumbers[lcio::ILDDetID::FTD] = 0;
- hitNumbers[lcio::ILDDetID::TPC] = 0;
- hitNumbers[lcio::ILDDetID::SET] = 0;
- hitNumbers[lcio::ILDDetID::ETD] = 0;
+ hitNumbers[UTIL::ILDDetID::VXD] = 0;
+ hitNumbers[UTIL::ILDDetID::SIT] = 0;
+ hitNumbers[UTIL::ILDDetID::FTD] = 0;
+ hitNumbers[UTIL::ILDDetID::TPC] = 0;
+ hitNumbers[UTIL::ILDDetID::SET] = 0;
+ hitNumbers[UTIL::ILDDetID::ETD] = 0;
for(unsigned int j=0; j<hit_list.size(); ++j) {
-
cellID_encoder.setValue(hit_list.at(j)->getCellID()) ;
int detID = cellID_encoder[UTIL::ILDCellID0::subdet];
++hitNumbers[detID];
- // streamlog_out( DEBUG1 ) << "Hit from Detector " << detID << std::endl;
+ //std::cout << "debug: " << "Hit from Detector " << detID << std::endl;
}
int offset = 2 ;
if ( hits_in_fit == false ) { // all hit atributed by patrec
offset = 1 ;
}
- track->addToSubDetectorHitNumbers(hitNumbers[lcio::ILDDetID::VXD]);
- track->addToSubDetectorHitNumbers(hitNumbers[lcio::ILDDetID::FTD]);
- track->addToSubDetectorHitNumbers(hitNumbers[lcio::ILDDetID::SIT]);
- track->addToSubDetectorHitNumbers(hitNumbers[lcio::ILDDetID::TPC]);
- track->addToSubDetectorHitNumbers(hitNumbers[lcio::ILDDetID::SET]);
- track->addToSubDetectorHitNumbers(hitNumbers[lcio::ILDDetID::ETD]);
+ track->addToSubDetectorHitNumbers(hitNumbers[UTIL::ILDDetID::VXD]);
+ track->addToSubDetectorHitNumbers(hitNumbers[UTIL::ILDDetID::FTD]);
+ track->addToSubDetectorHitNumbers(hitNumbers[UTIL::ILDDetID::SIT]);
+ track->addToSubDetectorHitNumbers(hitNumbers[UTIL::ILDDetID::TPC]);
+ track->addToSubDetectorHitNumbers(hitNumbers[UTIL::ILDDetID::SET]);
+ track->addToSubDetectorHitNumbers(hitNumbers[UTIL::ILDDetID::ETD]);
//track->subdetectorHitNumbers().resize(2 * lcio::ILDDetID::ETD);
//track->subdetectorHitNumbers()[ 2 * lcio::ILDDetID::VXD - offset ] = hitNumbers[lcio::ILDDetID::VXD];
//track->subdetectorHitNumbers()[ 2 * lcio::ILDDetID::FTD - offset ] = hitNumbers[lcio::ILDDetID::FTD];
@@ -730,9 +702,6 @@ namespace MarlinTrk {
//track->subdetectorHitNumbers()[ 2 * lcio::ILDDetID::TPC - offset ] = hitNumbers[lcio::ILDDetID::TPC];
//track->subdetectorHitNumbers()[ 2 * lcio::ILDDetID::SET - offset ] = hitNumbers[lcio::ILDDetID::SET];
//track->subdetectorHitNumbers()[ 2 * lcio::ILDDetID::ETD - offset ] = hitNumbers[lcio::ILDDetID::ETD];
-
-
-
}
void addHitNumbersToTrack(edm4hep::Track* track, std::vector<std::pair<edm4hep::TrackerHit* , double> >& hit_list, bool hits_in_fit, UTIL::BitField64& cellID_encoder){
@@ -741,7 +710,7 @@ namespace MarlinTrk {
// check inputs
///////////////////////////////////////////////////////
if( track == 0 ){
- throw EVENT::Exception( std::string("MarlinTrk::addHitsToTrack: TrackImpl == NULL ") ) ;
+ throw std::runtime_error( std::string("MarlinTrk::addHitsToTrack: TrackImpl == NULL ") ) ;
}
std::map<int, int> hitNumbers;
@@ -754,11 +723,10 @@ namespace MarlinTrk {
hitNumbers[lcio::ILDDetID::ETD] = 0;
for(unsigned int j=0; j<hit_list.size(); ++j) {
-
cellID_encoder.setValue(hit_list.at(j).first->getCellID()) ;
int detID = cellID_encoder[UTIL::ILDCellID0::subdet];
++hitNumbers[detID];
- // streamlog_out( DEBUG1 ) << "Hit from Detector " << detID << std::endl;
+ //std::cout << "debug: Hit from Detector " << detID << std::endl;
}
int offset = 2 ;
--
GitLab