ATLAS Offline Software
Loading...
Searching...
No Matches
MuonVertexValidationMacroUtils Namespace Reference

Functions

bool inBarrel (double eta)
bool inBarrel (const Amg::Vector3D &vtx)
bool inEndcaps (double eta)
bool inEndcaps (const Amg::Vector3D &vtx)
bool inDetectorRegion (const Amg::Vector3D &vtx)
int getNvtxBarrel (const std::vector< Amg::Vector3D > &vertices)
int getNvtxEndcaps (const std::vector< Amg::Vector3D > &vertices)
int getNvtxDetectorRegion (const std::vector< Amg::Vector3D > &vertices)
bool inFiducialVolBarrel (const Amg::Vector3D &vtx)
bool inFiducialVolEndcaps (const Amg::Vector3D &vtx)
bool inFiducialVol (const Amg::Vector3D &vtx)
int NvtxFiducialVol (const std::vector< Amg::Vector3D > &vertices)
bool isGoodVtx (const Amg::Vector3D &vtx)
double getMatchMetric (const Amg::Vector3D &vtx1, const Amg::Vector3D &vtx2)
Amg::Vector3D findBestMatch (const Amg::Vector3D &vtx, const std::vector< Amg::Vector3D > &candidates)
bool isValidMatch (const Amg::Vector3D &match_candidate)
bool hasMatch (const Amg::Vector3D &vtx1, const std::vector< Amg::Vector3D > &vtx2_vec)
std::vector< Amg::Vector3DgetVertexPos (const std::vector< double > &vtx_x, const std::vector< double > &vtx_y, const std::vector< double > &vtx_z)
std::vector< std::vector< Amg::Vector3D > > getConstituentPos (int Nvtx, const std::vector< int > &obj_vtx_link, const std::vector< double > &obj_x, const std::vector< double > &obj_y, const std::vector< double > &obj_z)

Variables

constexpr double fidVol_barrel_etaCut = 0.7
constexpr double fidVol_Lxy_low = 3*Gaudi::Units::m
constexpr double fidVol_Lxy_up = 8*Gaudi::Units::m
constexpr double fidVol_endcaps_etaCut_low = 1.3
constexpr double fidVol_endcaps_etaCut_up = 2.5
constexpr double fidVol_endcaps_Lxy_up = 10*Gaudi::Units::m
constexpr double fidVol_z_low = 5*Gaudi::Units::m
constexpr double fidVol_z_up = 15*Gaudi::Units::m
constexpr double match_max = 0.4

Function Documentation

◆ findBestMatch()

Amg::Vector3D MuonVertexValidationMacroUtils::findBestMatch ( const Amg::Vector3D & vtx,
const std::vector< Amg::Vector3D > & candidates )

Definition at line 102 of file MuonSpectrometer/MuonValidation/MuonVertexValidation/util/Utils.cxx.

102 {
103 // finds the best match for vtx among candidates
104 // The match metric is assumed to decreases for better matches. If no match could be found, a zero vector is returned
105 Amg::Vector3D best{Amg::Vector3D::Zero()}; // initialised with zeros
106 double aux(0), min(999);
107
108 for (unsigned int i=0; i<candidates.size(); ++i){
109 aux = getMatchMetric(vtx, candidates[i]);
110 if (aux > match_max) continue;
111 if (aux < min){
112 min = aux;
113 best = candidates[i];
114 }
115 }
116
117 return best;
118}
#define min(a, b)
Definition cfImp.cxx:40
Eigen::Matrix< double, 3, 1 > Vector3D
double getMatchMetric(const Amg::Vector3D &vtx1, const Amg::Vector3D &vtx2)
best(iterable, priorities=[3, 2, 1, -1, 0])

◆ getConstituentPos()

std::vector< std::vector< Amg::Vector3D > > MuonVertexValidationMacroUtils::getConstituentPos ( int Nvtx,
const std::vector< int > & obj_vtx_link,
const std::vector< double > & obj_x,
const std::vector< double > & obj_y,
const std::vector< double > & obj_z )

Definition at line 144 of file MuonSpectrometer/MuonValidation/MuonVertexValidation/util/Utils.cxx.

145 {
146 // returns the constituent position vectors for each reconstructed vertex
147 std::vector<std::vector<Amg::Vector3D>> consti_vec;
148 for (int j=0; j<Nvtx; ++j){
149 std::vector<Amg::Vector3D> consti;
150 for (unsigned int k=0; k<obj_x.size(); ++k){
151 if (obj_vtx_link[k] == j){consti.push_back(Amg::Vector3D(obj_x[k], obj_y[k], obj_z[k]));}
152 }
153 consti_vec.push_back(consti);
154 }
155
156 return consti_vec;
157}

◆ getMatchMetric()

double MuonVertexValidationMacroUtils::getMatchMetric ( const Amg::Vector3D & vtx1,
const Amg::Vector3D & vtx2 )

Definition at line 96 of file MuonSpectrometer/MuonValidation/MuonVertexValidation/util/Utils.cxx.

96 {
97 // Computes the metric on which matches are made between two vertices
98 return xAOD::P4Helpers::deltaR(vtx1.eta(), vtx1.phi(), vtx2.eta(), vtx2.phi());
99}
double deltaR(double rapidity1, double phi1, double rapidity2, double phi2)
from bare bare rapidity,phi

◆ getNvtxBarrel()

int MuonVertexValidationMacroUtils::getNvtxBarrel ( const std::vector< Amg::Vector3D > & vertices)

Definition at line 38 of file MuonSpectrometer/MuonValidation/MuonVertexValidation/util/Utils.cxx.

38 {
39 // find the number of vertices in the barrel region
40 return std::accumulate(vertices.begin(), vertices.end(), 0, [](unsigned int N, const Amg::Vector3D& vtx){return N+inBarrel(vtx);});
41}

◆ getNvtxDetectorRegion()

int MuonVertexValidationMacroUtils::getNvtxDetectorRegion ( const std::vector< Amg::Vector3D > & vertices)

Definition at line 50 of file MuonSpectrometer/MuonValidation/MuonVertexValidation/util/Utils.cxx.

50 {
51 // find the number of vertices in the either the barrel or endcaps
52 return getNvtxBarrel(vertices) + getNvtxEndcaps(vertices);
53}

◆ getNvtxEndcaps()

int MuonVertexValidationMacroUtils::getNvtxEndcaps ( const std::vector< Amg::Vector3D > & vertices)

Definition at line 44 of file MuonSpectrometer/MuonValidation/MuonVertexValidation/util/Utils.cxx.

44 {
45 // find the number of vertices in the endcap region
46 return std::accumulate(vertices.begin(), vertices.end(), 0, [](unsigned int N, const Amg::Vector3D& vtx){return N+inEndcaps(vtx);});
47}

◆ getVertexPos()

std::vector< Amg::Vector3D > MuonVertexValidationMacroUtils::getVertexPos ( const std::vector< double > & vtx_x,
const std::vector< double > & vtx_y,
const std::vector< double > & vtx_z )

Definition at line 133 of file MuonSpectrometer/MuonValidation/MuonVertexValidation/util/Utils.cxx.

133 {
134 std::vector<Amg::Vector3D> vertices;
135 for (unsigned int i=0; i<vtx_x.size(); ++i){
136 Amg::Vector3D vtx_pos(vtx_x[i], vtx_y[i], vtx_z[i]);
137 vertices.push_back(vtx_pos);
138 }
139
140 return vertices;
141}

◆ hasMatch()

bool MuonVertexValidationMacroUtils::hasMatch ( const Amg::Vector3D & vtx1,
const std::vector< Amg::Vector3D > & vtx2_vec )

Definition at line 126 of file MuonSpectrometer/MuonValidation/MuonVertexValidation/util/Utils.cxx.

126 {
127 // Checks if vtx1 has a match in vtx2_vec
128 Amg::Vector3D candidate = findBestMatch(vtx1, vtx2_vec);
129 return isValidMatch(candidate);
130}
Amg::Vector3D findBestMatch(const Amg::Vector3D &vtx, const std::vector< Amg::Vector3D > &candidates)

◆ inBarrel() [1/2]

bool MuonVertexValidationMacroUtils::inBarrel ( const Amg::Vector3D & vtx)

◆ inBarrel() [2/2]

bool MuonVertexValidationMacroUtils::inBarrel ( double eta)

Definition at line 13 of file MuonSpectrometer/MuonValidation/MuonVertexValidation/util/Utils.cxx.

13 {
14 return std::abs(eta) <= fidVol_barrel_etaCut;
15}
Scalar eta() const
pseudorapidity method

◆ inDetectorRegion()

bool MuonVertexValidationMacroUtils::inDetectorRegion ( const Amg::Vector3D & vtx)

◆ inEndcaps() [1/2]

bool MuonVertexValidationMacroUtils::inEndcaps ( const Amg::Vector3D & vtx)

Definition at line 28 of file MuonSpectrometer/MuonValidation/MuonVertexValidation/util/Utils.cxx.

28 {
29 return inEndcaps(vtx.eta());
30}

◆ inEndcaps() [2/2]

bool MuonVertexValidationMacroUtils::inEndcaps ( double eta)

◆ inFiducialVol()

bool MuonVertexValidationMacroUtils::inFiducialVol ( const Amg::Vector3D & vtx)

◆ inFiducialVolBarrel()

bool MuonVertexValidationMacroUtils::inFiducialVolBarrel ( const Amg::Vector3D & vtx)

◆ inFiducialVolEndcaps()

bool MuonVertexValidationMacroUtils::inFiducialVolEndcaps ( const Amg::Vector3D & vtx)

◆ isGoodVtx()

bool MuonVertexValidationMacroUtils::isGoodVtx ( const Amg::Vector3D & vtx)

Definition at line 86 of file MuonSpectrometer/MuonValidation/MuonVertexValidation/util/Utils.cxx.

86 {
87 // add further conditions for a vertex to be considered in efficiency calculations here
88 bool is_good = true;
89 if (!inFiducialVol(vtx)) is_good = false;
90 return is_good;
91}

◆ isValidMatch()

bool MuonVertexValidationMacroUtils::isValidMatch ( const Amg::Vector3D & match_candidate)

Definition at line 121 of file MuonSpectrometer/MuonValidation/MuonVertexValidation/util/Utils.cxx.

121 {
122 return match_candidate.mag() > 0.001;
123}

◆ NvtxFiducialVol()

int MuonVertexValidationMacroUtils::NvtxFiducialVol ( const std::vector< Amg::Vector3D > & vertices)

Definition at line 78 of file MuonSpectrometer/MuonValidation/MuonVertexValidation/util/Utils.cxx.

78 {
79 // find the number ofvertices in the fiducial volume
80 unsigned int N = 0;
81 for (const Amg::Vector3D &vtx : vertices) if (inFiducialVol(vtx)) ++N;
82 return N;
83}

Variable Documentation

◆ fidVol_barrel_etaCut

double MuonVertexValidationMacroUtils::fidVol_barrel_etaCut = 0.7
constexpr

◆ fidVol_endcaps_etaCut_low

double MuonVertexValidationMacroUtils::fidVol_endcaps_etaCut_low = 1.3
constexpr

◆ fidVol_endcaps_etaCut_up

double MuonVertexValidationMacroUtils::fidVol_endcaps_etaCut_up = 2.5
constexpr

◆ fidVol_endcaps_Lxy_up

double MuonVertexValidationMacroUtils::fidVol_endcaps_Lxy_up = 10*Gaudi::Units::m
constexpr

◆ fidVol_Lxy_low

double MuonVertexValidationMacroUtils::fidVol_Lxy_low = 3*Gaudi::Units::m
constexpr

◆ fidVol_Lxy_up

double MuonVertexValidationMacroUtils::fidVol_Lxy_up = 8*Gaudi::Units::m
constexpr

◆ fidVol_z_low

double MuonVertexValidationMacroUtils::fidVol_z_low = 5*Gaudi::Units::m
constexpr

◆ fidVol_z_up

double MuonVertexValidationMacroUtils::fidVol_z_up = 15*Gaudi::Units::m
constexpr

◆ match_max

double MuonVertexValidationMacroUtils::match_max = 0.4
constexpr