ATLAS Offline Software
Loading...
Searching...
No Matches
LArDigitRetriever.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#include "LArDigitRetriever.h"
6
7#include "CLHEP/Units/SystemOfUnits.h"
8
10
11#include "CaloDetDescr/CaloDetDescrElement.h"
19#include "GaudiKernel/ThreadLocalContext.h"
20#include <memory>
21#include <array>
22
23using CLHEP::GeV;
24
25namespace JiveXML {
26
33 LArDigitRetriever::LArDigitRetriever(const std::string& type,const std::string& name,const IInterface* parent):
34 AthAlgTool(type,name,parent),
35 m_calocell_id(nullptr)
36 {
37 //Only declare the interface
38 declareInterface<IDataRetriever>(this);
39
40 m_doDigit = false;
41 m_inputdpd = false;
42
43 declareProperty("CellThreshold", m_cellThreshold = 50.);
44 declareProperty("RetrieveLAr" , m_lar = true);
45 declareProperty("RetrieveHEC" , m_hec = true);
46 declareProperty("RetrieveFCAL" , m_fcal = true);
47 declareProperty("DoLArDigit", m_doLArDigit = false);
48 declareProperty("DoHECDigit", m_doHECDigit = false);
49 declareProperty("DoFCalDigit", m_doFCalDigit = false);
50 declareProperty("CellEnergyPrec", m_cellEnergyPrec = 3);
51 declareProperty("CellTimePrec", m_cellTimePrec = 3);
52
53// Check the cell conditions. Not present in MC data, so false by default. Switch to 'true'
54// for real (commissioning) cells.
55 declareProperty("CellConditionCut", m_cellConditionCut = false);
56 declareProperty("LArChannelsToIgnoreM5", m_LArChannelsToIgnoreM5);
57 declareProperty("DoMaskLArChannelsM5", m_doMaskLArChannelsM5 = false);
58 }
59
63
65
66 ATH_MSG_DEBUG( "Initialising Tool" );
67 ATH_CHECK( detStore()->retrieve (m_calocell_id, "CaloCell_ID") );
68 ATH_CHECK( m_adc2mevKey.initialize() );
69
70 ATH_CHECK( m_sgKey.initialize() );
71 ATH_CHECK( m_sgKeyLArDigit_raw.initialize() );
72 ATH_CHECK( m_sgKeyLArDigit_esd.initialize() );
73 ATH_CHECK( m_cablingKey.initialize() );
74
75 return StatusCode::SUCCESS;
76 }
77
81 StatusCode LArDigitRetriever::retrieve(ToolHandle<IFormatTool> &FormatTool) {
82
83 ATH_MSG_DEBUG( "in retrieve()" );
84
86 m_doDigit = true;
87 }
88
90 if (!cellContainer.isValid()){
91 ATH_MSG_WARNING( "Could not retrieve Calorimeter Cells " );
92 }
93 else{
94
95 if(m_fcal){
96 DataMap data = getLArDigitData(&(*cellContainer),"FCAL",CaloCell_ID::LARFCAL);
97 ATH_CHECK( FormatTool->AddToEvent("FCAL", m_sgKey.key(), &data) );
98 ATH_MSG_DEBUG( "FCAL retrieved" );
99 }
100
101 if(m_lar){
102 DataMap data = getLArDigitData(&(*cellContainer),"LAr",CaloCell_ID::LAREM);
103 ATH_CHECK( FormatTool->AddToEvent("LAr", m_sgKey.key(), &data) );
104 ATH_MSG_DEBUG( "LAr retrieved" );
105 }
106
107 if(m_hec){
108 DataMap data = getLArDigitData(&(*cellContainer),"HEC",CaloCell_ID::LARHEC);
109 ATH_CHECK( FormatTool->AddToEvent("HEC", m_sgKey.key(), &data) );
110 ATH_MSG_DEBUG( "HEC retrieved" );
111 }
112 }
113
114 //Tile cells retrieved okay
115 return StatusCode::SUCCESS;
116 }
117
118
124 const std::string& datatype,
125 CaloCell_ID::SUBCALO calotype)
126 {
127
128 //be verbose
129 ATH_MSG_DEBUG( "getLArDigitData()" );
130 const EventContext& ctx = Gaudi::Hive::currentContext();
131
132 char rndStr[30];
134
135 DataVect phi; phi.reserve(cellContainer->size());
136 DataVect eta; eta.reserve(cellContainer->size());
137 DataVect x; x.reserve(cellContainer->size());
138 DataVect y; y.reserve(cellContainer->size());
139 DataVect dx; dx.reserve(cellContainer->size());
140 DataVect dy; dy.reserve(cellContainer->size());
141
142 DataVect energy; energy.reserve(cellContainer->size());
143 DataVect idVec; idVec.reserve(cellContainer->size());
144 DataVect channel; channel.reserve(cellContainer->size());
145 DataVect feedThrough; feedThrough.reserve(cellContainer->size());
146 DataVect slotVec; slotVec.reserve(cellContainer->size());
147
148 DataVect cellTimeVec; cellTimeVec.reserve(cellContainer->size());
149 DataVect cellGain; cellGain.reserve(cellContainer->size());
150 DataVect cellPedestal; cellPedestal.reserve(cellContainer->size());
151 DataVect adc2Mev; adc2Mev.reserve(cellContainer->size());
152
153// m_sub; m_sub.reserve(cellContainer->size());
154 m_sub.clear(); // need to clear before each event, otherwise appended
155
156 DataVect LArSampleIndexVec; LArSampleIndexVec.reserve(cellContainer->size());
157
158 std::string LArSampleIndexStr="adcCounts multiple=\"0\"";
159
162 const LArDigitContainer* LArDigitCnt = nullptr;
163 if (LArDigitCnt_raw_handle.isValid()) {
164 LArDigitCnt = &(*LArDigitCnt_raw_handle);
165 }
166 else if (LArDigitCnt_esd_handle.isValid()) {
167 LArDigitCnt = &(*LArDigitCnt_esd_handle);
168 m_inputdpd = true;
169 }
170 else {
171 ATH_MSG_WARNING( "Could not retrieve LArDigits" );
172 }
173
174 const ILArPedestal* larPedestal = nullptr;
175 if ( detStore()->retrieve(larPedestal).isFailure()){
176 ATH_MSG_ERROR( "in getLArDigitData(), Could not retrieve LAr Pedestal" );
177 }
178
179 const LArOnlineID* onlineId = nullptr;
180 if ( detStore()->retrieve(onlineId, "LArOnlineID").isFailure()) {
181 ATH_MSG_ERROR( "in getLArDigitData(),Could not get LArOnlineID!" );
182 }
183
185
186//------------------------------------------------------
187//Loop over the digits and find Cell (LAr,HEC, FCAL)
188//------------------------------------------------------
189
191 const LArOnOffIdMapping* cabling{*cablingHdl};
192 if(!cabling) {
193 ATH_MSG_ERROR( "Do not have cabling mapping from key " << m_cablingKey.key() );
194 return DataMap;
195 }
196
197
198 if (LArDigitCnt && m_doDigit==true) {
199
200 auto pCellIndex = std::make_unique<std::array<int, 200000>>();
201 pCellIndex->fill(0);
202 int nLArSamples = 0;
203 HWIdentifier LArHardwareId;
204 Identifier LArId;
205
206 double energyGeV,cellTime;
207
208 for (const LArDigit* digit : *LArDigitCnt) {
209
210 LArHardwareId = digit->hardwareID();
211 if (!cabling->isOnlineConnected(LArHardwareId))continue;
212
213 LArId = cabling->cnvToIdentifier(LArHardwareId); //converter
214 const IdentifierHash cellhash=m_calocell_id->calo_cell_hash(LArId); //fast method to find cell
215
216 int Index = cellContainer->findIndex(cellhash); //find Cell Index
217 if (Index >= 0)
218 (*pCellIndex)[Index] = Index;
219
220 nLArSamples = digit->nsamples();
221 std::vector<short> LArSamples = digit->samples();
222 int largain = digit->gain();
223 int FT = onlineId->feedthrough(LArHardwareId);
224 int slot = onlineId->slot(LArHardwareId);
225 int larchan = onlineId->channel(LArHardwareId);
226 float pedestal=larPedestal->pedestal(LArHardwareId,largain);
227 float pedvalue=0;
228 if (pedestal >= (1.0+LArElecCalib::ERRORCODE)) pedvalue = pedestal;
229 else pedvalue = LArSamples[0];
230
231 LArVectorProxy polynom_adc2mev = adc2mev->ADC2MEV(LArId,largain);
232
233 if ( Index >= 0 ){ // can be -1
234 if ( (*cellContainer)[Index]->energy() >m_cellThreshold) {
235
236 if (((((*cellContainer)[Index]->provenance())&0xFF)!=0xA5)&&m_cellConditionCut) continue; // check full conition for LAr Cells
237
238 //ignore LAr, HEC, and/or FCAL cells that are to be masked
240 bool maskChannel = false;
241 for (size_t i = 0; i < m_LArChannelsToIgnoreM5.size(); i++){
242 if ( (*cellContainer)[Index]->ID().get_compact() == m_LArChannelsToIgnoreM5[i]){
243 maskChannel = true;
244 break; // exit loop over bad channels
245 }
246 }
247 if (maskChannel) continue; // continue loop over all channels
248 }
249
250
251 if (datatype == "LAr" && m_calocell_id->is_em(LArId)) {
252
253 calcEMLayerSub(LArId);
254
255 energyGeV = (*cellContainer)[Index]->energy() * (1./GeV);
256 energy.push_back(DataType( gcvt( energyGeV, m_cellEnergyPrec, rndStr) ));
257 idVec.push_back(DataType((*cellContainer)[Index]->ID().get_compact() ));
258 phi.push_back(DataType((*cellContainer)[Index]->phi()));
259 eta.push_back(DataType((*cellContainer)[Index]->eta()));
260
261 cellTime = (*cellContainer)[Index]->time();
262 cellTimeVec.push_back(DataType( gcvt( cellTime, m_cellTimePrec, rndStr) ) );
263 cellPedestal.push_back(DataType(pedvalue));
264 cellGain.push_back(DataType(largain));
265 channel.push_back(DataType(larchan));
266 feedThrough.push_back(DataType(FT));
267 slotVec.push_back(DataType(slot));
268 if (polynom_adc2mev.size()==0){
269 adc2Mev.push_back(DataType(-1));
270 }else{
271 adc2Mev.push_back(DataType(polynom_adc2mev[1]));
272 }
273 if ( m_doLArDigit ) {
274 LArSampleIndexStr="adcCounts multiple=\""+DataType(nLArSamples).toString()+"\"";
275 for(int i=0; i<nLArSamples; i++) LArSampleIndexVec.push_back(DataType(LArSamples[i]));
276 }
277 } //match datatype LAr
278
279 else if(datatype == "HEC" && m_calocell_id->is_hec(LArId)) {
280
281 calcHECLayerSub(LArId);
282 energyGeV = (*cellContainer)[Index]->energy() * (1./GeV);
283 energy.push_back(DataType( gcvt( energyGeV, m_cellEnergyPrec, rndStr) ));
284 idVec.push_back(DataType((*cellContainer)[Index]->ID().get_compact() ));
285 phi.push_back(DataType((*cellContainer)[Index]->phi()));
286 eta.push_back(DataType((*cellContainer)[Index]->eta()));
287
288
289 cellTime = (*cellContainer)[Index]->time();
290 cellTimeVec.push_back(DataType( gcvt( cellTime, m_cellTimePrec, rndStr) ) );
291 cellPedestal.push_back(DataType(pedvalue));
292 cellGain.push_back(DataType(largain));
293 channel.push_back(DataType(larchan));
294 feedThrough.push_back(DataType(FT));
295 slotVec.push_back(DataType(slot));
296 if (polynom_adc2mev.size()==0)
297 adc2Mev.push_back(DataType(-1));
298 else
299 adc2Mev.push_back(DataType(polynom_adc2mev[1]));
300 if (m_doHECDigit){
301 LArSampleIndexStr="adcCounts multiple=\""+DataType(nLArSamples).toString()+"\"";
302 for(int i=0; i<nLArSamples; i++) LArSampleIndexVec.push_back(DataType(LArSamples[i]));
303 }
304 }//match datatype HEC
305
306 else if(datatype == "FCAL" && m_calocell_id->is_fcal(LArId)) {
307
308 energyGeV = (*cellContainer)[Index]->energy() * (1./GeV);
309 energy.push_back(DataType( gcvt( energyGeV, m_cellEnergyPrec, rndStr) ));
310 idVec.push_back(DataType((*cellContainer)[Index]->ID().get_compact() ));
311 x.push_back(DataType((*cellContainer)[Index]->x() *0.1));
312 y.push_back(DataType((*cellContainer)[Index]->y() *0.1));
313
314 const CaloDetDescrElement* elt = (*cellContainer)[Index]->caloDDE();
315
316 dx.push_back(DataType(elt->dx()*0.1));
317 dy.push_back(DataType(elt->dy()*0.1));
318
319 if (m_calocell_id->pos_neg(LArId)==2) m_sub.push_back(DataType(1));
320 else m_sub.push_back(DataType(0));
321
322
323 cellTime = (*cellContainer)[Index]->time();
324 cellTimeVec.push_back(DataType( gcvt( cellTime, m_cellTimePrec, rndStr) ) );
325 cellPedestal.push_back(DataType(pedvalue));
326 cellGain.push_back(DataType(largain));
327 channel.push_back(DataType(larchan));
328 feedThrough.push_back(DataType(FT));
329 slotVec.push_back(DataType(slot));
330 if (polynom_adc2mev.size()==0)
331 adc2Mev.push_back(DataType(-1));
332 else
333 adc2Mev.push_back(DataType(polynom_adc2mev[1]));
334 if (m_doFCalDigit){
335 LArSampleIndexStr="adcCounts multiple=\""+DataType(nLArSamples).toString()+"\"";
336 for(int i=0; i<nLArSamples; i++) LArSampleIndexVec.push_back(DataType(LArSamples[i]));
337 }
338 }//match datatype FCAL
339 } //if cellThreshold
340 } // if Index >0
341 } //digit
342
343
344// If the digits are retrieved from DPD, retrieve the other cells which do not have the digits available
345
346 if (m_inputdpd) {
347
348 CaloCellContainer::const_iterator it1 = cellContainer->beginConstCalo(calotype);
349 CaloCellContainer::const_iterator it2 = cellContainer->endConstCalo(calotype);
350
351
352 for (;it1!=it2;++it1) {
353 //---------------------------------------------
354
355 if ( (*it1)->energy() < m_cellThreshold ) continue;
356
357 Identifier cellid = (*it1)->ID();
358
359 const IdentifierHash cellhash=m_calocell_id->calo_cell_hash(cellid); //fast method to find cell
360 int Index = cellContainer->findIndex(cellhash); //find Cell Index
361 if (Index >= 0 && (*pCellIndex)[Index] == Index)
362 continue; //test whether this cell was already retrieved
363
364 HWIdentifier LArhwid = cabling->createSignalChannelIDFromHash((*it1)->caloDDE()->calo_hash());
365
366 if (m_calocell_id->is_tile(cellid) ) continue;
367 if (((((*it1)->provenance())&0xFF)!=0xA5)&&m_cellConditionCut) continue; // check full condition for LAr cells
368
369 //ignore LAr, HEC, and/or FCAL cells that are to be masked
371 bool maskChannel = false;
372 for (size_t i = 0; i < m_LArChannelsToIgnoreM5.size(); i++){
373 if (cellid == m_LArChannelsToIgnoreM5[i]){
374 maskChannel = true;
375 break; // exit loop over bad channels
376 }
377 }
378 if (maskChannel) continue; // continue loop over all channels
379 }
380
381 energyGeV = (*it1)->energy()*(1./GeV);
382 energy.push_back(DataType( gcvt( energyGeV, m_cellEnergyPrec, rndStr) ));
383 idVec.push_back(DataType((*it1)->ID().get_compact() ));
384 channel.push_back(DataType(onlineId->channel(LArhwid)));
385 feedThrough.push_back(DataType(onlineId->feedthrough(LArhwid)));
386 slotVec.push_back(DataType(onlineId->slot(LArhwid)));
387
388 cellTime = (*it1)->time();
389 cellTimeVec.push_back(DataType( gcvt( cellTime, m_cellTimePrec, rndStr) ) );
390 cellGain.push_back(DataType((*it1)->gain()));
391 int largain = (*it1)->gain();
392 float pedvalue=0;
393 float pedestal=larPedestal->pedestal(LArhwid,largain);
394 if (pedestal >= (1.0+LArElecCalib::ERRORCODE)) pedvalue = pedestal;
395 else pedvalue = 0;
396 cellPedestal.push_back(DataType(pedvalue));
397
398 LArVectorProxy polynom_adc2mev = adc2mev->ADC2MEV(cellid,largain);
399 if (polynom_adc2mev.size()==0){
400 adc2Mev.push_back(DataType(-1));
401 }else{
402 adc2Mev.push_back(DataType(polynom_adc2mev[1]));
403 }
404
405
406 if (datatype == "LAr") {
407 calcEMLayerSub(cellid);
408 phi.push_back(DataType((*it1)->phi()));
409 eta.push_back(DataType((*it1)->eta()));
410 if ( m_doLArDigit ) {
411 LArSampleIndexStr="adcCounts multiple=\""+DataType(nLArSamples).toString()+"\"";
412 for(int i=0; i<nLArSamples; i++) LArSampleIndexVec.push_back(DataType(0)); }
413 }
414
415 else if (datatype == "HEC") {
416 calcHECLayerSub(cellid);
417 phi.push_back(DataType((*it1)->phi()));
418 eta.push_back(DataType((*it1)->eta()));
419 if (m_doHECDigit){
420 LArSampleIndexStr="adcCounts multiple=\""+DataType(nLArSamples).toString()+"\"";
421 for(int i=0; i<nLArSamples; i++) LArSampleIndexVec.push_back(DataType(0)); }
422 }
423
424 else if (datatype == "FCAL") {
425
426 x.push_back(DataType( (*it1)->x()*0.1 ));
427 y.push_back(DataType( (*it1)->y()*0.1 ));
428
429 const CaloDetDescrElement* elt = (*it1)->caloDDE();
430 dx.push_back(DataType( elt->dx()*0.1 ));
431 dy.push_back(DataType( elt->dy()*0.1 ));
432
433 if (m_calocell_id->pos_neg(cellid)==2) m_sub.push_back(DataType(1));
434 else m_sub.push_back(DataType(0));
435 if (m_doFCalDigit){
436 LArSampleIndexStr="adcCounts multiple=\""+DataType(nLArSamples).toString()+"\"";
437 for(int i=0; i<nLArSamples; i++) LArSampleIndexVec.push_back(DataType(0));
438 }
439 }
440
441 }//for(;it1!=it2;it1++)
442 }// if (m_inputdpd)
443 }
444
445//----------------above lines are trying to get the LAr,HEC and FCAL tags --------------
446
447 // write values into DataMap
448 const auto nEntries = phi.size();
449 if(!(datatype=="FCAL")){
450 DataMap["phi"] = std::move(phi);
451 DataMap["eta"] = std::move(eta);
452 } else {
453 DataMap["x"] = std::move(x);
454 DataMap["y"] = std::move(y);
455 DataMap["dx"] = std::move(dx);
456 DataMap["dy"] = std::move(dy);
457 }
458
459 DataMap["energy"] = std::move(energy);
460 DataMap["id"] = std::move(idVec);
461 DataMap["channel"] = std::move(channel);
462 DataMap["feedThrough"] = std::move(feedThrough);
463 DataMap["slot"] = std::move(slotVec);
464
465 // adc counts
466 DataMap["cellTime"] = std::move(cellTimeVec);
467 DataMap["cellGain"] = std::move(cellGain);
468 DataMap["cellPedestal"] = std::move(cellPedestal);
469 DataMap["adc2Mev"] = std::move(adc2Mev);
470
471 DataMap[LArSampleIndexStr] = std::move(LArSampleIndexVec); // adcCounts
472
473 //Be verbose
474 ATH_MSG_DEBUG( dataTypeName() << " retrieved with " << nEntries<< " entries" );
475
476 //All collections retrieved okay
477 return DataMap;
478
479 } // getTileData
480
481
482 //--------------------------------------------------------------------------
483
485 {
486 const auto posNeg = std::abs(m_calocell_id->pos_neg(cellid));
487 if (posNeg<1 or posNeg>3) return;
488 static constexpr std::array<int,3> datatypes{2,3,0};
489 m_sub.emplace_back(datatypes[posNeg-1]);
490 }
491
492 //--------------------------------------------------------------------------
493
495 {
496 if(m_calocell_id->pos_neg(cellid)==2)
497 m_sub.emplace_back(1);
498 else
499 m_sub.emplace_back(0);
500 }
501
502
503 //--------------------------------------------------------------------------
504
505} // JiveXML namespace
Scalar eta() const
pseudorapidity method
Scalar phi() const
phi method
#define ATH_CHECK
Evaluate an expression and check for errors.
#define ATH_MSG_ERROR(x)
#define ATH_MSG_WARNING(x)
#define ATH_MSG_DEBUG(x)
OFFLINE_FRAGMENTS_NAMESPACE::PointerType DataType
std::vector< Identifier > ID
char data[hepevt_bytes_allocation_ATLAS]
Definition HepEvt.cxx:11
IndexedConstituentUserInfo::Index Index
#define y
#define x
AthAlgTool(const std::string &type, const std::string &name, const IInterface *parent)
Constructor with parameters:
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T, V, H > &t)
const ServiceHandle< StoreGateSvc > & detStore() const
Container class for CaloCell.
CaloCellContainer::const_iterator beginConstCalo(CaloCell_ID::SUBCALO caloNum) const
get const iterators on cell of just one calo
CaloCellContainer::const_iterator endConstCalo(CaloCell_ID::SUBCALO caloNum) const
int findIndex(const IdentifierHash theHash) const
Return index of the cell with a given hash.
CaloCell_Base_ID::SUBCALO SUBCALO
Definition CaloCell_ID.h:50
This class groups all DetDescr information related to a CaloCell.
DataModel_detail::const_iterator< DataVector > const_iterator
Definition DataVector.h:838
size_type size() const noexcept
Returns the number of elements in the collection.
virtual float pedestal(const HWIdentifier &id, int gain) const =0
This is a "hash" representation of an Identifier.
Templated class to convert any object that is streamable in a ostringstream in a string.
Definition DataType.h:21
LArDigitRetriever(const std::string &type, const std::string &name, const IInterface *parent)
Standard Constructor.
SG::ReadCondHandleKey< LArADC2MeV > m_adc2mevKey
const CaloCell_ID * m_calocell_id
std::vector< Identifier::value_type > m_LArChannelsToIgnoreM5
SG::ReadCondHandleKey< LArOnOffIdMapping > m_cablingKey
virtual std::string dataTypeName() const
Return the name of the data type.
StatusCode initialize()
Default AthAlgTool methods.
const DataMap getLArDigitData(const CaloCellContainer *cellContainer, const std::string &datatype, CaloCell_ID::SUBCALO calotype)
Retrieve Tile cell location and details.
SG::ReadHandleKey< CaloCellContainer > m_sgKey
SG::ReadHandleKey< LArDigitContainer > m_sgKeyLArDigit_esd
virtual StatusCode retrieve(ToolHandle< IFormatTool > &FormatTool)
Retrieve all the data.
SG::ReadHandleKey< LArDigitContainer > m_sgKeyLArDigit_raw
Container class for LArDigit.
Liquid Argon digit base class.
Definition LArDigit.h:25
int feedthrough(const HWIdentifier id) const
Return the feedthrough of a hardware cell identifier : feedthrough = [0,31] Barrel - A/C side or H/...
int slot(const HWIdentifier id) const
Return the slot number of a hardware cell identifier: slot = [1,15] Slot-ID in top part of the crat...
int channel(const HWIdentifier id) const
Return the channel number of a hardware cell identifier channel = [0,127] in all FEB.
Proxy for accessing a range of float values like a vector.
virtual bool isValid() override final
Can the handle be successfully dereferenced?
This header is shared inbetween the C-style server thread and the C++ Athena ServerSvc.
std::map< std::string, DataVect > DataMap
Definition DataType.h:59
std::vector< DataType > DataVect
Defines a map with a key and a vector of DataType objects e.g.
Definition DataType.h:58