Loading [MathJax]/extensions/tex2jax.js
ATLAS Offline Software
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
StripClusterTool.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 "GaudiKernel/ConcurrencyFlags.h"
6 
8 
9 namespace NSWL1 {
10 
11  StripClusterTool::StripClusterTool( const std::string& type, const std::string& name, const IInterface* parent) :
13  {
14  declareInterface<NSWL1::IStripClusterTool>(this);
15  }
16 
18  ATH_MSG_DEBUG( "initializing " << name() );
19 
21 
22  // retrieve the MuonDetectormanager
24  ATH_CHECK(m_idHelperSvc.retrieve());
25  return StatusCode::SUCCESS;
26  }
27 
29  std::vector<std::unique_ptr<StripClusterData>>& clusters,
30  std::vector<std::shared_ptr<std::vector<std::unique_ptr<StripData> >> > &cluster_cache) const {
32  ATH_MSG_DEBUG("Cluster cache received " << cluster_cache.size());
33 
34  bool first_strip=true;
35  for(const auto &this_cl : cluster_cache){
36  float x_pos=0;
37  float y_pos=0;
38  float z_pos=0;
39 
40  float x_lpos=0;
41  float y_lpos=0;
42  float z_lpos=0;
43 
44  int charge=0;
45  int n_strip=0;
46 
47  first_strip=true;
48  float locx=-999999;
49  float locy=-999999;
50  if (this_cl->empty()){
51  ATH_MSG_WARNING("Zero size cluster!!");
52  continue;
53  }
54 
55  const MuonSimDataCollection* sdo_container = nullptr;
56  if(m_isMC) {
57  SG::ReadHandle<MuonSimDataCollection> readMuonSimDataCollection( m_sTgcSdoContainerKey, ctx );
58  if( !readMuonSimDataCollection.isValid() ){
59  ATH_MSG_WARNING("could not retrieve the sTGC SDO container: it will not be possible to associate the MC truth");
60  return StatusCode::FAILURE;
61  }
62  sdo_container = readMuonSimDataCollection.cptr();
63  }
64 
65  for(const auto &strip_cl : *this_cl){
66  n_strip++;
67  ATH_MSG_DEBUG("Start strip" << n_strip);
68  // Save truth deposits associated with cluster should be the same on for the whole strip, so take the first one need to work on this logic
69  if(m_isMC && first_strip) {
70  first_strip=false;
71  Identifier Id = strip_cl->Identity();
72  const MuonGM::sTgcReadoutElement* rdoEl = detManager->getsTgcReadoutElement(Id);
73  auto it = sdo_container->find(Id);
74  if(it == sdo_container->end()) continue;
75  const MuonSimData strip_sdo = it->second;
76  std::vector<MuonSimData::Deposit> deposits;
77  strip_sdo.deposits(deposits);
78  //retrieve the info of the first associated hit, i.e. the fastest in time
79  if (deposits.size()!=1) ATH_MSG_WARNING("Multiple cluster hits for strip!");
80  if (deposits.empty()){
81  ATH_MSG_WARNING("Empty hit here");
82  continue;
83  }
84 
85  int truth_barcode = deposits[0].first.barcode();
86  double truth_localPosX = deposits[0].second.firstEntry();
87  double truth_localPosY = deposits[0].second.secondEntry();
88  Amg::Vector3D hit_gpos(0.,0.,0.);
89  Amg::Vector2D lpos(truth_localPosX,truth_localPosY);
90  rdoEl->surface(Id).localToGlobal(lpos, hit_gpos,hit_gpos);
91  double truth_globalPosX = hit_gpos.x();
92  double truth_globalPosY = hit_gpos.y();
93  double truth_globalPosZ = hit_gpos.z();
94  float truth_energy = strip_sdo.word();
95 
96  if(std::abs(locx-lpos.x())>.001 || std::abs(locy - lpos.y())>.001){
97  ATH_MSG_DEBUG("OLD locx " << locx << " new locx " << lpos.x() << " b " << int(locx!=lpos.x()));
98  ATH_MSG_DEBUG("OLD locy " << locy << " new locy " << lpos.y() << " b " << int(locy!=lpos.y()));
99  ATH_MSG_DEBUG("Cluster hit, truth barcode = " << truth_barcode);
100  ATH_MSG_DEBUG("Cluster hit, truth globalPosX = " << truth_globalPosX
101  << ", truth globalPosY = " << truth_globalPosY
102  << ", truth globalPosZ = " << truth_globalPosZ
103  << ", truth enegy deposit = " << truth_energy);
104  ATH_MSG_DEBUG("Cluster hit, truth localPosX = " << lpos.x()
105  << ", truth localPosY = " << lpos.y()
106  << ", truth enegy deposit = " << truth_energy);
107  }
108  }
109 
110  float s_charge=strip_cl->strip_charge_6bit();
111  charge+=s_charge;
112  x_pos+=strip_cl->globX()*s_charge;
113  y_pos+=strip_cl->globY()*s_charge;
114  z_pos+=strip_cl->globZ()*s_charge;
115 
116  x_lpos+=(strip_cl->locX())*s_charge;
117  y_lpos+=(strip_cl->locY())*s_charge;
118  z_lpos+=(strip_cl->locZ())*s_charge;
119 
120 
121  ATH_MSG_DEBUG("Cluster ------------------------------------------" );
122  ATH_MSG_DEBUG("Cluster strip charge: " << s_charge);
123  ATH_MSG_DEBUG("Cluster strip loc X: " << strip_cl->locX());
124  ATH_MSG_DEBUG("Cluster strip loc Y: " << strip_cl->locY());
125  ATH_MSG_DEBUG("Cluster strip glob X: " << strip_cl->globX());
126  ATH_MSG_DEBUG("Cluster strip glob Y: " << strip_cl->globY());
127  ATH_MSG_DEBUG("Cluster strip glob Z: " << strip_cl->globZ());
128  ATH_MSG_DEBUG("Cluster strip locx dist: " << locx-strip_cl->locX());
129  ATH_MSG_DEBUG("Cluster strip charge o dist: " << s_charge/(locx-strip_cl->locX()));
130  ATH_MSG_DEBUG("Channel " << strip_cl->channelId());
131 
132  }//end of this_cl loop
133 
134  if (charge != 0){
135  if ( std::abs(x_pos/charge)<200. && std::abs(y_pos/charge)<200.){
136  ATH_MSG_WARNING("Cluster ------------------------------------------" );
137  ATH_MSG_WARNING("Cluster strip charge: " << charge );
138  ATH_MSG_WARNING("Cluster strip glob X: " << x_pos << x_pos/charge);
139  ATH_MSG_WARNING("Cluster strip glob Y: " << y_pos << y_pos/charge);
140  ATH_MSG_WARNING("Cluster strip glob Z: " << z_pos << z_pos/charge);
141  }
142  x_pos=x_pos/charge;
143  y_pos=y_pos/charge;
144  z_pos=z_pos/charge;
145  x_lpos=x_lpos/charge;
146  y_lpos=y_lpos/charge;
147  z_lpos=z_lpos/charge;
148  }
149  ATH_MSG_DEBUG("Cluster dump with X:" << x_pos << " Y: " << y_pos << " Z: " << z_pos << " cluster charge: " << charge);
150  ATH_MSG_DEBUG("Cluster dump with lX:" << x_lpos << " lY: " << y_lpos << " lZ: " << z_lpos << " cluster charge: " << charge);
151 
152  auto stripClOfflData=std::make_unique<StripClusterOfflineData>(
153  this_cl->at(0)->bandId(),
154  this_cl->at(0)->trig_BCID(),
155  this_cl->at(0)->sideId(),
156  this_cl->at(0)->phiId(),
157  this_cl->at(0)->isSmall(),
158  this_cl->at(0)->moduleId(),
159  this_cl->at(0)->sectorId(),
160  this_cl->at(0)->wedge(),
161  this_cl->at(0)->layer(),
162  n_strip,
163  charge,
164  x_pos,
165  y_pos,
166  z_pos);
167  clusters.push_back(std::move(stripClOfflData));
168  }
169  return StatusCode::SUCCESS;
170  }
171 
172 
173  StatusCode StripClusterTool::cluster_strip_data(const EventContext& ctx, std::vector<std::unique_ptr<StripData>>& strips, std::vector< std::unique_ptr<StripClusterData> >& clusters) const {
174 
175  if(strips.empty()){
176  ATH_MSG_WARNING("Received 0 strip hits... Skip event");
177  return StatusCode::SUCCESS;
178  }
179 
180  std::map<uint32_t, std::vector<std::unique_ptr<StripData>> > stripMap;
181  // Filter by layer:
182  for (auto& st : strips) {
183  // sideId [0, 1]
184  auto sideId = st->sideId();
185  // sectorType: small==0, large==1
186  auto sectorType = st->sectorType();
187  // sectorId [1,8]
188  auto sectorId = st->sectorId();
189  // etaId [1,3]
190  auto etaId = st->moduleId();
191  // wedgeId [1,2]
192  auto wedgeId = st->wedge();
193  // layer [1,4]
194  auto layerId = st->layer();
195 
196  // add 1 as the first placeholder, hence cache_hash always has 7 digits
197  std::string id_str = std::to_string(1)+std::to_string(sideId)+std::to_string(sectorType)+std::to_string(sectorId)+std::to_string(etaId)+std::to_string(wedgeId)+std::to_string(layerId);
198  const uint32_t cache_hash = atoi(id_str.c_str());
199 
200  auto it = stripMap.find(cache_hash);
201  if (it != stripMap.end()){
202  it->second.push_back(std::move(st));
203  }
204  else{
205  stripMap[cache_hash].push_back(std::move(st));
206  }
207 
208  }
209 
210  std::vector< std::shared_ptr<std::vector<std::unique_ptr<StripData> >> > cluster_cache;
211  for (auto &item : stripMap) {
212  std::vector<std::unique_ptr<StripData>>& stripList = item.second;
213 
214  //S.I sort strip w.r.t channelId in ascending order
215  std::sort(stripList.begin(), stripList.end(), [](const auto& s1,const auto& s2) { return s1->channelId()<s2->channelId(); });
216 
217  auto hit=stripList.begin();
218  ATH_MSG_DEBUG("Cluster Hits :" << (*hit)->channelId() << " " << m_idHelperSvc->stgcIdHelper().gasGap( (*hit)->Identity())
219  << " " << (*hit)->moduleId() << " " << (*hit)->sectorId() << " " << (*hit)->wedge() << " " << (*hit)->sideId() );
220  int first_ch=(*hit)->channelId();//channel id of the first strip
221  int prev_ch=-1;
222 
223  auto cr_cluster=std::make_shared< std::vector<std::unique_ptr<StripData>> >();
224 
225  for(auto& this_hit : stripList){
226  if(!(this_hit)->readStrip() ) continue;
227  if( ((this_hit)->bandId()==-1 || this_hit->phiId()==-1) ){
228  ATH_MSG_WARNING("Read Strip without BandId :" << (this_hit)->channelId() << " " << m_idHelperSvc->stgcIdHelper().gasGap( (this_hit)->Identity())
229  << " " << (this_hit)->moduleId() << " " << (this_hit)->sectorId() << " " << (this_hit)->wedge() << " " << (this_hit)->sideId() );
230  continue;
231  }
232 
233  if (prev_ch==-1){//for the first time...
234  prev_ch = first_ch;
235  cr_cluster->push_back(std::move(this_hit));
236  continue;
237  }
238 
239  int this_ch=(this_hit)->channelId();
240  if ( (this_ch < prev_ch)) {
241  ATH_MSG_ERROR("Hits Ordered incorrectly!!!" );
242  return StatusCode::FAILURE;
243  }
244 
245  if (this_ch == prev_ch || this_ch == prev_ch+1) cr_cluster->push_back(std::move(this_hit)); // form cluster with adjacent +-1 strips
246  else {
247  cluster_cache.push_back(std::move(cr_cluster)); //put the current cluster into the clusters buffer
248  cr_cluster=std::make_shared<std::vector<std::unique_ptr<StripData>>>(); //create a new empty cluster and assign this hit as the first hit
249  cr_cluster->push_back(std::move(this_hit));
250  }
251 
252  prev_ch=this_ch;
253  }
254 
255  if(!cr_cluster->empty()) cluster_cache.push_back(std::move(cr_cluster)); //don't forget the last cluster in the loop
256  }
257 
258  ATH_MSG_DEBUG("Found :" << cluster_cache.size() << " clusters");
259  ATH_CHECK(fill_strip_validation_id(ctx, clusters, cluster_cache));
260  cluster_cache.clear();
261 
262  return StatusCode::SUCCESS;
263  }
264 }
NSWL1::StripClusterTool::initialize
virtual StatusCode initialize() override
Definition: StripClusterTool.cxx:17
MuonSimData::word
int word() const
Definition: MuonSimData.h:89
SG::ReadCondHandle
Definition: ReadCondHandle.h:44
SG::ReadHandle::cptr
const_pointer_type cptr()
Dereference the pointer.
xAOD::uint32_t
setEventNumber uint32_t
Definition: EventInfo_v1.cxx:127
Amg::Vector2D
Eigen::Matrix< double, 2, 1 > Vector2D
Definition: GeoPrimitives.h:48
SG::ReadHandle
Definition: StoreGate/StoreGate/ReadHandle.h:67
NSWL1::StripClusterTool::m_sTgcSdoContainerKey
SG::ReadHandleKey< MuonSimDataCollection > m_sTgcSdoContainerKey
Definition: StripClusterTool.h:74
skel.it
it
Definition: skel.GENtoEVGEN.py:407
NSWL1::StripClusterTool::fill_strip_validation_id
StatusCode fill_strip_validation_id(const EventContext &ctx, std::vector< std::unique_ptr< StripClusterData >> &clusters, std::vector< std::shared_ptr< std::vector< std::unique_ptr< StripData > >> > &cluster_cache) const
Definition: StripClusterTool.cxx:28
MuonGM::MuonClusterReadoutElement::surface
virtual const Trk::PlaneSurface & surface() const override
access to chamber surface (phi orientation), uses the first gas gap
Definition: MuonClusterReadoutElement.h:123
python.CaloAddPedShiftConfig.type
type
Definition: CaloAddPedShiftConfig.py:42
NSWL1::StripClusterTool::m_detManagerKey
SG::ReadCondHandleKey< MuonGM::MuonDetectorManager > m_detManagerKey
Definition: StripClusterTool.h:70
MuonSimData::deposits
void deposits(std::vector< Deposit > &deposits) const
Definition: MuonSimData.h:99
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x)
Definition: AthMsgStreamMacros.h:33
EL::StatusCode
::StatusCode StatusCode
StatusCode definition for legacy code.
Definition: PhysicsAnalysis/D3PDTools/EventLoop/EventLoop/StatusCode.h:22
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition: AthMsgStreamMacros.h:29
MuonGM::sTgcReadoutElement
An sTgcReadoutElement corresponds to a single STGC module; therefore typicaly a barrel muon station c...
Definition: MuonDetDescr/MuonReadoutGeometry/MuonReadoutGeometry/sTgcReadoutElement.h:30
MuonSimDataCollection
Definition: MuonSimDataCollection.h:21
test_pyathena.parent
parent
Definition: test_pyathena.py:15
ATH_CHECK
#define ATH_CHECK
Definition: AthCheckMacros.h:40
StripClusterTool.h
SG::VarHandleKey::initialize
StatusCode initialize(bool used=true)
If this object is used as a property, then this should be called during the initialize phase.
Definition: AthToolSupport/AsgDataHandles/Root/VarHandleKey.cxx:103
CaloCondBlobAlgs_fillNoiseFromASCII.channelId
channelId
Definition: CaloCondBlobAlgs_fillNoiseFromASCII.py:122
SG::ReadHandle::isValid
virtual bool isValid() override final
Can the handle be successfully dereferenced?
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:240
ActsTrk::to_string
std::string to_string(const DetectorType &type)
Definition: GeometryDefs.h:34
charge
double charge(const T &p)
Definition: AtlasPID.h:931
NSWL1::StripClusterTool::cluster_strip_data
virtual StatusCode cluster_strip_data(const EventContext &ctx, std::vector< std::unique_ptr< StripData >> &strips, std::vector< std::unique_ptr< StripClusterData >> &clusters) const override
Definition: StripClusterTool.cxx:173
SG::CondHandleKey::initialize
StatusCode initialize(bool used=true)
MuonSimData
Definition: MuonSimData.h:62
item
Definition: ItemListSvc.h:43
Amg::Vector3D
Eigen::Matrix< double, 3, 1 > Vector3D
Definition: GeoPrimitives.h:47
NSWL1::StripClusterTool::StripClusterTool
StripClusterTool(const std::string &type, const std::string &name, const IInterface *parent)
Definition: StripClusterTool.cxx:11
ATH_MSG_WARNING
#define ATH_MSG_WARNING(x)
Definition: AthMsgStreamMacros.h:32
ReadCellNoiseFromCoolCompare.s2
s2
Definition: ReadCellNoiseFromCoolCompare.py:379
RunTileMonitoring.clusters
clusters
Definition: RunTileMonitoring.py:133
CxxUtils::atoi
int atoi(std::string_view str)
Helper functions to unpack numbers decoded in string into integers and doubles The strings are requir...
Definition: Control/CxxUtils/Root/StringUtils.cxx:85
AthAlgTool
Definition: AthAlgTool.h:26
Trk::PlaneSurface::localToGlobal
virtual void localToGlobal(const Amg::Vector2D &locp, const Amg::Vector3D &mom, Amg::Vector3D &glob) const override final
Specified for PlaneSurface: LocalToGlobal method without dynamic memory allocation.
Definition: PlaneSurface.cxx:204
NSWL1
A trigger trigger candidate for a stgc sector.
Definition: NSWL1Simulation.cxx:7
NSWL1::StripClusterTool::m_idHelperSvc
ServiceHandle< Muon::IMuonIdHelperSvc > m_idHelperSvc
Definition: StripClusterTool.h:71
NSWL1::StripClusterTool::m_isMC
Gaudi::Property< bool > m_isMC
Definition: StripClusterTool.h:73
Identifier
Definition: IdentifierFieldParser.cxx:14