gen/detail/MatScintSpecInserter.hh

0001 //------------------------------- -*- C++ -*- -------------------------------//
0002 // Copyright Celeritas contributors: see top-level COPYRIGHT file for details
0003 // SPDX-License-Identifier: (Apache-2.0 OR MIT)
0004 //---------------------------------------------------------------------------//
0005 //! \file celeritas/optical/gen/detail/MatScintSpecInserter.hh
0006 //---------------------------------------------------------------------------//
0007 #pragma once
0008
0009 #include "corecel/Assert.hh"
0010 #include "corecel/Types.hh"
0011 #include "corecel/data/CollectionBuilder.hh"
0012 #include "corecel/data/DedupeCollectionBuilder.hh"
0013 #include "corecel/math/SoftEqual.hh"
0014 #include "celeritas/io/ImportOpticalMaterial.hh"
0015
0016 #include "../ScintillationData.hh"
0017
0018 namespace celeritas
0019 {
0020 namespace detail
0021 {
0022 //---------------------------------------------------------------------------//
0023 /*!
0024  * Build scintillation spectrum data.
0025  */
0026 class MatScintSpecInserter
0027 {
0028   public:
0029     //!@{
0030     //! \name Type aliases
0031     using Data = HostVal<ScintillationData>;
0032     //!@}
0033
0034   public:
0035     // Construct with data to insert into
0036     explicit inline MatScintSpecInserter(Data* data);
0037
0038     // Add scintillation data for a single material
0039     auto operator()(ImportMaterialScintSpectrum const& mat);
0040
0041   private:
0042     using MatId = OptMatId;
0043
0044     CollectionBuilder<MatScintSpectrum, MemSpace::host, MatId> materials_;
0045     DedupeCollectionBuilder<real_type> reals_;
0046     CollectionBuilder<ScintRecord> scint_records_;
0047 };
0048
0049 //---------------------------------------------------------------------------//
0050 // INLINE DEFINITIONS
0051 //---------------------------------------------------------------------------//
0052 /*!
0053  * Construct with defaults.
0054  */
0055 MatScintSpecInserter::MatScintSpecInserter(Data* data)
0056     : materials_{&data->materials}
0057     , reals_{&data->reals}
0058     , scint_records_{&data->scint_records}
0059 {
0060     CELER_EXPECT(data);
0061 }
0062
0063 //---------------------------------------------------------------------------//
0064 /*!
0065  * Add scintillation data for a single material.
0066  */
0067 auto MatScintSpecInserter::operator()(ImportMaterialScintSpectrum const& mat)
0068 {
0069     CELER_EXPECT(mat);
0070
0071     CELER_VALIDATE(mat.yield_per_energy > 0,
0072                    << "invalid yield=" << mat.yield_per_energy
0073                    << " for scintillation (should be positive)");
0074
0075     double total_yield{0};
0076     std::vector<double> yield_pdf;
0077     auto const begin_components = scint_records_.size_id();
0078     for (ImportScintComponent const& comp : mat.components)
0079     {
0080         CELER_VALIDATE(comp.lambda_mean > 0,
0081                        << "invalid lambda_mean=" << comp.lambda_mean
0082                        << " for scintillation component (should be positive)");
0083         CELER_VALIDATE(comp.lambda_sigma > 0,
0084                        << "invalid lambda_sigma=" << comp.lambda_sigma
0085                        << " (should be positive)");
0086         CELER_VALIDATE(comp.rise_time >= 0,
0087                        << "invalid rise_time=" << comp.rise_time
0088                        << " (should be nonnegative)");
0089         CELER_VALIDATE(comp.fall_time > 0,
0090                        << "invalid fall_time=" << comp.fall_time
0091                        << " (should be positive)");
0092         CELER_VALIDATE(comp.yield_frac > 0,
0093                        << "invalid yield=" << comp.yield_frac);
0094
0095         ScintRecord scint;
0096         scint.lambda_mean = comp.lambda_mean;
0097         scint.lambda_sigma = comp.lambda_sigma;
0098         scint.rise_time = comp.rise_time;
0099         scint.fall_time = comp.fall_time;
0100         scint_records_.push_back(scint);
0101
0102         yield_pdf.push_back(comp.yield_frac);
0103         total_yield += comp.yield_frac;
0104     }
0105
0106     // Normalize yield PDF by total yield
0107     for (auto& y : yield_pdf)
0108     {
0109         y /= total_yield;
0110     }
0111
0112     MatScintSpectrum spectrum;
0113     spectrum.yield_per_energy = mat.yield_per_energy;
0114     spectrum.components = {begin_components, scint_records_.size_id()};
0115     spectrum.yield_pdf = reals_.insert_back(yield_pdf.begin(), yield_pdf.end());
0116
0117     CELER_ENSURE(spectrum.components.size() == mat.components.size());
0118     return materials_.push_back(std::move(spectrum));
0119 }
0120
0121 //---------------------------------------------------------------------------//
0122 }  // namespace detail
0123 }  // namespace celeritas