src/Material/BinnedSurfaceMaterialAccumulater.cpp

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 CERN for the benefit of the ACTS project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at https://mozilla.org/MPL/2.0/.
0008
0009 #include "Acts/Material/BinnedSurfaceMaterialAccumulater.hpp"
0010
0011 #include "Acts/Geometry/GeometryContext.hpp"
0012 #include "Acts/Material/BinnedSurfaceMaterial.hpp"
0013 #include "Acts/Material/ProtoSurfaceMaterial.hpp"
0014 #include "Acts/Surfaces/Surface.hpp"
0015 #include "Acts/Utilities/BinAdjustment.hpp"
0016 #include "Acts/Utilities/BinUtility.hpp"
0017
0018 Acts::BinnedSurfaceMaterialAccumulater::BinnedSurfaceMaterialAccumulater(
0019     const Config& cfg, std::unique_ptr<const Logger> mlogger)
0020     : m_cfg(cfg), m_logger(std::move(mlogger)) {}
0021
0022 std::unique_ptr<Acts::ISurfaceMaterialAccumulater::State>
0023 Acts::BinnedSurfaceMaterialAccumulater::createState() const {
0024   auto state = std::make_unique<State>();
0025
0026   /// Create the surface accumulation
0027   for (const auto& surface : m_cfg.materialSurfaces) {
0028     GeometryIdentifier geoID = surface->geometryId();
0029
0030     // Get the Surface Material
0031     const ISurfaceMaterial* surfaceMaterial = surface->surfaceMaterial();
0032     if (surfaceMaterial == nullptr) {
0033       throw std::invalid_argument(
0034           "Surface material is not set, inconsistent configuration.");
0035     }
0036
0037     // First attempt from ProtoSurfaceMaterial
0038     auto psm = dynamic_cast<const ProtoSurfaceMaterial*>(surfaceMaterial);
0039     if (psm != nullptr) {
0040       auto binUtility = psm->binning();
0041       // Screen output for Binned Surface material
0042       ACTS_DEBUG("       - (proto) binning from ProtoSurfaceMateria is "
0043                  << binUtility);
0044       // Now adjust to surface type
0045       binUtility = adjustBinUtility(binUtility, *surface, m_cfg.geoContext);
0046       // Screen output for Binned Surface material
0047       ACTS_DEBUG("       - adjusted binning is " << binUtility);
0048       state->accumulatedMaterial[geoID] =
0049           AccumulatedSurfaceMaterial(binUtility);
0050       // Material accumulation  is created for this
0051       continue;
0052     }
0053     // Second attempt from ProtoGridSurfaceMaterial
0054     auto psgm = dynamic_cast<const ProtoGridSurfaceMaterial*>(surfaceMaterial);
0055     if (psgm != nullptr) {
0056       auto binUtility = psgm->binning().toBinUtility();
0057       // Screen output for Binned Surface material
0058       ACTS_DEBUG("       - (proto) binning from ProtoGridSurfaceMaterial is "
0059                  << binUtility);
0060       // Now adjust to surface type
0061       binUtility = adjustBinUtility(binUtility, *surface, m_cfg.geoContext);
0062       // Screen output for Binned Surface material
0063       ACTS_DEBUG("       - adjusted binning is " << binUtility);
0064       state->accumulatedMaterial[geoID] =
0065           AccumulatedSurfaceMaterial(binUtility);
0066       // Material accumulation  is created for this
0067       continue;
0068     }
0069     // Third attempt: binned material
0070     auto bmp = dynamic_cast<const BinnedSurfaceMaterial*>(surfaceMaterial);
0071     if (bmp != nullptr) {
0072       // Screen output for Binned Surface material
0073       ACTS_DEBUG("       - binning from BinnedSurfaceMaterial is "
0074                  << bmp->binUtility());
0075       state->accumulatedMaterial[geoID] =
0076           AccumulatedSurfaceMaterial(bmp->binUtility());
0077       // Material accumulation  is created for this
0078       continue;
0079     }
0080     // Create a homogeneous type of material
0081     ACTS_DEBUG("       - this is homogeneous material.");
0082     state->accumulatedMaterial[geoID] = AccumulatedSurfaceMaterial();
0083   }
0084   return state;
0085 }
0086
0087 void Acts::BinnedSurfaceMaterialAccumulater::accumulate(
0088     ISurfaceMaterialAccumulater::State& state,
0089     const std::vector<MaterialInteraction>& interactions,
0090     const std::vector<IAssignmentFinder::SurfaceAssignment>&
0091         surfacesWithoutAssignment) const {
0092   // Cast into the right state object (guaranteed by upstream algorithm)
0093   State* cState = static_cast<State*>(&state);
0094   if (cState == nullptr) {
0095     throw std::invalid_argument(
0096         "Invalid state object provided, something is seriously wrong.");
0097   }
0098
0099   using MapBin =
0100       std::pair<AccumulatedSurfaceMaterial*, std::array<std::size_t, 3>>;
0101   std::map<AccumulatedSurfaceMaterial*, std::array<std::size_t, 3>>
0102       touchedMapBins;
0103
0104   // Assign the hits
0105   for (const auto& mi : interactions) {
0106     // Get the surface
0107     const Surface* surface = mi.surface;
0108     GeometryIdentifier geoID = surface->geometryId();
0109     // Get the accumulated material
0110     auto accMaterial = cState->accumulatedMaterial.find(geoID);
0111     if (accMaterial == cState->accumulatedMaterial.end()) {
0112       throw std::invalid_argument(
0113           "Surface material is not found, inconsistent configuration.");
0114     }
0115     // Accumulate the material - remember the touched bin
0116     auto tBin = accMaterial->second.accumulate(mi.intersection, mi.materialSlab,
0117                                                mi.pathCorrection);
0118     touchedMapBins.insert(MapBin(&(accMaterial->second), tBin));
0119   }
0120
0121   // After mapping this track, average the touched bins
0122   for (const auto& [key, value] : touchedMapBins) {
0123     std::vector<std::array<std::size_t, 3>> trackBins = {value};
0124     key->trackAverage(trackBins, true);
0125   }
0126
0127   // Empty bin correction
0128   if (m_cfg.emptyBinCorrection) {
0129     for (const auto& [surface, position, direction] :
0130          surfacesWithoutAssignment) {
0131       // Get the accumulated material
0132       auto missedMaterial =
0133           cState->accumulatedMaterial.find(surface->geometryId());
0134       if (missedMaterial == cState->accumulatedMaterial.end()) {
0135         throw std::invalid_argument(
0136             "Surface material is not found, inconsistent configuration.");
0137       }
0138       // Apply empty hit correction
0139       missedMaterial->second.trackAverage(position, true);
0140     }
0141   }
0142 }
0143
0144 std::map<Acts::GeometryIdentifier,
0145          std::shared_ptr<const Acts::ISurfaceMaterial>>
0146 Acts::BinnedSurfaceMaterialAccumulater::finalizeMaterial(
0147     ISurfaceMaterialAccumulater::State& state) const {
0148   std::map<GeometryIdentifier, std::shared_ptr<const ISurfaceMaterial>>
0149       sMaterials;
0150
0151   // Cast into the right state object (guaranteed by upstream algorithm)
0152   State* cState = static_cast<State*>(&state);
0153   if (cState == nullptr) {
0154     throw std::invalid_argument(
0155         "Invalid state object provided, something is seriously wrong.");
0156   }
0157
0158   // iterate over the map to call the total average
0159   for (auto& accMaterial : cState->accumulatedMaterial) {
0160     ACTS_DEBUG("Finalizing map for Surface " << accMaterial.first);
0161     auto sMaterial = accMaterial.second.totalAverage();
0162     sMaterials[accMaterial.first] = std::move(sMaterial);
0163   }
0164
0165   return sMaterials;
0166 }