Update AsgLeptonTrackSelectionAlg to write out d0 and z0sinTheta variables

Update AsgLeptonTrackSelectionAlg to write out d0 and z0sinTheta variables

PhysicsAnalysis/Algorithms/AnalysisAlgorithmsConfig/data/for_compare.yaml

@@ -30,6 +30,8 @@ Electrons:
         likelihoodWP: 'LooseBLayerLH'
         isolationWP: 'Loose_VarRad'
         recomputeLikelihood: False
+        writeTrackD0Z0: True
+
     PtEtaSelection:
         minPt: 10000.0
 
@@ -54,6 +56,7 @@ Muons:
         quality: 'Medium'
         isolation: 'Loose_VarRad'
         onlyRecoEffSF: True
+        writeTrackD0Z0: True
     PtEtaSelection: {}
 
 TauJets:

PhysicsAnalysis/Algorithms/AnalysisAlgorithmsConfig/python/ConfigText_unitTest.py

@@ -99,6 +99,7 @@ def compareTextBuilder(yamlPath='', *, checkOrder=False) :
     config.setOptions (likelihoodWP='LooseBLayerLH')
     config.setOptions (isolationWP='Loose_VarRad')
     config.setOptions (recomputeLikelihood=False)
+    config.setOptions (writeTrackD0Z0=True)
     # Electrons.PtEtaSelection
     config.addBlock ('Electrons.PtEtaSelection')
     config.setOptions (containerName='AnaElectrons')
@@ -133,6 +134,7 @@ def compareTextBuilder(yamlPath='', *, checkOrder=False) :
     config.setOptions (quality='Medium')
     config.setOptions (isolation='Loose_VarRad')
     config.setOptions (onlyRecoEffSF=True)
+    config.setOptions (writeTrackD0Z0=True)
     # Muons.PtEtaSelection
     config.addBlock ('Muons.PtEtaSelection')
     config.setOptions (containerName='AnaMuons')

PhysicsAnalysis/Algorithms/AnalysisAlgorithmsConfig/python/FullCPAlgorithmsTest.py

@@ -791,6 +791,7 @@ def makeSequenceBlocks (dataType, algSeq, forCompare, isPhyslite,
         configSeq.setOptionValue ('.likelihoodWP', 'LooseDNN')
     configSeq.setOptionValue ('.isolationWP', 'Loose_VarRad')
     configSeq.setOptionValue ('.recomputeLikelihood', recomputeLikelihood)
+    configSeq.setOptionValue ('.writeTrackD0Z0', True)
 
     configSeq += config.makeConfig ('Electrons.PtEtaSelection',
         containerName='AnaElectrons')
@@ -835,6 +836,8 @@ def makeSequenceBlocks (dataType, algSeq, forCompare, isPhyslite,
     configSeq.setOptionValue ('.isolation', 'Loose_VarRad')
     if forCompare :
         configSeq.setOptionValue ('.onlyRecoEffSF', True)
+    configSeq.setOptionValue ('.writeTrackD0Z0', True)
+
     # TODO: MCP should restore this when the recommendations for Tight WP exist in R23
     # configSeq += config.makeConfig ('Muons.Selection', 'AnaMuons.tight')
     # configSeq.setOptionValue ('.quality', 'Tight')
@@ -846,7 +849,6 @@ def makeSequenceBlocks (dataType, algSeq, forCompare, isPhyslite,
     configSeq.setOptionValue ('.minPt', muonMinPt)
     configSeq.setOptionValue ('.maxEta', muonMaxEta)
 
-
     # Include, and then set up the tau analysis algorithm sequence:
     configSeq += config.makeConfig ('TauJets',
         containerName='AnaTauJets')

PhysicsAnalysis/Algorithms/AsgAnalysisAlgorithms/AsgAnalysisAlgorithms/AsgLeptonTrackSelectionAlg.h

@@ -57,6 +57,9 @@ namespace CP
     Gaudi::Property<int> m_nMaxPixelHits {this, "nMaxPixelHits", -1, "maximum number of required Pixel hits (or -1 for no cut)"};
     Gaudi::Property<int> m_nMinSCTHits {this, "nMinSCTHits", -1, "minimum number of required SCT hits (or -1 for no cut)"};
     Gaudi::Property<int> m_nMaxSCTHits {this, "nMaxSCTHits", -1, "maximum number of required SCT hits (or -1 for no cut)"};
+    Gaudi::Property<bool> m_decorateTTVAVars{this, "decorateTTVAVars", false, "save the calculated d0sig and z0sinTheta variables"};
+    Gaudi::Property<std::string> m_d0sigDecoration {this, "d0sigDecoration", "", "the decoration name for d0 significance"};
+    Gaudi::Property<std::string> m_z0sinthetaDecoration {this, "z0sinthetaDecoration", "", "the decoration name for z0sintheta"};
 
     /// \}
 
@@ -90,6 +93,13 @@ namespace CP
   private:
     ServiceHandle<ISelectionNameSvc> m_nameSvc {"SelectionNameSvc", "AsgLeptonTrackSelectionAlg"};
 
+    /// \brief the name of the variable being decorated for d0significance
+  private:
+    std::unique_ptr<const SG::AuxElement::Decorator<float> > m_d0sigDecorator {};
+
+    /// \brief the name of the variable being decorated for z0sintheta
+  private:
+   std::unique_ptr<const SG::AuxElement::Decorator<float> > m_z0sinthetaDecorator {};
 
     /// \brief the \ref asg::AcceptInfo we are using
   private:

PhysicsAnalysis/Algorithms/AsgAnalysisAlgorithms/Root/AsgLeptonTrackSelectionAlg.cxx

@@ -48,6 +48,28 @@ namespace CP
     if (m_nMinSCTHits != -1 || m_nMaxSCTHits != -1)
       m_accept.addCut ("numSCTHits", "Minimum and/or maxiumum SCT hits");
 
+    if(m_decorateTTVAVars){
+      if (m_d0sigDecoration.empty()){
+        ANA_MSG_ERROR ("No d0significance decoration name set");
+        return StatusCode::FAILURE;
+      } else {
+        m_d0sigDecorator = std::make_unique<SG::AuxElement::Decorator<float> > (m_d0sigDecoration);
+      }
+      if (m_z0sinthetaDecoration.empty()){
+        ANA_MSG_ERROR ("No z0sintheta decoration name set");
+        return StatusCode::FAILURE;
+      } else {
+        m_z0sinthetaDecorator = std::make_unique<SG::AuxElement::Decorator<float> > (m_z0sinthetaDecoration);
+      }
+    } else{
+      if (!m_d0sigDecoration.empty()){
+        ANA_MSG_WARNING ("d0significance decoration name set, please set decorateTTVVars to True if you want this to be written out");
+      }
+      if (!m_z0sinthetaDecoration.empty()){
+        ANA_MSG_WARNING ("z0sintheta decoration name set, please set decorateTTVVars to True if you want this to be written out");
+      }
+    }
+
     ANA_CHECK (m_particlesHandle.initialize (m_systematicsList));
     ANA_CHECK (m_preselection.initialize (m_systematicsList, m_particlesHandle, SG::AllowEmpty));
     ANA_CHECK (m_selectionHandle.initialize (m_systematicsList, m_particlesHandle));
@@ -72,28 +94,25 @@ namespace CP
   execute ()
   {
     SG::ReadHandle<xAOD::EventInfo> eventInfo(m_eventInfoKey);
-
+    SG::ReadHandle<xAOD::VertexContainer> vertices(m_primaryVerticesKey);
     const xAOD::Vertex *primaryVertex {nullptr};
-    if (m_maxDeltaZ0SinTheta > 0)
+
+    for (const xAOD::Vertex *vertex : *vertices)
     {
-      SG::ReadHandle<xAOD::VertexContainer> vertices(m_primaryVerticesKey);
-      for (const xAOD::Vertex *vertex : *vertices)
+      if (vertex->vertexType() == xAOD::VxType::PriVtx)
       {
-        if (vertex->vertexType() == xAOD::VxType::PriVtx)
+        // The default "PrimaryVertex" container is ordered in
+        // sum-pt, and the tracking group recommends to pick the one
+        // with the maximum sum-pt, so this will do the right thing.
+        // If the user needs a different primary vertex, they need to
+        // provide a reordered primary vertex container and point
+        // this algorithm to it.  Currently there is no central
+        // algorithm to do that, so users will have to write their
+        // own (15 Aug 18).
+        if (primaryVertex == nullptr)
         {
-          // The default "PrimaryVertex" container is ordered in
-          // sum-pt, and the tracking group recommends to pick the one
-          // with the maximum sum-pt, so this will do the right thing.
-          // If the user needs a different primary vertex, he needs to
-          // provide a reordered primary vertex container and point
-          // this algorithm to it.  Currently there is no central
-          // algorithm to do that, so users will have to write their
-          // own (15 Aug 18).
-          if (primaryVertex == nullptr)
-          {
-            primaryVertex = vertex;
-            break;
-          }
+          primaryVertex = vertex;
+          break;
         }
       }
     }
@@ -105,44 +124,38 @@ namespace CP
       for (const xAOD::IParticle *particle : *particles)
       {
         asg::AcceptData acceptData (&m_accept);
+        float d0sig = -999;
+        float deltaZ0SinTheta = -999;
 
         if (m_preselection.getBool (*particle, sys))
         {
           std::size_t cutIndex {0};
-
+          
           const xAOD::TrackParticle *track {nullptr};
-          if (const xAOD::Muon *muon = dynamic_cast<const xAOD::Muon *>(particle))
+          if (const xAOD::Muon *muon = dynamic_cast<const xAOD::Muon *>(particle)){
             track = muon->primaryTrackParticle();
-          else if (const xAOD::Electron *electron = dynamic_cast<const xAOD::Electron *>(particle))
+          } else if (const xAOD::Electron *electron = dynamic_cast<const xAOD::Electron *>(particle)){
             track = electron->trackParticle();
-          else
-          {
+          } else {
             ANA_MSG_ERROR ("failed to cast input to electron or muon");
             return StatusCode::FAILURE;
           }
 
           acceptData.setCutResult (cutIndex ++, track != nullptr);
-          if (track != nullptr)
-          {
-            if (m_maxD0Significance > 0)
-            {
-              try
-              {
-                const float d0sig = xAOD::TrackingHelpers::d0significance
-                  (track, eventInfo->beamPosSigmaX(), eventInfo->beamPosSigmaY(),
-                   eventInfo->beamPosSigmaXY());
-                acceptData.setCutResult (cutIndex ++, fabs( d0sig ) < m_maxD0Significance);
-              } catch (const std::runtime_error &) {
-                acceptData.setCutResult (cutIndex ++, false);
-              }
-            }
-            if (m_maxDeltaZ0SinTheta > 0)
-            {
-              const double vertex_z = primaryVertex ? primaryVertex->z() : 0;
-              const float deltaZ0SinTheta
-                = (track->z0() + track->vz() - vertex_z) * sin (particle->p4().Theta());
-              acceptData.setCutResult (cutIndex ++, fabs (deltaZ0SinTheta) < m_maxDeltaZ0SinTheta);
+          
+          if (track != nullptr) {
+            try {
+              d0sig = xAOD::TrackingHelpers::d0significance(track, eventInfo->beamPosSigmaX(), eventInfo->beamPosSigmaY(), eventInfo->beamPosSigmaXY());
+              if (m_maxD0Significance > 0) acceptData.setCutResult (cutIndex ++, fabs( d0sig ) < m_maxD0Significance);
+            
+            } catch (const std::runtime_error &) {
+              acceptData.setCutResult (cutIndex ++, false);
             }
+            
+            const double vertex_z = primaryVertex ? primaryVertex->z() : 0;
+            deltaZ0SinTheta = (track->z0() + track->vz() - vertex_z) * sin (particle->p4().Theta());
+            if (m_maxDeltaZ0SinTheta > 0) acceptData.setCutResult (cutIndex ++, fabs (deltaZ0SinTheta) < m_maxDeltaZ0SinTheta);
+
             if (m_nMinPixelHits != -1 || m_nMaxPixelHits != -1) {
               uint8_t nPixelHits;
               track->summaryValue(nPixelHits, xAOD::numberOfPixelHits);
@@ -169,9 +182,11 @@ namespace CP
             }
           }
         }
-
-        m_selectionHandle.setBits
-          (*particle, selectionFromAccept (acceptData), sys);
+        if(m_decorateTTVAVars){
+          (*m_z0sinthetaDecorator)(*particle) = deltaZ0SinTheta;
+          (*m_d0sigDecorator)(*particle) = d0sig;
+        }
+        m_selectionHandle.setBits(*particle, selectionFromAccept (acceptData), sys);
       }
     }
 

PhysicsAnalysis/Algorithms/EgammaAnalysisAlgorithms/python/ElectronAnalysisConfig.py

@@ -207,6 +207,8 @@ class ElectronWorkingPointConfig (ConfigBlock) :
         self.addOption ('maxDeltaZ0SinTheta', 0.5, type=float,
             info="maximum z0sinTheta in mm used for the trackSelection"
             "The default is 0.5 mm")
+        self.addOption ('writeTrackD0Z0', False, type = bool,
+            info="save the d0 significance and z0sinTheta variables so they can be written out")
         self.addOption ('likelihoodWP', None, type=str,
             info="the ID WP (string) to use. Supported ID WPs: TightLH, "
             "MediumLH, LooseBLayerLH. ")
@@ -255,16 +257,23 @@ class ElectronWorkingPointConfig (ConfigBlock) :
             postfix = '_' + postfix
 
         # Set up the track selection algorithm:
-        if self.trackSelection :
+        if self.writeTrackD0Z0 or self.trackSelection :
             alg = config.createAlgorithm( 'CP::AsgLeptonTrackSelectionAlg',
                                         'ElectronTrackSelectionAlg' + postfix )
             alg.selectionDecoration = 'trackSelection' + postfix + ',as_bits'
             alg.maxD0Significance = self.maxD0Significance
-            alg.maxDeltaZ0SinTheta = self.maxDeltaZ0SinTheta
+            alg.maxDeltaZ0SinTheta = self.maxDeltaZ0SinTheta 
+            alg.decorateTTVAVars = self.writeTrackD0Z0
+            alg.d0sigDecoration = 'd0sig' + postfix
+            alg.z0sinthetaDecoration = 'z0sintheta' + postfix
             alg.particles = config.readName (self.containerName)
             alg.preselection = config.getPreselection (self.containerName, '')
-            config.addSelection (self.containerName, self.selectionName, alg.selectionDecoration)
-
+            if self.trackSelection : 
+                config.addSelection (self.containerName, self.selectionName, alg.selectionDecoration)
+            if self.writeTrackD0Z0 :
+                config.addOutputVar (self.containerName, alg.d0sigDecoration, alg.d0sigDecoration,noSys=True)
+                config.addOutputVar (self.containerName, alg.z0sinthetaDecoration, alg.z0sinthetaDecoration,noSys=True)
+                
         if 'LH' in self.likelihoodWP:
             # Set up the likelihood ID selection algorithm
             # It is safe to do this before calibration, as the cluster E is used

PhysicsAnalysis/Algorithms/MuonAnalysisAlgorithms/python/MuonAnalysisConfig.py

@@ -130,6 +130,8 @@ class MuonWorkingPointConfig (ConfigBlock) :
         self.addOption ('maxDeltaZ0SinTheta', 0.5, type=float,
             info="maximum Delta z0sinTheta in mm used for the trackSelection"
             "The default is 0.5 mm")
+        self.addOption ('writeTrackD0Z0', False, type = bool,
+            info="save the d0 significance and z0sinTheta variables so they can be written out")
         self.addOption ('quality', None, type=str,
             info="the ID WP (string) to use. Supported ID WPs: Tight, Medium, "
             "Loose, LowPt, HighPt.")
@@ -183,15 +185,22 @@ class MuonWorkingPointConfig (ConfigBlock) :
             postfix = '_' + postfix
 
         # Set up the track selection algorithm:
-        if self.trackSelection :
+        if self.writeTrackD0Z0 or self.trackSelection:
             alg = config.createAlgorithm( 'CP::AsgLeptonTrackSelectionAlg',
                                 'MuonTrackSelectionAlg' + postfix )
             alg.selectionDecoration = 'trackSelection' + postfix + ',as_bits'
+            alg.decorateTTVAVars = self.writeTrackD0Z0
+            alg.d0sigDecoration = 'd0sig' + postfix
+            alg.z0sinthetaDecoration = 'z0sintheta' + postfix
             alg.maxD0Significance = self.maxD0Significance
             alg.maxDeltaZ0SinTheta = self.maxDeltaZ0SinTheta
             alg.particles = config.readName (self.containerName)
             alg.preselection = config.getPreselection (self.containerName, '')
-            config.addSelection (self.containerName, self.selectionName, alg.selectionDecoration, preselection=self.qualitySelectionOutput)
+            if self.trackSelection :
+                config.addSelection (self.containerName, self.selectionName, alg.selectionDecoration, preselection=self.qualitySelectionOutput)
+            if self.writeTrackD0Z0 :
+                config.addOutputVar (self.containerName, alg.d0sigDecoration, alg.d0sigDecoration,noSys=True)
+                config.addOutputVar (self.containerName, alg.z0sinthetaDecoration, alg.z0sinthetaDecoration,noSys=True)
 
         # Setup the muon quality selection
         alg = config.createAlgorithm( 'CP::MuonSelectionAlgV2',