Improve interface to speckled lasers

docs/source/api/profiles/index.rst

@@ -10,4 +10,4 @@ Laser Profiles
    combined_profile
    longitudinal/index
    transverse/index
-   speckled
+   speckled/index

docs/source/api/profiles/speckled/fm_ssd.rst

@@ -0,0 +1,5 @@
+Frequency-modulated Smoothing by Spectral Dispersion (FM-SSD) Laser Profile
+===========================================================================
+
+.. autoclass:: lasy.profiles.speckled.FM_SSD_Profile
+    :members: evaluate

docs/source/api/profiles/speckled/gp_isi.rst

@@ -0,0 +1,5 @@
+Smoothing by Induced Spatial Incoherence (ISI) Laser Profile
+============================================================
+
+.. autoclass:: lasy.profiles.speckled.GP_ISI_Profile
+    :members: evaluate

docs/source/api/profiles/speckled/gp_rpm_ssd.rst

@@ -0,0 +1,5 @@
+Random Phase Modulated (RPM) Smoothing by Spectral Dispersion (RPM-SSD) Laser Profile
+=====================================================================================
+
+.. autoclass:: lasy.profiles.speckled.GP_RPM_SSD_Profile
+    :members: evaluate

docs/source/api/profiles/speckled/index.rst

@@ -0,0 +1,11 @@
+Speckled Laser Profiles
+=======================
+
+.. toctree::
+   :maxdepth: 4
+   :hidden:
+
+   rpp_cpp
+   fm_ssd
+   gp_rpm_ssd
+   gp_isi

docs/source/api/profiles/speckled/rpp_cpp.rst

@@ -0,0 +1,5 @@
+RPP/CPP only Laser Profile
+==========================
+
+.. autoclass:: lasy.profiles.speckled.PhasePlateProfile
+    :members: evaluate

lasy/profiles/__init__.py

@@ -3,13 +3,11 @@
 from .from_array_profile import FromArrayProfile
 from .from_openpmd_profile import FromOpenPMDProfile
 from .from_insight_file import FromInsightFile
-from .speckle_profile import SpeckleProfile
 
 __all__ = [
     "CombinedLongitudinalTransverseProfile",
     "GaussianProfile",
     "FromArrayProfile",
     "FromOpenPMDProfile",
     "FromInsightFile",
-    "SpeckleProfile",
 ]

lasy/profiles/speckled/__init__.py

@@ -0,0 +1,13 @@
+from .speckle_profile import SpeckleProfile
+from .rpp_cpp import PhasePlateProfile
+from .fm_ssd import FM_SSD_Profile
+from .gp_rpm_ssd import GP_RPM_SSD_Profile
+from .gp_isi import GP_ISI_Profile
+
+__all__ = [
+    "SpeckleProfile",
+    "PhasePlateProfile",
+    "FM_SSD_Profile",
+    "GP_RPM_SSD_Profile",
+    "GP_ISI_Profile",
+]

lasy/profiles/speckled/documentation_splice.py

@@ -0,0 +1,41 @@
+class _DocumentedMetaClass(type):
+    """Metaclass that combines the __doc__ of the SpeckleProfile base and of the implementation."""
+
+    def __new__(cls, name, bases, attrs):
+        # "if bases" skips this for the _ClassWithInit (which has no bases)
+        # "if bases[0].__doc__ is not None" skips this for the picmistandard classes since their bases[0] (i.e. _ClassWithInit)
+        # has no __doc__.
+        if bases and bases[0].__doc__ is not None:
+            implementation_doc = attrs.get("__doc__", "")
+            # print('implementation doc',implementation_doc)
+            base_doc = bases[0].__doc__
+            param_delimiter = "Parameters\n    ----------\n"
+            opt_param_delimiter = "    do_include_transverse_envelope"
+
+            if implementation_doc:
+                # The format of the added string is intentional.
+                # The double return "\n\n" is needed to start a new section in the documentation.
+                # Then the four spaces matches the standard level of indentation for doc strings
+                # (assuming PEP8 formatting).
+                # The final return "\n" assumes that the implementation doc string begins with a return,
+                # i.e. a line with only three quotes, """.
+                implementation_notes, implementation_params = implementation_doc.split(
+                    param_delimiter
+                )
+                base_doc_notes, base_doc_params = base_doc.split(param_delimiter)
+                base_doc_needed_params, base_doc_opt_params = base_doc_params.split(
+                    opt_param_delimiter
+                )
+                attrs["__doc__"] = (
+                    base_doc_notes
+                    + implementation_notes
+                    + param_delimiter
+                    + base_doc_needed_params
+                    + implementation_params[:-5]
+                    + "\n\n"
+                    + opt_param_delimiter
+                    + base_doc_opt_params
+                )
+            else:
+                attrs["__doc__"] = base_doc
+        return super(_DocumentedMetaClass, cls).__new__(cls, name, bases, attrs)

lasy/profiles/speckled/fm_ssd.py

@@ -0,0 +1,144 @@
+import numpy as np
+from .speckle_profile import SpeckleProfile
+from .documentation_splice import _DocumentedMetaClass
+
+
+class FM_SSD_Profile(SpeckleProfile, metaclass=_DocumentedMetaClass):
+    r"""Speckled laser profile information specific to smoothing by frequency modulated (FM) spectral dispersion (SSD).
+
+    In frequency-modulated smoothing by spectral dispersion, or FM-SSD, the amplitude of the beamlets is always :math:`A_{ml}(t)=1`.
+    There are two contributions to the phase :math:`\phi_{ml}` of each beamlet:
+
+    .. math::
+
+        \phi_{ml}(t)=\phi_{PP,ml}+\phi_{SSD,ml}.
+
+    The phase plate part :math:`\phi_{PP,ml}` is the initial phase delay from the randomly sized phase plate sections,
+    drawn from uniform distribution on the interval :math:`[0,2\pi]`.
+    The temporal smoothing is from the SSD term:
+
+    .. math::
+
+        \begin{aligned}
+        \phi_{SSD,ml}(t)&=\delta_{x} \sin\left(\omega_{x} t + 2\pi\frac{mN_{cc,x}}{N_{bx}}\right)\\
+        &+\delta_{y} \sin\left(\omega_{y} t + 2\pi\frac{lN_{cc,y}}{N_{by}}\right).
+        \end{aligned}
+
+    The modulation frequencies :math:`\omega_x,\omega_y` are determined by the
+    laser bandwidth and modulation amplitudes according to the relation
+
+    .. math::
+
+        \omega_x = \frac{\Delta_\nu r_x }{2\delta_x},
+        \omega_y = \frac{\Delta_\nu r_y }{2\delta_y},
+
+    where :math:`\Delta_\nu` is the relative bandwidth of the laser pulse
+    and :math:`r_x, r_y` are additional rotation factors supplied by the user
+    in the `transverse_bandwidth_distribution` parameter that determine
+    how much of the modulation is in x and how much is in y. [Michel, Eqn. 9.69]
+
+    Parameters
+    ----------
+    relative_laser_bandwidth : float
+        Resulting bandwidth :math:`\Delta_\nu` of the laser pulse, relative to central frequency, due to the frequency modulation.
+
+    phase_modulation_amplitude : list of 2 floats
+        Amplitudes :math:`\delta_{x},\delta_{y}` of phase modulation in each transverse direction.
+
+    number_color_cycles : list of 2 floats
+        Number of color cycles :math:`N_{cc,x},N_{cc,y}` of SSD spectrum to include in modulation
+
+    transverse_bandwidth_distribution : list of 2 floats
+        Determines how much SSD is distributed in the :math:`x` and :math:`y` directions.
+        if `transverse_bandwidth_distribution=[a,b]`, then the SSD frequency modulation is :math:`r_x=a/\sqrt{a^2+b^2}` in :math:`x` and :math:`r_y=b/\sqrt{a^2+b^2}` in :math:`y`.
+    """
+
+    def __init__(
+        self,
+        wavelength,
+        pol,
+        laser_energy,
+        focal_length,
+        beam_aperture,
+        n_beamlets,
+        relative_laser_bandwidth,
+        phase_modulation_amplitude,
+        number_color_cycles,
+        transverse_bandwidth_distribution,
+        do_include_transverse_envelope=True,
+        long_profile=None,
+    ):
+        super().__init__(
+            wavelength,
+            pol,
+            laser_energy,
+            focal_length,
+            beam_aperture,
+            n_beamlets,
+            do_include_transverse_envelope,
+            long_profile,
+        )
+        self.laser_bandwidth = relative_laser_bandwidth
+        # the amplitude of phase along each direction
+        self.phase_modulation_amplitude = phase_modulation_amplitude
+        # number of color cycles
+        self.number_color_cycles = number_color_cycles
+        # bandwidth distributed with respect to the two transverse direction
+        self.transverse_bandwidth_distribution = transverse_bandwidth_distribution
+        normalization = np.sqrt(
+            self.transverse_bandwidth_distribution[0] ** 2
+            + self.transverse_bandwidth_distribution[1] ** 2
+        )
+        frac = [
+            self.transverse_bandwidth_distribution[0] / normalization,
+            self.transverse_bandwidth_distribution[1] / normalization,
+        ]
+        self.phase_modulation_frequency = [
+            self.laser_bandwidth * sf * 0.5 / pma
+            for sf, pma in zip(frac, self.phase_modulation_amplitude)
+        ]
+        self.time_delay = (
+            (
+                self.number_color_cycles[0] / self.phase_modulation_frequency[0]
+                if self.phase_modulation_frequency[0] > 0
+                else 0
+            ),
+            (
+                self.number_color_cycles[1] / self.phase_modulation_frequency[1]
+                if self.phase_modulation_frequency[1] > 0
+                else 0
+            ),
+        )
+        self.x_y_dephasing = np.random.standard_normal(2) * np.pi
+        self.phase_plate = np.random.uniform(
+            -np.pi, np.pi, size=self.n_beamlets[0] * self.n_beamlets[1]
+        ).reshape(self.n_beamlets)
+
+    def beamlets_complex_amplitude(
+        self,
+        t_now,
+    ):
+        """Calculate complex amplitude of the beamlets in the near-field, before propagating to the focal plane.
+
+        Parameters
+        ----------
+        t_now: float, time at which to evaluate complex amplitude
+
+        Returns
+        -------
+        array of complex numbers giving beamlet amplitude and phases in the near-field
+        """
+        phase_t = self.phase_modulation_amplitude[0] * np.sin(
+            self.x_y_dephasing[0]
+            + 2
+            * np.pi
+            * self.phase_modulation_frequency[0]
+            * (t_now - self.X_lens_matrix * self.time_delay[0] / self.n_beamlets[0])
+        ) + self.phase_modulation_amplitude[1] * np.sin(
+            self.x_y_dephasing[1]
+            + 2
+            * np.pi
+            * self.phase_modulation_frequency[1]
+            * (t_now - self.Y_lens_matrix * self.time_delay[1] / self.n_beamlets[1])
+        )
+        return np.exp(1j * (self.phase_plate + phase_t))

lasy/profiles/speckled/gp_isi.py

@@ -0,0 +1,115 @@
+import numpy as np
+from .speckle_profile import SpeckleProfile
+from .stochastic_process_utilities import gen_gaussian_time_series
+from .documentation_splice import _DocumentedMetaClass
+
+
+class GP_ISI_Profile(SpeckleProfile, metaclass=_DocumentedMetaClass):
+    r"""Speckled laser profile information specific to smoothing inspired by Induced Spatial Incoherence (ISI).
+
+    This is a smoothing technique with temporal stochastic variation in the beamlet phases and amplitudes
+    to simulate the random phase differences and amplitudes the beamlets pick up when passing through ISI echelons.
+    In this case, :math:`\phi_{ml}(t)` and :math:`A_{ml}(t)` are chosen randomly.
+    Practically, this is done by drawing the complex amplitudes :math:\tilde A_{ml}(t)`
+    from a stochastic process with Guassian power spectral density having mean of 1 and FWHM of twice the laser bandwidth.
+
+    Parameters
+    ----------
+    relative_laser_bandwidth : float
+        Bandwidth :math:`\Delta_\nu` of the incoming laser pulse, relative to the central frequency.
+
+    """
+
+    def __init__(
+        self,
+        wavelength,
+        pol,
+        laser_energy,
+        focal_length,
+        beam_aperture,
+        n_beamlets,
+        relative_laser_bandwidth,
+        do_include_transverse_envelope=True,
+        long_profile=None,
+    ):
+        super().__init__(
+            wavelength,
+            pol,
+            laser_energy,
+            focal_length,
+            beam_aperture,
+            n_beamlets,
+            do_include_transverse_envelope,
+            long_profile,
+        )
+        self.laser_bandwidth = relative_laser_bandwidth
+        self.dt_update = 1 / self.laser_bandwidth / 50
+        return
+
+    def init_gaussian_time_series(
+        self,
+        series_time,
+    ):
+        r"""Initialize a time series sampled from a Gaussian process.
+
+        At every time specified by the input `series_time`, calculate the random phases and/or amplitudes.
+
+        * This function returns a time series with random phase offsets in x and y at each time.
+            The phase offsets are real-valued and centered around the user supplied ``phase_modulation_amplitude``
+            :math:`\delta_{x},\delta_{y}`, with distribution FWHM ``phase_modulation_frequency``.
+
+        Parameters
+        ----------
+        series_time: array of times at which to sample from Gaussian process
+
+        Returns
+        -------
+        array-like, the supplied `series_time`
+        array-like, either with 2 random numbers at every time
+        """
+        complex_amp = np.stack(
+            [
+                np.stack(
+                    [
+                        gen_gaussian_time_series(
+                            series_time.size,
+                            self.dt_update,
+                            2 * self.laser_bandwidth,
+                            1,
+                        )
+                        for _i in range(self.n_beamlets[1])
+                    ]
+                )
+                for _j in range(self.n_beamlets[0])
+            ]
+        )
+        return complex_amp
+
+    def setup_for_evaluation(self, t_norm):
+        """Create or update data used in evaluation."""
+        self.x_y_dephasing = np.random.standard_normal(2) * np.pi
+        self.phase_plate = np.random.uniform(
+            -np.pi, np.pi, size=self.n_beamlets[0] * self.n_beamlets[1]
+        ).reshape(self.n_beamlets)
+
+        t_max = t_norm[-1]
+        series_time = np.arange(0, t_max + self.dt_update, self.dt_update)
+
+        self.time_series = self.init_gaussian_time_series(series_time)
+        return
+
+    def beamlets_complex_amplitude(
+        self,
+        t_now,
+    ):
+        """Calculate complex amplitude of the beamlets in the near-field, before propagating to the focal plane.
+
+        Parameters
+        ----------
+        t_now: float, time at which to evaluate complex amplitude
+
+        Returns
+        -------
+        array of complex numbers giving beamlet amplitude and phases in the near-field
+        """
+        return self.time_series[:, :, int(round(t_now / self.dt_update))]