11#include "CLHEP/Random/RandomEngine.h"
12#include "CLHEP/Random/RandGaussZiggurat.h"
14#include "Identifier/Identifier.h"
50 return StatusCode::SUCCESS;
60 ATH_MSG_ERROR(
"Could not get background CscRawDataContainer called " << bkgContainer.
name() <<
" from store " << bkgContainer.
store());
61 return StatusCode::FAILURE;
63 ATH_MSG_DEBUG(
"Found background CscRawDataContainer called " << bkgContainer.
name() <<
" in store " << bkgContainer.
store());
66 if(!signalContainer.
isValid()) {
67 ATH_MSG_ERROR(
"Could not get signal CscRawOverlayContainer called " << signalContainer.
name() <<
" from store " << signalContainer.
store());
68 return StatusCode::FAILURE;
70 ATH_MSG_DEBUG(
"Found signal CscRawOverlayContainer called " << signalContainer.
name() <<
" in store " << signalContainer.
store());
73 ATH_CHECK(outputContainer.
record(std::make_unique<CscRawDataContainer>(bkgContainer->size())));
74 if (!outputContainer.
isValid()) {
75 ATH_MSG_ERROR(
"Could not record output CscRawOverlayContainer called " << outputContainer.
name() <<
" to store " << outputContainer.
store());
76 return StatusCode::FAILURE;
78 ATH_MSG_DEBUG(
"Recorded output CscRawOverlayContainer called " << outputContainer.
name() <<
" in store " << outputContainer.
store());
85 return StatusCode::SUCCESS;
98 std::vector < std::pair<IdentifierHash, bool> > overlapMap;
101 overlapMap.emplace_back(hashId,
false);
106 auto search = std::lower_bound( overlapMap.begin(), overlapMap.end(), hashId,
107 [](
const std::pair<IdentifierHash, bool> &lhs,
IdentifierHash rhs) ->
bool { return lhs.first < rhs; } );
108 if (
search == overlapMap.end() ||
search->first != hashId) {
110 std::unique_ptr<CscRawDataCollection> bkgCollection =
copyCollection(ptr);
112 if (outputContainer->
addCollection(bkgCollection.get(), hashId).isFailure()) {
113 ATH_MSG_ERROR(
"Adding background Collection with hashId " << hashId <<
" failed");
114 return StatusCode::FAILURE;
116 (void)bkgCollection.release();
126 rngWrapper->
setSeed( name(), Gaudi::Hive::currentContext() );
127 CLHEP::HepRandomEngine *rndmEngine(*rngWrapper);
131 for (
const auto &[hashId, overlap] : overlapMap) {
136 std::unique_ptr<CscRawDataCollection> outputCollection{};
147 mergeCollections(bkgCollection, signalCollection, outputCollection.get(), rndmEngine);
152 uint16_t numSamples = signalCollection->
numSamples();
161 std::vector<uint16_t> samples;
164 uint32_t hashOffset =
data->hashId();
167 for (uint16_t j = 0; j <
width; ++j) {
168 uint32_t stripHash = hashOffset + j;
171 std::vector<uint16_t> stripSamples;
172 bool extractSamplesStatus =
data->samples(j, numSamples, stripSamples);
173 if (!extractSamplesStatus) {
175 <<
" Online Cluster width = " <<
width
176 <<
" for number of Samples = " << numSamples
177 <<
" continuing ...");
179 for (uint16_t sample : stripSamples) {
180 double sampleNoise = CLHEP::RandGaussZiggurat::shoot(rndmEngine, 0.0, stripNoise);
181 float adcCount = sample + sampleNoise;
183 ATH_MSG_DEBUG(
"value out of range (copying over signal): " << adcCount <<
" "
184 <<
" Setting it to max value = " <<
MAX_AMPL
185 <<
" IdentifierHash is " << stripHash);
188 samples.push_back( (uint16_t) std::rint(adcCount) );
194 auto rdo = std::make_unique<CscRawData>(samples,
data->address(),
data->identify(),
data->time(),
data->rpuID(),
data->width(),
data->hashId());
198 for (uint16_t j = 0; j <
width; ++j) {
207 outputCollection->push_back(rdo.release());
212 if (outputContainer->
addCollection(outputCollection.get(), hashId).isFailure()) {
213 ATH_MSG_ERROR(
"Adding overlaid Collection with hashId " << hashId <<
" failed");
214 return StatusCode::FAILURE;
218 (void)outputCollection.release();
223 return StatusCode::SUCCESS;
227std::unique_ptr<CscRawDataCollection>
229 bool propertiesOnly)
const
231 auto outputCollection = std::make_unique<CscRawDataCollection>(collection->
identify());
234 outputCollection->setIdentifyHash(collection->
identifyHash());
235 outputCollection->set_eventType(collection->
eventType());
236 outputCollection->setRodId(collection->
rodId());
237 outputCollection->setSubDetectorId(collection->
subDetectorId());
238 outputCollection->set_samplingPhase(collection->
samplingPhase());
239 outputCollection->set_triggerType(collection->
triggerType());
241 outputCollection->set_scaAddress(collection->
scaAddress());
242 for (uint8_t dataType : collection->
dataType()) {
243 outputCollection->addDataType(dataType);
247 if (propertiesOnly) {
248 return outputCollection;
251 for (
const CscRawData *existingDatum : *collection) {
253 auto *datumCopy =
new CscRawData(*existingDatum);
254 outputCollection->push_back(datumCopy);
257 return outputCollection;
262 data.clear();
if ( !coll )
return;
265 for ( ; idata != edata; ++idata ) {
266 if ( (*idata)->rpuID() == spuID )
data.push_back( *idata );
268 ATH_MSG_DEBUG(
"spuData(): made data vector of size "<<
data.size()<<
" for SPU "<<spuID);
273 int measuresPhi = ( (address & 0x00000100) >> 8);
274 if (address<2147483640 && measuresPhi) {
275 int stationEta = ( ((address & 0x00001000) >> 12 ) == 0x0) ? -1 : 1;
276 if (stationEta>0) {
return true; }
286 CLHEP::HepRandomEngine *rndmEngine)
const
291 unsigned int nSigSamples = bkgCollection->
numSamples();
292 unsigned int nOvlSamples = signalCollection->
numSamples();
295 unsigned int dataSamplingTime = bkgCollection->
rate();
296 unsigned int ovlSamplingTime = signalCollection->
rate();
298 if ( dataSamplingTime != ovlSamplingTime ) {
299 ATH_MSG_ERROR(
"Overlay of inconsistent data - sampling times not the same "
300 << dataSamplingTime <<
" ns " << ovlSamplingTime <<
" ns");
301 throw std::runtime_error(
"mergeCollections(): sampling time mismatch");
304 if ( nSigSamples != nOvlSamples ) {
305 ATH_MSG_ERROR(
"Overlay of inconsistent data - number of samples not the same "
306 << nSigSamples <<
" " << nOvlSamples);
307 throw std::runtime_error(
"mergeCollections(): number of samples mismatch");
312 uint16_t clusterCounts[] = {0,0,0,0,0,0,0,0,0,0};
313 uint16_t rpuCount[] = {0,0};
314 for ( uint16_t spuID=0; spuID<10; ++spuID) {
316 std::vector<const CscRawData*> sigData;
317 this->
spuData(bkgCollection, spuID, sigData);
319 std::vector<const CscRawData*> ovlData;
320 this->
spuData(signalCollection, spuID, ovlData);
325 if ( spuID == 4 || spuID == 9 ) layer=4;
326 for (
int j=0; j<=layer; ++j ) {
327 std::map< int,std::vector<uint16_t> > sigSamples;
328 std::map< int,std::vector<uint16_t> > ovlSamples;
332 uint32_t ovlAddress = this->
stripData( ovlData, nOvlSamples, ovlSamples, ovlHash, spuID, j ,
false);
333 if (sigSamples.size()==0 && ovlSamples.size()==0)
continue;
335 uint32_t hash = std::min( sigHash, ovlHash );
336 uint32_t address = std::min( sigAddress, ovlAddress );
337 if (sigSamples.size()!=0 && ovlSamples.size()!=0 &&
needtoflip(address)){
338 ATH_MSG_DEBUG(
"Looking for overlap of hashes and addresses within widths because needtoflip");
339 msg() << MSG::VERBOSE ;
340 std::set<int> sig;
int lastindex=-1;
341 for (std::map<
int,std::vector<uint16_t> >
::const_iterator si=sigSamples.begin(); si!=sigSamples.end(); ++si) {
342 if (si!=sigSamples.begin() && si->first-lastindex!=1)
break;
344 sig.insert(si->first);
msg() << si->first <<
" ";
348 msg() <<MSG::VERBOSE ;
349 for (std::map<
int,std::vector<uint16_t> >
::const_iterator so=ovlSamples.begin(); so!=ovlSamples.end(); ++so) {
351 msg() << (so->first)-1 <<
" ";
352 if (sig.find((so->first)-1)!=sig.end()) {overlap=
true;
msg() <<
"!!";}
353 msg() << (so->first) <<
" ";
354 if (sig.find((so->first))!=sig.end()) {overlap=
true;
msg() <<
"!!";}
355 msg() << (so->first)+1 <<
" ";
356 if (sig.find((so->first)+1)!=sig.end()) {overlap=
true;
msg() <<
"!!";}
360 ATH_MSG_DEBUG(
"Taking max of hashes and addresses because needtoflip and no overlap");
361 hash = std::max( sigHash, ovlHash );
362 address = std::max( sigAddress, ovlAddress );
367 std::set<int> insertedstrips, readstrips;
368 for (std::map<
int,std::vector<uint16_t> >
::const_iterator s=sigSamples.begin(); s!=sigSamples.end(); ++s){readstrips.insert(s->first);}
369 for (std::map<
int,std::vector<uint16_t> >
::const_iterator si=ovlSamples.begin(); si!=ovlSamples.end(); ++si){readstrips.insert(si->first);}
371 std::vector<CscRawData*> datums = this->
overlay(sigSamples, ovlSamples,address, spuID, outputCollection->identify(), hash, rndmEngine);
373 for (
unsigned int di=0; di<datums.size(); ++di){
377 int stripstart = ( address & 0x000000FF) + 1 + 0;
378 ATH_MSG_DEBUG(
"Datum in layer="<<j<<
" has hash="<<hash<<
" address="<<address<<
" stripstart="<<stripstart<<
", "<< *datum );
379 if (datum->
width()==0) {
385 int stationName = ( ( address & 0x00010000) >> 16 ) + 50;
386 int stationEta = ( ((address & 0x00001000) >> 12 ) == 0x0) ? -1 : 1;
387 int stationPhi = ( ( address & 0x0000E000) >> 13 ) + 1;
389 ATH_MSG_VERBOSE(
"stationName,Eta,Phi="<<stationName<<
","<<stationEta<<
","<<stationPhi<<
" - me="<<me);
391 for (
unsigned int j=0; j<datum->
width(); ++j) {
392 int chamberLayer = ( (address & 0x00000800) >> 11) + 1;
393 int wireLayer = ( (address & 0x00000600) >> 9) + 1;
394 int measuresPhi = ( (address & 0x00000100) >> 8);
395 int strip = ( address & 0x000000FF) + 1 + j;
396 ATH_MSG_VERBOSE(
"chamberlayer,wirelayer,measuresphi,strip="<<chamberLayer<<
","<<wireLayer<<
","<<measuresPhi<<
","<<
strip);
399 int stationEta = ( ((address & 0x00001000) >> 12 ) == 0x0) ? -1 : 1;
405 insertedstrips.insert(
strip);
410 ATH_MSG_WARNING(
"Invalid CSC Identifier in merge! - skipping " << channelId );
413 else{
ATH_MSG_DEBUG(
"Valid CSC Identifier in merge " << channelId);}
415 if (good){ outputCollection->push_back(datum); }
419 if (spuID <10) clusterCounts[spuID] += 1;
420 if ( spuID <= 4 ) rpuCount[0] = 5;
421 else if ( spuID > 4 && spuID <= 9 ) rpuCount[1] = 11;
425 if (readstrips!=insertedstrips){
427 std::ostringstream readstream;
428 for (std::set<int>::const_iterator i = readstrips.begin(); i!=readstrips.end(); ++i){readstream<<*i<<
" ";}
430 std::ostringstream insertstream;
431 for (std::set<int>::const_iterator i = insertedstrips.begin(); i!=insertedstrips.end(); ++i){insertstream<<*i<<
" ";}
437 for (
unsigned int i=0; i<10; ++i) outputCollection->set_spuCount(i,clusterCounts[i]);
438 for (
unsigned int i=0; i<2; ++i) {
if (rpuCount[i] != 0) outputCollection->addRPU(rpuCount[i]); }
444 const unsigned int numSamples,
445 std::map<
int,std::vector<uint16_t> >& samples,
447 const uint16_t spuID,
448 const int gasLayer,
bool isdata)
const
450 ATH_MSG_DEBUG(
"stripData<>() begin: gasLayer="<<gasLayer<<
" spuID="<<spuID<<
" isdata="<<isdata);
455 uint32_t maxInt = 2147483640;
456 uint32_t address = maxInt;
460 std::vector<const CscRawData*>::const_iterator idata =
data.begin();
461 std::vector<const CscRawData*>::const_iterator edata =
data.end();
462 for ( ; idata != edata; ++idata ) {
464 uint32_t hashOffset = datum->
hashId();
468 m_idHelperSvc->cscIdHelper().get_id(hashOffset, stripId, &context);
470 int layer =
m_idHelperSvc->cscIdHelper().wireLayer( stripId );
476 bool non_precision = (gasLayer==layer) && (spuID==4 || spuID==9);
477 bool precision = (gasLayer==0) && (!(spuID==4 || spuID==9));
478 bool check = precision || non_precision;
485 unsigned int newaddress = datum->
address();
488 ATH_MSG_VERBOSE(
"needtoflip in stripdata, newaddress was = "<<newaddress<<
", strip was = "<<
strip);
494 uint32_t oldFirstStrip = uint32_t (newaddress & 0x000000FF);
495 uint32_t newFirstStrip = uint32_t (47-oldFirstStrip) -
width +1;
496 newaddress=newaddress - oldFirstStrip + newFirstStrip;
498 uint32_t oldStrip = uint32_t (
strip & 0x000000FF);
499 uint32_t newStrip = uint32_t (49-oldStrip);
502 ATH_MSG_VERBOSE(
"needtoflip in stripdata, newaddress now = "<<newaddress<<
", strip now = "<<
strip);
506 if (hash == maxInt) hash=0;
507 if (address == maxInt) address=0;
508 if ( hashOffset > hash ) hash = hashOffset;
509 if ( newaddress > address ) address = newaddress;
512 if ( hashOffset < hash ) hash = hashOffset;
513 if ( newaddress < address ) address = newaddress;
516 ATH_MSG_DEBUG(
"stripData(): width="<<
width<<
" hashOffset="<<hashOffset<<
" datumaddress="<<datum->
address()<<
" layer="<<layer<<
" strip="<<
strip<<
", hash="<<hash<<
" address="<<address);
518 for (
unsigned int j=0; j<
width; ++j) {
519 std::vector<uint16_t> adcs;
520 bool extractSamples = datum->
samples(j, numSamples, adcs);
521 if ( !extractSamples ) {
523 <<
" Online Cluster width = " <<
width <<
" for number of Samples = " << numSamples);
526 int newstrip = (
strip+j);
531 samples.insert ( std::make_pair( newstrip, adcs) );
536 ATH_MSG_DEBUG(
"stripData<>() end: hash=" << hash <<
" address=" << address);
541 const std::map<
int,std::vector<uint16_t> >& ovlSamples,
542 const uint32_t address,
543 const uint16_t spuID,
544 const uint16_t collId,
546 CLHEP::HepRandomEngine *rndmEngine)
const
548 ATH_MSG_DEBUG(
"overlay<>() begin: hash="<<hash<<
" address="<<address);
549 std::vector<CscRawData*> datas;
552 if ( spuID == 4 || spuID == 9 )
max = 48;
553 std::vector<uint16_t> samples;
557 int myhash=hash;
int myaddress=address;
560 for (
int i=1; i<=
max; ++i) {
561 sig = sigSamples.find(i);
562 ovl = ovlSamples.find(i);
565 if ( sig != sigSamples.end() && ovl == ovlSamples.end() ) {
567 for (
unsigned int j=0; j<(*sig).second.size(); ++j ) {
568 samples.push_back( (*sig).second.at(j) );
569 assert((*sig).second.at(j)<=
MAX_AMPL);
573 else if ( sig == sigSamples.end() && ovl != ovlSamples.end() ) {
577 for (
unsigned int j=0; j<(*ovl).second.size(); ++j ) {
578 double theNoise = CLHEP::RandGaussZiggurat::shoot(rndmEngine, 0.0, noise);
579 float adcCount = (*ovl).second.at(j) + theNoise ;
581 ATH_MSG_DEBUG(
"value out of range (adding noise): " << adcCount <<
" "
582 <<
" Setting it to max value = " <<
MAX_AMPL
583 <<
" IdentifierHash is " << (myhashw));
586 samples.push_back( (uint16_t) rint(adcCount) );
590 else if ( sig != sigSamples.end() && ovl != ovlSamples.end() ) {
593 double pedestal =
m_cscCalibTool->stripPedestal( (myhashw),
false );
594 for (
unsigned int j=0; j<(*sig).second.size(); ++j ) {
595 float adcCount = (*sig).second.at(j) + (*ovl).second.at(j) - pedestal ;
597 ATH_MSG_DEBUG(
"value out of range (adding data+MC samples - pedestal): " << adcCount <<
" "
598 <<
" Setting it to max value = " <<
MAX_AMPL
599 <<
" IdentifierHash is " << (myhashw));
602 samples.push_back( (uint16_t) rint(adcCount) );
607 if (
used==
false && datas.size()>0 ){
608 if (
needtoflip(myaddress)) {myhash-=1; myaddress-=1;}
609 else {myhash+=1; myaddress+=1;}
614 if ( (
used==
false||i==
max) && samples.size()>0){
619 datas.push_back(rawData);
620 ATH_MSG_DEBUG(
"overlay<>() add datum: hash="<<myhash<<
" address="<<myaddress<<
" width="<<
width);
627 ATH_MSG_DEBUG(
"overlay<>() end: CscRawDatas size="<<datas.size());
#define ATH_CHECK
Evaluate an expression and check for errors.
#define ATH_MSG_VERBOSE(x)
#define ATH_MSG_WARNING(x)
char data[hepevt_bytes_allocation_ATLAS]
Handle class for reading from StoreGate.
Handle class for recording to StoreGate.
A wrapper class for event-slot-local random engines.
void setSeed(const std::string &algName, const EventContext &ctx)
Set the random seed using a string (e.g.
An algorithm that can be simultaneously executed in multiple threads.
void spuData(const CscRawDataCollection *coll, const uint16_t spuID, std::vector< const CscRawData * > &data) const
get the data in one SPU of a chamber
StatusCode overlayContainer(const CscRawDataContainer *bkgContainer, const CscRawDataContainer *signalContainer, CscRawDataContainer *outputContainer) const
Overlay signal on the background container and record to the output one.
ToolHandle< ICscCalibTool > m_cscCalibTool
ToolHandle< Muon::ICSC_RDO_Decoder > m_cscRdoDecoderTool
std::unique_ptr< CscRawDataCollection > copyCollection(const CscRawDataCollection *collection, bool propertiesOnly=false) const
Copy CscRawDataCollection, optionally only copy properties.
CscOverlay(const std::string &name, ISvcLocator *pSvcLocator)
std::vector< CscRawData * > overlay(const std::map< int, std::vector< uint16_t > > &sigSamples, const std::map< int, std::vector< uint16_t > > &ovlSamples, const uint32_t address, const uint16_t spuID, const uint16_t collId, const uint32_t hash, CLHEP::HepRandomEngine *rndmEngine) const
do the overlay - summing the ADC samples on one plane if there is overlap between zero bias data and ...
ServiceHandle< IAthRNGSvc > m_rndmSvc
SG::WriteHandleKey< CscRawDataContainer > m_outputKey
void mergeCollections(const CscRawDataCollection *bkgCollection, const CscRawDataCollection *signalCollection, CscRawDataCollection *outputCollection, CLHEP::HepRandomEngine *rndmEngine) const
In case of overlap merge signal and background collections.
virtual StatusCode execute(const EventContext &ctx) const override final
uint32_t stripData(const std::vector< const CscRawData * > &data, const unsigned int numSamples, std::map< int, std::vector< uint16_t > > &samples, uint32_t &hash, const uint16_t spuID, const int gasLayer, bool isdata) const
data in one gas lauer
ServiceHandle< Muon::IMuonIdHelperSvc > m_idHelperSvc
virtual StatusCode initialize() override final
bool needtoflip(const int address) const
SG::ReadHandleKey< CscRawDataContainer > m_bkgInputKey
Gaudi::Property< bool > m_isDataOverlay
SG::ReadHandleKey< CscRawDataContainer > m_signalInputKey
Collection of CSC Raw Hits, arranged according to CSC Detector Elements Author: Ketevi A.
uint32_t scaAddress() const
bool samplingPhase() const
IdentifierHash identifyHash() const
Returns the OFFLINE identifier hash for this collection.
uint8_t firstBitSummary() const
uint16_t identify() const
access methods
uint16_t subDetectorId() const
uint16_t numSamples() const
uint32_t eventType() const
const std::vector< uint8_t > & dataType() const
uint8_t rate() const
the rate could be 25 or 50 ns
This container provides access to collections of CSC RDOs and a mechanism for recording them.
Class to hold the electronic output for a single CSC readout channel: n sampling ADC data + the addre...
void setHashID(uint32_t hash)
const std::vector< uint16_t > & samples() const
void setTime(uint16_t time)
DataModel_detail::const_iterator< DataVector > const_iterator
const_iterator end() const noexcept
Return a const_iterator pointing past the end of the collection.
const_iterator begin() const noexcept
Return a const_iterator pointing at the beginning of the collection.
This class saves the "context" of an expanded identifier (ExpandedIdentifier) for compact or hash ver...
virtual size_t numberOfCollections() const override final
return number of collections
auto GetAllHashPtrPair() const
virtual const T * indexFindPtr(IdentifierHash hashId) const override final
return pointer on the found entry or null if out of range using hashed index - fast version,...
virtual StatusCode addCollection(const T *coll, IdentifierHash hashId) override final
insert collection into container with id hash if IDC should not take ownership of collection,...
This is a "hash" representation of an Identifier.
virtual bool isValid() override final
Can the handle be successfully dereferenced?
const_pointer_type cptr()
Dereference the pointer.
std::string store() const
Return the name of the store holding the object we are proxying.
const std::string & name() const
Return the StoreGate ID for the referenced object.
StatusCode record(std::unique_ptr< T > data)
Record a const object to the store.
virtual bool isValid() override final
Can the handle be successfully dereferenced?
pointer_type ptr()
Dereference the pointer.
holding In fact this class is here in order to allow STL container for all features This class is sho...
void search(TDirectory *td, const std::string &s, std::string cwd, node *n)
recursive directory search for TH1 and TH2 and TProfiles