Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AccumulatedSurfaceMaterial.cpp
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file AccumulatedSurfaceMaterial.cpp
1 // This file is part of the Acts project.
2 //
3 // Copyright (C) 2018-2020 CERN for the benefit of the Acts project
4 //
5 // This Source Code Form is subject to the terms of the Mozilla Public
6 // License, v. 2.0. If a copy of the MPL was not distributed with this
7 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
8 
10 
13 
14 #include <utility>
15 
16 // Default Constructor - for homogeneous material
18  : m_splitFactor(splitFactor) {
20  m_accumulatedMaterial = {{accMat}};
21 }
22 
23 // Binned Material constructor with split factor
25  const BinUtility& binUtility, double splitFactor)
26  : m_binUtility(binUtility), m_splitFactor(splitFactor) {
27  size_t bins0 = m_binUtility.bins(0);
28  size_t bins1 = m_binUtility.bins(1);
31 }
32 
33 // Assign a material properties object
35  const Vector2& lp, const MaterialSlab& mp, double pathCorrection) {
36  if (m_binUtility.dimensions() == 0) {
37  m_accumulatedMaterial[0][0].accumulate(mp, pathCorrection);
38  return {0, 0, 0};
39  }
40  size_t bin0 = m_binUtility.bin(lp, 0);
41  size_t bin1 = m_binUtility.bin(lp, 1);
42  m_accumulatedMaterial[bin1][bin0].accumulate(mp, pathCorrection);
43  return {bin0, bin1, 0};
44 }
45 
46 // Assign a material properties object
48  const Vector3& gp, const MaterialSlab& mp, double pathCorrection) {
49  if (m_binUtility.dimensions() == 0) {
50  m_accumulatedMaterial[0][0].accumulate(mp, pathCorrection);
51  return {0, 0, 0};
52  }
53  std::array<size_t, 3> bTriple = m_binUtility.binTriple(gp);
54  m_accumulatedMaterial[bTriple[1]][bTriple[0]].accumulate(mp, pathCorrection);
55  return bTriple;
56 }
57 
58 // Void average for vacuum assignment
60  MaterialSlab slabReference,
61  bool emptyHit) {
62  if (m_binUtility.dimensions() == 0) {
63  m_accumulatedMaterial[0][0].trackVariance(slabReference, emptyHit);
64  return;
65  }
66  std::array<size_t, 3> bTriple = m_binUtility.binTriple(gp);
67  std::vector<std::array<size_t, 3>> trackBins = {bTriple};
68  trackVariance(trackBins, slabReference);
69 }
70 
71 // Average the information accumulated during one event
73  const std::vector<std::array<size_t, 3>>& trackBins,
74  MaterialSlab slabReference, bool emptyHit) {
75  // the homogeneous material case
76  if (m_binUtility.dimensions() == 0) {
77  m_accumulatedMaterial[0][0].trackVariance(slabReference, emptyHit);
78  return;
79  }
80  // The touched bins are known, so you can access them directly
81  if (not trackBins.empty()) {
82  for (auto bin : trackBins) {
83  m_accumulatedMaterial[bin[1]][bin[0]].trackVariance(slabReference);
84  }
85  } else {
86  // Touched bins are not known: Run over all bins
87  for (auto& matVec : m_accumulatedMaterial) {
88  for (auto& mat : matVec) {
89  mat.trackVariance(slabReference);
90  }
91  }
92  }
93 }
94 
95 // Void average for vacuum assignment
97  bool emptyHit) {
98  if (m_binUtility.dimensions() == 0) {
99  m_accumulatedMaterial[0][0].trackAverage(emptyHit);
100  return;
101  }
102 
103  std::array<size_t, 3> bTriple = m_binUtility.binTriple(gp);
104  std::vector<std::array<size_t, 3>> trackBins = {bTriple};
105  trackAverage(trackBins, emptyHit);
106 }
107 
108 // Average the information accumulated during one event
110  const std::vector<std::array<size_t, 3>>& trackBins, bool emptyHit) {
111  // the homogeneous material case
112  if (m_binUtility.dimensions() == 0) {
113  m_accumulatedMaterial[0][0].trackAverage(emptyHit);
114  return;
115  }
116 
117  // The touched bins are known, so you can access them directly
118  if (not trackBins.empty()) {
119  for (auto bin : trackBins) {
120  m_accumulatedMaterial[bin[1]][bin[0]].trackAverage(emptyHit);
121  }
122  } else {
123  // Touched bins are not known: Run over all bins
124  for (auto& matVec : m_accumulatedMaterial) {
125  for (auto& mat : matVec) {
126  mat.trackAverage(emptyHit);
127  }
128  }
129  }
130 }
131 
133 std::unique_ptr<const Acts::ISurfaceMaterial>
135  if (m_binUtility.bins() == 1) {
136  // Return HomogeneousSurfaceMaterial
137  return std::make_unique<HomogeneousSurfaceMaterial>(
138  m_accumulatedMaterial[0][0].totalAverage().first, m_splitFactor);
139  }
140  // Create the properties matrix
141  MaterialSlabMatrix mpMatrix(
142  m_binUtility.bins(1),
143  MaterialSlabVector(m_binUtility.bins(0), MaterialSlab()));
144  // Loop over and fill
145  for (size_t ib1 = 0; ib1 < m_binUtility.bins(1); ++ib1) {
146  for (size_t ib0 = 0; ib0 < m_binUtility.bins(0); ++ib0) {
147  mpMatrix[ib1][ib0] = m_accumulatedMaterial[ib1][ib0].totalAverage().first;
148  }
149  }
150  // Now return the BinnedSurfaceMaterial
151  return std::make_unique<const BinnedSurfaceMaterial>(
152  m_binUtility, std::move(mpMatrix), m_splitFactor);
153 }