ATLAS Offline Software
FPGATrackSimSectorSlice.cxx
Go to the documentation of this file.
1 // Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
2 
11 
12 
13 #include <TFile.h>
14 #include <TTree.h>
15 #include <TBits.h>
16 #include <TClonesArray.h>
18 
19 
20 using namespace asg::msgUserCode;
21 
22 
23 void largest_region_wrap(std::vector<bool> const & good_bin, int & i_lower, int & i_upper);
24 
25 namespace{
26  // Max is exclusive
27  int clamp(int bin, int min, int max)
28  {
29  if (bin < min) return min;
30  if (bin >= max) return max - 1;
31  return bin;
32  }
33 }
34 
36 // Constructor/Initialization
37 
38 
41  :
42  m_nBins(nBins),
43  m_min(min),
44  m_max(max)
45 {
46  m_bits_phi = new TClonesArray("TBits", m_nBins.phi);
47  m_bits_c = new TClonesArray("TBits", m_nBins.qOverPt);
48  m_bits_d0 = new TClonesArray("TBits", m_nBins.d0);
49  m_bits_z0 = new TClonesArray("TBits", m_nBins.z0);
50  m_bits_eta = new TClonesArray("TBits", m_nBins.eta);
51 
52  for (int i = 0; i < m_nBins.phi; i++) new ((*m_bits_phi)[i]) TBits(nSectors);
53  for (int i = 0; i < m_nBins.qOverPt; i++) new ((*m_bits_c)[i]) TBits(nSectors);
54  for (int i = 0; i < m_nBins.d0; i++) new ((*m_bits_d0)[i]) TBits(nSectors);
55  for (int i = 0; i < m_nBins.z0; i++) new ((*m_bits_z0)[i]) TBits(nSectors);
56  for (int i = 0; i < m_nBins.eta; i++) new ((*m_bits_eta)[i]) TBits(nSectors);
57 
59 }
60 
61 
63 {
64  ANA_MSG_INFO("Reading " << filepath);
65  TFile *slice_file = TFile::Open(filepath.c_str());
66 
67  m_bits_phi = dynamic_cast<TClonesArray*>(slice_file->Get("c_bits_phi"));
68  m_bits_c = dynamic_cast<TClonesArray*>(slice_file->Get("c_bits_c"));
69  m_bits_d0 = dynamic_cast<TClonesArray*>(slice_file->Get("c_bits_d"));
70  m_bits_z0 = dynamic_cast<TClonesArray*>(slice_file->Get("c_bits_z0"));
71 
72  // kludge for backward compatability
73  if ((slice_file->Get("c_bits_eta")) != nullptr)
74  m_bits_eta = dynamic_cast<TClonesArray*>(slice_file->Get("c_bits_eta"));
75  else
76  m_bits_eta = dynamic_cast<TClonesArray*>(slice_file->Get("c_bits_ctheta"));
77 
78  if (not m_bits_eta){
79  ANA_MSG_ERROR("m_bits_eta is null in FPGATrackSimSectorSlice constructor");
80  delete slice_file;
81  return;
82  }
83  TTree *slice_tree = dynamic_cast<TTree*>(slice_file->Get("bin_info"));
84  if (not slice_tree){
85  ANA_MSG_ERROR("slice_tree is null in FPGATrackSimSectorSlice constructor");
86  delete slice_file;
87  return;
88  }
89  slice_tree->SetBranchAddress("qOverPt_max", &m_max.qOverPt);
90  slice_tree->SetBranchAddress("d0_max", &m_max.d0);
91  slice_tree->SetBranchAddress("phi_max", &m_max.phi);
92  slice_tree->SetBranchAddress("z0_max", &m_max.z0);
93  slice_tree->SetBranchAddress("eta_max", &m_max.eta);
94 
95  slice_tree->SetBranchAddress("qOverPt_min", &m_min.qOverPt);
96  slice_tree->SetBranchAddress("d0_min", &m_min.d0);
97  slice_tree->SetBranchAddress("phi_min", &m_min.phi);
98  slice_tree->SetBranchAddress("z0_min", &m_min.z0);
99  slice_tree->SetBranchAddress("eta_min", &m_min.eta);
100 
101  slice_tree->SetBranchAddress("qOverPt_bins",&m_nBins.qOverPt);
102  slice_tree->SetBranchAddress("d0_bins", &m_nBins.d0);
103  slice_tree->SetBranchAddress("phi_bins", &m_nBins.phi);
104  slice_tree->SetBranchAddress("z0_bins", &m_nBins.z0);
105  slice_tree->SetBranchAddress("eta_bins", &m_nBins.eta);
106 
107  // For backwards compatability
108  if (slice_tree->GetBranch("qOverPt_max") == nullptr) slice_tree->SetBranchAddress("halfInvPt_max", &m_max.qOverPt);
109  if (slice_tree->GetBranch("qOverPt_min") == nullptr) slice_tree->SetBranchAddress("halfInvPt_min", &m_min.qOverPt);
110  if (slice_tree->GetBranch("qOverPt_bins") == nullptr) slice_tree->SetBranchAddress("halfInvPt_bins", &m_nBins.qOverPt);
111 
112  slice_tree->GetEntry(0);
113 
115 
116  delete slice_file;
117 }
118 
119 
120 
122 {
123 
126  m_step.d0 = (m_max.d0 - m_min.d0) / m_nBins.d0;
127  m_step.z0 = (m_max.z0 - m_min.z0) /m_nBins.z0;
129 }
130 
131 
133 // Slice writing functions
134 
136 {
137  int bin;
138 
139  bin = clamp(bins.phi, 0, m_nBins.phi);
140  static_cast<TBits*>(m_bits_phi->UncheckedAt(bin))->SetBitNumber(sector);
141 
142  bin = clamp(bins.qOverPt, 0, m_nBins.qOverPt);
143  static_cast<TBits*>(m_bits_c->UncheckedAt(bin))->SetBitNumber(sector);
144 
145  bin = clamp(bins.d0, 0, m_nBins.d0);
146  static_cast<TBits*>(m_bits_d0->UncheckedAt(bin))->SetBitNumber(sector);
147 
148  bin = clamp(bins.z0, 0, m_nBins.z0);
149  static_cast<TBits*>(m_bits_z0->UncheckedAt(bin))->SetBitNumber(sector);
150 
151  bin = clamp(bins.eta, 0, m_nBins.eta);
152  static_cast<TBits*>(m_bits_eta->UncheckedAt(bin))->SetBitNumber(sector);
153 }
154 
155 
156 void FPGATrackSimSectorSlice::saveSlices(const std::string & filepath)
157 {
158  TFile ofile(filepath.c_str(), "recreate");
159 
160  m_bits_phi->Write("c_bits_phi", TObject::kSingleKey);
161  m_bits_c->Write("c_bits_c", TObject::kSingleKey);
162  m_bits_d0->Write("c_bits_d", TObject::kSingleKey);
163  m_bits_z0->Write("c_bits_z0", TObject::kSingleKey);
164  m_bits_eta->Write("c_bits_eta", TObject::kSingleKey);
165 
166  TTree *slice_tree = new TTree("bin_info", "Slice binning information");
167 
168  slice_tree->Branch("qOverPt_max", &m_max.qOverPt);
169  slice_tree->Branch("d0_max", &m_max.d0);
170  slice_tree->Branch("phi_max", &m_max.phi);
171  slice_tree->Branch("z0_max", &m_max.z0);
172  slice_tree->Branch("eta_max", &m_max.eta);
173 
174  slice_tree->Branch("qOverPt_min", &m_min.qOverPt);
175  slice_tree->Branch("d0_min", &m_min.d0);
176  slice_tree->Branch("phi_min", &m_min.phi);
177  slice_tree->Branch("z0_min", &m_min.z0);
178  slice_tree->Branch("eta_min", &m_min.eta);
179 
180  slice_tree->Branch("qOverPt_bins",&m_nBins.qOverPt);
181  slice_tree->Branch("d0_bins", &m_nBins.d0);
182  slice_tree->Branch("phi_bins", &m_nBins.phi);
183  slice_tree->Branch("z0_bins", &m_nBins.z0);
184  slice_tree->Branch("eta_bins", &m_nBins.eta);
185 
186  slice_tree->Fill();
187 
188  ofile.Write();
189  ofile.Close();
190 }
191 
192 
194 // searchSectors and helpers
195 
196 
198 {
199  if (pars.phi < m_min.phi || pars.phi > m_max.phi) return false;
200  if (pars.qOverPt < m_min.qOverPt || pars.qOverPt > m_max.qOverPt) return false;
201  if (pars.d0 < m_min.d0 || pars.d0 > m_max.d0) return false;
202  if (pars.z0 < m_min.z0 || pars.z0 > m_max.z0) return false;
203  if (pars.eta < m_min.eta || pars.eta > m_max.eta) return false;
204  return true;
205 }
206 
207 
208 // If the bin is underrun or overrun set to 0 or clamp - 1
209 int getBin(double x, double min, double step, int clamp_max)
210 {
211  int bin = static_cast<int>((x - min) / step);
212  return clamp(bin, 0, clamp_max);
213 }
214 
215 
217 {
218  std::vector<sector_t> sectors;
219  if (!checkTrackPars(pars))
220  {
221  ANA_MSG_ERROR("searchSectors() Bad track parameters: " << pars);
222  return sectors;
223  }
224 
225  // Special case for phi
226  int bin = static_cast<int>((pars.phi - m_min.phi) / m_step.phi);
227  bin = (bin + m_nBins.phi) % m_nBins.phi; // modulo operator: phi wraps around
228  TBits result_bits = *(TBits*)(m_bits_phi->UncheckedAt(bin));
229 
230  result_bits &= *(TBits*)(m_bits_c->UncheckedAt(getBin(pars.qOverPt, m_min.qOverPt, m_step.qOverPt, m_nBins.qOverPt)));
231  result_bits &= *(TBits*)(m_bits_d0->UncheckedAt(getBin(pars.d0, m_min.d0, m_step.d0, m_nBins.d0)));
232  result_bits &= *(TBits*)(m_bits_z0->UncheckedAt(getBin(pars.z0, m_min.z0, m_step.z0, m_nBins.z0)));
233  result_bits &= *(TBits*)(m_bits_eta->UncheckedAt(getBin(pars.eta, m_min.eta, m_step.eta, m_nBins.eta)));
234 
235  unsigned int curPos = result_bits.FirstSetBit(0);
236  while (curPos != result_bits.GetNbits())
237  {
238  sectors.push_back(static_cast<sector_t>(curPos));
239  curPos = result_bits.FirstSetBit(curPos + 1);
240  }
241 
242  return sectors;
243 }
244 
245 
246 
248 // getBoundaries and helpers
249 
250 // Calculates the parameter boundaries for where bitmasks are non-empty.
251 // Returns (min, max).
252 std::pair<FPGATrackSimTrackPars, FPGATrackSimTrackPars> FPGATrackSimSectorSlice::getBoundaries() const
253 {
256 
257  getBoundary(m_bits_phi, m_min.phi, m_max.phi, min.phi, max.phi, true, "phi");
258  getBoundary(m_bits_c, m_min.qOverPt, m_max.qOverPt, min.qOverPt, max.qOverPt, false, "c");
259  getBoundary(m_bits_d0, m_min.d0, m_max.d0, min.d0, max.d0, false, "d0");
260  getBoundary(m_bits_z0, m_min.z0, m_max.z0, min.z0, max.z0, false, "z0");
261  getBoundary(m_bits_eta, m_min.eta, m_max.eta, min.eta, max.eta, false, "eta");
262 
263  ANA_MSG_INFO("getBoundaries() returning " << min << "; " << max);
264  return { min, max };
265 }
266 
267 
276 void FPGATrackSimSectorSlice::getBoundary(const TClonesArray *bitmasks, double x_min, double x_max,
277  double &bound_min, double &bound_max, bool wraps, const char *debug) const
278 {
279  int nbin = bitmasks->GetEntries();
280 
281  // Find the bounds of where there are bits, i_lower and i_upper are inclusive.
282  int i_lower, i_upper;
283  if (!wraps)
284  {
285  for (i_lower = 0; i_lower < nbin-1; i_lower++)
286  if (dynamic_cast<TBits const *>(bitmasks->UncheckedAt(i_lower))->CountBits() > 0) break;
287  for (i_upper = nbin-1; i_upper > 0; i_upper--)
288  if (dynamic_cast<TBits const *>(bitmasks->UncheckedAt(i_upper))->CountBits() > 0) break;
289  }
290  else // special case for cyclic values that wrap around, like phi
291  {
292  // Mark bins with nonzero bits, with a safety margin
293  std::vector<bool> good_bins(nbin);
294  for (int j = 0; j < nbin; j++)
295  {
296  if (dynamic_cast<TBits const *>(bitmasks->UncheckedAt(j))->CountBits() > 0)
297  {
298  // add safety margin: +/- 1 bin
299  // if the good bins are not actually continuous, the region found below will leave out some bins
300  good_bins[(j+nbin-1) % nbin] = true; // mod to handle wrap-around
301  good_bins[j] = true;
302  good_bins[(j+1) % nbin] = true;
303  }
304  }
305 
306  // Locate largest continuous region
307  largest_region_wrap(good_bins, i_lower, i_upper);
308  }
309 
310  if (i_lower >= i_upper) ANA_MSG_WARNING("getBoundary() Completely empty bitmask array");
311  bound_min = x_min + ((x_max-x_min) * i_lower) / nbin;
312  bound_max = x_min + ((x_max-x_min) * (i_upper+1)) / nbin;
313 
314  ANA_MSG_INFO("getBoundary() " << debug << " nbin=" << nbin
315  << " range=[" << x_min << ", " << x_max
316  << "], index=[" << i_lower << ", " << i_upper
317  << "], result=[" << bound_min << ", " << bound_max << "]");
318 }
319 
320 
321 void largest_region_wrap(std::vector<bool> const & good_bins, int & i_lower, int & i_upper)
322 {
323  int start = -1; // start of current region
324  int best_start = -1; // start of best region
325  int best_size = 0; // size of best region
326  for (int j = 0; j < 2*(int)good_bins.size(); j++) // j < 2*size to get regions across j=0
327  {
328  if (good_bins[j % good_bins.size()])
329  {
330  if (start < 0) start = j; // start region
331  int size = j + 1 - start;
332  if (size > best_size)
333  {
334  best_start = start;
335  best_size = size;
336  if (best_size == (int)good_bins.size()) break; // All bins are filled
337  }
338  }
339  else // end region
340  {
341  start = -1;
342  }
343  }
344  i_lower = best_start;
345  i_upper = best_start + best_size - 1; // inclusive
346  // May be larger than nbin, which causes max phi to be larger than 2pi, e.g.
347 }
348 
349 
make_hlt_rep.pars
pars
Definition: make_hlt_rep.py:90
FPGATrackSimTrackPars::phi
double phi
Definition: FPGATrackSimTrackPars.h:24
max
#define max(a, b)
Definition: cfImp.cxx:41
CaloCellPos2Ntuple.int
int
Definition: CaloCellPos2Ntuple.py:24
FPGATrackSimTrackParsI::eta
int eta
Definition: FPGATrackSimTrackPars.h:62
python.App.bins
bins
Definition: App.py:410
FPGATrackSimTrackPars
Definition: FPGATrackSimTrackPars.h:22
FPGATrackSimTrackPars::qOverPt
double qOverPt
Definition: FPGATrackSimTrackPars.h:25
FPGATrackSimSectorSlice::FPGATrackSimSectorSlice
FPGATrackSimSectorSlice(size_t nSectors, FPGATrackSimTrackParsI const &nBins, FPGATrackSimTrackPars const &min, FPGATrackSimTrackPars const &max)
Definition: FPGATrackSimSectorSlice.cxx:39
FPGATrackSimSectorSlice::getBoundary
void getBoundary(const TClonesArray *bitmasks, double x_min, double x_max, double &autoMin, double &autoMax, bool wraps, const char *debug) const
Given a range [x_min, x_max] split in bins corresponding to the size of bitmasks, finds a tighter ran...
Definition: FPGATrackSimSectorSlice.cxx:276
FPGATrackSimSectorSlice::calcDependentVals
void calcDependentVals()
Definition: FPGATrackSimSectorSlice.cxx:121
mergePhysValFiles.start
start
Definition: DataQuality/DataQualityUtils/scripts/mergePhysValFiles.py:14
FPGATrackSimSectorSlice::m_bits_d0
TClonesArray * m_bits_d0
Definition: FPGATrackSimSectorSlice.h:63
FPGATrackSimSectorSlice::m_bits_c
TClonesArray * m_bits_c
Definition: FPGATrackSimSectorSlice.h:62
validation.ofile
ofile
Definition: validation.py:96
FPGATrackSimTrackPars::d0
double d0
Definition: FPGATrackSimTrackPars.h:26
FPGATrackSimSectorSlice.h
Stores the range of eta/phi/etc. of each sector.
FPGATrackSimSectorSlice::m_step
FPGATrackSimTrackPars m_step
Definition: FPGATrackSimSectorSlice.h:75
bin
Definition: BinsDiffFromStripMedian.h:43
FPGATrackSimSectorSlice::m_nBins
FPGATrackSimTrackParsI m_nBins
Definition: FPGATrackSimSectorSlice.h:68
getBin
int getBin(double x, double min, double step, int clamp_max)
Definition: FPGATrackSimSectorSlice.cxx:209
ANA_MSG_ERROR
#define ANA_MSG_ERROR(xmsg)
Macro printing error messages.
Definition: Control/AthToolSupport/AsgMessaging/AsgMessaging/MessageCheck.h:294
FPGATrackSimTrackParsI::d0
int d0
Definition: FPGATrackSimTrackPars.h:60
FPGATrackSimSectorSlice::getBoundaries
std::pair< FPGATrackSimTrackPars, FPGATrackSimTrackPars > getBoundaries() const
Definition: FPGATrackSimSectorSlice.cxx:252
FPGATrackSimSectorSlice::m_min
FPGATrackSimTrackPars m_min
Definition: FPGATrackSimSectorSlice.h:71
FPGATrackSimTrackPars::eta
double eta
Definition: FPGATrackSimTrackPars.h:28
x
#define x
python.setupRTTAlg.size
int size
Definition: setupRTTAlg.py:39
lumiFormat.i
int i
Definition: lumiFormat.py:92
FPGATrackSimSectorSlice::checkTrackPars
bool checkTrackPars(FPGATrackSimTrackPars const &pars) const
Definition: FPGATrackSimSectorSlice.cxx:197
ANA_MSG_INFO
#define ANA_MSG_INFO(xmsg)
Macro printing info messages.
Definition: Control/AthToolSupport/AsgMessaging/AsgMessaging/MessageCheck.h:290
MessageCheck.h
macros for messaging and checking status codes
ANA_MSG_WARNING
#define ANA_MSG_WARNING(xmsg)
Macro printing warning messages.
Definition: Control/AthToolSupport/AsgMessaging/AsgMessaging/MessageCheck.h:292
FPGATrackSimSectorSlice::m_max
FPGATrackSimTrackPars m_max
Definition: FPGATrackSimSectorSlice.h:72
min
#define min(a, b)
Definition: cfImp.cxx:40
sector_t
int32_t sector_t
Definition: FPGATrackSimTypes.h:21
dumpTgcDigiJitter.nBins
list nBins
Definition: dumpTgcDigiJitter.py:29
FPGATrackSimSectorSlice::m_bits_eta
TClonesArray * m_bits_eta
Definition: FPGATrackSimSectorSlice.h:65
debug
const bool debug
Definition: MakeUncertaintyPlots.cxx:53
plotBeamSpotVxVal.bin
int bin
Definition: plotBeamSpotVxVal.py:83
FPGATrackSimSectorSlice::m_bits_phi
TClonesArray * m_bits_phi
Definition: FPGATrackSimSectorSlice.h:61
FPGATrackSimTrackParsI::phi
int phi
Definition: FPGATrackSimTrackPars.h:58
FPGATrackSimTrackPars::z0
double z0
Definition: FPGATrackSimTrackPars.h:27
FPGATrackSimTrackParsI::z0
int z0
Definition: FPGATrackSimTrackPars.h:61
LArCellBinning.step
step
Definition: LArCellBinning.py:158
largest_region_wrap
void largest_region_wrap(std::vector< bool > const &good_bin, int &i_lower, int &i_upper)
Definition: FPGATrackSimSectorSlice.cxx:321
FPGATrackSimSectorSlice::addSectorToSlice
void addSectorToSlice(sector_t sector, FPGATrackSimTrackParsI const &bins)
Definition: FPGATrackSimSectorSlice.cxx:135
FPGATrackSimSectorSlice::searchSectors
std::vector< sector_t > searchSectors(FPGATrackSimTrackPars const &pars) const
Definition: FPGATrackSimSectorSlice.cxx:216
FPGATrackSimSectorSlice::m_bits_z0
TClonesArray * m_bits_z0
Definition: FPGATrackSimSectorSlice.h:64
FPGATrackSimTrackParsI::qOverPt
int qOverPt
Definition: FPGATrackSimTrackPars.h:59
FPGATrackSimTrackParsI
Definition: FPGATrackSimTrackPars.h:56
FPGATrackSimSectorSlice::saveSlices
void saveSlices(const std::string &filepath)
Definition: FPGATrackSimSectorSlice.cxx:156