ATLAS Offline Software
Loading...
Searching...
No Matches
gFEXJwoJAlgo.cxx
Go to the documentation of this file.
1/*
2 Copyright (C) 2002-2026 CERN for the benefit of the ATLAS collaboration
3*/
4//***************************************************************************
5// gFEXJwoJAlg - Jets without jets algorithm for gFEX
6// -------------------
7// begin : 21 12 2023
8// email : cecilia.tosciri@cern.ch
9//***************************************************************************
10
11
12
13#include "gFEXJwoJAlgo.h"
14#include "L1CaloFEXSim/gFEXJwoJTOB.h"
15
16#include <cmath> //for std::sqrt
17#include <vector>
18
19namespace LVL1 {
20
21 // default constructor for persistency
22gFEXJwoJAlgo::gFEXJwoJAlgo(const std::string& type, const std::string& name, const IInterface* parent):
23 AthAlgTool(type, name, parent)
24 {
25 declareInterface<IgFEXJwoJAlgo>(this);
26 }
27
28
29StatusCode gFEXJwoJAlgo::initialize(){
30
31 ATH_CHECK(m_DBToolKey.initialize());
32 ATH_CHECK(m_l1MenuKey.initialize());
33
34 return StatusCode::SUCCESS;
35}
36
37
38std::vector<std::unique_ptr<gFEXJwoJTOB>> gFEXJwoJAlgo::jwojAlgo(const EventContext& ctx,
39 const gTowersType& Atwr, int pucA_JWJ,
40 const gTowersType& Btwr, int pucB_JWJ,
41 const gTowersType& Ctwr, int pucC_JWJ,
42 std::array<int32_t, 4> & outTOB) const {
43
44
45 SG::ReadCondHandle<gFEXDBCondData> myDBTool = SG::ReadCondHandle<gFEXDBCondData>(m_DBToolKey, ctx);
46 if (!myDBTool.isValid()) {
47 ATH_MSG_ERROR("Could not retrieve DB tool " << m_DBToolKey);
48 throw std::runtime_error("Could not retrieve DB tool");
49 }
50
51 std::string fwVersion = myDBTool->get_FWVersion();
52 int major = std::stoi(fwVersion);
53 bool SumETfast = (major >= 1);
54 bool metRho = (major >= 2);
55
56 // Retrieve the L1 menu configuration
57 SG::ReadHandle<TrigConf::L1Menu> l1Menu (m_l1MenuKey, ctx);
58 if (!l1Menu.isValid()) {
59 ATH_MSG_ERROR("Could not retrieve L1Menu " << m_l1MenuKey);
60 throw std::runtime_error("Could not retrieve L1Menu");
61 }
62
63 //Parameters related to gXE (MET objects, both JwoJ and alternative MET calculation)
64 const auto & thr_gXE = l1Menu->thrExtraInfo().gXE();
65 int gBlockthresholdA = thr_gXE.seedThr('A'); //defined in counts
66 int gBlockthresholdB = thr_gXE.seedThr('B'); //defined in counts
67 int gBlockthresholdC = thr_gXE.seedThr('C'); //defined in counts
68 int aFPGA_A = thr_gXE.JWOJ_param('A','a');// 1003
69 int bFPGA_A = thr_gXE.JWOJ_param('A','b');// 409
70 int aFPGA_B = thr_gXE.JWOJ_param('B','a');// 1003
71 int bFPGA_B = thr_gXE.JWOJ_param('B','b');// 409
72 int aFPGA_C = thr_gXE.JWOJ_param('C','a');// 1003
73 int bFPGA_C = thr_gXE.JWOJ_param('C','b');// 409
74
75 // find gBlocks
76 gTowersType AgBlk;
77 gTowersType Ascaled;
78
79 gTowersType BgBlk;
80 gTowersType Bscaled;
81
82 gTowersType CgBlk;
83 gTowersType Cscaled;
84
85 gTowersType hasSeed;
86
87 gBlockAB(Atwr, AgBlk, hasSeed, gBlockthresholdA);
88 gBlockAB(Btwr, BgBlk, hasSeed, gBlockthresholdB);
89 gBlockAB(Ctwr, CgBlk, hasSeed, gBlockthresholdC);
90
91
92 // switch to 10 bit number
93 // DMS -- do we eventaully need to check for overflows here?
94
95 for(int irow = 0; irow<FEXAlgoSpaceDefs::ABCrows; irow++){
96 for(int jcolumn = 0; jcolumn<FEXAlgoSpaceDefs::ABcolumns; jcolumn++){
97 Ascaled[irow][jcolumn] = Atwr[irow][jcolumn] >> 2;
98 AgBlk[irow][jcolumn] = AgBlk[irow][jcolumn] >> 2;
99
100 Bscaled[irow][jcolumn] = Btwr[irow][jcolumn] >> 2;
101 BgBlk[irow][jcolumn] = BgBlk[irow][jcolumn] >> 2;
102
103 Cscaled[irow][jcolumn] = Ctwr[irow][jcolumn] >> 2;
104 CgBlk[irow][jcolumn] = CgBlk[irow][jcolumn] >> 2;
105
106 }
107 }
108
109
110 //FPGA A observables
111 int A_MHT_x = 0x0;
112 int A_MHT_y = 0x0;
113 int A_MST_x = 0x0;
114 int A_MST_y = 0x0;
115 int A_MET_x = 0x0;
116 int A_MET_y = 0x0;
117
118 int A_eth = 0x0;
119 int A_ets = 0x0;
120 int A_etw = 0x0;
121
122 //FPGA B observables
123 int B_MHT_x = 0x0;
124 int B_MHT_y = 0x0;
125 int B_MST_x = 0x0;
126 int B_MST_y = 0x0;
127 int B_MET_x = 0x0;
128 int B_MET_y = 0x0;
129
130 int B_eth = 0x0;
131 int B_ets = 0x0;
132 int B_etw = 0x0;
133
134 //FPGA C observables
135 int C_MHT_x = 0x0;
136 int C_MHT_y = 0x0;
137 int C_MST_x = 0x0;
138 int C_MST_y = 0x0;
139 int C_MET_x = 0x0;
140 int C_MET_y = 0x0;
141
142 int C_eth = 0x0;
143 int C_ets = 0x0;
144 int C_etw = 0x0;
145
146 //Global observables
147 int MHT_x = 0x0;
148 int MHT_y = 0x0;
149 int MST_x = 0x0;
150 int MST_y = 0x0;
151 int MET_x = 0x0;
152 int MET_y = 0x0;
153
154 int ETH = 0x0;
155 int ETS = 0x0;
156 int ETW = 0x0;
157
158 int total_sumEt = 0x0;
159 int MET = 0x0;
160
161 // will need to hard code etFPGA ,a's and b's
162 int etBprime = 0;
163
164 if (metRho) metFPGA_rho(0, Ascaled, pucA_JWJ, AgBlk, gBlockthresholdA, aFPGA_A, bFPGA_A, A_MHT_x, A_MHT_y, A_MST_x, A_MST_y, A_MET_x, A_MET_y);
165 else metFPGA(0, Ascaled, AgBlk, gBlockthresholdA, aFPGA_A, bFPGA_A, A_MHT_x, A_MHT_y, A_MST_x, A_MST_y, A_MET_x, A_MET_y);
166 if (SumETfast) etFastFPGA(0, Ascaled, AgBlk, gBlockthresholdA, aFPGA_A, etBprime, A_eth, A_ets, A_etw);
167 else etFPGA(0, Ascaled, AgBlk, gBlockthresholdA, aFPGA_A, etBprime, A_eth, A_ets, A_etw);
168
169 if (metRho) metFPGA_rho(1, Bscaled, pucB_JWJ, BgBlk, gBlockthresholdB, aFPGA_B, bFPGA_B, B_MHT_x, B_MHT_y, B_MST_x, B_MST_y, B_MET_x, B_MET_y);
170 else metFPGA(1, Bscaled, BgBlk, gBlockthresholdB, aFPGA_B, bFPGA_B, B_MHT_x, B_MHT_y, B_MST_x, B_MST_y, B_MET_x, B_MET_y);
171 if (SumETfast) etFastFPGA(1, Bscaled, BgBlk, gBlockthresholdB, aFPGA_B, etBprime, B_eth, B_ets, B_etw);
172 else etFPGA(1, Bscaled, BgBlk, gBlockthresholdB, aFPGA_B, etBprime, B_eth, B_ets, B_etw);
173
174 if (metRho) metFPGA_rho(2, Cscaled, pucC_JWJ, CgBlk, gBlockthresholdC, aFPGA_C, bFPGA_C, C_MHT_x, C_MHT_y, C_MST_x, C_MST_y, C_MET_x, C_MET_y);
175 else metFPGA(2, Cscaled, CgBlk, gBlockthresholdC, aFPGA_C, bFPGA_C, C_MHT_x, C_MHT_y, C_MST_x, C_MST_y, C_MET_x, C_MET_y);
176 if (SumETfast) etFastFPGA(2, Cscaled, CgBlk, gBlockthresholdC, aFPGA_C, etBprime, C_eth, C_ets, C_etw);
177 else etFPGA(2, Cscaled, CgBlk, gBlockthresholdC, aFPGA_C, etBprime, C_eth, C_ets, C_etw);
178
179 metTotal(A_MHT_x, A_MHT_y, B_MHT_x, B_MHT_y, C_MHT_x, C_MHT_y, MHT_x, MHT_y);
180 metTotal(A_MST_x, A_MST_y, B_MST_x, B_MST_y, C_MST_x, C_MST_y, MST_x, MST_y);
181 metTotal(A_MET_x, A_MET_y, B_MET_x, B_MET_y, C_MET_x, C_MET_y, MET_x, MET_y);
182
183 etTotal(A_eth, B_eth, C_eth, ETH);
184 etTotal(A_ets, B_ets, C_ets, ETS);
185 etTotal(A_etw, B_etw, C_etw, ETW);
186 total_sumEt = ETW;
187
188 // components should all be less than 12 bits at this point with 200 MeV LSB
189 int MET2 = MET_x * MET_x + MET_y * MET_y;
190
191 if (MET2 > 0x0FFFFFF) {
192 MET = 0x000FFF;
193 } else {
194 // repeat the byte stream converter calculation here -- not what the hardware actually does
195 MET = std::sqrt(MET2);
196
197
198 // best guess at current hardware. Note that this is computed in the bytestream converter
199 // take most signficant 12 bits
200 //int MET12 = MET2 >> 12;
201 // simulate the look up -- only 6 most signficant bits currently set -- to be checked
202 //MET = ( (int)(std::sqrt(MET12)) << 6) & 0x00000FC0 ;
203 }
204
205
206 //Define a vector to be filled with all the TOBs of one event
207 std::vector<std::unique_ptr<gFEXJwoJTOB>> tobs_v;
208 tobs_v.resize(4);
209
210
211 // fill in TOBs
212 // The order of the TOBs is given according to the TOB ID (TODO: check how it's done in fw)
213
214 // First TOB is (MET, SumEt)
215 outTOB[0] = (total_sumEt & 0x00000FFF) << 0; //set the Quantity2 to the corresponding slot (LSB)
216 outTOB[0] = outTOB[0] | (MET & 0x00000FFF) << 12;//Quantity 1 (in bit number 12)
217 if (total_sumEt != 0) outTOB[0] = outTOB[0] | 0x00000001 << 24;//Status bit for Quantity 2 (0 if quantity is null)
218 if (MET != 0) outTOB[0] = outTOB[0] | 0x00000001 << 25;//Status bit for Quantity 1 (0 if quantity is null)
219 outTOB[0] = outTOB[0] | (1 & 0x0000001F) << 26;//TOB ID is 1 for scalar values (5 bits starting at 26)
220
221 // std::cout << "DMS MET " << std::hex << MET << " total_sumEt " << total_sumEt << std::endl << std::dec;
222
223// Second TOB is (MET_x, MET_y)
224 outTOB[1] = (MET_y & 0x00000FFF) << 0; //set the Quantity2 to the corresponding slot (LSB)
225 outTOB[1] = outTOB[1] | (MET_x & 0x00000FFF) << 12;//Quantity 1 (in bit number 12)
226 if (MET_y != 0) outTOB[1] = outTOB[1] | 0x00000001 << 24;//Status bit for Quantity 2 (0 if quantity is null)
227 if (MET_x != 0) outTOB[1] = outTOB[1] | 0x00000001 << 25;//Status bit for Quantity 1 (0 if quantity is null)
228 outTOB[1] = outTOB[1] | (2 & 0x0000001F) << 26;//TOB ID is 2 for MET_x, MET_y (5 bits starting at 26)
229
230// Third TOB is hard components (MHT_x, MHT_y)
231 outTOB[2] = (MHT_y & 0x00000FFF) << 0; //set the Quantity2 to the corresponding slot (LSB)
232 outTOB[2] = outTOB[2] | (MHT_x & 0x00000FFF) << 12;//Quantity 1 (in bit number 12)
233 if (MHT_y != 0) outTOB[2] = outTOB[2] | 0x00000001 << 24;//Status bit for Quantity 2 (0 if quantity is null)
234 if (MHT_x != 0) outTOB[2] = outTOB[2] | 0x00000001 << 25;//Status bit for Quantity 1 (0 if quantity is null)
235 outTOB[2] = outTOB[2] | (3 & 0x0000001F) << 26;//TOB ID is 3 for hard components (5 bits starting at 26)
236
237 // Fourth TOB is hard components (MST_x, MST_y)
238 outTOB[3] = (MST_y & 0x00000FFF) << 0; //set the Quantity2 to the corresponding slot (LSB)
239 outTOB[3] = outTOB[3] | (MST_x & 0x00000FFF) << 12;//Quantity 1 (in bit number 12)
240 if (MST_y != 0) outTOB[3] = outTOB[3] | 0x00000001 << 24;//Status bit for Quantity 2 (0 if quantity is null)
241 if (MST_x != 0) outTOB[3] = outTOB[3] | 0x00000001 << 25;//Status bit for Quantity 1 (0 if quantity is null)
242 outTOB[3] = outTOB[3] | (4 & 0x0000001F) << 26;//TOB ID is 4 for soft components (5 bits starting at 26)
243
244
245 tobs_v[0] = std::make_unique<gFEXJwoJTOB>();
246 tobs_v[0]->setWord(outTOB[0]);
247 tobs_v[0]->setQuantity1(MET);
248 tobs_v[0]->setQuantity2(total_sumEt);
249 tobs_v[0]->setSaturation(0); //Always 0 for now, need a threshold?
250 tobs_v[0]->setTobID(1);
251 if( MET != 0 ) tobs_v[0]->setStatus1(1);
252 else tobs_v[0]->setStatus1(0);
253 if(total_sumEt!= 0) tobs_v[0]->setStatus2(1);
254 else tobs_v[0]->setStatus2(0);
255
256 tobs_v[1] = std::make_unique<gFEXJwoJTOB>();
257 tobs_v[1]->setWord(outTOB[1]);
258 tobs_v[1]->setQuantity1(MET_x);
259 tobs_v[1]->setQuantity2(MET_y);
260 tobs_v[1]->setSaturation(0); //Always 0 for now, need a threshold?
261 tobs_v[1]->setTobID(2);
262 if( MET_x != 0 ) tobs_v[1]->setStatus1(1);
263 else tobs_v[1]->setStatus1(0);
264 if(MET_y!= 0) tobs_v[1]->setStatus2(1);
265 else tobs_v[1]->setStatus2(0);
266
267 tobs_v[2] = std::make_unique<gFEXJwoJTOB>();
268 tobs_v[2]->setWord(outTOB[2]);
269 tobs_v[2]->setQuantity1(MHT_x);
270 tobs_v[2]->setQuantity2(MHT_y);
271 tobs_v[2]->setSaturation(0); //Always 0 for now, need a threshold?
272 tobs_v[2]->setTobID(3);
273 if( MHT_x != 0 ) tobs_v[2]->setStatus1(1);
274 else tobs_v[2]->setStatus1(0);
275 if(MHT_y!= 0) tobs_v[2]->setStatus2(1);
276 else tobs_v[2]->setStatus2(0);
277
278 tobs_v[3] = std::make_unique<gFEXJwoJTOB>();
279 tobs_v[3]->setWord(outTOB[3]);
280 tobs_v[3]->setQuantity1(MST_x);
281 tobs_v[3]->setQuantity2(MST_y);
282 tobs_v[3]->setSaturation(0); //Always 0 for now, need a threshold?
283 tobs_v[3]->setTobID(4);
284 if( MST_x != 0 ) tobs_v[3]->setStatus1(1);
285 else tobs_v[2]->setStatus1(0);
286 if(MST_y!= 0) tobs_v[3]->setStatus2(1);
287 else tobs_v[3]->setStatus2(0);
288
289
290 return tobs_v;
291
292}
293
294
295
296void gFEXJwoJAlgo::gBlockAB(const gTowersType & twrs, gTowersType & gBlkSum, gTowersType & hasSeed, int seedThreshold) const {
297
298 const int rows = twrs.size();
299 const int cols = twrs[0].size();
300 for( int irow = 0; irow < rows; irow++ ){
301 for(int jcolumn = 0; jcolumn<cols; jcolumn++){
302 // zero jet sum here
303 gBlkSum[irow][jcolumn] = 0;
304 int krowUp = (irow + 1)%32;
305 int krowDn = (irow - 1 +32)%32;
306 if( (jcolumn == 0) || (jcolumn == 6) ) {
307 //left edge case
308 gBlkSum[irow][jcolumn] =
309 twrs[irow][jcolumn] + twrs[krowUp][jcolumn] + twrs[krowDn][jcolumn] +
310 twrs[irow][jcolumn+1] + twrs[krowUp][jcolumn+1] + twrs[krowDn][jcolumn+1];
311 } else if( (jcolumn == 5) || (jcolumn == 11)) {
312 // right edge case
313 gBlkSum[irow][jcolumn] =
314 twrs[irow][jcolumn] + twrs[krowUp][jcolumn] + twrs[krowDn][jcolumn] +
315 twrs[irow][jcolumn-1] + twrs[krowUp][jcolumn-1] + twrs[krowDn][jcolumn-1];
316 } else{
317 // normal case; jcolumn is not 11 so does not overrun
318 //coverity[OVERRUN:FALSE]
319 gBlkSum[irow][jcolumn] =
320 twrs[irow][jcolumn] + twrs[krowUp][jcolumn] + twrs[krowDn][jcolumn] +
321 twrs[irow][jcolumn-1] + twrs[krowUp][jcolumn-1] + twrs[krowDn][jcolumn-1] +
322 twrs[irow][jcolumn+1] + twrs[krowUp][jcolumn+1] + twrs[krowDn][jcolumn+1];
323 }
324
325 if( gBlkSum[irow][jcolumn] > seedThreshold) {
326 hasSeed[irow][jcolumn] = 1;
327 } else {
328 hasSeed[irow][jcolumn] = 0;
329 }
330
331 if ( gBlkSum[irow][jcolumn] < 0 )
332 gBlkSum[irow][jcolumn] = 0;
333
334 // was bits 11+3 downto 3, now is 11 downto 0
335 if ( gBlkSum[irow][jcolumn] > FEXAlgoSpaceDefs::gBlockMax ) {
336 gBlkSum[irow][jcolumn] = FEXAlgoSpaceDefs::gBlockMax;
337 }
338 }
339 }
340}
341void gFEXJwoJAlgo::metFPGA_rho(int FPGAnum, const gTowersType& twrs, int puc_jwj,
342 const gTowersType & gBlkSum, int gBlockthreshold,
343 int aFPGA, int bFPGA,
344 int & MHT_x, int & MHT_y,
345 int & MST_x, int & MST_y,
346 int & MET_x, int & MET_y) const {
347 gBlockthreshold = gBlockthreshold * 200 / 800;
348
349 int64_t h_tx_hi = 0;
350 int64_t h_ty_hi = 0;
351 int64_t h_tx_lw = 0;
352 int64_t h_ty_lw = 0;
353
354 int64_t e_tx_hi = 0;
355 int64_t e_ty_hi = 0;
356 int64_t e_tx_lw = 0;
357 int64_t e_ty_lw = 0;
358
359
360 int64_t RHO_SUM_OF_COS_h_tx_hi = 0;
361 int64_t RHO_SUM_OF_SIN_h_ty_hi = 0;
362 int64_t RHO_SUM_OF_COS_h_tx_lw = 0;
363 int64_t RHO_SUM_OF_SIN_h_ty_lw = 0;
364
365 int64_t RHO_SUM_OF_COS_e_tx_hi = 0;
366 int64_t RHO_SUM_OF_SIN_e_ty_hi = 0;
367 int64_t RHO_SUM_OF_COS_e_tx_lw = 0;
368 int64_t RHO_SUM_OF_SIN_e_ty_lw = 0;
369
370
371 for( int irow = 0; irow < FEXAlgoSpaceDefs::ABCrows; irow++ ){
372 for(int jcolumn = 6; jcolumn<12; jcolumn++){
373 if( FPGAnum == 2){
374 int frow = 2*(irow/2) + 1;
375
376 if(gBlkSum[irow][jcolumn] > gBlockthreshold){
377 h_tx_hi += (twrs[irow][jcolumn])*(cosLUT(frow, 5));
378 h_ty_hi += (twrs[irow][jcolumn])*(sinLUT(frow, 5));
379 RHO_SUM_OF_COS_h_tx_hi += (cosLUT(frow, 5));
380 RHO_SUM_OF_SIN_h_ty_hi += (sinLUT(frow, 5));
381
382 } else {
383 e_tx_hi += (twrs[irow][jcolumn])*(cosLUT(frow, 5));
384 e_ty_hi += (twrs[irow][jcolumn])*(sinLUT(frow, 5));
385 RHO_SUM_OF_COS_e_tx_hi += (cosLUT(frow, 5));
386 RHO_SUM_OF_SIN_e_ty_hi += (sinLUT(frow, 5));
387 }
388
389 } else {
390
391 if(gBlkSum[irow][jcolumn] > gBlockthreshold){
392 h_tx_hi += (twrs[irow][jcolumn])*(cosLUT(irow, 5));
393 h_ty_hi += (twrs[irow][jcolumn])*(sinLUT(irow, 5));
394 RHO_SUM_OF_COS_h_tx_hi += (cosLUT(irow, 5));
395 RHO_SUM_OF_SIN_h_ty_hi += (sinLUT(irow, 5));
396 } else {
397 e_tx_hi += (twrs[irow][jcolumn])*(cosLUT(irow, 5));
398 e_ty_hi += (twrs[irow][jcolumn])*(sinLUT(irow, 5));
399 RHO_SUM_OF_COS_e_tx_hi += (cosLUT(irow, 5));
400 RHO_SUM_OF_SIN_e_ty_hi += (sinLUT(irow, 5));
401 }
402 }
403 }
404
405 for(int jcolumn = 0; jcolumn<6; jcolumn++){
406 if( FPGAnum == 2){
407 int frow = 2*(irow/2) + 1;
408
409 if(gBlkSum[irow][jcolumn] > gBlockthreshold){
410 h_tx_lw += (twrs[irow][jcolumn])*(cosLUT(frow, 5));
411 h_ty_lw += (twrs[irow][jcolumn])*(sinLUT(frow, 5));
412 RHO_SUM_OF_COS_h_tx_lw += (cosLUT(frow, 5));
413 RHO_SUM_OF_SIN_h_ty_lw += (sinLUT(frow, 5));
414 } else{
415 e_tx_lw += (twrs[irow][jcolumn])*(cosLUT(frow, 5));
416 e_ty_lw += (twrs[irow][jcolumn])*(sinLUT(frow, 5));
417 RHO_SUM_OF_COS_e_tx_lw += (cosLUT(frow, 5));
418 RHO_SUM_OF_SIN_e_ty_lw += (sinLUT(frow, 5));
419 }
420 } else {
421
422 if(gBlkSum[irow][jcolumn] > gBlockthreshold){
423 h_tx_lw += (twrs[irow][jcolumn])*(cosLUT(irow, 5));
424 h_ty_lw += (twrs[irow][jcolumn])*(sinLUT(irow, 5));
425 RHO_SUM_OF_COS_h_tx_lw += (cosLUT(irow, 5));
426 RHO_SUM_OF_SIN_h_ty_lw += (sinLUT(irow, 5));
427 } else {
428 e_tx_lw += (twrs[irow][jcolumn])*(cosLUT(irow, 5));
429 e_ty_lw += (twrs[irow][jcolumn])*(sinLUT(irow, 5));
430 RHO_SUM_OF_COS_e_tx_lw += (cosLUT(irow, 5));
431 RHO_SUM_OF_SIN_e_ty_lw += (sinLUT(irow, 5));
432 }
433 }
434 }
435 }
436
437 // REMEMBER TO DO BIT ADJUSTMENTS FOR SUBTRACTION
438
439 long int fMHT_x = (h_tx_hi + h_tx_lw) ;
440 long int fMHT_y = (h_ty_hi + h_ty_lw) ;
441 long int fMST_x = (e_tx_hi + e_tx_lw) ;
442 long int fMST_y = (e_ty_hi + e_ty_lw) ;
443
444 long int RHO_MULTIPLIED_BY_SUM_OF_COS_HARD_RESULT_hi = ( puc_jwj * (RHO_SUM_OF_COS_h_tx_hi) ) >> 4 ; // could be >> 2, or >> 14
445 long int RHO_MULTIPLIED_BY_SUM_OF_SIN_HARD_RESULT_hi = ( puc_jwj * (RHO_SUM_OF_SIN_h_ty_hi) ) >> 4 ;
446 long int RHO_MULTIPLIED_BY_SUM_OF_COS_SOFT_RESULT_hi = ( puc_jwj * (RHO_SUM_OF_COS_e_tx_hi) ) >> 4 ;
447 long int RHO_MULTIPLIED_BY_SUM_OF_SIN_SOFT_RESULT_hi = ( puc_jwj * (RHO_SUM_OF_SIN_e_ty_hi) ) >> 4 ;
448
449
450 long int RHO_MULTIPLIED_BY_SUM_OF_COS_HARD_RESULT_lw = ( puc_jwj * (RHO_SUM_OF_COS_h_tx_lw) ) >> 4 ;
451 long int RHO_MULTIPLIED_BY_SUM_OF_SIN_HARD_RESULT_lw = ( puc_jwj * (RHO_SUM_OF_SIN_h_ty_lw) ) >> 4 ;
452 long int RHO_MULTIPLIED_BY_SUM_OF_COS_SOFT_RESULT_lw = ( puc_jwj * (RHO_SUM_OF_COS_e_tx_lw) ) >> 4 ;
453 long int RHO_MULTIPLIED_BY_SUM_OF_SIN_SOFT_RESULT_lw = ( puc_jwj * (RHO_SUM_OF_SIN_e_ty_lw) ) >> 4 ;
454
455 long int RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_h_tx_hi = (h_tx_hi - RHO_MULTIPLIED_BY_SUM_OF_COS_HARD_RESULT_hi) ;
456 long int RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_h_ty_hi = (h_ty_hi - RHO_MULTIPLIED_BY_SUM_OF_SIN_HARD_RESULT_hi) ;
457 long int RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_e_tx_hi = (e_tx_hi - RHO_MULTIPLIED_BY_SUM_OF_COS_SOFT_RESULT_hi) ;
458 long int RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_e_ty_hi = (e_ty_hi - RHO_MULTIPLIED_BY_SUM_OF_SIN_SOFT_RESULT_hi) ;
459
460 long int RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_h_tx_lw = (h_tx_lw - RHO_MULTIPLIED_BY_SUM_OF_COS_HARD_RESULT_lw) ;
461 long int RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_h_ty_lw = (h_ty_lw - RHO_MULTIPLIED_BY_SUM_OF_SIN_HARD_RESULT_lw) ;
462 long int RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_e_tx_lw = (e_tx_lw - RHO_MULTIPLIED_BY_SUM_OF_COS_SOFT_RESULT_lw) ;
463 long int RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_e_ty_lw = (e_ty_lw - RHO_MULTIPLIED_BY_SUM_OF_SIN_SOFT_RESULT_lw) ;
464
465 MHT_x = (RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_h_tx_hi + RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_h_tx_lw) >> 3;
466 MHT_y = (RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_h_ty_hi + RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_h_ty_lw) >> 3;
467 MST_x = (RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_e_tx_hi + RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_e_tx_lw) >> 3;
468 MST_y = (RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_e_ty_hi + RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_e_ty_lw) >> 3;
469
470 fMHT_x = (RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_h_tx_hi + RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_h_tx_lw) ;
471 fMHT_y = (RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_h_ty_hi + RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_h_ty_lw) ;
472 fMST_x = (RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_e_tx_hi + RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_e_tx_lw) ;
473 fMST_y = (RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_e_ty_hi + RHO_SUBTRACTED_BEFORE_FINAL_MULTIPLY_e_ty_lw) ;
474
475 long int fMET_x = ( aFPGA * (fMHT_x) + bFPGA * (fMST_x) ) >> 13 ;
476 long int fMET_y = ( aFPGA * (fMHT_y) + bFPGA * (fMST_y) ) >> 13 ;
477
478 MET_x = fMET_x;
479 MET_y = fMET_y;
480
481}
482
483
484void gFEXJwoJAlgo::metFPGA(int FPGAnum, const gTowersType& twrs,
485 const gTowersType & gBlkSum, int gBlockthreshold,
486 int aFPGA, int bFPGA,
487 int & MHT_x, int & MHT_y,
488 int & MST_x, int & MST_y,
489 int & MET_x, int & MET_y) const {
490
491 gBlockthreshold = gBlockthreshold * 200 / 800; //gBlockthreshold is provided in counts with a resolution of 200 MeV, but here needs to be applied with a resolution of 800 GeV
492 // in the RTL code these are 19+ 5 = 24 bits
493 int64_t h_tx_hi = 0;
494 int64_t h_ty_hi = 0;
495 int64_t h_tx_lw = 0;
496 int64_t h_ty_lw = 0;
497
498 int64_t e_tx_hi = 0;
499 int64_t e_ty_hi = 0;
500 int64_t e_tx_lw = 0;
501 int64_t e_ty_lw = 0;
502
503 for( int irow = 0; irow < FEXAlgoSpaceDefs::ABCrows; irow++ ){
504 for(int jcolumn = 6; jcolumn<12; jcolumn++){
505 if( FPGAnum == 2){
506 int frow = 2*(irow/2) + 1;
507
508 if(gBlkSum[irow][jcolumn] > gBlockthreshold){
509 h_tx_hi += (twrs[irow][jcolumn])*(cosLUT(frow, 5));
510 h_ty_hi += (twrs[irow][jcolumn])*(sinLUT(frow, 5));
511 } else {
512 e_tx_hi += (twrs[irow][jcolumn])*(cosLUT(frow, 5));
513 e_ty_hi += (twrs[irow][jcolumn])*(sinLUT(frow, 5));
514 }
515
516 } else {
517
518 if(gBlkSum[irow][jcolumn] > gBlockthreshold){
519 h_tx_hi += (twrs[irow][jcolumn])*(cosLUT(irow, 5));
520 h_ty_hi += (twrs[irow][jcolumn])*(sinLUT(irow, 5));
521 } else {
522 e_tx_hi += (twrs[irow][jcolumn])*(cosLUT(irow, 5));
523 e_ty_hi += (twrs[irow][jcolumn])*(sinLUT(irow, 5));
524 }
525 }
526 }
527
528 for(int jcolumn = 0; jcolumn<6; jcolumn++){
529 if( FPGAnum == 2){
530 int frow = 2*(irow/2) + 1;
531
532 if(gBlkSum[irow][jcolumn] > gBlockthreshold){
533 h_tx_lw += (twrs[irow][jcolumn])*(cosLUT(frow, 5));
534 h_ty_lw += (twrs[irow][jcolumn])*(sinLUT(frow, 5));
535 } else{
536 e_tx_lw += (twrs[irow][jcolumn])*(cosLUT(frow, 5));
537 e_ty_lw += (twrs[irow][jcolumn])*(sinLUT(frow, 5));
538 }
539 } else {
540
541 if(gBlkSum[irow][jcolumn] > gBlockthreshold){
542 h_tx_lw += (twrs[irow][jcolumn])*(cosLUT(irow, 5));
543 h_ty_lw += (twrs[irow][jcolumn])*(sinLUT(irow, 5));
544 } else {
545 e_tx_lw += (twrs[irow][jcolumn])*(cosLUT(irow, 5));
546 e_ty_lw += (twrs[irow][jcolumn])*(sinLUT(irow, 5));
547 }
548 }
549 }
550 }
551
552 // According to https://gitlab.cern.ch/atlas-l1calo/gfex/firmware/-/issues/406#note_5662344
553 // there is no division -- so LSB is indeed 25 MeV
554
555 //but later changed to 200 MeV so divide by 8
556
557 long int fMHT_x = (h_tx_hi + h_tx_lw) ;
558 long int fMHT_y = (h_ty_hi + h_ty_lw) ;
559 long int fMST_x = (e_tx_hi + e_tx_lw) ;
560 long int fMST_y = (e_ty_hi + e_ty_lw) ;
561
562 MHT_x = (h_tx_hi + h_tx_lw) >> 3;
563 MHT_y = (h_ty_hi + h_ty_lw) >> 3;
564 MST_x = (e_tx_hi + e_tx_lw) >> 3;
565 MST_y = (e_ty_hi + e_ty_lw) >> 3;
566
567 // a and b coffecients are 10 bits
568 // multiplication has an addtional 2^10
569 // constant JWJ_OW : integer := 35;--Out width
570 // values are 35 bits long and top 16 bits are taken -- so divide by 2^19
571 // 2^10/2^19 = 1/2^9 = 1/512
572
573 long int fMET_x = ( aFPGA * (fMHT_x) + bFPGA * (fMST_x) ) >> 13 ;
574 long int fMET_y = ( aFPGA * (fMHT_y) + bFPGA * (fMST_y) ) >> 13 ;
575
576 MET_x = fMET_x;
577 MET_y = fMET_y;
578
579}
580
581void gFEXJwoJAlgo::etFPGA(int FPGAnum, const gTowersType& twrs, gTowersType &gBlkSum,
582 int gBlockthreshold, int A, int B, int &eth, int &ets, int &etw) const {
583
584 gBlockthreshold = gBlockthreshold * 200 / 800; //gBlockthreshold is provided in counts with a resolution of 200 MeV, but here needs to be applied with a resolution of 800 GeV
585
586 int64_t ethard_hi = 0;
587 int64_t etsoft_hi = 0;
588 int64_t ethard_lo = 0;
589 int64_t etsoft_lo = 0;
590
591 int64_t ethard = 0.0;
592 int64_t etsoft = 0.0;
593
594 int multiplicitiveFactor = 0;
595
596 if(FPGAnum < 2 ) {
597 multiplicitiveFactor = cosLUT(0, 5);
598 } else{
599 multiplicitiveFactor = cosLUT(1, 5);
600 }
601
602// firmware treats upper and lower columns differnetly
603
604 for( int irow = 0; irow < FEXAlgoSpaceDefs::ABCrows; irow++ ){
605 for(int jcolumn = 0; jcolumn<6; jcolumn++){
606 if(gBlkSum[irow][jcolumn] > gBlockthreshold){
607 ethard_lo = ethard_lo + twrs[irow][jcolumn]*multiplicitiveFactor;
608 } else {
609 etsoft_lo = etsoft_lo + twrs[irow][jcolumn]*multiplicitiveFactor;
610 }
611 }
612 }
613
614 for( int irow = 0; irow < FEXAlgoSpaceDefs::ABCrows; irow++ ){
615 for(int jcolumn = 6; jcolumn<12; jcolumn++){
616 if(gBlkSum[irow][jcolumn] > gBlockthreshold){
617 ethard_hi = ethard_hi + twrs[irow][jcolumn]*multiplicitiveFactor;
618 } else {
619 etsoft_hi = etsoft_hi + twrs[irow][jcolumn]*multiplicitiveFactor;
620 }
621 }
622 }
623
624 ethard = ethard_hi + ethard_lo;
625 etsoft = etsoft_hi + etsoft_lo;
626
627
628 int64_t etsum_hi = ethard_hi*A + etsoft_hi*B ;
629 if ( etsum_hi < 0 ) etsum_hi = 0;
630
631 int64_t etsum_lo = ethard_lo*A + etsoft_lo*B ;
632 if ( etsum_lo < 0 ) etsum_lo = 0;
633
634 int64_t etsum = etsum_hi + etsum_lo;
635
636
637 // convert 200 MeV LSB here
638 eth = ethard>>3;
639 ets = etsoft>>3;
640 etw = (etsum >>13 ) ;
641
642 if( etw < 0 ) etw = 0;
643 // max value is 15 bits with 800 MeV LSB -- so 17 bits here
644 if( etw > 0X001FFFF ) etw = 0X001FFFF ;
645
646
647 if(msgLvl(MSG::DEBUG)) {
648 std::cout << "DMS FPGA gTEJWOJ " << std::hex << FPGAnum << "et sum hard " << eth << "etsum soft" << ets << " A " << A << " B " << B << " weighted term " << etw << std::endl << std::dec;
649 }
650}
651
652void gFEXJwoJAlgo::etFastFPGA(int FPGAnum, const gTowersType& twrs, gTowersType &gBlkSum,
653 int gBlockthreshold, int A, int B, int &eth, int &ets, int &etw) const {
654
655 gBlockthreshold = gBlockthreshold * 200 / 800; //gBlockthreshold is provided in counts with a resolution of 200 MeV, but here needs to be applied with a resolution of 800 GeV
656 int64_t ethard_hi = 0;
657 int64_t etsoft_hi = 0;
658 int64_t ethard_lo = 0;
659 int64_t etsoft_lo = 0;
660
661 int64_t ethard = 0.0;
662 int64_t etsoft = 0.0;
663
664 // firmware treats upper and lower columns differently
665 for( int irow = 0; irow < FEXAlgoSpaceDefs::ABCrows; irow++ ){
666 for(int jcolumn = 0; jcolumn<6; jcolumn++){
667 if(gBlkSum[irow][jcolumn] > gBlockthreshold){
668 ethard_lo = ethard_lo + twrs[irow][jcolumn];
669 }
670 else {
671 etsoft_lo = etsoft_lo + twrs[irow][jcolumn];
672 }
673 }
674 }
675
676 for( int irow = 0; irow < FEXAlgoSpaceDefs::ABCrows; irow++ ){
677 for(int jcolumn = 6; jcolumn<12; jcolumn++){
678 if(gBlkSum[irow][jcolumn] > gBlockthreshold){
679 ethard_hi = ethard_hi + twrs[irow][jcolumn];
680 }
681 else {
682 etsoft_hi = etsoft_hi + twrs[irow][jcolumn];
683 }
684 }
685 }
686
687 ethard = ethard_hi + ethard_lo;
688 etsoft = etsoft_hi + etsoft_lo;
689
690 // convert 200 MeV LSB here
691 eth = ethard;
692 ets = etsoft; // Keep for reference -- not used in etFast
693 etw = ethard; // For EtFast the weighted term is just the hard term
694
695 // 16 bits signed, set max and min
696 if( etw < -32768 ) etw = -32768;
697 if( etw > 32767 ) etw = 32767;
698
699 if(msgLvl(MSG::DEBUG)) {
700 std::cout << "DMS FPGA gTEJWOJ " << std::hex << FPGAnum << "et sum hard " << eth << "etsum soft" << ets << " A " << A << " B " << B << " weighted term " << etw << std::endl << std::dec;
701 }
702}
703
704void gFEXJwoJAlgo::metTotal(int A_MET_x, int A_MET_y,
705 int B_MET_x, int B_MET_y,
706 int C_MET_x, int C_MET_y,
707 int & MET_x, int & MET_y) const {
708
709
710 MET_x = A_MET_x + B_MET_x + C_MET_x;
711 MET_y = A_MET_y + B_MET_y+ C_MET_y;
712
713 // Truncation of the result, as the individual quantities are 16 bits, while the TOB field is 12 bits
714 // MET_x = MET_x >> 4;
715 // MET_y = MET_y >> 4;
716
717 if (MET_x < -0x000800) MET_x = -0x000800; //-2048
718 if (MET_y < -0x000800) MET_y = -0x000800; //-2048
719
720 if (MET_x > 0x0007FF) MET_x = 0x0007FF; //2047
721 if (MET_y > 0x0007FF) MET_y = 0x0007FF; //2047
722
723}
724
725void gFEXJwoJAlgo::etTotal(int A_ET,
726 int B_ET,
727 int C_ET,
728 int & ET ) const {
729
730 // leave at 200 MeV scale
731 if (A_ET < 0 ) A_ET = 0;
732 if (B_ET < 0 ) B_ET = 0;
733 if (C_ET < 0 ) C_ET = 0;
734
735 ET = (A_ET + B_ET + C_ET);
736
737 // main value of ET is always positive
738 if( ET > 0x0000FFF) ET = 0x0000FFF;
739
740}
741
742//----------------------------------------------------------------------------------
743// bitwise simulation of sine LUT in firmware
744//----------------------------------------------------------------------------------
745float gFEXJwoJAlgo::sinLUT(unsigned int phiIDX, unsigned int aw) const
746{
747 float c = static_cast<float>(phiIDX)/std::pow(2,aw);
748 float rad = (2*M_PI) *c;
749 float rsin = std::sin(rad);
750 int isin = std::round(rsin*(std::pow(2,aw) - 1));
751 return isin;
752
753}
754
755//----------------------------------------------------------------------------------
756// bitwise simulation cosine LUT in firmware
757//----------------------------------------------------------------------------------
758float gFEXJwoJAlgo::cosLUT(unsigned int phiIDX, unsigned int aw) const
759{
760 float c = static_cast<float>(phiIDX)/std::pow(2,aw);
761 float rad = (2*M_PI) *c;
762 float rcos = std::cos(rad);
763 int icos = std::round(rcos*(std::pow(2,aw) - 1));
764 return icos;
765}
766
767} // namespace LVL1
M_PI
#define M_PI
Definition ActiveFraction.h:14
ATH_CHECK
#define ATH_CHECK
Evaluate an expression and check for errors.
Definition AthCheckMacros.h:40
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x)
Definition AthMsgStreamMacros.h:33
AthAlgTool::AthAlgTool
AthAlgTool(const std::string &type, const std::string &name, const IInterface *parent)
Constructor with parameters:
Definition AthAlgTool.cxx:16
AthCommonMsg< AlgTool >::msgLvl
bool msgLvl(const MSG::Level lvl) const
Definition AthCommonMsg.h:30
LVL1::FEXAlgoSpaceDefs::ABCrows
static constexpr int ABCrows
Definition FEXAlgoSpaceDefs.h:157
LVL1::FEXAlgoSpaceDefs::gBlockMax
static constexpr int gBlockMax
Definition FEXAlgoSpaceDefs.h:161
LVL1::FEXAlgoSpaceDefs::ABcolumns
static constexpr int ABcolumns
Definition FEXAlgoSpaceDefs.h:158
LVL1::gFEXJwoJAlgo::etFastFPGA
void etFastFPGA(int FPGAnum, const gTowersType &twrs, gTowersType &gBlkSum, int gBlockthreshold, int A, int B, int &eth, int &ets, int &etw) const
Definition gFEXJwoJAlgo.cxx:652
LVL1::gFEXJwoJAlgo::metFPGA
void metFPGA(int FPGAnum, const gTowersType &twrs, const gTowersType &gBlkSum, int gBlockthreshold, int aFPGA, int bFPGA, int &MHT_x, int &MHT_y, int &MST_x, int &MST_y, int &MET_x, int &MET_y) const
Definition gFEXJwoJAlgo.cxx:484
LVL1::gFEXJwoJAlgo::m_DBToolKey
SG::ReadCondHandleKey< gFEXDBCondData > m_DBToolKey
Definition gFEXJwoJAlgo.h:48
LVL1::gFEXJwoJAlgo::etFPGA
void etFPGA(int FPGAnum, const gTowersType &twrs, gTowersType &gBlkSum, int gBlockthreshold, int A, int B, int &eth, int &ets, int &etw) const
Definition gFEXJwoJAlgo.cxx:581
LVL1::gFEXJwoJAlgo::initialize
virtual StatusCode initialize() override
standard Athena-Algorithm method
Definition gFEXJwoJAlgo.cxx:29
LVL1::gFEXJwoJAlgo::cosLUT
float cosLUT(unsigned int phiIDX, unsigned int aw) const
Definition gFEXJwoJAlgo.cxx:758
LVL1::gFEXJwoJAlgo::metFPGA_rho
void metFPGA_rho(int FPGAnum, const gTowersType &twrs, int puc_jwj, const gTowersType &gBlkSum, int gBlockthreshold, int aFPGA, int bFPGA, int &MHT_x, int &MHT_y, int &MST_x, int &MST_y, int &MET_x, int &MET_y) const
Definition gFEXJwoJAlgo.cxx:341
LVL1::gFEXJwoJAlgo::m_l1MenuKey
SG::ReadHandleKey< TrigConf::L1Menu > m_l1MenuKey
Definition gFEXJwoJAlgo.h:50
LVL1::gFEXJwoJAlgo::metTotal
void metTotal(int A_MET_x, int A_MET_y, int B_MET_x, int B_MET_y, int C_MET_x, int C_MET_y, int &MET_x, int &MET_y) const
Definition gFEXJwoJAlgo.cxx:704
LVL1::gFEXJwoJAlgo::sinLUT
float sinLUT(unsigned int phiIDX, unsigned int aw) const
Definition gFEXJwoJAlgo.cxx:745
LVL1::gFEXJwoJAlgo::gBlockAB
void gBlockAB(const gTowersType &twrs, gTowersType &gBlkSum, gTowersType &hasSeed, int seedThreshold) const
Definition gFEXJwoJAlgo.cxx:296
LVL1::gFEXJwoJAlgo::jwojAlgo
virtual std::vector< std::unique_ptr< gFEXJwoJTOB > > jwojAlgo(const EventContext &ctx, const gTowersType &Atwr, int pucA_JWJ, const gTowersType &Btwr, int pucB_JWJ, const gTowersType &Ctwr, int pucC_JWJ, std::array< int32_t, 4 > &outTOB) const override
Definition gFEXJwoJAlgo.cxx:38
LVL1::gFEXJwoJAlgo::gFEXJwoJAlgo
gFEXJwoJAlgo(const std::string &type, const std::string &name, const IInterface *parent)
Constructors.
Definition gFEXJwoJAlgo.cxx:22
LVL1::gFEXJwoJAlgo::etTotal
void etTotal(int A_ET, int B_ET, int C_ET, int &ET) const
Definition gFEXJwoJAlgo.cxx:725
SG::ReadHandle::isValid
virtual bool isValid() override final
Can the handle be successfully dereferenced?
LVL1
eFexTowerBuilder creates xAOD::eFexTowerContainer from supercells (LATOME) and triggerTowers (TREX) i...
Definition IZdcDataAccess.h:13
LVL1::gTowersType
std::array< std::array< int, 12 >, 32 > gTowersType
Definition IgFEXFPGA.h:25
MET
Definition MET.py:1
A
hold the test vectors and ease the comparison
type
Generated on for ATLAS Offline Software by doxygen 1.17.0