Loading [MathJax]/extensions/tex2jax.js
ATLAS Offline Software
•All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
AnomDetVAE.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2025 CERN for the benefit of the ATLAS collaboration
3 */
4 /*********************************
5  * AnomDetVAE.cxx
6  * Created by Sagar Addepalli on 25/11/2024.
7  *
8  * @brief algorithm uses a variational auto-encoder based anomaly detection network
9  * Uses 6 jets, 4 muons, 4 taus, and MET to calculate the anomaly score
10  * based on the KL-divergence in a lower dimension latent space
11  *
12  * @param AnomalyScoreThresh
13 **********************************/
14 
15 #include <cmath>
16 
17 #include "L1TopoAlgorithms/AnomDetVAE.h"
18 #include "L1TopoCommon/Exception.h"
19 #include "L1TopoInterfaces/Decision.h"
20 #include <VAENetwork.h>
21 
22 REGISTER_ALG_TCS(ADVAE_2A)
23 
24 
25 TCS::ADVAE_2A::ADVAE_2A(const std::string & name) : DecisionAlg(name)
26 {
27  defineParameter("InputWidth1", 6);
28  defineParameter("InputWidth2", 6);
29  defineParameter("InputWidth3", 6);
30  defineParameter("InputWidth4", 1);
31  defineParameter("MaxTob1", 6);
32  defineParameter("MaxTob2", 4);
33  defineParameter("MaxTob3", 4);
34  defineParameter("MaxTob4", 1);
35  defineParameter("NumResultBits", 2);
36 
37  // Version parameter, used for L1TopoFW book-keeping, no practical application
38  defineParameter("ADVAEVersion", 1);
39  //minEt cuts, one per input list (TOBs failing these are set to ET = eta = phi = 0)
40  defineParameter("MinET1",0);
41  defineParameter("MinET2",0);
42  defineParameter("MinET3",0);
43  defineParameter("MinET4",0);
44  //The value used is sum of squares of the NN result vector elements, bit-shifted by 8 to get all decimal bits
45  defineParameter("AnomalyScoreThresh", 1000000, 0);
46  defineParameter("AnomalyScoreThresh", 1000000, 1);
47 
48  setNumberOutputBits(2);
49 }
50 
51 TCS::ADVAE_2A::~ADVAE_2A(){}
52 
53 
54 TCS::StatusCode
55 TCS::ADVAE_2A::initialize() {
56  p_NumberLeading1 = parameter("InputWidth1").value();
57  p_NumberLeading2 = parameter("InputWidth2").value();
58  p_NumberLeading3 = parameter("InputWidth3").value();
59  p_NumberLeading4 = parameter("InputWidth4").value();
60 
61  if(parameter("MaxTob1").value() > 0) p_NumberLeading1 = parameter("MaxTob1").value();
62  if(parameter("MaxTob2").value() > 0) p_NumberLeading2 = parameter("MaxTob2").value();
63  if(parameter("MaxTob3").value() > 0) p_NumberLeading3 = parameter("MaxTob3").value();
64  if(parameter("MaxTob4").value() > 0) p_NumberLeading4 = parameter("MaxTob4").value();
65 
66  p_minEt1 = parameter("MinET1").value();
67  p_minEt2 = parameter("MinET2").value();
68  p_minEt3 = parameter("MinET3").value();
69  p_minEt4 = parameter("MinET4").value();
70 
71  for(unsigned int i=0; i<numberOutputBits(); ++i) {
72  p_AnomalyScoreThresh[i] = parameter("AnomalyScoreThresh", i).value();
73  }
74 
75  TRG_MSG_INFO("number output : " << numberOutputBits());
76 
77  // book histograms
78  for(unsigned int i=0; i<numberOutputBits(); ++i) {
79  std::string hname_accept = "hAnomalyScore_accept_bit"+std::to_string((int)i);
80  std::string hname_reject = "hAnomalyScore_reject_bit"+std::to_string((int)i);
81  // score
82  bookHist(m_histAccept, hname_accept, "ADScore", 2000, 0, 2000000);
83  bookHist(m_histReject, hname_reject, "ADScore", 2000, 0, 2000000);
84  }
85 
86  return StatusCode::SUCCESS;
87 }
88 
89 
90 
91 TCS::StatusCode
92 TCS::ADVAE_2A::processBitCorrect( const std::vector<TCS::TOBArray const *> & input,
93  const std::vector<TCS::TOBArray *> & output,
94  Decision & decision )
95 {
96 
97 
98  if( input.size() == 4) {
99 
100  TCS::TOBArray const* jets = input[0];
101  TCS::TOBArray const* taus = input[1];
102  TCS::TOBArray const* mus = input[2];
103  TCS::TOBArray const* met = input[3];
104  TRG_MSG_DEBUG("Number of jets are " << (*jets).size());
105  TRG_MSG_DEBUG("Number of taus are " << (*taus).size());
106  TRG_MSG_DEBUG("Number of mus are " << (*mus).size());
107  TRG_MSG_DEBUG("Number of met are " << (*met).size());
108 
109  std::vector<u_int> jet_pt(6,0), tau_pt(4,0), mu_pt(4,0), met_pt(1,0);
110  std::vector<int> jet_eta(6,0), tau_eta(4,0), mu_eta(4,0); //no met_eta
111  std::vector<int> jet_phi(6,0), tau_phi(4,0), mu_phi(4,0), met_phi(1,0);
112 
113  for (u_int i = 0; i<(*jets).size() && i<6; ++i) {
114  if ( parType_t( (*jets)[i].Et() ) <= p_minEt1 ) continue; //ET cut, leave NN inputs at default values (0)
115  jet_pt[i] = (*jets)[i].Et();
116  jet_eta[i] = (*jets)[i].eta();
117  jet_phi[i] = (*jets)[i].phi();
118  }
119  for (u_int i = 0; i < (*taus).size() && i<4; ++i) {
120  if ( parType_t( (*taus)[i].Et() ) <= p_minEt2 ) continue; //ET cut, leave NN inputs at default values (0)
121  tau_pt[i] = (*taus)[i].Et();
122  tau_eta[i] = (*taus)[i].eta();
123  tau_phi[i] = (*taus)[i].phi();
124  }
125  for (u_int i = 0; i < (*mus).size() && i<4; ++i) {
126  if ( parType_t( (*mus)[i].Et() ) <= p_minEt3 ) continue; //ET cut, leave NN inputs at default values (0)
127  mu_pt[i] = (*mus)[i].Et();
128  mu_eta[i] = (*mus)[i].eta();
129  mu_phi[i] = (*mus)[i].phi();
130  }
131  for (u_int i = 0; i < (*met).size() && i<1; ++i) {
132  if ( parType_t( (*met)[i].Et() ) <= p_minEt4 ) continue; //ET cut, leave NN inputs at default values (0)
133  met_pt[i] = (*met)[i].Et();
134  met_phi[i] = (*met)[i].phi();
135  }
136 
137 
138  TRG_MSG_DEBUG("Jet0: " << jet_pt[0] << ", " << jet_eta[0] << ", " << jet_phi[0] );
139  TRG_MSG_DEBUG("Jet1: " << jet_pt[1] << ", " << jet_eta[1] << ", " << jet_phi[1] );
140  TRG_MSG_DEBUG("Jet2: " << jet_pt[2] << ", " << jet_eta[2] << ", " << jet_phi[2] );
141  TRG_MSG_DEBUG("Jet3: " << jet_pt[3] << ", " << jet_eta[3] << ", " << jet_phi[3] );
142  TRG_MSG_DEBUG("Jet4: " << jet_pt[4] << ", " << jet_eta[4] << ", " << jet_phi[4] );
143  TRG_MSG_DEBUG("Jet5: " << jet_pt[5] << ", " << jet_eta[5] << ", " << jet_phi[5] );
144 
145  TRG_MSG_DEBUG("Tau0: " << tau_pt[0] << ", " << tau_eta[0] << ", " << tau_phi[0] );
146  TRG_MSG_DEBUG("Tau1: " << tau_pt[1] << ", " << tau_eta[1] << ", " << tau_phi[1] );
147  TRG_MSG_DEBUG("Tau2: " << tau_pt[2] << ", " << tau_eta[2] << ", " << tau_phi[2] );
148  TRG_MSG_DEBUG("Tau3: " << tau_pt[3] << ", " << tau_eta[3] << ", " << tau_phi[3] );
149 
150  TRG_MSG_DEBUG("Mu0: " << mu_pt[0] << ", " << mu_eta[0] << ", " << mu_phi[0] );
151  TRG_MSG_DEBUG("Mu1: " << mu_pt[1] << ", " << mu_eta[1] << ", " << mu_phi[1] );
152  TRG_MSG_DEBUG("Mu2: " << mu_pt[2] << ", " << mu_eta[2] << ", " << mu_phi[2] );
153  TRG_MSG_DEBUG("Mu3: " << mu_pt[3] << ", " << mu_eta[3] << ", " << mu_phi[3] );
154 
155  TRG_MSG_DEBUG("MET: " << met_pt[0] << ", " << met_phi[0] << std::endl);
156 
157  ADVAE2A::VAENetwork AD_Network( jet_pt[0], jet_eta[0], jet_phi[0],
158  jet_pt[1], jet_eta[1], jet_phi[1],
159  jet_pt[2], jet_eta[2], jet_phi[2],
160  jet_pt[3], jet_eta[3], jet_phi[3],
161  jet_pt[4], jet_eta[4], jet_phi[4],
162  jet_pt[5], jet_eta[5], jet_phi[5],
163  tau_pt[0], tau_eta[0], tau_phi[0],
164  tau_pt[1], tau_eta[1], tau_phi[1],
165  tau_pt[2], tau_eta[2], tau_phi[2],
166  tau_pt[3], tau_eta[3], tau_phi[3],
167  mu_pt [0], mu_eta [0], mu_phi [0],
168  mu_pt [1], mu_eta [1], mu_phi [1],
169  mu_pt [2], mu_eta [2], mu_phi [2],
170  mu_pt [3], mu_eta [3], mu_phi [3],
171  met_pt[0], met_phi[0] );
172  int64_t anomScoreInt64 = AD_Network.getAnomalyScoreInt64();
173 
174  for(u_int i=0; i<numberOutputBits(); ++i) {
175  bool accept = false;
176  int32_t threshold = int32_t ( p_AnomalyScoreThresh[i] );
177  if ( anomScoreInt64 > threshold ) {
178  accept = true;
179  decision.setBit(i, true);
180  for ( u_int j = 0; j<6 && j<(*jets).size(); ++j ) output[i]->push_back((*jets)[j]);
181  for ( u_int j = 0; j<4 && j<(*taus).size(); ++j ) output[i]->push_back((*taus)[j]);
182  for ( u_int j = 0; j<4 && j<(*mus).size() ; ++j ) output[i]->push_back((*mus) [j]);
183  output[i]->push_back((*met)[0]);
184  }
185 
186  if(fillHistos() and accept) {
187  fillHist1D(m_histAccept[i],anomScoreInt64);
188  } else if(fillHistos() && !accept) {
189  fillHist1D(m_histReject[i],anomScoreInt64);
190  }
191 
192  TRG_MSG_DEBUG("Decision for bit" << i << ": " << (accept?"pass":"fail") << " anomaly score = " << anomScoreInt64 << std::endl);
193  }
194  } else {
195  TCS_EXCEPTION("ADVAE_2A alg must have 4 inputs, but got " << input.size());
196  }
197 
198  return TCS::StatusCode::SUCCESS;
199 }
200 
201 TCS::StatusCode
202 TCS::ADVAE_2A::process( const std::vector<TCS::TOBArray const *> & input,
203  const std::vector<TCS::TOBArray *> & output,
204  Decision & decision )
205 {
206 
207 
208  if( input.size() == 4) {
209 
210  TCS::TOBArray const* jets = input[0];
211  TCS::TOBArray const* taus = input[1];
212  TCS::TOBArray const* mus = input[2];
213  TCS::TOBArray const* met = input[3];
214  TRG_MSG_DEBUG("Number of jets are " << (*jets).size());
215  TRG_MSG_DEBUG("Number of taus are " << (*taus).size());
216  TRG_MSG_DEBUG("Number of mus are " << (*mus).size());
217  TRG_MSG_DEBUG("Number of met are " << (*met).size());
218 
219  std::vector<u_int> jet_pt(6,0), tau_pt(4,0), mu_pt(4,0), met_pt(1,0);
220  std::vector<int> jet_eta(6,0), tau_eta(4,0), mu_eta(4,0); //no met_eta
221  std::vector<int> jet_phi(6,0), tau_phi(4,0), mu_phi(4,0), met_phi(1,0);
222 
223  for (u_int i = 0; i<(*jets).size() && i<6; ++i) {
224  if ( parType_t( (*jets)[i].Et() ) <= p_minEt1 ) continue; //ET cut, leave NN inputs at default values (0)
225  jet_pt[i] = (*jets)[i].Et();
226  jet_eta[i] = (*jets)[i].eta();
227  jet_phi[i] = (*jets)[i].phi();
228  }
229  for (u_int i = 0; i < (*taus).size() && i<4; ++i) {
230  if ( parType_t( (*taus)[i].Et() ) <= p_minEt2 ) continue; //ET cut, leave NN inputs at default values (0)
231  tau_pt[i] = (*taus)[i].Et();
232  tau_eta[i] = (*taus)[i].eta();
233  tau_phi[i] = (*taus)[i].phi();
234  }
235  for (u_int i = 0; i < (*mus).size() && i<4; ++i) {
236  if ( parType_t( (*mus)[i].Et() ) <= p_minEt3 ) continue; //ET cut, leave NN inputs at default values (0)
237  mu_pt[i] = (*mus)[i].Et();
238  mu_eta[i] = (*mus)[i].eta();
239  mu_phi[i] = (*mus)[i].phi();
240  }
241  for (u_int i = 0; i < (*met).size() && i<1; ++i) {
242  if ( parType_t( (*met)[i].Et() ) <= p_minEt4 ) continue; //ET cut, leave NN inputs at default values (0)
243  met_pt[i] = (*met)[i].Et();
244  met_phi[i] = (*met)[i].phi();
245  }
246 
247  ADVAE2A::VAENetwork AD_Network( jet_pt[0], jet_eta[0], jet_phi[0],
248  jet_pt[1], jet_eta[1], jet_phi[1],
249  jet_pt[2], jet_eta[2], jet_phi[2],
250  jet_pt[3], jet_eta[3], jet_phi[3],
251  jet_pt[4], jet_eta[4], jet_phi[4],
252  jet_pt[5], jet_eta[5], jet_phi[5],
253  tau_pt[0], tau_eta[0], tau_phi[0],
254  tau_pt[1], tau_eta[1], tau_phi[1],
255  tau_pt[2], tau_eta[2], tau_phi[2],
256  tau_pt[3], tau_eta[3], tau_phi[3],
257  mu_pt [0], mu_eta [0], mu_phi [0],
258  mu_pt [1], mu_eta [1], mu_phi [1],
259  mu_pt [2], mu_eta [2], mu_phi [2],
260  mu_pt [3], mu_eta [3], mu_phi [3],
261  met_pt[0], met_phi[0] );
262  int64_t anomScoreInt64 = AD_Network.getAnomalyScoreInt64();
263 
264  for(u_int i=0; i<numberOutputBits(); ++i) {
265  bool accept = false;
266  int32_t threshold = int32_t ( p_AnomalyScoreThresh[i] );
267  if ( anomScoreInt64 > threshold ) {
268  accept = true;
269  decision.setBit(i, true);
270  for ( u_int j = 0; j<6 && j<(*jets).size(); ++j ) output[i]->push_back((*jets)[j]);
271  for ( u_int j = 0; j<4 && j<(*taus).size(); ++j ) output[i]->push_back((*taus)[j]);
272  for ( u_int j = 0; j<4 && j<(*mus).size() ; ++j ) output[i]->push_back((*mus) [j]);
273  output[i]->push_back((*met)[0]);
274  }
275 
276  if(fillHistos() and accept) {
277  fillHist1D(m_histAccept[i],anomScoreInt64);
278  } else if(fillHistos() && !accept) {
279  fillHist1D(m_histReject[i],anomScoreInt64);
280  }
281 
282  TRG_MSG_DEBUG("Decision for bit" << i << ": " << (accept?"pass":"fail") << " anomaly score = " << anomScoreInt64 << std::endl);
283  }
284  } else {
285  TCS_EXCEPTION("ADVAE_2A alg must have 4 inputs, but got " << input.size());
286  }
287 
288  return TCS::StatusCode::SUCCESS;
289 }
TCS::ADVAE_2A::initialize
virtual StatusCode initialize()
Definition: AnomDetVAE.cxx:55
TCS::StatusCode::SUCCESS
@ SUCCESS
Definition: Trigger/TrigT1/L1Topo/L1TopoCommon/L1TopoCommon/StatusCode.h:17
TCS::parType_t
uint32_t parType_t
Definition: Parameter.h:22
CutsMETMaker::accept
StatusCode accept(const xAOD::Muon *mu)
Definition: CutsMETMaker.cxx:18
defineDB.jets
jets
Definition: JetTagCalibration/share/defineDB.py:24
athena.value
value
Definition: athena.py:124
TCS::ADVAE_2A
Definition: AnomDetVAE.h:16
Decision.h
const
bool const RAWDATA *ch2 const
Definition: LArRodBlockPhysicsV0.cxx:560
TCS::DecisionAlg
Definition: Trigger/TrigT1/L1Topo/L1TopoInterfaces/L1TopoInterfaces/DecisionAlg.h:22
TCS::Decision::setBit
void setBit(unsigned int index, bool value)
Definition: L1Topo/L1TopoInterfaces/Root/Decision.cxx:12
met
Definition: IMETSignificance.h:24
lumiFormat.i
int i
Definition: lumiFormat.py:85
TCS_EXCEPTION
#define TCS_EXCEPTION(MSG)
Definition: Trigger/TrigT1/L1Topo/L1TopoCommon/L1TopoCommon/Exception.h:14
TCS::TOBArray
Definition: TOBArray.h:24
plotIsoValidation.mu_phi
mu_phi
Definition: plotIsoValidation.py:152
TCS::Decision
Definition: L1Topo/L1TopoInterfaces/L1TopoInterfaces/Decision.h:19
TRG_MSG_INFO
#define TRG_MSG_INFO(x)
Definition: Trigger/TrigConfiguration/TrigConfBase/TrigConfBase/MsgStreamMacros.h:27
TCS::ADVAE_2A::processBitCorrect
virtual StatusCode processBitCorrect(const std::vector< TCS::TOBArray const * > &input, const std::vector< TCS::TOBArray * > &output, Decision &decison)
Definition: AnomDetVAE.cxx:92
TCS::ADVAE_2A::process
virtual StatusCode process(const std::vector< TCS::TOBArray const * > &input, const std::vector< TCS::TOBArray * > &output, Decision &decison)
Definition: AnomDetVAE.cxx:202
REGISTER_ALG_TCS
#define REGISTER_ALG_TCS(CLASS)
Definition: AlgFactory.h:62
AnomDetVAE.h
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:240
ActsTrk::to_string
std::string to_string(const DetectorType &type)
Definition: GeometryDefs.h:34
threshold
Definition: chainparser.cxx:74
TCS
Definition: Global/GlobalSimulation/src/IO/Decision.h:18
plotIsoValidation.mu_eta
mu_eta
Definition: plotIsoValidation.py:151
Trk::jet_phi
@ jet_phi
Definition: JetVtxParamDefs.h:28
Exception.h
TCS::ADVAE_2A::~ADVAE_2A
virtual ~ADVAE_2A()
Definition: AnomDetVAE.cxx:51
TRG_MSG_DEBUG
#define TRG_MSG_DEBUG(x)
Definition: Trigger/TrigConfiguration/TrigConfBase/TrigConfBase/MsgStreamMacros.h:25
TCS::StatusCode
Definition: Trigger/TrigT1/L1Topo/L1TopoCommon/L1TopoCommon/StatusCode.h:15
plotIsoValidation.mu_pt
mu_pt
Definition: plotIsoValidation.py:150
Generated on Thu May 8 2025 21:06:09 for ATLAS Offline Software by   doxygen 1.8.18