From 42a7984a0f6211d4700f70680972e5f23be4a04d Mon Sep 17 00:00:00 2001
From: Vadim Kostyukhin <vkost@cern.ch>
Date: Tue, 13 Dec 2016 17:01:00 +0100
Subject: [PATCH] Increase high pt jet SV efficiency
(InDetVKalVxInJetTool-00-08-13)
* Increase high pt jet SV efficiency
* Retune pseudo-vertex finder
* InDetVKalVxInJetTool-00-08-13
Former-commit-id: e1de95e31856e0841684466a04d6a33497234802
---
.../InDetVKalVxInJetTool.h | 6 +-
.../InDetVKalVxInJetTool/src/BTagVrtSec.cxx | 62 ++++++++++++-------
.../src/InDetVKalVxInJetTool.cxx | 4 +-
.../InDetVKalVxInJetTool/src/Utilities.cxx | 47 +++++++++++---
4 files changed, 86 insertions(+), 33 deletions(-)
diff --git a/InnerDetector/InDetRecTools/InDetVKalVxInJetTool/InDetVKalVxInJetTool/InDetVKalVxInJetTool.h b/InnerDetector/InDetRecTools/InDetVKalVxInJetTool/InDetVKalVxInJetTool/InDetVKalVxInJetTool.h
index b4aae8b21ad..47033a0102b 100755
--- a/InnerDetector/InDetRecTools/InDetVKalVxInJetTool/InDetVKalVxInJetTool/InDetVKalVxInJetTool.h
+++ b/InnerDetector/InDetRecTools/InDetVKalVxInJetTool/InDetVKalVxInJetTool/InDetVKalVxInJetTool.h
@@ -214,6 +214,7 @@ namespace InDet {
double m_AntiFake2trVrtCut;
double m_JetPtFractionCut;
int m_TrackInJetNumberLimit;
+ double m_MaterialPtCut;
double m_pseudoSigCut;
bool m_FillHist;
@@ -355,7 +356,8 @@ namespace InDet {
TLorentzVector GetBDir( const xAOD::TrackParticle* trk1,
const xAOD::TrackParticle* trk2,
- const xAOD::Vertex & PrimVrt) const;
+ const xAOD::Vertex & PrimVrt,
+ Amg::Vector3D &V1, Amg::Vector3D &V2) const;
int nTrkCommon( std::vector<WrkVrt> *WrkVrtSet, int V1, int V2) const;
double minVrtVrtDist( std::vector<WrkVrt> *WrkVrtSet, int & V1, int & V2) const;
@@ -469,6 +471,8 @@ namespace InDet {
void getPixelLayers(const xAOD::TrackParticle* Part, int &blHit, int &l1Hit, int &l2Hit, int &nLay) const;
void getPixelDiscs(const xAOD::TrackParticle* Part, int &d0Hit, int &d1Hit, int &d2Hit) const;
void getPixelDiscs(const Rec::TrackParticle* Part, int &d0Hit, int &d1Hit, int &d2Hit) const;
+ void getPixelProblems(const xAOD::TrackParticle* Part, int &splshIBL, int &splshBL ) const;
+ void getPixelProblems(const Rec::TrackParticle* Part, int &splshIBL, int &splshBL ) const;
StatusCode VKalVrtFitBase(const std::vector<const Rec::TrackParticle*> & listPart,
diff --git a/InnerDetector/InDetRecTools/InDetVKalVxInJetTool/src/BTagVrtSec.cxx b/InnerDetector/InDetRecTools/InDetVKalVxInJetTool/src/BTagVrtSec.cxx
index 56fc331375c..ae63a4d56f9 100755
--- a/InnerDetector/InDetRecTools/InDetVKalVxInJetTool/src/BTagVrtSec.cxx
+++ b/InnerDetector/InDetRecTools/InDetVKalVxInJetTool/src/BTagVrtSec.cxx
@@ -4,6 +4,7 @@
// Header include
#include "InDetVKalVxInJetTool/InDetVKalVxInJetTool.h"
+#include "GeoPrimitives/GeoPrimitivesHelpers.h"
//-------------------------------------------------
// Other stuff
#include "AnalysisUtils/AnalysisMisc.h"
@@ -299,7 +300,7 @@ namespace InDet{
long int tChrg=0;
double Chi2=0.;
std::vector<const xAOD::TrackParticle*> tTrkForFit(2,0);
- std::vector<float> tmpCov(15,0.); tmpCov[0]=1e-4; tmpCov[2]=4.e-4; tmpCov[5]=1.e-4; tmpCov[9]=1.e-4; tmpCov[14]=1.e-12;
+ std::vector<float> tmpCov(15,0.); tmpCov[0]=1e-4; tmpCov[2]=4.e-4; tmpCov[5]=4.e-4; tmpCov[9]=4.e-4; tmpCov[14]=1.e-10;
StatusCode scode; scode.setChecked();
std::vector<const xAOD::TrackParticle*> reducedTrkSet(SelectedTracks.begin(),SelectedTracks.begin()+nTrkLead);
double maxImp=RemoveNegImpact(reducedTrkSet,PrimVrt,JetDir,m_pseudoSigCut);
@@ -318,9 +319,12 @@ namespace InDet{
tTrkForFit[1]=tmpBTP;
TLorentzVector sumB(0.,0.,0.,0.);
int nvFitted=0;
- for(int it=0; it<(int)reducedTrkSet.size(); it++){
+ int nPRT=reducedTrkSet.size();
+ for(int it=0; it<nPRT; it++){
tTrkForFit[0]=reducedTrkSet[it];
if(tTrkForFit[0]->pt()<2000.)continue;
+ if(nPRT==1 && tTrkForFit[0]->pt()<300.*log(JetDir.Pt()))continue;
+ m_fitSvc->setApproximateVertex(0., 0., 0.);
scode=VKalVrtFitBase( tTrkForFit, FitVertex, Momentum, tChrg, ErrorMatrix, Chi2PerTrk, TrkAtVrt, Chi2); nvFitted++;
if(scode.isFailure() || Chi2>6.)continue;
if(FitVertex.perp()>reducedTrkSet[it]->radiusOfFirstHit())continue;
@@ -328,11 +332,11 @@ namespace InDet{
double Signif3D; VrtVrtDist(PrimVrt, FitVertex, ErrorMatrix, Signif3D);
if(FitVertex.perp()>m_Rbeampipe && Signif3D<2.) continue; // Cleaning of material interactions
if(m_FillHist)m_hb_DST_JetTrkSV->Fill(Signif3D,1.);
- if(tTrkForFit[0]->pt()<5.e4){
+ if(tTrkForFit[0]->pt()<2.e4){
if(!Check2TrVertexInPixel(tTrkForFit[0],tTrkForFit[0],FitVertex,Signif3D)) continue;
}
- sumSelTrk += MomAtVrt(TrkAtVrt[0]) ; sumB += MomAtVrt(TrkAtVrt[1]) ;
- if(Signif3D>selDST){ sel1T=it; selDST=Signif3D;}
+ //sumSelTrk += MomAtVrt(TrkAtVrt[0]) ; sumB += MomAtVrt(TrkAtVrt[1]) ;
+ if(Signif3D>selDST){ sel1T=it; selDST=Signif3D; sumSelTrk=MomAtVrt(TrkAtVrt[0]); sumB=MomAtVrt(TrkAtVrt[1]);}
}
if(sumSelTrk.Pt()>0. && sel1T>=0 ){
Results.resize(4);
@@ -352,22 +356,24 @@ namespace InDet{
//
//------ Plane-Plane crossing doesn't provide good vertices for the moment
// int sel2TI=-1,sel2TJ=-1;
-// if(reducedTrkSet.size()>1 && sel1T<0){
-// int nPRT=reducedTrkSet.size();
+// if(reducedTrkSet.size()<1 && sel1T<0){
+// int nPRT=reducedTrkSet.size(); Amg::Vector3D VB1,VB2; float distMax=1.e10;
// for(int it=0; it<nPRT-1; it++){ for(int jt=it+1; jt<nPRT; jt++){
-// TLorentzVector pseudoB=GetBDir( reducedTrkSet[it],reducedTrkSet[jt],PrimVrt);
-// if(pseudoB.DeltaR(JetDir)>0.3)continue;
-// sel2TI=it; sel2TJ=jt;
+// TLorentzVector pseudoB=GetBDir( reducedTrkSet[it], reducedTrkSet[jt], PrimVrt, VB1, VB2);
+// if( VB1.dot(VB2)==0. || pseudoB.DeltaR(JetDir)>0.3 ) continue;
+// if(Amg::distance(VB1,VB2)<distMax){ sel2TI=it; sel2TJ=jt; distMax=Amg::distance(VB1,VB2);}
// } }
-// if(sel1T<0 && sel2TI>=0){
+// if(sel2TI>=0){
// tTrkForFit[0]=reducedTrkSet[sel2TI]; tTrkForFit[1]=reducedTrkSet[sel2TJ];
+// //float RMHIT=TMath::Min(tTrkForFit[0]->radiusOfFirstHit(),tTrkForFit[1]->radiusOfFirstHit()); // Closest hit on both tracks
+// m_fitSvc->setApproximateVertex(0., 0., 0.);
// scode=VKalVrtFitBase( tTrkForFit, FitVertex, Momentum, tChrg, ErrorMatrix, Chi2PerTrk, TrkAtVrt, Chi2);
-// if(scode.isSuccess() && ProjPosT(FitVertex-PrimVrt.position(),JetDir) > 0 ){
+// if(scode.isSuccess() && ProjPosT(FitVertex-PrimVrt.position(),JetDir)>0 && FitVertex.perp()<180.){
// sumSelTrk += MomAtVrt(TrkAtVrt[0]) ; sumSelTrk += MomAtVrt(TrkAtVrt[1]) ;
-// Results.resize(4);
-// Results[0]=sumSelTrk.M(); //Invariant mass
-// Results[1]=sumSelTrk.Pt()/TrkJet.Pt(); //Ratio
-// Results[2]=0.; //Should be
+// Results.resize(4);
+// Results[0]=Momentum.M(); //Invariant mass
+// Results[1]=Momentum.Pt()/TrkJet.Pt(); //Ratio
+// Results[2]=0.; //Should be
// Results[3]=nPRT; //Found leading tracks with high positive impact
// } }
// }
@@ -878,10 +884,10 @@ namespace InDet{
ImpactSignif = sqrt( SignifR*SignifR + SignifZ*SignifZ);
}else if(m_getNegativeTag){
ImpactSignif = sqrt( (SignifR-0.6)*(SignifR-0.6)
- + (SignifZ-0.6)*(SignifZ-0.6) );
+ + (SignifZ-0.5)*(SignifZ-0.5) );
}else{
ImpactSignif = sqrt( (SignifR+0.6)*(SignifR+0.6)
- + (SignifZ+0.6)*(SignifZ+0.6) );
+ + (SignifZ+0.5)*(SignifZ+0.5) );
}
if(fabs(SignifR) < m_AntiPileupSigRCut) { // cut against tracks from pileup vertices
if(SignifZ > 1.+m_AntiPileupSigZCut ) ImpactSignif=0.;
@@ -1050,7 +1056,7 @@ namespace InDet{
//
float xvt=FitVertex.x(); float yvt=FitVertex.y();
if(m_FillHist){m_hb_r2d->Fill( Dist2D, m_w_1);}
- if(m_useMaterialRejection){
+ if(m_useMaterialRejection && (SelectedTracks[i]->pt()<m_MaterialPtCut || SelectedTracks[j]->pt()<m_MaterialPtCut) ){
//if(m_materialMap){
// if(m_materialMap->inMaterial(FitVertex)) BadTracks=4;
// if(msgLvl(MSG::DEBUG))msg(MSG::DEBUG)<<" MaterialMap test="<< BadTracks<<endreq;
@@ -1161,11 +1167,15 @@ namespace InDet{
const
{
int blTrk[2]={0,0};
+ int blP[2]={0,0};
int l1Trk[2]={0,0};
+ int l1P[2]={0,0};
int l2Trk[2]={0,0};
int nLays[2]={0,0};
getPixelLayers( p1, blTrk[0] , l1Trk[0], l2Trk[0], nLays[0] );
getPixelLayers( p2, blTrk[1] , l1Trk[1], l2Trk[1], nLays[1] ); // Very close to PV. Both b-layer hits are mandatory.
+ getPixelProblems(p1, blP[0], l1P[0] );
+ getPixelProblems(p2, blP[1], l1P[1] );
if( Signif3D<15. && FitVertex.perp()<15. ){
if( blTrk[0]<1 && l1Trk[0]<1 ) return false;
if( blTrk[1]<1 && l1Trk[1]<1 ) return false;
@@ -1179,8 +1189,9 @@ namespace InDet{
if( nLays[1] <2 ) return false; // Less than 2 layers on track 1
return true;
}else if(Dist2DBL > m_RlayerB+m_SVResolutionR){ //----------------------------------------- Outside b-layer
- if( blTrk[0]>0 || blTrk[1]>0 ) return false; // Hit in b-layer is present
- }
+ if( blTrk[0]>0 && blP[0]==0 ) return false; // Good hit in b-layer is present
+ if( blTrk[1]>0 && blP[1]==0 ) return false; // Good hit in b-layer is present
+ }
//
// L1 and L2 are considered only if vertex is in acceptance
//
@@ -1195,7 +1206,8 @@ namespace InDet{
if( l1Trk[1]<1 && l2Trk[1]<1 ) return false; // Less than 1 hits on track 1
return true;
}else if(Dist2DL1 > m_Rlayer1+m_SVResolutionR) { //------------------------------------------- Outside 1st-layer
- if( l1Trk[0]>0 || l1Trk[1]>0 ) return false; // L1 hits are present
+ if( l1Trk[0]>0 && l1P[0]==0 ) return false; // Good L1 hit is present
+ if( l1Trk[1]>0 && l1P[1]==0 ) return false; // Good L1 hit is present
}
if (Dist2DL2 < m_Rlayer2-m_SVResolutionR) { //------------------------------------------- Inside 2nd-layer
@@ -1238,7 +1250,7 @@ namespace InDet{
std::vector<double> ImpR(NTracks);
std::vector<double> Impact,ImpactError;
AmgVector(5) tmpPerigee; tmpPerigee<<0.,0.,0.,0.,0.;
- double maxImp=-1.e10;
+ double maxImp=-1.e10, zImp=0.;
for (i=0; i<NTracks; i++) {
m_fitSvc->VKalGetImpact(inTracks[i], PrimVrt.position(), 1, Impact, ImpactError);
tmpPerigee = GetPerigee(inTracks[i])->parameters();
@@ -1250,7 +1262,8 @@ namespace InDet{
if(fabs(SignifR) < m_AntiPileupSigRCut) { // cut against tracks from pileup vertices
if( fabs(Impact[1])/sqrt(ImpactError[2]) > m_AntiPileupSigZCut ) ImpR[i]=-9999.;
}
- if(ImpR[i]>maxImp)maxImp=ImpR[i];
+ if(fabs(Impact[1])/sqrt(ImpactError[2])>fabs(ImpR[i]) && Impact[1]<0 && ImpR[i]>0) ImpR[i]*=-1.;
+ if(ImpR[i]>maxImp){maxImp=ImpR[i]; zImp=Impact[1]/sqrt(ImpactError[2]);}
}
if(maxImp<Limit){ inTracks.clear(); return maxImp;}
double rmin=1.e6;
@@ -1258,6 +1271,7 @@ namespace InDet{
for(i=0; i<(int)ImpR.size(); i++){ if(rmin>ImpR[i]){rmin=ImpR[i]; jpm=i;}}; if(rmin>Limit)continue;
ImpR.erase(ImpR.begin()+jpm); inTracks.erase(inTracks.begin()+jpm);
}while(rmin<=Limit);
+ if(inTracks.size()==1 && zImp<1.){ inTracks.clear(); maxImp=0.;}
return maxImp;
}
diff --git a/InnerDetector/InDetRecTools/InDetVKalVxInJetTool/src/InDetVKalVxInJetTool.cxx b/InnerDetector/InDetRecTools/InDetVKalVxInJetTool/src/InDetVKalVxInJetTool.cxx
index 74a9f0254c5..7e07688536e 100755
--- a/InnerDetector/InDetRecTools/InDetVKalVxInJetTool/src/InDetVKalVxInJetTool.cxx
+++ b/InnerDetector/InDetRecTools/InDetVKalVxInJetTool/src/InDetVKalVxInJetTool.cxx
@@ -52,6 +52,7 @@ InDetVKalVxInJetTool::InDetVKalVxInJetTool(const std::string& type,
m_AntiFake2trVrtCut(0.5),
m_JetPtFractionCut(0.01),
m_TrackInJetNumberLimit(25),
+ m_MaterialPtCut(5.e3),
m_pseudoSigCut(3.),
m_FillHist(false),
m_existIBL(true),
@@ -68,7 +69,7 @@ InDetVKalVxInJetTool::InDetVKalVxInJetTool(const std::string& type,
m_RlayerB (0.), // in jobO or initialize()
m_Rlayer1 (0.),
m_Rlayer2 (0.),
- m_SVResolutionR(0.),
+ m_SVResolutionR(3.),
m_useMaterialRejection(true),
m_useVertexCleaning(true),
m_MassType (1),
@@ -118,6 +119,7 @@ InDetVKalVxInJetTool::InDetVKalVxInJetTool(const std::string& type,
declareProperty("AntiFake2trVrtCut", m_AntiFake2trVrtCut, "Cut to reduce fake 2-track vertices contribution.Single Vertex Finder only" );
declareProperty("JetPtFractionCut", m_JetPtFractionCut, "Reduce high Pt fakes. Jet HLV input is mandatory, direction is not enough. Multi and single vertex versions are affected" );
declareProperty("TrackInJetNumberLimit", m_TrackInJetNumberLimit, " Use only limited number of highest pT tracks in jet for vertex search" );
+ declareProperty("MaterialPtCut", m_MaterialPtCut, " To be a material vertex at least one track pt should be below this cut" );
declareProperty("PseudoSigCut", m_pseudoSigCut, " Cut on track impact significance for pseudo-vertex search" );
declareProperty("FillHist", m_FillHist, "Fill technical histograms" );
diff --git a/InnerDetector/InDetRecTools/InDetVKalVxInJetTool/src/Utilities.cxx b/InnerDetector/InDetRecTools/InDetVKalVxInJetTool/src/Utilities.cxx
index 27b7f6c7cf8..85722183d0b 100755
--- a/InnerDetector/InDetRecTools/InDetVKalVxInJetTool/src/Utilities.cxx
+++ b/InnerDetector/InDetRecTools/InDetVKalVxInJetTool/src/Utilities.cxx
@@ -16,7 +16,8 @@ namespace InDet{
TLorentzVector InDetVKalVxInJetTool::GetBDir( const xAOD::TrackParticle* trk1,
const xAOD::TrackParticle* trk2,
- const xAOD::Vertex & PrimVrt)
+ const xAOD::Vertex & PrimVrt,
+ Amg::Vector3D &V1, Amg::Vector3D &V2)
const
{ // B hadron flight direction based on 2 separate tracks and PV. Calculated via plane-plane crossing
Amg::Vector3D PVRT(PrimVrt.x(),PrimVrt.y(),PrimVrt.z());
@@ -26,18 +27,33 @@ namespace InDet{
////Amg::Vector3D new1(trk1->vx()-d01*sin(phi1),trk1->vy()+d01*cos(phi1),trk1->vz()+z01); // Wrong sign!!!!!
////Amg::Vector3D new2(trk2->vx()-d02*sin(phi2),trk2->vy()+d02*cos(phi2),trk2->vz()+z02);
//----------------------------------------------------------------------------
- Amg::Vector3D pnt1=trk1->perigeeParameters().position();
+ Amg::Vector3D pnt1=trk1->perigeeParameters().position()-PVRT;
Amg::Vector3D mom1((trk1->p4()).Px(),(trk1->p4()).Py(),(trk1->p4()).Pz());
- Amg::Vector3D pnt2=trk2->perigeeParameters().position();
+ Amg::Vector3D pnt2=trk2->perigeeParameters().position()-PVRT;
Amg::Vector3D mom2((trk2->p4()).Px(),(trk2->p4()).Py(),(trk2->p4()).Pz());
- pnt1 -= PVRT; pnt2 -= PVRT;
pnt1.normalize(); pnt2.normalize(); mom1.normalize(); mom2.normalize();
//------------------------------------------------------------------------
+ const float dRLim=0.4;
Amg::Vector3D norm1=pnt1.cross(mom1);
Amg::Vector3D norm2=pnt2.cross(mom2);
- Amg::Vector3D t=norm1.cross(norm2); t.normalize();
+ Amg::Vector3D t=norm1.cross(norm2); t.normalize(); if(t.dot(mom1+mom2)<0.) t*=-1.;
double aveP=(trk1->p4()+trk2->p4()).P()/2.;
- TLorentzVector tl; tl.SetXYZM(t.x()*aveP,t.y()*aveP,t.z()*aveP,139.57);
+ TLorentzVector tl; tl.SetXYZM(t.x()*aveP,t.y()*aveP,t.z()*aveP,139.57); //Crossing line of 2 planes
+ if( tl.DeltaR(trk1->p4()) >dRLim || tl.DeltaR(trk2->p4()) >dRLim ) {V1*=0.; V2*=0.; return tl;}//Too big dR between tracks and found "B line"
+//------------------------------------------------------------------------
+ double X;
+ pnt1=trk1->perigeeParameters().position()-PVRT;
+ pnt2=trk2->perigeeParameters().position()-PVRT;
+ fabs(mom1[1]*t[2]-mom1[2]*t[1])>fabs(mom1[0]*t[2]-mom1[2]*t[0]) ? X=(t[1]*pnt1[2]-t[2]*pnt1[1])/(mom1[1]*t[2]-mom1[2]*t[1])
+ : X=(t[0]*pnt1[2]-t[2]*pnt1[0])/(mom1[0]*t[2]-mom1[2]*t[0]);
+ V1=pnt1+mom1*X; // First particle vertex
+ fabs(mom2[1]*t[2]-mom2[2]*t[1])>fabs(mom2[0]*t[2]-mom2[2]*t[0]) ? X=(t[1]*pnt2[2]-t[2]*pnt2[1])/(mom2[1]*t[2]-mom2[2]*t[1])
+ : X=(t[0]*pnt2[2]-t[2]*pnt2[0])/(mom2[0]*t[2]-mom2[2]*t[0]);
+ V2=pnt2+mom2*X; // Second particle vertex
+//------------------------------------------------------------------------
+ if(V1.dot(t)<0. && V2.dot(t)<0.) {V1*=0.;V2*=0.;} // Check correctness of topology
+ else {V1+=PVRT; V2+=PVRT;} // Transform to detector frame
+//------------------------------------------------------------------------
return tl;
}
@@ -571,7 +587,24 @@ namespace InDet{
}
}
-
+ void InDetVKalVxInJetTool::getPixelProblems(const xAOD::TrackParticle* Part, int &splshIBL, int &splshBL ) const
+ {
+ splshIBL=splshBL=0;
+ if(m_existIBL){ // 4-layer pixel detector
+ uint8_t share,split;
+ //if(!Part->summaryValue( IBLout, xAOD::numberOfInnermostPixelLayerOutliers ) ) IBLout = 0;
+ if(!Part->summaryValue( share, xAOD::numberOfInnermostPixelLayerSharedHits ) ) share = 0;
+ if(!Part->summaryValue( split, xAOD::numberOfInnermostPixelLayerSplitHits ) ) split = 0;
+ splshIBL=share+split;
+ if(!Part->summaryValue( share, xAOD::numberOfNextToInnermostPixelLayerSharedHits ) ) share = 0;
+ if(!Part->summaryValue( split, xAOD::numberOfNextToInnermostPixelLayerSplitHits ) ) split = 0;
+ splshBL=share+split;
+ }
+ }
+ void InDetVKalVxInJetTool::getPixelProblems(const Rec::TrackParticle* , int &splshIBL, int &splshBL ) const
+ {
+ splshIBL=splshBL=0; // Temporary implementation
+ }
void InDetVKalVxInJetTool::getPixelDiscs(const xAOD::TrackParticle* Part, int &d0Hit, int &d1Hit, int &d2Hit) const
{
uint32_t HitPattern=Part->hitPattern();
--
GitLab