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