Vertex Info Added:

Addition of th e secondary vertex information to the output stream and a tidy/comment of code

Reconstruction/VKalVrt/GNNVertexConstructor/GNNVertexConstructor/GNNVertexConstructorTool.h

@@ -24,6 +24,7 @@ namespace Trk {
 class TrkVKalVrtFitter;
 class IVertexFitter;
 class IVKalState;
+class VxSecVKalVertexInfo;
 } // namespace Trk
 
 namespace Rec {
@@ -31,7 +32,6 @@ namespace Rec {
 struct workVectorArrxAOD {
   std::vector<const xAOD::TrackParticle *> listSelTracks; // Selected tracks after quality cuts
   std::vector<const xAOD::TrackParticle *> tmpListTracks;
-  std::vector<const xAOD::TrackParticle *> inpTrk; // All tracks provided to tool
   double beamX = 0.;
   double beamY = 0.;
   double beamZ = 0.;
@@ -50,7 +50,10 @@ public:
   StatusCode finalize();
 
   virtual StatusCode decorateJets(const xAOD::JetContainer*) const;
-  virtual StatusCode performVertexFit(const xAOD::JetContainer*, xAOD::VertexContainer*, const EventContext&) const;
+  virtual StatusCode performVertexFit(const xAOD::JetContainer*, 
+                                      xAOD::VertexContainer*, 
+                                      const xAOD::Vertex & primaryVertex, 
+                                      const EventContext&) const;
 
   // Tools
   ToolHandle<FlavorTagDiscriminants::GNNTool> m_gnn_Tool{this, "gnn_Tool", "",
@@ -72,10 +75,12 @@ public:
   // Conditions
   SG::ReadCondHandleKey<InDet::BeamSpotData> m_beamSpotKey{this, "BeamSpotKey", "BeamSpotData",
                                                            "SG key for beam spot"};
-
+  
+  //Access the Primary Vertex Info
+  const xAOD::Vertex* m_thePV;
+  
 private:
-  double m_w_1{};
-  float m_chiScale[11]{};
+  
   std::string m_jetCollection;
 
   struct WrkVrt {
@@ -88,10 +93,7 @@ private:
     std::vector<double> chi2PerTrk;
     std::vector<std::vector<double>> trkAtVrt;
     double chi2{};
-    double projectedVrt = 0.;
-    int detachedTrack = -1;
-    double BDT = 1.1;
-  };
+    };//end WrkVrt
 };
 } // namespace Rec
 

Reconstruction/VKalVrt/GNNVertexConstructor/GNNVertexConstructor/IGNNVertexConstructorInterface.h

@@ -14,7 +14,10 @@ namespace Rec {
     public:
       static const InterfaceID& interfaceID() { return IID_IGNNVertexConstructorInterface;}
 
-      virtual StatusCode performVertexFit( const xAOD::JetContainer* jetCont, xAOD::VertexContainer* vertexCont, const EventContext& ctx ) const = 0;
+      virtual StatusCode performVertexFit( const xAOD::JetContainer* jetCont, 
+                                           xAOD::VertexContainer* vertexCont, 
+                                           const xAOD::Vertex & primaryVertex, 
+                                           const EventContext& ctx ) const = 0;
       
   };
 

Reconstruction/VKalVrt/GNNVertexConstructor/Root/GNNVertexConstructorTool.cxx

@@ -18,7 +18,8 @@ GNNVertexConstructorTool::GNNVertexConstructorTool(const std::string &type, cons
                                                    const IInterface *parent)
     : AthAlgTool(type, name, parent),
       m_vertexFitterTool("Trk::TrkVKalVrtFitter/VertexFitterTool", this),
-      m_jetCollection("AntiKt4EMPFlowJets") {
+      m_jetCollection("AntiKt4EMPFlowJets")
+      /*m_thePV(nullptr)*/{
   declareInterface<IGNNVertexConstructorInterface>(this);
 
   declareProperty("JetTrackLinks", m_trackLinksKey = "BTagging_AntiKt4EMPFlowAuxDyn.TrackLinks");
@@ -30,7 +31,6 @@ GNNVertexConstructorTool::GNNVertexConstructorTool(const std::string &type, cons
 }
 
 /* Destructor */
-
 GNNVertexConstructorTool::~GNNVertexConstructorTool() {
   ATH_MSG_DEBUG("GNNVertexConstructorTool destructor called");
 }
@@ -50,16 +50,15 @@ StatusCode GNNVertexConstructorTool::initialize() {
   ATH_CHECK(m_gnn_Tool.retrieve());
   ATH_CHECK(m_vertexFitterTool.retrieve());
 
+  //Additional Info for Vertex Fit
   ATH_CHECK(m_beamSpotKey.initialize());
-
   ATH_CHECK(m_eventInfoKey.initialize());
 
-  // @Luke TODO: what is this parameter?
-  m_w_1 = 1.;
-
   return StatusCode::SUCCESS;
 }
 
+//Use GNN to decorate the tracks
+//May be removed in future
 StatusCode GNNVertexConstructorTool::decorateJets(const xAOD::JetContainer *jetCont) const {
 
   for (auto jet : *jetCont) {
@@ -68,55 +67,67 @@ StatusCode GNNVertexConstructorTool::decorateJets(const xAOD::JetContainer *jetC
   return StatusCode::SUCCESS;
 }
 
+//Perform Vertex fit using the jet decorations of the GNN
 StatusCode GNNVertexConstructorTool::performVertexFit(const xAOD::JetContainer *inJetContainer,
                             xAOD::VertexContainer *outVertexContainer,
+                            const xAOD::Vertex & primVrt,
                             const EventContext &ctx) const {
 
   using TLC = std::vector<ElementLink<DataVector<xAOD::TrackParticle_v1>>>;
   using TL  = ElementLink<DataVector<xAOD::TrackParticle_v1>>;
 
+  //Read Decor Handle for Track links and Vertex links
   SG::ReadDecorHandle<xAOD::JetContainer, TLC> trackLinksHandle(m_trackLinksKey, ctx);
   SG::ReadDecorHandle<xAOD::JetContainer, std::vector<char, std::allocator<char>>>
       vertexLinksHandle(m_vertexLinksKey, ctx);
   SG::WriteDecorHandle< xAOD::JetContainer, std::vector<ElementLink<xAOD::VertexContainer>>> 
       jetWriteDecorHandleVertexLink (m_jetWriteDecorKeyVertexLink, ctx);
+  
   // Vertex decorators
   SG::AuxElement::Decorator<float> decor_mass("mass");
   SG::AuxElement::Decorator<float> decor_pT("pt");
   SG::AuxElement::Decorator<float> decor_charge("charge");
   SG::AuxElement::Decorator<float> decor_vPos("vPos");
-
+  SG::AuxElement::Decorator<float> decor_Lxy("Lxy");
+  SG::AuxElement::Decorator<float> decor_significance3d("significance3d");
+  SG::AuxElement::Decorator<float> decor_deltaR("deltaR");
+  SG::AuxElement::Decorator<float> decor_NGTinSvx("NGTinSvx");
+  SG::AuxElement::Decorator<float> decor_L3D("L3d");
+  SG::AuxElement::Decorator<float> decor_N2Tpair("N2Tpair");
+  //SG::AuxElement::Decorator<float> decor_efracsv("efracsv");
+  
   // Create a map of track links and track vertexing values (Using mutlimap)
   std::multimap<int, TL> vertexMap;
 
   // Loop over the jets
   for (const auto &jet : *inJetContainer) {
-
+    
+    //Ensure map is empty from previous iterations
     vertexMap.clear();
 
+    //Retrieve the Vertex and Track Collections
     auto vertexCollection = vertexLinksHandle(*jet);
     auto trackCollection = trackLinksHandle(*jet);
 
+    //Fill the map
     int i = 0;
-
     for (auto v : vertexCollection) {
       vertexMap.insert(std::pair<int, TL>(v, (trackCollection[i])));
       i++;
     }
     std::multimap<int, TL>::iterator itr;
 
+    //Working xAOD
     workVectorArrxAOD *xAODwrk = new workVectorArrxAOD();
     SG::ReadCondHandle<InDet::BeamSpotData> beamSpotHandle{m_beamSpotKey, ctx};
 
+    //Beam Conditions
     xAODwrk->beamX = beamSpotHandle->beamPos().x();
     xAODwrk->beamY = beamSpotHandle->beamPos().y();
     xAODwrk->beamZ = beamSpotHandle->beamPos().z();
-
     xAODwrk->tanBeamTiltX = tan(beamSpotHandle->beamTilt(0));
     xAODwrk->tanBeamTiltY = tan(beamSpotHandle->beamTilt(1));
 
-    std::vector<xAOD::Vertex *> finalVertices(0);
-
     std::unique_ptr<std::vector<WrkVrt>> wrkVrtSet = std::make_unique<std::vector<WrkVrt>>();
     WrkVrt newvrt;
     newvrt.Good = true;
@@ -133,26 +144,41 @@ StatusCode GNNVertexConstructorTool::performVertexFit(const xAOD::JetContainer *
         // Need at least 2 tracks to perform a fit
         auto elements = vertexMap.equal_range(k);
 
+        //Retrieve the tracks and push to working xAOD
         for (auto i = elements.first; i != elements.second; ++i) {
           xAODwrk->listSelTracks.push_back(*(i->second));
         }
-
+        
+        //Perform the Vertex Fit
         StatusCode sc = (m_vertexFitterTool->VKalVrtFit(
             xAODwrk->listSelTracks, neutralPartDummy, newvrt.vertex, newvrt.vertexMom,
             newvrt.vertexCharge, newvrt.vertexCov, newvrt.chi2PerTrk, newvrt.trkAtVrt, newvrt.chi2,
             *state, false));
         if (sc.isFailure())
           continue;
-
+        
+        //Chi2 Cut
+        if (newvrt.chi2<20){
         ATH_MSG_DEBUG("Found IniVertex=" << newvrt.vertex[0] << ", " << newvrt.vertex[1] << ", "
-                                         << newvrt.vertex[2]);
-
-        Amg::Vector3D vDist = newvrt.vertex; // - m_thePV->position();
-
-        double vPos = (vDist.x() * newvrt.vertexMom.Px() + vDist.y() * newvrt.vertexMom.Py() +
-                       vDist.z() * newvrt.vertexMom.Pz()) /
+                                         << newvrt.vertex[2] << " trks " << newvrt.trkAtVrt.size());
+
+        Amg::Vector3D vDir = newvrt.vertex - primVrt.position();  //Vertex Dirction in relation to Primary
+        
+        Amg::Vector3D jetDir(jet->p4().Px(),jet->p4().Py(),jet->p4().Pz()); //Jet Direction
+        
+        
+        double vPos = (vDir.x() * newvrt.vertexMom.Px() + vDir.y() * newvrt.vertexMom.Py() +
+                       vDir.z() * newvrt.vertexMom.Pz()) /
                       newvrt.vertexMom.Rho();
-
+        
+        double L3D =sqrt(vDir[0]*vDir[0]+vDir[1]*vDir[1]+vDir[2]*vDir[2]);
+        ATH_MSG_DEBUG("L3D" << L3D);
+        
+        
+        double drJPVSV = Amg::deltaR(jetDir,vDir); //DeltaR
+        
+        int NGTatVtx=newvrt.trkAtVrt.size();
+        
         xAOD::Vertex *GNNvertex = new xAOD::Vertex;
         outVertexContainer->emplace_back(GNNvertex);
 
@@ -164,21 +190,33 @@ StatusCode GNNVertexConstructorTool::performVertexFit(const xAOD::JetContainer *
           // Register the link to the vertex
           GNNvertex->addTrackAtVertex( link_trk, 1. );
         }
-
+       
+        //Add Vertex Info into Container        
         GNNvertex->setVertexType(xAOD::VxType::SecVtx);
         GNNvertex->setPosition(newvrt.vertex);
         GNNvertex->setFitQuality(newvrt.chi2, 1);
 
-        decor_mass(*GNNvertex) = newvrt.vertexMom.M();
-        decor_pT(*GNNvertex) = newvrt.vertexMom.Perp();
-        decor_charge(*GNNvertex) = newvrt.vertexCharge;
-        decor_vPos(*GNNvertex) = vPos;
-
+        decor_mass(*GNNvertex)            = newvrt.vertexMom.M();
+        decor_pT(*GNNvertex)              = newvrt.vertexMom.Perp();
+        decor_charge(*GNNvertex)          = newvrt.vertexCharge;
+        decor_vPos(*GNNvertex)            = vPos;
+        decor_Lxy(*GNNvertex)             = sqrt(vDir[0]*vDir[0]+vDir[1]*vDir[1]);
+        decor_L3D(*GNNvertex)             = L3D;
+        decor_significance3d(*GNNvertex)  = L3D/newvrt.chi2;
+        decor_NGTinSvx(*GNNvertex)        = NGTatVtx;
+        decor_deltaR(*GNNvertex)          = drJPVSV;
+        
+        if (newvrt.trkAtVrt.size()==2){
+          decor_N2Tpair(*GNNvertex)=newvrt.trkAtVrt.size();
+        }
+        
         ElementLink< xAOD::VertexContainer> linkVertex;
         linkVertex.setElement(GNNvertex);
         linkVertex.setStorableObject(*outVertexContainer);
         jetWriteDecorHandleVertexLink(*jet).push_back(linkVertex);
-      }
+      
+      }//end of Chi2 cut
+      }//end of 2 Track requirement
     }
     delete xAODwrk;
   } // end loop over jets

Reconstruction/VKalVrt/GNNVertexConstructor/src/GNNVertexConstructorAlg.cxx

@@ -19,7 +19,8 @@ StatusCode GNNVertexConstructorAlg::initialize() {
   // Initializing Keys
   ATH_CHECK(m_inJetsKey.initialize());
   ATH_CHECK(m_outVertexKey.initialize());
-
+  ATH_CHECK( m_pvContainerKey.initialize() );
+  
   return StatusCode::SUCCESS;
 }
 
@@ -27,12 +28,14 @@ StatusCode GNNVertexConstructorAlg::execute(const EventContext &ctx) const {
 
   ATH_MSG_DEBUG("In GNNVertexConstructorAlg::execute()");
 
+  //Extract Jets
   SG::ReadHandle<xAOD::JetContainer> inJetContainer(m_inJetsKey, ctx);
   if (!inJetContainer.isValid()) {
     ATH_MSG_WARNING("No xAOD::JetContainer named " << m_inJetsKey.key() << " found in StoreGate");
     return StatusCode::FAILURE;
   }
-
+  
+  //Write new GNN Vertice Container
   SG::WriteHandle<xAOD::VertexContainer> outVertexContainer(m_outVertexKey, ctx);
   if (outVertexContainer
           .record(std::make_unique<xAOD::VertexContainer>(),
@@ -42,8 +45,24 @@ StatusCode GNNVertexConstructorAlg::execute(const EventContext &ctx) const {
     return StatusCode::FAILURE;
   }
 
+  const xAOD::Vertex* pv = nullptr;
+  //-- Extract Primary Vertices
+   SG::ReadHandle<xAOD::VertexContainer> pv_cont(m_pvContainerKey, ctx);
+   if ( !pv_cont.isValid() ) {
+     ATH_MSG_WARNING( "No Primary Vertices container found in TDS" );
+   }else{
+     //-- Extract PV itself
+     for ( auto v : *pv_cont ) {
+       if (v->vertexType()==xAOD::VxType::PriVtx) {    pv = v;   break; }
+     }
+   }
+
+  //Decorates the Jets using the GNN model
+  //May b removed in future with development of GNN model
   ATH_CHECK(m_VtxTool->decorateJets(inJetContainer.ptr()));
-  ATH_CHECK(m_VtxTool->performVertexFit(inJetContainer.ptr(), outVertexContainer.ptr(), ctx));
+  
+  //Perform a Vertex fit
+  ATH_CHECK(m_VtxTool->performVertexFit(inJetContainer.ptr(), outVertexContainer.ptr(), *pv,  ctx));
 
   return StatusCode::SUCCESS;
 }

Reconstruction/VKalVrt/GNNVertexConstructor/src/GNNVertexConstructorAlg.h

@@ -17,13 +17,18 @@ public:
 
 private:
   // Tools
-  ToolHandle<Rec::GNNVertexConstructorTool> m_VtxTool;
+  ToolHandle<Rec::GNNVertexConstructorTool>   m_VtxTool;
   // Input jets
-  SG::ReadHandleKey<xAOD::JetContainer> m_inJetsKey{
+  SG::ReadHandleKey<xAOD::JetContainer>       m_inJetsKey{
       this, "inputJetContainer", "AntiKt4EMPFlowJets", "Input jet container"};
   // Output vertices
-  SG::WriteHandleKey<xAOD::VertexContainer> m_outVertexKey{
+  SG::WriteHandleKey<xAOD::VertexContainer>   m_outVertexKey{
       this, "outputVertexContainer", "GNNVertices", "Output vertex container"};
+  
+  //Input Primary Vertices
+  SG::ReadHandleKey<xAOD::VertexContainer>    m_pvContainerKey{
+      this,"PrimaryVertexContainer","PrimaryVertices","Read PrimaryVertices container"};
+
 };
 
 } // namespace Rec