optical/gen/WavelengthShiftGenerator.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/WavelengthShiftGenerator.hh
0006 //---------------------------------------------------------------------------//
0007 #pragma once
0008
0009 #include "corecel/Macros.hh"
0010 #include "corecel/Types.hh"
0011 #include "corecel/data/StackAllocator.hh"
0012 #include "corecel/random/distribution/ExponentialDistribution.hh"
0013 #include "corecel/random/distribution/IsotropicDistribution.hh"
0014 #include "corecel/random/distribution/UniformRealDistribution.hh"
0015 #include "celeritas/Types.hh"
0016 #include "celeritas/grid/NonuniformGridCalculator.hh"
0017 #include "celeritas/optical/Interaction.hh"
0018 #include "celeritas/optical/ParticleTrackView.hh"
0019 #include "celeritas/optical/TrackInitializer.hh"
0020 #include "celeritas/optical/WavelengthShiftData.hh"
0021 #include "celeritas/phys/InteractionUtils.hh"
0022
0023 namespace celeritas
0024 {
0025 namespace optical
0026 {
0027 //---------------------------------------------------------------------------//
0028 /*!
0029  * Sample optical photons from the wavelength shift process.
0030  *
0031  * A wavelength shifter absorbs an incident light and reemits secondary lights
0032  * isotropically at longer wavelengths. It usually shifts the ultraviolet
0033  * region of the radiation spectrum to the visible region, which enhances the
0034  * light collection or reduces the self-absorption of the optical production.
0035  * The number of the reemitted photons follows the Poisson distribution with
0036  * the mean number of the characteristic light production, which depends on the
0037  * optical property of wavelength shifters. The polarization of the reemitted
0038  * lights is assumed to be incoherent with respect to the polarization of the
0039  * primary optical photon.
0040  *
0041  * \note This performs the same sampling routine as in the G4OpWLS class of
0042  * the Geant4 release 11.2.
0043  */
0044 class WavelengthShiftGenerator
0045 {
0046   public:
0047     // Construct with shared and state data
0048     inline CELER_FUNCTION WavelengthShiftGenerator(
0049         NativeCRef<WavelengthShiftData> const& shared,
0050         WlsDistributionData const& distribution);
0051
0052     // Sample a WLS photon
0053     template<class Engine>
0054     inline CELER_FUNCTION TrackInitializer operator()(Engine& rng);
0055
0056   private:
0057     //// TYPES ////
0058
0059     using Energy = units::MevEnergy;
0060
0061     //// DATA ////
0062
0063     WlsDistributionData const& distribution_;
0064     real_type time_constant_;
0065     WlsDistribution time_profile_;
0066     NonuniformGridCalculator calc_cdf_;
0067 };
0068
0069 //---------------------------------------------------------------------------//
0070 // INLINE DEFINITIONS
0071 //---------------------------------------------------------------------------//
0072 /*!
0073  * Construct with shared and state data.
0074  */
0075 CELER_FUNCTION WavelengthShiftGenerator::WavelengthShiftGenerator(
0076     NativeCRef<WavelengthShiftData> const& shared,
0077     WlsDistributionData const& distribution)
0078     : distribution_(distribution)
0079     , time_constant_(shared.wls_record[distribution_.material].time_constant)
0080     , time_profile_(shared.time_profile)
0081     , calc_cdf_(shared.energy_cdf[distribution_.material], shared.reals)
0082 {
0083     CELER_EXPECT(distribution_);
0084     CELER_EXPECT(distribution_.energy.value() > calc_cdf_.grid().front());
0085 }
0086
0087 //---------------------------------------------------------------------------//
0088 // INLINE DEFINITIONS
0089 //---------------------------------------------------------------------------//
0090 /*!
0091  * Sampling the wavelength shift (WLS) photons.
0092  */
0093 template<class Engine>
0094 CELER_FUNCTION TrackInitializer WavelengthShiftGenerator::operator()(
0095     Engine& rng)
0096 {
0097     using namespace celeritas::literals;
0098
0099     // Sample wavelength shifted optical photon
0100     TrackInitializer result;
0101
0102     // Sample the emitted energy from the inverse cumulative distribution
0103     // TODO: add CDF sampler; see
0104     // https://github.com/celeritas-project/celeritas/pull/1507/files#r1844973621
0105     NonuniformGridCalculator calc_energy = calc_cdf_.make_inverse();
0106     real_type energy = calc_energy(generate_canonical(rng));
0107     if (CELER_UNLIKELY(energy >= value_as<Energy>(distribution_.energy)))
0108     {
0109         // Sample a restricted energy below the incident photon energy
0110         real_type cdf_max = calc_cdf_(value_as<Energy>(distribution_.energy));
0111         UniformRealDistribution<real_type> sample_cdf(0, cdf_max);
0112         energy = calc_energy(sample_cdf(rng));
0113     }
0114     CELER_ENSURE(energy < value_as<Energy>(distribution_.energy));
0115     result.energy = Energy{energy};
0116
0117     // Use the post-step position
0118     result.position = distribution_.position;
0119
0120     // Sample the emitted photon (incoherent) direction and polarization
0121     result.direction = IsotropicDistribution{}(rng);
0122     result.polarization = TransversePolarizationSampler{result.direction}(rng);
0123
0124     // Sample the delta time (based on the exponential relaxation)
0125     result.time = distribution_.time
0126                   + (time_profile_ == WlsDistribution::delta
0127                          ? time_constant_
0128                          : ExponentialDistribution(1_r / time_constant_)(rng));
0129     result.primary = distribution_.primary;
0130
0131     return result;
0132 }
0133
0134 //---------------------------------------------------------------------------//
0135 }  // namespace optical
0136 }  // namespace celeritas