ATLAS Offline Software
Loading...
Searching...
No Matches
TrigTRTHTHCounter.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#include "TrigTRTHTHCounter.h"
8
9#include <cmath>
10#include <memory>
11
12//Function to calculate distance for road algorithm
13float dist2COR(float R, float phi1, float phi2){
14 float PHI= std::abs(phi1-phi2);
15 return std::abs(R*std::sin(PHI));
16}
17
18//TRT hit struct used for convenience
19struct TRT_hit {
20 float phi;
21 float R;
22 bool isHT;
23};
24
25//Constructor of above struct
26TRT_hit make_hit(float phi, float R, bool isHT){
27 TRT_hit my_hit={phi,R,isHT};
28 return my_hit;
29}
30
31//hth_eta_match is used for matching hits to RoI eta
32bool hth_eta_match(float caleta, float trteta, float etaWindow){
33 return std::abs(caleta)<etaWindow or caleta*trteta>0.;
34}
35//---------------------------------------------------------------------------------
36
37
38TrigTRTHTHCounter::TrigTRTHTHCounter(const std::string & name, ISvcLocator* pSvcLocator)
39 : AthReentrantAlgorithm(name, pSvcLocator)
40{}
41
43{
44 ATH_MSG_DEBUG ( "Initialising this TrigTRTHTHCounter: " << name());
45
46 // Get a TRT identifier helper
47 ATH_CHECK( detStore()->retrieve(m_trtHelper, "TRT_ID") );
48
49 if (!m_doFullScan){
50 ATH_MSG_INFO ( "PhiHalfWidth: " << m_phiHalfWidth << " EtaHalfWidth: "<< m_etaHalfWidth);
51 } else {
52 ATH_MSG_INFO ( "FullScan mode");
53 }
54
55 ATH_CHECK( m_roiCollectionKey.initialize() );
56 ATH_CHECK( m_trtDCContainerKey.initialize() );
57 ATH_CHECK( m_trigRNNOutputKey.initialize() );
58
59 return StatusCode::SUCCESS;
60}
61
62
63StatusCode TrigTRTHTHCounter::execute(const EventContext& ctx) const {
64
65 ATH_MSG_DEBUG( "Executing " <<name());
66
67 auto trigRNNOutputColl = SG::makeHandle (m_trigRNNOutputKey, ctx);
68 ATH_CHECK(trigRNNOutputColl.record (std::make_unique<xAOD::TrigRNNOutputContainer>(),
69 std::make_unique<xAOD::TrigRNNOutputAuxContainer>()));
70
71 ATH_MSG_DEBUG( "Made WriteHandle " << m_trigRNNOutputKey );
72
73 auto roiCollection = SG::makeHandle(m_roiCollectionKey, ctx);
74 ATH_MSG_DEBUG( "Made handle " << m_roiCollectionKey );
75
76 if (roiCollection->size()==0) {
77 ATH_MSG_DEBUG(" RoI collection size = 0");
78 return StatusCode::SUCCESS;
79 }
80
81 const TrigRoiDescriptor* roiDescriptor = *(roiCollection->begin());
82
83 ATH_MSG_DEBUG(" RoI ID = " << (roiDescriptor)->roiId()
84 << ": Eta = " << (roiDescriptor)->eta()
85 << ", Phi = " << (roiDescriptor)->phi());
86
87
88 std::vector<float> trththits{0,-999,0,-999,-999,-999};
89
90 //Vectors to hold the count of total and HT TRT hits in the coarse bins
91 std::vector<int> count_httrt_c(m_nBinCoarse);
92 std::vector<int> count_tottrt_c(m_nBinCoarse);
93
94 //Vectors to hold the count of total and HT TRT hits in the fine bins
95 std::vector<int> count_httrt(3*m_nBinFine);
96 std::vector<int> count_tottrt(3*m_nBinFine);
97
98 //Vector to hold TRT hits that are within RoI
99 std::vector<TRT_hit> hit;
100
101 //Sanity check of the ROI size
102 double deltaEta= std::abs(roiDescriptor->etaPlus()-roiDescriptor->etaMinus());
103 double deltaPhi= std::abs(CxxUtils::deltaPhi(roiDescriptor->phiPlus(),roiDescriptor->phiMinus()));
104
105 ATH_MSG_DEBUG( "roiDescriptor->etaPlus() in TrigTRTHTHCounter:"<<roiDescriptor->etaPlus());
106 ATH_MSG_DEBUG( "roiDescriptor->etaMinus() in TrigTRTHTHCounter:"<<roiDescriptor->etaMinus());
107 ATH_MSG_DEBUG( "deltaEta in TrigTRTHTHCounter:"<<deltaEta);
108
109 float phiTolerance = 0.001;
110 float etaTolerance = 0.001;
111
112 if((m_etaHalfWidth - deltaEta/2.) > etaTolerance)
113 ATH_MSG_WARNING ( "ROI eta range too small : " << deltaEta);
114
115 if((m_phiHalfWidth - deltaPhi/2.) > phiTolerance)
116 ATH_MSG_WARNING ( "ROI phi range too small : " << deltaPhi);
117
118 float coarseWedgeHalfWidth = m_phiHalfWidth/m_nBinCoarse;
119 float fineWedgeHalfWidth = coarseWedgeHalfWidth/m_nBinFine;
120
121 //Code will only proceed if the RoI eta is not too large; used to limit rate from endcap
122 if ( std::abs(roiDescriptor->eta())<=m_maxCaloEta ){
123
125 ATH_MSG_DEBUG( "Made handle " << m_trtDCContainerKey );
126
127 if (trtDCC->size() == 0){
128 return StatusCode::SUCCESS; // Exit early if there are no tracks
129 }
130
131 InDet::TRT_DriftCircleContainer::const_iterator trtdriftContainerItr = trtDCC->begin();
132 InDet::TRT_DriftCircleContainer::const_iterator trtdriftContainerItrE = trtDCC->end();
133
134 for (; trtdriftContainerItr != trtdriftContainerItrE; ++trtdriftContainerItr) {
135
136 InDet::TRT_DriftCircleCollection::const_iterator trtItr = (*trtdriftContainerItr)->begin();
137 InDet::TRT_DriftCircleCollection::const_iterator trtEnd = (*trtdriftContainerItr)->end();
138
139 for(; trtItr!=trtEnd; ++trtItr){
140
141 // find out which detector element the hit belongs to
142 const InDetDD::TRT_BaseElement *det = (*trtItr)->detectorElement();
143 Identifier ID = (*trtItr)->identify();
144 const Amg::Vector3D& strawcenter = det->strawCenter(m_trtHelper->straw(ID));
145
146 bool hth = false;
147 float hphi = strawcenter.phi();
148 float heta = strawcenter.eta();
149 float R = strawcenter.perp();
150
151 if ((*trtItr)->highLevel()) hth = true;
152
153 //hit needs to be stored
154 if(hth_eta_match((roiDescriptor)->eta(), heta, m_etaHalfWidth))
155 hit.push_back(make_hit(hphi,R,hth));
156
157
158 //First, define coarse wedges in phi, and count the TRT hits in these wedges
159 int countbin=0;
160 if(std::abs(CxxUtils::deltaPhi(hphi, static_cast<float>(roiDescriptor->phi()))) < 0.1){
161 float startValue = roiDescriptor->phi() - m_phiHalfWidth + coarseWedgeHalfWidth;
162 float endValue = roiDescriptor->phi() + m_phiHalfWidth;
163 float increment = 2*coarseWedgeHalfWidth;
164 for(float roibincenter = startValue; roibincenter < endValue; roibincenter += increment){
165 if (std::abs(CxxUtils::deltaPhi(hphi,roibincenter))<=coarseWedgeHalfWidth && hth_eta_match((roiDescriptor)->eta(), heta, m_etaHalfWidth)) {
166 if(hth) count_httrt_c.at(countbin) += 1.;
167 count_tottrt_c.at(countbin) += 1.;
168 break; //the hit has been assigned to one of the coarse wedges, so no need to continue the for loop
169 }
170 countbin++;
171 }
172 }
173 ATH_MSG_VERBOSE ( "timeOverThreshold=" << (*trtItr)->timeOverThreshold()
174 << " highLevel=" << (*trtItr)->highLevel()
175 << " rawDriftTime=" << (*trtItr)->rawDriftTime()
176 << " barrel_ec=" << m_trtHelper->barrel_ec(ID)
177 << " phi_module=" << m_trtHelper->phi_module(ID)
178 << " layer_or_wheel=" << m_trtHelper->layer_or_wheel(ID)
179 << " straw_layer=" << m_trtHelper->straw_layer(ID)
180 << " straw=" << m_trtHelper->straw(ID)
181 << " scR=" << det->strawCenter(m_trtHelper->straw(ID)).perp()
182 << " scPhi=" << hphi
183 << " scEta=" << heta);
184 } // end of loop over TRT drift circle collection
185 } //end of loop over TRT drift circle container
186 } //end of check to see if RoI eta is not too large
187
188 //Now figure out which of the coarse bins in phi has the max number of HT TRT hits
189 int maxHits = 0; //used to keep track of the max number of HT TRT hits in a coarse bin
190 int dist = 0; //used to keep track of which coarse bin has the max number of HT TRT hits
191
192 for (size_t iw = 0; iw < count_httrt_c.size(); iw++){
193 if(maxHits <= count_httrt_c[iw]){
194 maxHits = count_httrt_c[iw];
195 dist = iw;
196 }
197 }
198
199 //Value of dist can range between 0 and (m_nBinCoarse-1)
200 float center_pos_phi=roiDescriptor->phi()+(2*dist+1-m_nBinCoarse)*coarseWedgeHalfWidth;
201
202 //Now, define fine wedges in phi, centered around the best coarse wedge, and count the TRT hits in these fine wedges
203 for(size_t v=0;v<hit.size();v++){
204 int countbin=0;
205 if(std::abs(CxxUtils::deltaPhi(hit[v].phi, center_pos_phi)) < 0.01){
206 float startValue = center_pos_phi - 3*coarseWedgeHalfWidth + fineWedgeHalfWidth;
207 float endValue = center_pos_phi + 3*coarseWedgeHalfWidth;
208 float increment = 2*fineWedgeHalfWidth;
209 for(float roibincenter = startValue; roibincenter < endValue; roibincenter += increment){
210 if (std::abs(CxxUtils::deltaPhi(hit[v].phi,roibincenter))<=fineWedgeHalfWidth) {
211 if(hit[v].isHT) count_httrt.at(countbin) += 1.;
212 count_tottrt.at(countbin) += 1.;
213 break; //the hit has been assigned to one of the fine wedges, so no need to continue the for loop
214 }
215 countbin++;
216 }
217 }
218 }
219
220 //Now figure out which of the fine bins in phi has the max number of HT TRT hits
221 maxHits = 0; //used to keep track of the max number of HT TRT hits in a fine bin
222 dist = 0; //used to keep track of which fine bin has the max number of HT TRT hits
223
224 for (size_t iw = 0; iw < count_httrt.size(); iw++){
225 if(maxHits <= count_httrt[iw]){
226 maxHits = count_httrt[iw];
227 dist = iw;
228 }
229 }
230
231 //Value of dist can range between 0 and (3*m_nBinFine-1)
232 center_pos_phi+=(2*dist+1-3*m_nBinFine)*fineWedgeHalfWidth;
233
234 //Count the number of total and HT TRT hits for the road algorithm
235 int trthit=0, trthit_ht=0;
236 for(size_t v=0;v<hit.size();v++){
237 if (dist2COR(hit[v].R,hit[v].phi,center_pos_phi)<=m_roadWidth){
238 if(hit[v].isHT) trthit_ht+=1;
239 trthit+=1;
240 }
241 }
242
243 if (trthit!=0&&(std::abs(roiDescriptor->eta())<=m_roadMaxEta)){
244 trththits[0] = trthit_ht;
245 trththits[1] = (float)trthit_ht/trthit;
246 }
247
248 //Count the number of total and HT TRT hits for the wedge algorithm
249 trthit = 0;
250 trthit_ht = 0;
251
252 for (int k=0;k<(2*m_wedgeNBin+1);k++){
253 int binNumber = dist+k-m_wedgeNBin;
254 //Even though we are supposed to use 2*m_wedgeNBin+1 fine bins,
255 //need to make sure that binNumber is sensible
256 if(binNumber >= 0 && binNumber < (int)count_httrt.size()){
257 trthit += count_tottrt.at(binNumber);
258 trthit_ht += count_httrt.at(binNumber);
259 }
260 }
261
262 if (trthit!=0&&(std::abs(roiDescriptor->eta())>=m_wedgeMinEta)){
263 trththits[2] = trthit_ht;
264 trththits[3] = (float)trthit_ht/trthit;
265 }
266
267 trththits[4]=roiDescriptor->eta();
268 trththits[5]=roiDescriptor->phi();
269
270 ATH_MSG_VERBOSE ( "trthits with road algorithm : " << trththits[0]);
271 ATH_MSG_VERBOSE ( "fHT with road algorithm : " << trththits[1]);
272 ATH_MSG_VERBOSE ( "trthits with wedge algorithm : " << trththits[2]);
273 ATH_MSG_VERBOSE ( "fHT with wedge algorithm : " << trththits[3]);
274 ATH_MSG_VERBOSE ( "ROI eta : " << trththits[4]);
275
276 //Writing to xAOD
277 xAOD::TrigRNNOutput* rnnOutput = new xAOD::TrigRNNOutput();
278 trigRNNOutputColl->push_back(rnnOutput);
279 rnnOutput->setRnnDecision(trththits);
280
281 ATH_MSG_DEBUG("REGTEST: returning an xAOD::TrigRNNOutputContainer with size "<< trigRNNOutputColl->size() << ".");
282
283 return StatusCode::SUCCESS;
284}
285
286
Scalar eta() const
pseudorapidity method
Scalar deltaPhi(const MatrixBase< Derived > &vec) const
Scalar phi() const
phi method
#define ATH_CHECK
Evaluate an expression and check for errors.
#define ATH_MSG_INFO(x)
#define ATH_MSG_VERBOSE(x)
#define ATH_MSG_WARNING(x)
#define ATH_MSG_DEBUG(x)
Header file to be included by clients of the Monitored infrastructure.
float dist2COR(float R, float phi1, float phi2)
bool hth_eta_match(float caleta, float trteta, float etaWindow)
TRT_hit make_hit(float phi, float R, bool isHT)
const ServiceHandle< StoreGateSvc > & detStore() const
An algorithm that can be simultaneously executed in multiple threads.
Virtual base class of TRT readout elements.
virtual double etaMinus() const override final
gets eta at zMinus
virtual double etaPlus() const override final
gets eta at zedPlus
virtual double phi() const override final
Methods to retrieve data members.
virtual double phiMinus() const override final
gets phiMinus
virtual double eta() const override final
virtual double phiPlus() const override final
gets phiPlus
nope - should be used for standalone also, perhaps need to protect the class def bits ifndef XAOD_ANA...
virtual StatusCode initialize() override
Gaudi::Property< bool > m_doFullScan
Gaudi::Property< int > m_wedgeNBin
Gaudi::Property< float > m_etaHalfWidth
Gaudi::Property< int > m_nBinCoarse
Gaudi::Property< float > m_roadWidth
virtual StatusCode execute(const EventContext &ctx) const override
SG::ReadHandleKey< TrigRoiDescriptorCollection > m_roiCollectionKey
const TRT_ID * m_trtHelper
TRT ID helper.
Gaudi::Property< float > m_phiHalfWidth
SG::WriteHandleKey< xAOD::TrigRNNOutputContainer > m_trigRNNOutputKey
Gaudi::Property< float > m_roadMaxEta
Gaudi::Property< float > m_wedgeMinEta
SG::ReadHandleKey< InDet::TRT_DriftCircleContainer > m_trtDCContainerKey
Gaudi::Property< int > m_nBinFine
TrigTRTHTHCounter(const std::string &name, ISvcLocator *pSvcLocator)
void setRnnDecision(const std::vector< float > &d)
Eigen::Matrix< double, 3, 1 > Vector3D
T deltaPhi(T phiA, T phiB)
Return difference phiA - phiB in range [-pi, pi].
Definition phihelper.h:42
SG::ReadCondHandle< T > makeHandle(const SG::ReadCondHandleKey< T > &key, const EventContext &ctx=Gaudi::Hive::currentContext())
TrigRNNOutput_v2 TrigRNNOutput
Define the latest version of the RingerRings class.