ATLAS Offline Software
L2OverlapRemoverMon.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 #include "L2OverlapRemoverMon.h"
6 
7 #include "MuonMatchingTool.h"
8 
9 L2OverlapRemoverMon :: L2OverlapRemoverMon(const std::string& name, ISvcLocator* pSvcLocator )
10  : TrigMuonMonitorAlgorithm(name, pSvcLocator)
11 {}
12 
13 
14 StatusCode L2OverlapRemoverMon :: fillVariablesPerChain(const EventContext &, const std::string &chain) const {
15 
16  ATH_MSG_DEBUG ("Filling histograms for " << name() << "...");
17  ATH_CHECK( fillVariablesOverlapRemoverPlots<xAOD::L2StandAloneMuon>(chain, "L2SA") );
18  ATH_CHECK( fillVariablesOverlapRemoverPlots<xAOD::L2CombinedMuon>(chain, "L2CB") );
19 
20  return StatusCode::SUCCESS;
21 }
22 
23 
25 
26  const float ZERO_LIMIT = 0.00001;
27 
28  auto dR = Monitored::Scalar<float>("L2SA_"+chain+"_dR", -999.);
29  auto invMass = Monitored::Scalar<float>("L2SA_"+chain+"_invMass", -999.);
30  auto dRLog10 = Monitored::Scalar<float>("L2SA_"+chain+"_dRLog10", -999.);
31  auto invMassLog10 = Monitored::Scalar<float>("L2SA_"+chain+"_invMassLog10", -999.);
32 
33  const auto [mu1Pt, mu1Eta, mu1Phi] = L2ORPosForMatchFunc((*muEL1));
34  const auto [mu2Pt, mu2Eta, mu2Phi] = L2ORPosForMatchFunc((*muEL2));
35 
36  if( ( (std::abs(mu1Eta) < ZERO_LIMIT) && (std::abs(mu1Phi) < ZERO_LIMIT) ) ||
37  ( (std::abs(mu2Eta) < ZERO_LIMIT) && (std::abs(mu2Phi) < ZERO_LIMIT) ) ||
38  (std::abs(mu1Pt) < ZERO_LIMIT) || (std::abs(mu2Pt) < ZERO_LIMIT) ) return false;
39 
40  double dRThres = 9999.;
41  double massThres = 9999.;
42  bool isBarrel1 = (*muEL1)->sAddress() != -1;
43  bool isBarrel2 = (*muEL2)->sAddress() != -1;
44 
45  if( isBarrel1 && isBarrel2 ) { // BB
46  dRThres =m_dRSAThresBB;
47  massThres=m_massSAThresBB;
48  }
49  else if( (isBarrel1 && ! isBarrel2) || (!isBarrel1 && isBarrel2) ) { // BE
50  dRThres =m_dRSAThresBE;
51  massThres=m_massSAThresBE;
52  }
53  else { // EE
54  double absEta = (std::abs(mu1Pt) > std::abs(mu2Pt)) ? std::abs(mu1Eta) : std::abs(mu2Eta);
55  unsigned int iThres=0;
56  for(unsigned int i=0; i<(m_etaBins.size()-1); i++) {
57  if ( m_etaBins[i] <= absEta && absEta < m_etaBins[i+1] ) iThres = i;
58  }
59  dRThres = m_dRSAThresEC[iThres];
60  massThres = m_massSAThresEC[iThres];
61  }
62 
63 
64  // same sign cut
65  bool sameSign = mu1Pt*mu2Pt > 0;
66 
67  // dR cut
68  float deta = mu1Eta - mu2Eta;
69  float dphi = xAOD::P4Helpers::deltaPhi(mu1Phi, mu2Phi);
70  dR = sqrt(deta*deta + dphi*dphi);
71  bool dRisClose = dR < dRThres;
72 
73  // mass cut
74  invMass = calcinvMass(0., mu1Pt, mu1Eta, mu1Phi, 0., mu2Pt, mu2Eta, mu2Phi);
75  bool massIsClose = invMass < massThres;
76 
77  // for monitorinng log10 plot
78  const float monitor_limit = 1e-4;
79  dRLog10 = ( dR >= monitor_limit ) ? log10(dR) : log10(monitor_limit);
80  invMassLog10 = ( invMass >= monitor_limit ) ? log10(invMass) : log10(monitor_limit);
81 
82 
83  // total judge
84  bool overlap = sameSign && dRisClose && massIsClose;
85  ATH_MSG_DEBUG( " ...=> isOverlap=" << overlap );
86 
87  fill(m_group+"_"+chain, dR, invMass, dRLog10, invMassLog10);
88 
89  return overlap;
90 }
91 
92 
94 
95  const float ZERO_LIMIT = 0.00001;
96 
97  auto dR = Monitored::Scalar<float>("L2CB_"+chain+"_dR", -999.);
98  auto invMass = Monitored::Scalar<float>("L2CB_"+chain+"_invMass", -999.);
99  auto dRLog10 = Monitored::Scalar<float>("L2CB_"+chain+"_dRLog10", -999.);
100  auto invMassLog10 = Monitored::Scalar<float>("L2CB_"+chain+"_invMassLog10", -999.);
101 
102  const auto [mu1Pt, mu1Eta, mu1Phi] = L2ORPosForMatchFunc((*muEL1));
103  const auto [mu2Pt, mu2Eta, mu2Phi] = L2ORPosForMatchFunc((*muEL2));
104 
105  if( ( (std::abs(mu1Eta) < ZERO_LIMIT) && (std::abs(mu1Phi) < ZERO_LIMIT) ) ||
106  ( (std::abs(mu2Eta) < ZERO_LIMIT) && (std::abs(mu2Phi) < ZERO_LIMIT) ) ||
107  (std::abs(mu1Pt) < ZERO_LIMIT) || (std::abs(mu2Pt) < ZERO_LIMIT) ) return false;
108 
109  double absEta = (std::abs(mu1Pt) > std::abs(mu2Pt)) ? std::abs(mu1Eta) : std::abs(mu2Eta);
110  unsigned int iThres=0;
111  for(unsigned int i=0; i<(m_etaBins.size()-1); i++) {
112  if ( m_etaBins[i] <= absEta && absEta < m_etaBins[i+1] ) iThres = i;
113  }
114  float dRThres = m_dRCBThres[iThres];
115  float dRbySAThres = m_dRbySAThres[iThres];
116  float massThres = m_massCBThres[iThres];
117 
118 
119  // same sign cut
120  bool sameSign = mu1Pt*mu2Pt > 0;
121 
122  // dR cut
123  float deta = mu1Eta - mu2Eta;
124  float dphi = xAOD::P4Helpers::deltaPhi(mu1Phi, mu2Phi);
125  dR = sqrt(deta*deta + dphi*dphi);
126  bool dRisClose = dR < dRThres;
127 
128  // dR(by L2SA) cut
129  bool dRbySAisClose = false;
130  const xAOD::L2StandAloneMuon* muSA1 = (*muEL1)->muSATrack();
131  const xAOD::L2StandAloneMuon* muSA2 = (*muEL2)->muSATrack();
132  if( muSA1 == nullptr || muSA2 == nullptr ) return false;
133  else {
134  float deta = muSA1->etaMS() - muSA2->etaMS();
135  float dphi = xAOD::P4Helpers::deltaPhi(muSA1->phiMS(), muSA2->phiMS());
136  float dRBySA = sqrt(deta*deta + dphi*dphi);
137  if( dRBySA < dRbySAThres ) dRbySAisClose = true;
138  }
139 
140  // mass cut
141  invMass = calcinvMass(0., mu1Pt, mu1Eta, mu1Phi, 0., mu2Pt, mu2Eta, mu2Phi);
142  bool massIsClose = invMass < massThres;
143 
144  // for monitorinng log10 plot
145  const float monitor_limit = 1e-4;
146  dRLog10 = ( dR >= monitor_limit ) ? log10(dR) : log10(monitor_limit);
147  invMassLog10 = ( invMass >= monitor_limit ) ? log10(invMass) : log10(monitor_limit);
148 
149 
150  // total judge
151  bool overlap = sameSign && dRisClose && massIsClose && dRbySAisClose;
152  ATH_MSG_DEBUG( " ...=> isOverlap=" << overlap );
153 
154  fill(m_group+"_"+chain, dR, invMass, dRLog10, invMassLog10);
155 
156  return overlap;
157 }
158 
159 
160 StatusCode L2OverlapRemoverMon::chooseBestMuon(const std::string &chain, const std::vector< TrigCompositeUtils::LinkInfo<xAOD::L2StandAloneMuonContainer> >& featureCont, const std::vector<unsigned int>& muResult) const
161 {
162  unsigned int i,j,k;
163 
164  std::vector<float> vec_RemovedEta, vec_RemovedPhi, vec_RemovedPt;
165  auto RemovedEta = Monitored::Collection("L2SA_"+chain+"_RemovedEta", vec_RemovedEta);
166  auto RemovedPhi = Monitored::Collection("L2SA_"+chain+"_RemovedPhi", vec_RemovedPhi);
167  auto RemovedPt = Monitored::Collection("L2SA_"+chain+"_RemovedPt", vec_RemovedPt);
168 
169  for(i=0; i<featureCont.size(); i++) {
170  ATH_MSG_DEBUG( "++ i=" << i << ": result=" << muResult[i] );
171  if( muResult[i] != i ) {
172  ATH_MSG_DEBUG( " overlap to some one. skip." );
173  continue;
174  }
175  std::vector<unsigned int> others;
176  for(j=0; j<featureCont.size(); j++) {
177  if( muResult[j] == muResult[i] ) others.emplace_back(j);
178  }
179  if( others.size() == 1 ) {
180  ATH_MSG_DEBUG( " unique object. keep it active." );
181  continue;
182  }
183  else { // must choose one best
184  ATH_MSG_DEBUG( " overlapped objects among: " << others );
185  unsigned int BestMuon = 0;
186  float maxPt = 0.;
187  float maxPtRoI = 0.;
188  for(k=0; k<others.size(); k++) {
189  j=others[k];
190  const ElementLink<xAOD::L2StandAloneMuonContainer> muEL = featureCont[j].link;
191  float pt = std::abs((*muEL)->pt());
192  float ptRoI = (*muEL)->roiThreshold();
193  if( (ptRoI-maxPtRoI) > 0.1 ) {
194  maxPtRoI = ptRoI;
195  maxPt = pt;
196  BestMuon = j;
197  }
198  else if( std::abs(ptRoI-maxPtRoI) < 0.1 ) {
199  if( pt > maxPt ) {
200  maxPtRoI = ptRoI;
201  maxPt = pt;
202  BestMuon = j;
203  }
204  }
205  }
206 
207  for(k=0; k<others.size(); k++) {
208  j=others[k];
209  if( j != BestMuon ) { // removed muon
210  const ElementLink<xAOD::L2StandAloneMuonContainer> muEL = featureCont[j].link;
211  vec_RemovedPt.push_back( (*muEL)->pt() );
212  vec_RemovedEta.push_back( (*muEL)->etaMS() );
213  vec_RemovedPhi.push_back( (*muEL)->phiMS() );
214  }
215  }
216  }
217  }
218 
219  fill(m_group+"_"+chain, RemovedEta, RemovedPhi, RemovedPt);
220 
221  return StatusCode::SUCCESS;
222 }
223 
224 StatusCode L2OverlapRemoverMon::chooseBestMuon(const std::string &chain, const std::vector< TrigCompositeUtils::LinkInfo<xAOD::L2CombinedMuonContainer> >& featureCont, const std::vector<unsigned int>& muResult) const
225 {
226  unsigned int i,j,k;
227 
228  std::vector<float> vec_RemovedEta, vec_RemovedPhi, vec_RemovedPt;
229  auto RemovedEta = Monitored::Collection("L2CB_"+chain+"_RemovedEta", vec_RemovedEta);
230  auto RemovedPhi = Monitored::Collection("L2CB_"+chain+"_RemovedPhi", vec_RemovedPhi);
231  auto RemovedPt = Monitored::Collection("L2CB_"+chain+"_RemovedPt", vec_RemovedPt);
232 
233  for(i=0; i<featureCont.size(); i++) {
234  ATH_MSG_DEBUG( "++ i=" << i << ": result=" << muResult[i] );
235  if( muResult[i] != i ) {
236  ATH_MSG_DEBUG( " overlap to some one. skip." );
237  continue;
238  }
239  std::vector<unsigned int> others;
240  for(j=0; j<featureCont.size(); j++) {
241  if( muResult[j] == muResult[i] ) others.emplace_back(j);
242  }
243  if( others.size() == 1 ) {
244  ATH_MSG_DEBUG( " unique object. keep it active." );
245  continue;
246  }
247  else { // must choose one best
248  ATH_MSG_DEBUG( " overlapped objects among: " << others );
249  unsigned int BestMuon = 0;
250  float maxPt = 0.;
251  for(k=0; k<others.size(); k++) {
252  j=others[k];
253  const ElementLink<xAOD::L2CombinedMuonContainer> muEL = featureCont[j].link;
254  float pt = std::abs((*muEL)->pt()/1e3);
255  if( pt > maxPt ) {
256  maxPt = pt;
257  BestMuon = j;
258  }
259  }
260 
261  for(k=0; k<others.size(); k++) {
262  j=others[k];
263  if( j != BestMuon ) { // removed muon
264  const ElementLink<xAOD::L2CombinedMuonContainer> muEL = featureCont[j].link;
265  vec_RemovedPt.push_back( (*muEL)->pt()/1e3 * (*muEL)->charge() );
266  vec_RemovedEta.push_back( (*muEL)->eta() );
267  vec_RemovedPhi.push_back( (*muEL)->phi() );
268  }
269  }
270  }
271  }
272 
273  fill(m_group+"_"+chain, RemovedEta, RemovedPhi, RemovedPt);
274 
275  return StatusCode::SUCCESS;
276 }
277 
278 
279 
280 float L2OverlapRemoverMon::calcinvMass(double m1, double pt1, double eta1, double phi1,
281  double m2, double pt2, double eta2, double phi2)
282 {
283  const double ZERO_LIMIT = 1e-12;
284 
285  double theta1 = 2*atan2((double)exp(-eta1),1.);
286  double theta2 = 2*atan2((double)exp(-eta2),1.);
287 
288  double fpt1 = fabs(pt1);
289  double fpt2 = fabs(pt2);
290 
291  double px1 = fpt1*cos(phi1);
292  double py1 = fpt1*sin(phi1);
293  double pz1 = fpt1/tan(theta1);
294  double e1 = sqrt(px1*px1+py1*py1+pz1*pz1+m1*m1);
295 
296  double px2 = fpt2*cos(phi2);
297  double py2 = fpt2*sin(phi2);
298  double pz2 = fpt2/tan(theta2);
299  double e2 = sqrt(px2*px2+py2*py2+pz2*pz2+m2*m2);
300 
301  double pxsum = px1 + px2;
302  double pysum = py1 + py2;
303  double pzsum = pz1 + pz2;
304  double esum = e1 + e2;
305 
306  double mass = 0;
307  double mass2 = esum*esum - pxsum*pxsum - pysum*pysum - pzsum*pzsum;
308  if( mass2 > ZERO_LIMIT ) mass = sqrt(mass2);
309 
310  return mass;
311 }
312 
313 
314 std::tuple<float,float,float> L2OverlapRemoverMon :: L2ORPosForMatchFunc(const xAOD::L2StandAloneMuon *trig){
315  return std::forward_as_tuple(trig->pt(), trig->etaMS(), trig->phiMS());
316 }
317 
318 std::tuple<float,float,float> L2OverlapRemoverMon :: L2ORPosForMatchFunc(const xAOD::L2CombinedMuon *trig){
319  return std::forward_as_tuple( (trig->pt()/1e3 * trig->charge() ), trig->eta(), trig->phi());
320 }
xAOD::L2CombinedMuon_v1::phi
virtual double phi() const
The azimuthal angle ( ) of the particle.
AllowedVariables::e
e
Definition: AsgElectronSelectorTool.cxx:37
L2OverlapRemoverMon::isOverlap
bool isOverlap(const std::string &chain, const ElementLink< xAOD::L2StandAloneMuonContainer > &muEL1, const ElementLink< xAOD::L2StandAloneMuonContainer > &muEL2) const
Definition: L2OverlapRemoverMon.cxx:24
TrigMuonMonitorAlgorithm::m_group
Gaudi::Property< std::string > m_group
Name of monitored group.
Definition: TrigMuonMonitorAlgorithm.h:141
python.SystemOfUnits.m2
int m2
Definition: SystemOfUnits.py:92
L2OverlapRemoverMon::L2OverlapRemoverMon
L2OverlapRemoverMon(const std::string &name, ISvcLocator *pSvcLocator)
Definition: L2OverlapRemoverMon.cxx:9
runLayerRecalibration.chain
chain
Definition: runLayerRecalibration.py:175
ParticleGun_SamplingFraction.eta2
eta2
Definition: ParticleGun_SamplingFraction.py:96
P4Helpers::invMass
double invMass(const I4Momentum &pA, const I4Momentum &pB)
invariant mass from two I4momentum references
Definition: P4Helpers.h:252
xAOD::L2StandAloneMuon_v2::etaMS
float etaMS() const
Get the eta at muon spectrometer.
Base_Fragment.mass
mass
Definition: Sherpa_i/share/common/Base_Fragment.py:59
egammaEnergyPositionAllSamples::e1
double e1(const xAOD::CaloCluster &cluster)
return the uncorrected cluster energy in 1st sampling
xAOD::L2StandAloneMuon_v2
Class describing standalone muons reconstructed in the LVL2 trigger.
Definition: L2StandAloneMuon_v2.h:36
L2OverlapRemoverMon::L2ORPosForMatchFunc
static std::tuple< float, float, float > L2ORPosForMatchFunc(const xAOD::L2StandAloneMuon *trig)
Definition: L2OverlapRemoverMon.cxx:314
xAOD::L2CombinedMuon_v1::charge
float charge() const
get seeding muon charge
test_pyathena.pt
pt
Definition: test_pyathena.py:11
L2OverlapRemoverMon::m_massSAThresBB
Gaudi::Property< float > m_massSAThresBB
Definition: L2OverlapRemoverMon.h:51
xAOD::P4Helpers::deltaPhi
double deltaPhi(double phiA, double phiB)
delta Phi in range [-pi,pi[
Definition: xAODP4Helpers.h:69
xAOD::eta1
setEt setPhi setE277 setWeta2 eta1
Definition: TrigEMCluster_v1.cxx:41
xAOD::L2StandAloneMuon_v2::phiMS
float phiMS() const
Get the phi at muon spectrometer.
drawFromPickle.cos
cos
Definition: drawFromPickle.py:36
drawFromPickle.exp
exp
Definition: drawFromPickle.py:36
L2OverlapRemoverMon::calcinvMass
static float calcinvMass(double m1, double pt1, double eta1, double phi1, double m2, double pt2, double eta2, double phi2)
Definition: L2OverlapRemoverMon.cxx:280
Monitored::Collection
ValuesCollection< T > Collection(std::string name, const T &collection)
Declare a monitored (double-convertible) collection.
Definition: MonitoredCollection.h:38
L2OverlapRemoverMon::m_etaBins
Gaudi::Property< std::vector< float > > m_etaBins
Definition: L2OverlapRemoverMon.h:54
L2OverlapRemoverMon::fillVariablesPerChain
virtual StatusCode fillVariablesPerChain(const EventContext &ctx, const std::string &chain) const override
Function that fills variables of trigger objects associated to specified trigger chains.
Definition: L2OverlapRemoverMon.cxx:14
xAOD::L2CombinedMuon_v1
Class describing combined muon reconstructed in the LVL2 trigger.
Definition: L2CombinedMuon_v1.h:41
xAOD::L2CombinedMuon_v1::eta
virtual double eta() const
The pseudorapidity ( ) of the particle.
python.changerun.m1
m1
Definition: changerun.py:32
CheckAppliedSFs.e3
e3
Definition: CheckAppliedSFs.py:264
MuonMatchingTool.h
lumiFormat.i
int i
Definition: lumiFormat.py:85
EL::StatusCode
::StatusCode StatusCode
StatusCode definition for legacy code.
Definition: PhysicsAnalysis/D3PDTools/EventLoop/EventLoop/StatusCode.h:22
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition: AthMsgStreamMacros.h:29
L2OverlapRemoverMon.h
xAOD::L2CombinedMuon_v1::pt
virtual double pt() const
The transverse momentum ( ) of the particle.
ATH_CHECK
#define ATH_CHECK
Definition: AthCheckMacros.h:40
TrigMuonMonitorAlgorithm
Base class from which analyzers can define a derived class to do specific analysis.
Definition: TrigMuonMonitorAlgorithm.h:22
drawFromPickle.tan
tan
Definition: drawFromPickle.py:36
L2OverlapRemoverMon::m_massSAThresBE
Gaudi::Property< float > m_massSAThresBE
Definition: L2OverlapRemoverMon.h:53
AthMonitorAlgorithm::fill
void fill(const ToolHandle< GenericMonitoringTool > &groupHandle, std::vector< std::reference_wrapper< Monitored::IMonitoredVariable >> &&variables) const
Fills a vector of variables to a group by reference.
L2OverlapRemoverMon::chooseBestMuon
StatusCode chooseBestMuon(const std::string &chain, const std::vector< TrigCompositeUtils::LinkInfo< xAOD::L2StandAloneMuonContainer > > &featureCont, const std::vector< unsigned int > &muResult) const
Function that choose best muon.
Definition: L2OverlapRemoverMon.cxx:160
L2OverlapRemoverMon::m_dRSAThresBE
Gaudi::Property< float > m_dRSAThresBE
Definition: L2OverlapRemoverMon.h:52
ZERO_LIMIT
const float ZERO_LIMIT
Definition: VP1TriggerHandleL2.cxx:37
L2OverlapRemoverMon::m_dRSAThresEC
Gaudi::Property< std::vector< float > > m_dRSAThresEC
Definition: L2OverlapRemoverMon.h:55
L2OverlapRemoverMon::m_dRCBThres
Gaudi::Property< std::vector< float > > m_dRCBThres
Definition: L2OverlapRemoverMon.h:57
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:221
TrigCompositeUtils::LinkInfo
Helper to keep a Decision object, ElementLink and ActiveState (with respect to some requested ChainGr...
Definition: LinkInfo.h:28
L2OverlapRemoverMon::m_dRSAThresBB
Gaudi::Property< float > m_dRSAThresBB
Definition: L2OverlapRemoverMon.h:50
egammaEnergyPositionAllSamples::e2
double e2(const xAOD::CaloCluster &cluster)
return the uncorrected cluster energy in 2nd sampling
L2OverlapRemoverMon::m_massCBThres
Gaudi::Property< std::vector< float > > m_massCBThres
Definition: L2OverlapRemoverMon.h:59
TauGNNUtils::Variables::absEta
bool absEta(const xAOD::TauJet &tau, double &out)
Definition: TauGNNUtils.cxx:234
xAOD::L2StandAloneMuon_v2::pt
virtual double pt() const
The transverse momentum ( ) of the particle.
drawFromPickle.sin
sin
Definition: drawFromPickle.py:36
Monitored::Scalar
Declare a monitored scalar variable.
Definition: MonitoredScalar.h:34
L2OverlapRemoverMon::m_massSAThresEC
Gaudi::Property< std::vector< float > > m_massSAThresEC
Definition: L2OverlapRemoverMon.h:56
L2OverlapRemoverMon::m_dRbySAThres
Gaudi::Property< std::vector< float > > m_dRbySAThres
Definition: L2OverlapRemoverMon.h:58
fitman.k
k
Definition: fitman.py:528