From b9e2d1f8d6ce8257db522b6b95e302f131ff55ca Mon Sep 17 00:00:00 2001 From: Alex Kastanas <Alex.Kastanas@cern.ch> Date: Wed, 5 Oct 2016 11:22:28 +0200 Subject: [PATCH] Added SDO based SiHit association and tagging broken clusters (InDetPrepRawDataToxAOD-00-01-44) * Added a SiHit-to-cluster matching scheme based on SDO barcodes, rather than geometric matching * Added a flag for clusters with truth particles contributing that also contribute to other clusters * Tagging as InDetPrepRawDataToxAOD-00-01-39-01 2016-08-22 Alex Alonso * New sequencer for the ID xAOD * Tagging as InDetPrepRawDataToxAOD-00-01-43 2016-08-19 Alex Alonso * Change the Split track collection name, as it overlapes with muon track collection * Tag as InDetPrepRawDataToxAOD-00-01-42 2016-07-14 Leigh Schaefer <leigh.schaefer@cern.ch> * Remove initialization of TRT dE/dx tool from share/InDetDxAOD.py since this tool is no longer included in DerivationFrameworkInDet * Change TRT_DriftCircles.localXError to store uncertainty rather than covariance in src/TRT_PrepDataToxAOD.cxx * Tag as InDetPrepRawDataToxAOD-00-01-40 Former-commit-id: f45bd7628085e5224aefc00dcc303c50886755d7 --- .../share/InDetDxAOD.py | 44 +++--- .../src/PixelPrepDataToxAOD.cxx | 149 +++++++++++++----- .../src/PixelPrepDataToxAOD.h | 14 +- .../src/TRT_PrepDataToxAOD.cxx | 4 +- 4 files changed, 138 insertions(+), 73 deletions(-) diff --git a/InnerDetector/InDetEventCnv/InDetPrepRawDataToxAOD/share/InDetDxAOD.py b/InnerDetector/InDetEventCnv/InDetPrepRawDataToxAOD/share/InDetDxAOD.py index 1288ddfbebd..6b91c6e6860 100644 --- a/InnerDetector/InDetEventCnv/InDetPrepRawDataToxAOD/share/InDetDxAOD.py +++ b/InnerDetector/InDetEventCnv/InDetPrepRawDataToxAOD/share/InDetDxAOD.py @@ -62,6 +62,9 @@ prefixName = "" ################# ### Setup tools ################# +from AthenaCommon import CfgMgr +IDDerivationSequence = CfgMgr.AthSequencer("InDetDxAOD_seq") + if dumpTrtInfo: from TRT_ConditionsServices.TRT_ConditionsServicesConf import TRT_StrawNeighbourSvc TRTStrawNeighbourSvc=TRT_StrawNeighbourSvc() @@ -71,17 +74,11 @@ if dumpTrtInfo: TRTCalibDBSvc=TRT_CalDbSvc() ServiceMgr += TRTCalibDBSvc - from TRT_ToT_Tools.TRT_ToT_ToolsConf import TRT_ToT_dEdx - TRT_dEdx_Tool = TRT_ToT_dEdx(name="TRT_ToT_dEdx") - ToolSvc += TRT_dEdx_Tool - #Setup charge->ToT back-conversion to restore ToT info as well if dumpPixInfo: - from AthenaCommon.AlgSequence import AlgSequence - topSequence = AlgSequence() from PixelCalibAlgs.PixelCalibAlgsConf import PixelChargeToTConversion PixelChargeToTConversionSetter = PixelChargeToTConversion(name = "PixelChargeToTConversionSetter") - topSequence += PixelChargeToTConversionSetter + IDDerivationSequence += PixelChargeToTConversionSetter if (printIdTrkDxAODConf): print PixelChargeToTConversionSetter print PixelChargeToTConversionSetter.properties() @@ -105,8 +102,8 @@ if makeSplitTracks: splittercomb=InDet__InDetSplittedTracksCreator(name='CombinedTrackSplitter', TrackSplitterTool = splittertoolcomb, TrackCollection = "Tracks", - OutputTrackCollection = "Tracks_split") - topSequence+=splittercomb + OutputTrackCollection = "Tracks_splitID") + IDDerivationSequence +=splittercomb if (printIdTrkDxAODConf): print splittercomb print splittercomb.properties() @@ -122,20 +119,20 @@ if makeSplitTracks: # The following adds truth information, but needs further testing #include ("InDetRecExample/ConfiguredInDetTrackTruth.py") #if isIdTrkDxAODSimulation: - # InDetSplitTracksTruth = ConfiguredInDetTrackTruth("Tracks_split",'SplitTrackDetailedTruth','SplitTrackTruth') + # InDetSplitTracksTruth = ConfiguredInDetTrackTruth("Tracks_splitID",'SplitTrackDetailedTruth','SplitTrackTruth') xAODSplitTrackParticleCnvAlg = xAODMaker__TrackParticleCnvAlg('InDetSplitTrackParticles') xAODSplitTrackParticleCnvAlg.xAODContainerName = 'InDetSplitTrackParticles' xAODSplitTrackParticleCnvAlg.xAODTrackParticlesFromTracksContainerName = 'InDetSplitTrackParticles' xAODSplitTrackParticleCnvAlg.TrackParticleCreator = InDetxAODSplitParticleCreatorTool - xAODSplitTrackParticleCnvAlg.TrackContainerName = 'Tracks_split' + xAODSplitTrackParticleCnvAlg.TrackContainerName = 'Tracks_splitID' xAODSplitTrackParticleCnvAlg.ConvertTrackParticles = False xAODSplitTrackParticleCnvAlg.ConvertTracks = True xAODSplitTrackParticleCnvAlg.AddTruthLink = False #isIdTrkDxAODSimulation if (isIdTrkDxAODSimulation): xAODSplitTrackParticleCnvAlg.TrackTruthContainerName = 'SplitTrackTruth' xAODSplitTrackParticleCnvAlg.PrintIDSummaryInfo = True - topSequence += xAODSplitTrackParticleCnvAlg + IDDerivationSequence += xAODSplitTrackParticleCnvAlg ################# @@ -150,7 +147,7 @@ if dumpTrtInfo: xAOD_TRT_PrepDataToxAOD.UseTruthInfo = dumpTruthInfo #xAOD_TRT_PrepDataToxAOD.WriteSDOs = True - topSequence += xAOD_TRT_PrepDataToxAOD + IDDerivationSequence += xAOD_TRT_PrepDataToxAOD if (printIdTrkDxAODConf): print xAOD_TRT_PrepDataToxAOD print xAOD_TRT_PrepDataToxAOD.properties() @@ -165,7 +162,7 @@ if dumpSctInfo: #xAOD_SCT_PrepDataToxAOD.WriteSDOs = True #xAOD_SCT_PrepDataToxAOD.WriteSiHits = True # if available - topSequence += xAOD_SCT_PrepDataToxAOD + IDDerivationSequence += xAOD_SCT_PrepDataToxAOD if (printIdTrkDxAODConf): print xAOD_SCT_PrepDataToxAOD print xAOD_SCT_PrepDataToxAOD.properties() @@ -183,7 +180,7 @@ if dumpPixInfo: if InDetFlags.doSLHC(): xAOD_PixelPrepDataToxAOD.WriteNNinformation=False - topSequence += xAOD_PixelPrepDataToxAOD + IDDerivationSequence += xAOD_PixelPrepDataToxAOD if (printIdTrkDxAODConf): print xAOD_PixelPrepDataToxAOD print xAOD_PixelPrepDataToxAOD.properties() @@ -213,9 +210,6 @@ DFTSOS = DerivationFramework__TrackStateOnSurfaceDecorator(name = "DFTrackStateO StorePixel = dumpPixInfo, IsSimulation = isIdTrkDxAODSimulation, OutputLevel = INFO) -if dumpTrtInfo: - #Add tool to calculate TRT-based dEdx - DFTSOS.TRT_ToT_dEdx = TRT_dEdx_Tool ToolSvc += DFTSOS augmentationTools+=[DFTSOS] @@ -283,7 +277,7 @@ if dumpLArCollisionTime: from RecExConfig.ObjKeyStore import cfgKeyStore # We can only do this if we have the cell container. if cfgKeyStore.isInInput ('CaloCellContainer', 'AllCalo'): - LArCollisionTimeGetter (topSequence) + LArCollisionTimeGetter (IDDerivationSequence) from DerivationFrameworkInDet.DerivationFrameworkInDetConf import DerivationFramework__LArCollisionTimeDecorator lArCollisionTimeDecorator = DerivationFramework__LArCollisionTimeDecorator (name ='lArCollisionTimeDecorator', @@ -340,18 +334,16 @@ thinningTools.append(IDTRKThinningTool) #==================================================================== # Add the derivation job to the top AthAlgSeqeuence # DerivationJob is COMMON TO ALL DERIVATIONS -DerivationFrameworkJob = CfgMgr.AthSequencer("MySeq2") -from DerivationFrameworkCore.DerivationFrameworkCoreConf import DerivationFramework__CommonAugmentation -DerivationFrameworkJob += CfgMgr.DerivationFramework__DerivationKernel("DFTSOS_KERN", +IDDerivationSequence += CfgMgr.DerivationFramework__DerivationKernel("DFTSOS_KERN", AugmentationTools = augmentationTools, SkimmingTools = skimmingTools, ThinningTools = thinningTools, OutputLevel = INFO) -topSequence += DerivationFrameworkJob +topSequence += IDDerivationSequence if (printIdTrkDxAODConf): - print DerivationFrameworkJob - print DerivationFrameworkJob.properties() + print IDDerivationSequence + print IDDerivationSequence.properties() ################# ### Steer output file content @@ -418,6 +410,8 @@ if dumpTruthInfo: IDTRKVALIDStream.AddItem("xAOD::TruthVertexAuxContainer#*") IDTRKVALIDStream.AddItem("xAOD::TruthEventContainer#*") IDTRKVALIDStream.AddItem("xAOD::TruthEventAuxContainer#*") + IDTRKVALIDStream.AddItem("xAOD::TruthPileupEventContainer#*") + IDTRKVALIDStream.AddItem("xAOD::TruthPileupEventAuxContainer#*") # add pseudo tracking in xAOD IDTRKVALIDStream.AddItem("xAOD::TrackParticleContainer#InDetPseudoTrackParticles") IDTRKVALIDStream.AddItem("xAOD::TrackParticleAuxContainer#InDetPseudoTrackParticlesAux."+excludedAuxData) diff --git a/InnerDetector/InDetEventCnv/InDetPrepRawDataToxAOD/src/PixelPrepDataToxAOD.cxx b/InnerDetector/InDetEventCnv/InDetPrepRawDataToxAOD/src/PixelPrepDataToxAOD.cxx index d1113a1fcaa..c457a04b441 100644 --- a/InnerDetector/InDetEventCnv/InDetPrepRawDataToxAOD/src/PixelPrepDataToxAOD.cxx +++ b/InnerDetector/InDetEventCnv/InDetPrepRawDataToxAOD/src/PixelPrepDataToxAOD.cxx @@ -47,11 +47,13 @@ PixelPrepDataToxAOD::PixelPrepDataToxAOD(const std::string &name, ISvcLocator *p m_pixelDCSSvc("PixelDCSSvc", name), m_pixelBSErrorsSvc("PixelByteStreamErrorsSvc", name), m_PixelHelper(0), - m_firstEventWarnings(true) + m_firstEventWarnings(true), + m_useSiHitsGeometryMatching(true) { // --- Steering and configuration flags declareProperty("UseTruthInfo", m_useTruthInfo=false); + declareProperty("UseSiHitsGeometryMatching", m_useSiHitsGeometryMatching=true); declareProperty("WriteSDOs", m_writeSDOs = true); declareProperty("WriteSiHits", m_writeSiHits = true); declareProperty("WriteNNinformation", m_writeNNinformation = true); @@ -213,7 +215,7 @@ StatusCode PixelPrepDataToxAOD::execute() //May want to addinformation about the individual hits here } xprd->setRdoIdentifierList(rdoIdentifierList); - + //Add pixel cluster properties xprd->auxdata<int>("bec") = m_PixelHelper->barrel_ec(clusterId) ; xprd->auxdata<int>("layer") = m_PixelHelper->layer_disk(clusterId) ; @@ -280,23 +282,75 @@ StatusCode PixelPrepDataToxAOD::execute() } xprd->auxdata< std::vector<int> >("truth_barcode") = barcodes; } - + + std::vector< std::vector< int > > sdo_tracks; // Use the SDO Collection to get a list of all true particle contributing to the cluster per readout element // Also get the energy deposited by each true particle per readout element if( sdoCollection ){ - addSDOInformation( xprd, prd,sdoCollection); + sdo_tracks = addSDOInformation( xprd, prd,sdoCollection); } // Now Get the most detailed truth from the SiHits // Note that this could get really slow if there are a lot of hits and clusters if( sihitCollection ){ - if(m_writeSiHits) - addSiHitInformation( xprd, prd, sihitCollection); - if(m_writeNNinformation) - addNNTruthInfo( xprd, prd, sihitCollection ); + const std::vector<SiHit> matched_hits = findAllHitsCompatibleWithCluster( prd, sihitCollection, sdo_tracks ); + + if(m_writeSiHits) + { + if ( ! sdoCollection ) + ATH_MSG_WARNING("Si hit truth information requested, but SDO collection not available!" ); + addSiHitInformation( xprd, prd, matched_hits ); + } + if(m_writeNNinformation) + { + if ( ! sdoCollection ) + ATH_MSG_WARNING("Si hit truth information requested, but SDO collection not available!" ); + addNNTruthInfo( xprd, prd, matched_hits ); + } } } } + + for ( auto clusItr = xaod->begin(); clusItr != xaod->end(); clusItr++ ) + (*clusItr)->auxdata<char>("broken") = false; + + for ( auto clusItr = xaod->begin(); clusItr != xaod->end(); clusItr++ ) + { + auto pixelCluster = *clusItr; + int layer = pixelCluster->auxdata< int >("layer"); + std::vector<int> barcodes = pixelCluster->auxdata< std::vector< int > >("sihit_barcode"); + + for ( auto clusItr2 = clusItr + 1; clusItr2 != xaod->end(); clusItr2++ ) + { + auto pixelCluster2 = *clusItr2; + if ( pixelCluster2->auxdata< int >("layer") != layer ) + continue; + if ( pixelCluster->auxdata< int >("eta_module") != pixelCluster2->auxdata< int >("eta_module") ) + continue; + if ( pixelCluster->auxdata< int >("phi_module") != pixelCluster2->auxdata< int >("phi_module") ) + continue; + + std::vector<int> barcodes2 = pixelCluster2->auxdata< std::vector< int > >("sihit_barcode"); + + bool broken = false; + for ( auto bc : barcodes ) + { + for ( auto bc2 : barcodes2 ) + { + if ( bc2 == bc ) + { + broken = true; + pixelCluster->auxdata<char>("broken") = true; + pixelCluster2->auxdata<char>("broken") = true; + break; + } + } + if ( broken ) + break; + } + } + } + ATH_MSG_DEBUG( " recorded PixelPrepData objects: size " << xaod->size() ); m_firstEventWarnings = false; @@ -305,9 +359,9 @@ StatusCode PixelPrepDataToxAOD::execute() } -void PixelPrepDataToxAOD::addSDOInformation( xAOD::TrackMeasurementValidation* xprd, - const InDet::PixelCluster* prd, - const InDetSimDataCollection* sdoCollection ) const +std::vector< std::vector< int > > PixelPrepDataToxAOD::addSDOInformation( xAOD::TrackMeasurementValidation* xprd, + const InDet::PixelCluster* prd, + const InDetSimDataCollection* sdoCollection ) const { std::vector<int> sdo_word; std::vector< std::vector< int > > sdo_depositsBarcode; @@ -336,19 +390,17 @@ void PixelPrepDataToxAOD::addSDOInformation( xAOD::TrackMeasurementValidation* x xprd->auxdata< std::vector<int> >("sdo_words") = sdo_word; xprd->auxdata< std::vector< std::vector<int> > >("sdo_depositsBarcode") = sdo_depositsBarcode; xprd->auxdata< std::vector< std::vector<float> > >("sdo_depositsEnergy") = sdo_depositsEnergy; + + return sdo_depositsBarcode; } void PixelPrepDataToxAOD::addSiHitInformation( xAOD::TrackMeasurementValidation* xprd, - const InDet::PixelCluster* prd, - const SiHitCollection* sihitCollection) const - - + const InDet::PixelCluster* prd, + const std::vector<SiHit> & matchingHits ) const { - std::vector<SiHit> matchingHits = findAllHitsCompatibleWithCluster( prd, sihitCollection ); - int numHits = matchingHits.size(); std::vector<float> sihit_energyDeposit(numHits,0); @@ -415,7 +467,8 @@ void PixelPrepDataToxAOD::addSiHitInformation( xAOD::TrackMeasurementValidation std::vector<SiHit> PixelPrepDataToxAOD::findAllHitsCompatibleWithCluster( const InDet::PixelCluster* prd, - const SiHitCollection* collection) const + const SiHitCollection* collection, + std::vector< std::vector< int > > & trkBCs ) const { ATH_MSG_VERBOSE( "Got " << collection->size() << " SiHits to look through" ); std::vector<SiHit> matchingHits; @@ -431,7 +484,7 @@ std::vector<SiHit> PixelPrepDataToxAOD::findAllHitsCompatibleWithCluster( const // Check if it is a Pixel hit if( !siHit.isPixel() ) continue; - + //Check if it is on the correct module Identifier clusterId = prd->identify(); @@ -445,26 +498,43 @@ std::vector<SiHit> PixelPrepDataToxAOD::findAllHitsCompatibleWithCluster( const // Must be within +/- 1 hits of any hit in the cluster to be included ATH_MSG_DEBUG("Hit is on the same module"); - HepGeom::Point3D<double> averagePosition = siHit.localStartPosition() + siHit.localEndPosition(); - averagePosition *= 0.5; - Amg::Vector2D pos = de->hitLocalToLocal( averagePosition.z(), averagePosition.y() ); - InDetDD::SiCellId diode = de->cellIdOfPosition(pos); - - for( const auto &hitIdentifier : prd->rdoList() ){ - - ATH_MSG_DEBUG("Truth Phi " << diode.phiIndex() << " Cluster Phi " << m_PixelHelper->phi_index( hitIdentifier ) ); - ATH_MSG_DEBUG("Truth Eta " << diode.etaIndex() << " Cluster Eta " << m_PixelHelper->eta_index( hitIdentifier ) ); - - if( abs( int(diode.etaIndex()) - m_PixelHelper->eta_index( hitIdentifier ) ) <=1 - && abs( int(diode.phiIndex()) - m_PixelHelper->phi_index( hitIdentifier ) ) <=1 ) - { - multiMatchingHits.push_back(&siHit); - break; - } + if ( m_useSiHitsGeometryMatching ) + { + HepGeom::Point3D<double> averagePosition = siHit.localStartPosition() + siHit.localEndPosition(); + averagePosition *= 0.5; + Amg::Vector2D pos = de->hitLocalToLocal( averagePosition.z(), averagePosition.y() ); + InDetDD::SiCellId diode = de->cellIdOfPosition(pos); + + for( const auto &hitIdentifier : prd->rdoList() ){ + ATH_MSG_DEBUG("Truth Phi " << diode.phiIndex() << " Cluster Phi " << m_PixelHelper->phi_index( hitIdentifier ) ); + ATH_MSG_DEBUG("Truth Eta " << diode.etaIndex() << " Cluster Eta " << m_PixelHelper->eta_index( hitIdentifier ) ); + if( abs( int(diode.etaIndex()) - m_PixelHelper->eta_index( hitIdentifier ) ) <=1 + && abs( int(diode.phiIndex()) - m_PixelHelper->phi_index( hitIdentifier ) ) <=1 ) + { + multiMatchingHits.push_back(&siHit); + break; + } + } + } + else + { + bool foundHit = false; + for ( const auto barcodeSDOColl : trkBCs ) + { + for ( const auto barcode : barcodeSDOColl ) + { + if ( siHit.particleLink().barcode() == barcode ) + { + multiMatchingHits.push_back(&siHit); + foundHit = true; + break; + } + } + if ( foundHit ) + break; + } } } - - //Now we will now make 1 SiHit for each true particle if the SiHits "touch" other std::vector<const SiHit* >::iterator siHitIter = multiMatchingHits.begin(); std::vector<const SiHit* >::iterator siHitIter2 = multiMatchingHits.begin(); @@ -823,10 +893,9 @@ void PixelPrepDataToxAOD::addNNInformation(xAOD::TrackMeasurementValidation* xpr void PixelPrepDataToxAOD::addNNTruthInfo( xAOD::TrackMeasurementValidation* xprd, const InDet::PixelCluster* pixelCluster, - const SiHitCollection* sihitCollection ) const + const std::vector<SiHit> & matchingHits ) const { - std::vector<SiHit> matchingHits = findAllHitsCompatibleWithCluster( pixelCluster, sihitCollection ); unsigned int numberOfSiHits = matchingHits.size(); @@ -880,7 +949,7 @@ void PixelPrepDataToxAOD::addNNTruthInfo( xAOD::TrackMeasurementValidation* xp if (fabs(YposC)>design->width()/2 && fabs(averagePosition.y())<design->width()/2) - { + { if (YposC>design->width()/2) { YposC=design->width()/2-1e-6; diff --git a/InnerDetector/InDetEventCnv/InDetPrepRawDataToxAOD/src/PixelPrepDataToxAOD.h b/InnerDetector/InDetEventCnv/InDetPrepRawDataToxAOD/src/PixelPrepDataToxAOD.h index 5b9c85e487a..6c124e2aae5 100644 --- a/InnerDetector/InDetEventCnv/InDetPrepRawDataToxAOD/src/PixelPrepDataToxAOD.h +++ b/InnerDetector/InDetEventCnv/InDetPrepRawDataToxAOD/src/PixelPrepDataToxAOD.h @@ -54,22 +54,23 @@ public: private: - void addSDOInformation( xAOD::TrackMeasurementValidation* xprd, - const InDet::PixelCluster* prd, - const InDetSimDataCollection* sdoCollection ) const; + std::vector< std::vector< int > > addSDOInformation( xAOD::TrackMeasurementValidation* xprd, + const InDet::PixelCluster* prd, + const InDetSimDataCollection* sdoCollection ) const; void addSiHitInformation( xAOD::TrackMeasurementValidation* xprd, const InDet::PixelCluster* prd, - const SiHitCollection* collection) const; + const std::vector<SiHit> & matchingHits ) const; std::vector<SiHit> findAllHitsCompatibleWithCluster(const InDet::PixelCluster* prd, - const SiHitCollection* collection) const; + const SiHitCollection* collection, + std::vector< std::vector< int > > & trkBCs) const; void addNNTruthInfo( xAOD::TrackMeasurementValidation* xprd, const InDet::PixelCluster* prd, - const SiHitCollection* collection ) const; + const std::vector<SiHit> & matchingHits ) const; void addNNInformation( xAOD::TrackMeasurementValidation* xprd, const InDet::PixelCluster* pixelCluster, @@ -99,6 +100,7 @@ private: bool m_writeSiHits; bool m_writeNNinformation; bool m_writeRDOinformation; + bool m_useSiHitsGeometryMatching; ServiceHandle<IPixelCalibSvc> m_calibSvc; ServiceHandle<IPixelDCSSvc> m_pixelDCSSvc; diff --git a/InnerDetector/InDetEventCnv/InDetPrepRawDataToxAOD/src/TRT_PrepDataToxAOD.cxx b/InnerDetector/InDetEventCnv/InDetPrepRawDataToxAOD/src/TRT_PrepDataToxAOD.cxx index 0587f9de85f..96600635617 100644 --- a/InnerDetector/InDetEventCnv/InDetPrepRawDataToxAOD/src/TRT_PrepDataToxAOD.cxx +++ b/InnerDetector/InDetEventCnv/InDetPrepRawDataToxAOD/src/TRT_PrepDataToxAOD.cxx @@ -208,9 +208,9 @@ StatusCode TRT_PrepDataToxAOD::execute() const Amg::MatrixX& localCov = prd->localCovariance(); if(localCov.size() == 1){ - xprd->setLocalPositionError( localCov(0,0), 0., 0. ); + xprd->setLocalPositionError( sqrt(localCov(0,0)), 0., 0. ); } else if(localCov.size() == 2){ - xprd->setLocalPositionError( localCov(0,0), localCov(1,1), localCov(0,1) ); + xprd->setLocalPositionError( sqrt(localCov(0,0)), sqrt(localCov(1,1)), sqrt(localCov(0,1)) ); } else { xprd->setLocalPositionError(0.,0.,0.); } -- GitLab