Wavefront

Wavefront `dataclass`

Wavefront(Ex: ndarray | None = None, Ey: ndarray | None = None, dx: float = 1.0, dy: float = 1.0, dz: float = 1.0, wavelength: float = 1.0, attrs: WavefrontAttrs = WavefrontAttrs(), xmid: float = 0.0, ymid: float = 0.0, zmid: float = 0.0, s_position: float = 0.0)

Bases: WavefrontBase


              flowchart TD
              beamphysics.Wavefront[Wavefront]
              beamphysics.wavefront.wavefront.WavefrontBase[WavefrontBase]

                              beamphysics.wavefront.wavefront.WavefrontBase --> beamphysics.Wavefront
                


              click beamphysics.Wavefront href "" "beamphysics.Wavefront"
              click beamphysics.wavefront.wavefront.WavefrontBase href "" "beamphysics.wavefront.wavefront.WavefrontBase"

Real-space representation of electromagnetic wavefront fields.

This class represents electromagnetic fields in real space (position domain) for optical wavefront propagation and analysis. It stores complex electric field components Ex and Ey on a 3D Cartesian grid with uniform spacing.

The wavefront propagates in the +z direction, with the pulse head at max(z). Field intensities, energies, and other derived quantities follow standard electromagnetic definitions using SI units.

Parameters:

Ex (ndarray, default: None ) –

x-polarized electric field component in V/m, shape (nx, ny, nz)
Ey (ndarray, default: None ) –

y-polarized electric field component in V/m, shape (nx, ny, nz)
dx (float, default: 1.0 ) –

Grid spacing in x direction (m)
dy (float, default: 1.0 ) –

Grid spacing in y direction (m)
dz (float, default: 1.0 ) –

Grid spacing in z direction (m)
wavelength (float, default: 1.0 ) –

Central wavelength (m)
attrs (WavefrontAttrs or mapping, default: WavefrontAttrs() ) –

openPMD EXT_Wavefront attributes carried through file I/O, such as beamline. Not interpreted by this class. A mapping keyed by either openPMD or field names is coerced to WavefrontAttrs.
xmid (float, default: 0.0 ) –

Midpoint of the grid in m. The default of 0.0 gives a grid centered on the origin, so that xvec runs from -(nx-1)*dx/2 to +(nx-1)*dx/2. Note that z here is the intra-pulse coordinate, not the beamline position: see the note on frames below.
ymid (float, default: 0.0 ) –

Midpoint of the grid in m. The default of 0.0 gives a grid centered on the origin, so that xvec runs from -(nx-1)*dx/2 to +(nx-1)*dx/2. Note that z here is the intra-pulse coordinate, not the beamline position: see the note on frames below.
zmid (float, default: 0.0 ) –

Midpoint of the grid in m. The default of 0.0 gives a grid centered on the origin, so that xvec runs from -(nx-1)*dx/2 to +(nx-1)*dx/2. Note that z here is the intra-pulse coordinate, not the beamline position: see the note on frames below.
s_position (float, default: 0.0 ) –

Position of this plane along the beamline in m, in whatever frame the caller picked. Advanced by drift. See the note on frames below.

Attributes:

spatial_domain (SpatialDomain) –

Always returns SpatialDomain.R (real space)
intensity (ndarray) –

Total field intensity c ε₀/2 (|Ex|² + |Ey|²) in W/m²
energy (float) –

Total field energy ε₀/2 ∫∫∫ |E|² dx dy dz in J
fluence (ndarray) –

Fluence F(x,y) = ε₀/2 ∫ |E(x,y,z)|² dz in J/m²
power (ndarray) –

Power profile P(z) = ∫∫ I(x,y,z) dx dy in W
k0 (float) –

Central wavenumber 2π/λ in rad/m
photon_energy (float) –

Central photon energy ℏck₀ in eV

Methods:

to_kspace –

Transform to k-space representation (WavefrontK)
drift –

Propagate wavefront by distance z in free space
pad –

Zero-pad the field arrays
crop –

Crop the field arrays by removing elements from edges
auto_crop –

Determine and optionally apply symmetric crop based on intensity thresholds
estimate_curvature –

Estimate wavefront radius of curvature
plot_power –

Plot longitudinal power profile
plot_fluence –

Plot transverse fluence distribution
from_genesis4 –

Create Wavefront from Genesis4 HDF5 field file
write_genesis4 –

Write Wavefront to Genesis4 HDF5 format
from_openpmd –

Create Wavefront from an openPMD EXT_Wavefront file
write_openpmd –

Write Wavefront to an openPMD EXT_Wavefront file

Notes

At least one of Ex or Ey must be provided
All spacing parameters (dx, dy, dz, wavelength) must be positive
Field arrays must be complex dtype
Fourier transforms use norm='ortho' to preserve Plancherel's theorem: ∫∫∫ |E(x,y,z)|² dx dy dz = ∫∫∫ |Ẽ(kx,ky,kz)|² dkx dky dkz
Statistical moments (mean, sigma) are intensity-weighted

Frames. There are two longitudinal quantities and they belong to different frames. The mesh's own z axis, positioned by zmid, is the intra-pulse coordinate: it is co-moving, and drift leaves it alone because the paraxial kernel is transverse, applied slice by slice with no piston term. s_position is the position of this plane along the beamline, and drift does advance it, by the propagation distance.

That split follows ParticleGroup, where drift moves the coordinates that exist rather than a separate register: there, z records that the particles propagated. Here the mesh is co-moving and cannot record it, so s_position is the only place the information can go. Operations that do not move the wavefront along the beamline -- crop, pad, to_kspace, an applied lens phase -- leave it alone, just as they leave ParticleGroup.z alone.

s_position is measured in whatever frame the caller picked, so the default of 0.0 is the origin of that frame rather than an unknown. It is written as the EXT_Wavefront zCoordinate. Genesis4's refposition is not this quantity: it is the coordinate of the first slice, and corresponds to zmin.

Returned by:

API Wavefront
- Wavefront
- WavefrontK
  - crop
  - to_rspace

Used by:

API Wavefront openpmd write_wavefront_openpmd

dkx `property`

dkx: float

Transverse wavevector spacing in kx direction in rad/m.

dkx = 2π / (nx * dx)

dky `property`

dky

Transverse wavevector spacing in ky direction in rad/m.

dky = 2π / (ny * dy)

dkz `property`

dkz

Longitudinal wavevector spacing in kz direction in rad/m.

dkz = 2π / (nz * dz)

dthetax `property`

dthetax

Angular spacing in thetax direction in rad.

dthetax = wavelength / (nx * dx) = dkx / k0

dthetay `property`

dthetay

Angular spacing in thetay direction in rad.

dthetay = wavelength / (ny * dy) = dky / k0

energy `property`

energy

Total (time-averaged) energy in J

energy = ϵ0/2 ∫∫∫ |E|^2 dx dy dz

   = ϵ0/2  ∫∫∫ |Ẽ|^2 dkx dky dkz

energy_density `property`

energy_density

Energy density ϵ0/2 |E|^2 in J/m^3

3D real array

fluence `property`

fluence

Transverse fluence profile in J/m².

F(x,y) = ϵ0/2 ∫ |E(x,y,z)|² dz

2D real array with shape (nx, ny).

fluence_profile_x `property`

fluence_profile_x

1D fluence profile along x in J/m.

F_x(x) = ϵ0/2 ∫∫ |E(x,y,z)|² dy dz

Integrated over y and z directions.

fluence_profile_y `property`

fluence_profile_y

1D fluence profile along y in J/m.

F_y(y) = ϵ0/2 ∫∫ |E(x,y,z)|² dx dz

Integrated over x and z directions.

in_kspace `property`

in_kspace

True if the wavefront is in k-space (Fourier domain).

in_rspace `property`

in_rspace

True if the wavefront is in real space (position domain).

intensity `property`

intensity: ndarray | float

total field intensity in W/m^2 I = c ϵ0/2 (|Ex|^2 + |Ey|^2)

intensity_x `property`

intensity_x: ndarray | float

x polarization field intensity in W/m^2 Intensity = c ϵ0/2 |Ex|^2

intensity_y `property`

intensity_y: ndarray | float

y polarization field intensity in W/m^2 Intensity = c ϵ0/2 |Ey|^2

k0 `property`

k0

Central wavenumber in rad/m.

k0 = 2π / wavelength

kxmax `property`

kxmax

Maximum kx value in rad/m.

kxmin `property`

kxmin

Minimum kx value in rad/m.

kxvec `property`

kxvec

Array of kx (transverse wavenumber) values in rad/m.

kymax `property`

kymax

Maximum ky value in rad/m.

kymin `property`

kymin

Minimum ky value in rad/m.

kyvec `property`

kyvec

Array of ky (transverse wavenumber) values in rad/m.

kzmax `property`

kzmax

Maximum kz value in rad/m.

kzmin `property`

kzmin

Minimum kz value in rad/m.

kzvec `property`

kzvec

Array of kz (longitudinal wavenumber) values in rad/m.

mean_x `property`

mean_x

Intensity-weighted mean x position in m.

= ∫ x I(x,y,z) dx dy dz / ∫ I(x,y,z) dx dy dz

mean_y `property`

mean_y

Intensity-weighted mean y position in m.

= ∫ y I(x,y,z) dx dy dz / ∫ I(x,y,z) dx dy dz

mean_z `property`

mean_z

Intensity-weighted mean z position in m.

= ∫ z I(x,y,z) dx dy dz / ∫ I(x,y,z) dx dy dz

nx `property`

nx

Number of grid points in the x direction.

ny `property`

ny

Number of grid points in the y direction.

nz `property`

nz

Number of grid points in the z direction.

photon_energy `property`

photon_energy: float

Central photon energy in eV

power `property`

power

Longitudinal power profile along z in W.

P(z) = ∫∫ I(x, y, z) dx dy = ∫∫ (c ϵ0/2) |E(x, y, z)|² dx dy

1D real array with shape (nz,).

shape `property`

shape

Shape of the field arrays as (nx, ny, nz).

sigma_x `property`

sigma_x

RMS beam size in x in m.

σ_x = sqrt( - ²)

sigma_y `property`

sigma_y

RMS beam size in y in m.

σ_y = sqrt( - ²)

sigma_z `property`

sigma_z

RMS pulse length in z in m.

σ_z = sqrt( - ²)

thetaxmax `property`

thetaxmax

Maximum thetax value in rad.

thetaxmin `property`

thetaxmin

Minimum thetax value in rad.

thetaxvec `property`

thetaxvec

Array of thetax (angular deviation) values in rad.

thetax = kx / k0

thetaymax `property`

thetaymax

Maximum thetay value in rad.

thetaymin `property`

thetaymin

Minimum thetay value in rad.

thetayvec `property`

thetayvec

Array of thetay (angular deviation) values in rad.

thetay = ky / k0

xmax `property`

xmax

Maximum x coordinate in m.

xmax = xmid + (nx-1) * dx / 2

xmin `property`

xmin

Minimum x coordinate in m.

xmin = xmid - (nx-1) * dx / 2

xvec `property`

xvec

Array of x coordinates in m.

ymax `property`

ymax

Maximum y coordinate in m.

ymax = ymid + (ny-1) * dy / 2

ymin `property`

ymin

Minimum y coordinate in m.

ymin = ymid - (ny-1) * dy / 2

yvec `property`

yvec

Array of y coordinates in m.

zmax `property`

zmax

Maximum z coordinate in m.

zmax = zmid + (nz-1) * dz / 2

zmin `property`

zmin

Minimum z coordinate in m.

zmin = zmid - (nz-1) * dz / 2

zvec `property`

zvec

Array of z coordinates in m.

__post_init__

__post_init__()

Validate inputs after dataclass initialization

Source code in beamphysics/wavefront/wavefront.py

def __post_init__(self):
    """
    Validate inputs after dataclass initialization
    """
    # Validate field shapes match
    if self.Ex is not None and self.Ey is not None:
        if self.Ex.shape != self.Ey.shape:
            raise ValueError(
                f"Ex shape {self.Ex.shape} != Ey shape {self.Ey.shape}"
            )

    # Validate at least one field exists
    if self.Ex is None and self.Ey is None:
        raise ValueError("At least one of Ex or Ey must be provided")

    # Validate positive spacing
    for attr in ["dx", "dy", "dz", "wavelength"]:
        val = getattr(self, attr)
        if val <= 0:
            raise ValueError(f"{attr} must be positive, got {val}")

setattr

__setattr__(name, value)

Keep attrs a WavefrontAttrs however it is set.

A mapping keyed by either openPMD or field names is accepted for convenience, but the attribute itself is always the typed object. Coercing only in __post_init__ would leave a plain assignment such as w.attrs = {...} holding a raw dict, so that w.attrs.beamline raises AttributeError and the mapping's contents go unvalidated until write time.

Source code in beamphysics/wavefront/wavefront.py

def __setattr__(self, name, value):
    """
    Keep `attrs` a `WavefrontAttrs` however it is set.

    A mapping keyed by either openPMD or field names is accepted for
    convenience, but the attribute itself is always the typed object. Coercing
    only in `__post_init__` would leave a plain assignment such as
    ``w.attrs = {...}`` holding a raw dict, so that `w.attrs.beamline`
    raises `AttributeError` and the mapping's contents go unvalidated until
    write time.
    """
    if name == "attrs" and not isinstance(value, WavefrontAttrs):
        value = WavefrontAttrs.from_pmd(value)
    object.__setattr__(self, name, value)

auto_crop

auto_crop(threshold: float = 1e-06, apply: bool = True)

Determine symmetric crop values based on intensity profile thresholds.

Analyzes the 1D intensity projections along each axis to find regions where the signal falls below a relative threshold, then computes symmetric crop amounts for each axis.

Parameters:

threshold (float, default: 1e-6 ) –

Relative threshold (fraction of maximum) below which data is considered negligible and can be cropped.
apply (bool, default: True ) –

If True, apply the crop and return the cropped wavefront. If False, return only the crop amounts as a dictionary.

Returns:

Wavefront, WavefrontK, or dict –

If apply=True: New wavefront with cropped field arrays. If apply=False: Dictionary with keys 'nx', 'ny', 'nz' containing the symmetric crop amounts for each axis.

Examples:

>>> from beamphysics.wavefront import Wavefront
>>> w = Wavefront.from_gaussian(
...     shape=(64, 64, 16), dx=1e-6, dy=1e-6, dz=1e-6,
...     wavelength=1e-9, sigma0=3e-6, sigma_z=2e-6,
... )

Get the crop amounts without applying them:

>>> crop_info = w.auto_crop(threshold=1e-4, apply=False)
>>> crop_info
{'nx': 19, 'ny': 19, 'nz': 0}

Apply the crop directly, removing that many elements from each end:

>>> w.auto_crop(threshold=1e-4).shape
(26, 26, 16)

The two steps are equivalent, so the amounts can be inspected or adjusted first:

>>> w.crop(**crop_info).shape
(26, 26, 16)

Source code in beamphysics/wavefront/wavefront.py

def auto_crop(self, threshold: float = 1e-6, apply: bool = True):
    """
    Determine symmetric crop values based on intensity profile thresholds.

    Analyzes the 1D intensity projections along each axis to find
    regions where the signal falls below a relative threshold, then
    computes symmetric crop amounts for each axis.

    Parameters
    ----------
    threshold : float, default=1e-6
        Relative threshold (fraction of maximum) below which data
        is considered negligible and can be cropped.
    apply : bool, default=True
        If True, apply the crop and return the cropped wavefront.
        If False, return only the crop amounts as a dictionary.

    Returns
    -------
    Wavefront, WavefrontK, or dict
        If apply=True: New wavefront with cropped field arrays.
        If apply=False: Dictionary with keys 'nx', 'ny', 'nz' containing
        the symmetric crop amounts for each axis.

    Examples
    --------
    >>> from beamphysics.wavefront import Wavefront
    >>> w = Wavefront.from_gaussian(
    ...     shape=(64, 64, 16), dx=1e-6, dy=1e-6, dz=1e-6,
    ...     wavelength=1e-9, sigma0=3e-6, sigma_z=2e-6,
    ... )

    Get the crop amounts without applying them:

    >>> crop_info = w.auto_crop(threshold=1e-4, apply=False)
    >>> crop_info
    {'nx': 19, 'ny': 19, 'nz': 0}

    Apply the crop directly, removing that many elements from each end:

    >>> w.auto_crop(threshold=1e-4).shape
    (26, 26, 16)

    The two steps are equivalent, so the amounts can be inspected or adjusted
    first:

    >>> w.crop(**crop_info).shape
    (26, 26, 16)
    """

    def find_symmetric_crop(profile: np.ndarray, threshold: float) -> int:
        """Find symmetric crop amount for a 1D profile."""
        if profile.size == 0:
            return 0

        max_val = np.max(profile)
        if max_val == 0:
            return 0

        # Normalize and find where signal is above threshold
        normalized = profile / max_val
        above_threshold = normalized > threshold

        if not np.any(above_threshold):
            return 0

        # Find first and last indices above threshold
        indices = np.where(above_threshold)[0]
        first_idx = indices[0]
        last_idx = indices[-1]

        # Calculate symmetric crop (use minimum of both sides)
        crop_start = first_idx
        crop_end = len(profile) - 1 - last_idx

        # A plain int, so that the returned dict reprs and serializes cleanly.
        return int(min(crop_start, crop_end))

    # Get 1D intensity projections along each axis
    intensity = self.intensity
    profile_x = np.sum(intensity, axis=_axis_for_sum[self.axis_labels[0]])
    profile_y = np.sum(intensity, axis=_axis_for_sum[self.axis_labels[1]])
    profile_z = np.sum(intensity, axis=_axis_for_sum[self.axis_labels[2]])

    # Analyze each axis
    nx = find_symmetric_crop(profile_x, threshold)
    ny = find_symmetric_crop(profile_y, threshold)
    nz = find_symmetric_crop(profile_z, threshold)

    crop_amounts = {"nx": nx, "ny": ny, "nz": nz}

    if apply:
        return self.crop(**crop_amounts)
    else:
        return crop_amounts

axis_index

axis_index(key)

Returns axis index for a named axis label key.

Examples:

.axis_labels == ('x', 'y', 'z')
.axis_index('z') returns 2

Source code in beamphysics/wavefront/wavefront.py

def axis_index(self, key):
    """
    Returns axis index for a named axis label key.

    Examples:

    - `.axis_labels == ('x', 'y', 'z')`
    - `.axis_index('z')` returns `2`
    """
    return self.axis_labels.index(key)

copy

copy()

Returns a deep copy

Source code in beamphysics/wavefront/wavefront.py

def copy(self):
    """Returns a deep copy"""
    return deepcopy(self)

crop

crop(nx=(0, 0), ny=(0, 0), nz=(0, 0))

Crop the field arrays by removing elements from edges.

Parameters:

nx (int or tuple of (int, int), default: (0, 0) ) –

Number of elements to remove from (start, end) of first axis. If a single int, removes symmetrically from both sides.
ny (int or tuple of (int, int), default: (0, 0) ) –

Number of elements to remove from (start, end) of second axis. If a single int, removes symmetrically from both sides.
nz (int or tuple of (int, int), default: (0, 0) ) –

Number of elements to remove from (start, end) of third axis. If a single int, removes symmetrically from both sides.

Returns:

Wavefront or WavefrontK –

New wavefront with cropped field arrays.

Raises:

ValueError –

If crop amounts exceed array dimensions.

Notes

The surviving samples keep their physical coordinates: an asymmetric crop moves the grid origin rather than re-centering the grid.

Source code in beamphysics/wavefront/wavefront.py

def crop(self, nx=(0, 0), ny=(0, 0), nz=(0, 0)):
    """
    Crop the field arrays by removing elements from edges.

    Parameters
    ----------
    nx : int or tuple of (int, int), default=(0, 0)
        Number of elements to remove from (start, end) of first axis.
        If a single int, removes symmetrically from both sides.
    ny : int or tuple of (int, int), default=(0, 0)
        Number of elements to remove from (start, end) of second axis.
        If a single int, removes symmetrically from both sides.
    nz : int or tuple of (int, int), default=(0, 0)
        Number of elements to remove from (start, end) of third axis.
        If a single int, removes symmetrically from both sides.

    Returns
    -------
    Wavefront or WavefrontK
        New wavefront with cropped field arrays.

    Raises
    ------
    ValueError
        If crop amounts exceed array dimensions.

    Notes
    -----
    The surviving samples keep their physical coordinates: an asymmetric crop
    moves the grid origin rather than re-centering the grid.
    """
    nx = (nx, nx) if np.isscalar(nx) else nx
    ny = (ny, ny) if np.isscalar(ny) else ny
    nz = (nz, nz) if np.isscalar(nz) else nz

    # Calculate slice indices
    x_start = nx[0]
    x_end = self.shape[0] - nx[1] if nx[1] > 0 else None
    y_start = ny[0]
    y_end = self.shape[1] - ny[1] if ny[1] > 0 else None
    z_start = nz[0]
    z_end = self.shape[2] - nz[1] if nz[1] > 0 else None

    # Validate crop amounts
    if nx[0] + nx[1] >= self.shape[0]:
        raise ValueError(
            f"Crop in x ({nx}) exceeds array dimension ({self.shape[0]})"
        )
    if ny[0] + ny[1] >= self.shape[1]:
        raise ValueError(
            f"Crop in y ({ny}) exceeds array dimension ({self.shape[1]})"
        )
    if nz[0] + nz[1] >= self.shape[2]:
        raise ValueError(
            f"Crop in z ({nz}) exceeds array dimension ({self.shape[2]})"
        )

    Ex = (
        self.Ex[x_start:x_end, y_start:y_end, z_start:z_end]
        if self.Ex is not None
        else None
    )
    Ey = (
        self.Ey[x_start:x_end, y_start:y_end, z_start:z_end]
        if self.Ey is not None
        else None
    )

    return replace(
        self,
        Ex=Ex,
        Ey=Ey,
        attrs=self.attrs.copy(),
        **self._resized_mids(nx, ny, nz),
    )

estimate_curvature

estimate_curvature(axis: str = 'x', polarization: str | None = None, ix: int | None = None, iy: int | None = None, iz: int | None = None, plot: bool = False, cutoff: float = 0.0001) -> float | tuple[float, Figure]

Estimates the wavefront curvature by looking at the phase across an axis.

The method fits a parabola to the unwrapped phase: phase = k r^2 / 2R + phi0 where R is the radius of curvature and k is the wavenumber.

Parameters:

axis (str, default: 'x' ) –

Axis along which to estimate curvature ('x' or 'y')
polarization (str, default: None ) –

Polarization to analyze ('x' or 'y'). If None, uses first available.
ix (int, default: None ) –

Indices for slicing. iz defaults to peak intensity slice.
iy (int, default: None ) –

Indices for slicing. iz defaults to peak intensity slice.
iz (int, default: None ) –

Indices for slicing. iz defaults to peak intensity slice.
plot (bool, default: False ) –

Whether to create a diagnostic plot
cutoff (float, default: 1e-4 ) –

Relative intensity cutoff for data selection

Returns:

float –

Curvature in 1/m.

Raises:

ValueError –

If insufficient data points or invalid parameters

Source code in beamphysics/wavefront/wavefront.py

def estimate_curvature(
    self,
    axis: str = "x",
    polarization: str | None = None,
    ix: int | None = None,
    iy: int | None = None,
    iz: int | None = None,
    plot: bool = False,
    cutoff: float = 1e-4,
) -> float | tuple[float, plt.Figure]:
    """
    Estimates the wavefront curvature by looking at the phase across an axis.

    The method fits a parabola to the unwrapped phase: phase = k r^2 / 2R + phi0
    where R is the radius of curvature and k is the wavenumber.

    Parameters
    ----------
    axis : str, default='x'
        Axis along which to estimate curvature ('x' or 'y')
    polarization : str, optional
        Polarization to analyze ('x' or 'y'). If None, uses first available.
    ix, iy, iz : int, optional
        Indices for slicing. iz defaults to peak intensity slice.
    plot : bool, default=False
        Whether to create a diagnostic plot
    cutoff : float, default=1e-4
        Relative intensity cutoff for data selection

    Returns
    -------
    float
        Curvature in 1/m.

    Raises
    ------
    ValueError
        If insufficient data points or invalid parameters
    """
    w = self

    if polarization is None:
        polarization = "x" if w.Ex is not None else "y"
    field = getattr(w, "E" + polarization)
    if field is None:
        raise ValueError(
            f"E{polarization} field is None, cannot use polarization='{polarization}'"
        )

    if iz is None:
        iz = np.sum(np.abs(field) ** 2, axis=(0, 1)).argmax()

    field_xy = field[:, :, iz]
    field_xy2 = np.abs(field_xy) ** 2

    if axis == "x":
        if iy is None:
            iy = np.sum(field_xy2, axis=1).argmax()

        if ix is not None:
            raise ValueError(f"{ix=} must be None when requesting axis='x'")

        weights = field_xy2[:, iy]

        phase = np.unwrap(np.angle(field_xy[:, iy]))

    elif axis == "y":
        if ix is None:
            ix = np.sum(field_xy2, axis=0).argmax()

        if iy is not None:
            raise ValueError(f"{iy=} must be None when requesting axis='y'")

        weights = field_xy2[ix, :]

        phase = np.unwrap(np.angle(field_xy[ix, :]))
    else:
        raise ValueError(f"{axis=} must be 'x' or 'y'")

    x = getattr(w, axis + "vec")
    y = phase - phase.min()

    weights = weights / weights.max()
    mask = weights > cutoff

    x = x[mask]
    y = y[mask]
    weights = weights[mask]

    if len(x) < 3:
        raise ValueError(
            f"Insufficient points ({len(x)}) after cutoff filtering for polynomial fit"
        )

    a, b, c = np.polyfit(x, y, 2, w=weights)

    # curvature = 1/R = a λ / π

    curvature = a * w.wavelength / np.pi

    if plot:
        fig, ax = plt.subplots()
        norm = 360 / (2 * np.pi)
        ax.plot(x, y * norm)

        radius = 1 / curvature if abs(curvature) > 1e-10 else np.inf
        label = f"curvature_{axis} = {curvature:0.6f} 1/m \n radius_{axis} = {radius:0.1f} m"
        ax.plot(x, (a * x**2 + b * x + c) * norm, "--", label=label)
        ax.set_xlabel(f"{axis} (m)")

        ax.set_ylabel("phase (deg)")
        # ax.set_ylim(0, None)
        ax.legend()

        ax2 = ax.twinx()
        ax2.fill_between(x, 0, weights / weights.max(), alpha=0.5, color="purple")
        ax2.set_ylim(0, 2)
        ax2.set_ylabel("weight")

    return curvature

from_gaussian `classmethod`

from_gaussian(shape: tuple[int, int, int], dx: float = 1.0, dy: float = 1.0, dz: float = 1.0, wavelength: float = 1.0, sigma0: float | None = None, z0: float = 0.0, x0: float = 0.0, y0: float = 0.0, sigma_z: float | None = None, energy: float = 1.0, phase: float = 0.0, polarization: str = 'x') -> 'Wavefront'

Create a Wavefront with a Gaussian beam.

Parameters:

shape (tuple[int, int, int]) –

Grid shape (nx, ny, nz)
dx (float, default: 1.0 ) –

Grid spacing in x direction (m)
dy (float, default: 1.0 ) –

Grid spacing in y direction (m)
dz (float, default: 1.0 ) –

Grid spacing in z direction (m)
wavelength (float, default: 1.0 ) –

Central wavelength (m)
sigma0 (float, default: None ) –

RMS transverse beam size at the waist in meters (round beam). Related to waist size w0 via: w0 = 2·σ₀ Related to Rayleigh length via: z_R = 4π·σ₀² / λ = π·w0² / λ
z0 (float, default: 0.0 ) –

Propagation distance from the waist (m), the same quantity drift advances. Sets the transverse size, curvature and Gouy phase. It does not shift the pulse along the mesh z axis (set zmid), nor set the beamline position (set s_position).
x0 (float, default: 0.0 ) –

Beam center position in x (m)
y0 (float, default: 0.0 ) –

Beam center position in y (m)
sigma_z (float or None, default: None ) –

Longitudinal Gaussian width (m). If None, use constant longitudinal profile.
energy (float, default: 1.0 ) –

Total beam energy in Joules
phase (float, default: 0.0 ) –

Global phase offset (radians)
polarization (str, default: 'x' ) –

Polarization direction ('x' or 'y')

Returns:

Wavefront –

New Wavefront instance with Gaussian beam

Raises:

ValueError –

If sigma0 is not specified If polarization is not 'x' or 'y' If sigma_z is negative

Notes

The Gaussian beam uses the complex beam parameter formalism in the exp(-i·ω·t + i·k·z) convention shared by drift and openPMD EXT_Wavefront:

q(z) = z - i·z_R where z_R = 4π·σ₀² / λ = π·w0² / λ

u_xy(x,y,z) = (1/q) · exp(i·k·(x² + y²)/(2q))

For sigma_z > 0: u_z(z) = exp(-z²/(4σ_z²)) / √(√(2π)·σ_z) For sigma_z = None: u_z(z) = 1 (constant profile, normalized over grid)

The field is normalized such that: ∫∫∫ |u|² dx dy dz = energy / (2·Z₀·c)

The relationship between sigma0 and common Gaussian beam parameters: - Waist size (1/e² intensity radius): w0 = 2·σ₀ - Rayleigh length: z_R = π·w0²/λ = 4π·σ₀²/λ

References

Wikipedia: Gaussian beam, Complex beam parameter
Siegman "Lasers" 1986, Chapter 16.3

Source code in beamphysics/wavefront/wavefront.py

@classmethod
def from_gaussian(
    cls,
    shape: tuple[int, int, int],
    dx: float = 1.0,
    dy: float = 1.0,
    dz: float = 1.0,
    wavelength: float = 1.0,
    sigma0: float | None = None,
    z0: float = 0.0,
    x0: float = 0.0,
    y0: float = 0.0,
    sigma_z: float | None = None,
    energy: float = 1.0,
    phase: float = 0.0,
    polarization: str = "x",
) -> "Wavefront":
    """
    Create a Wavefront with a Gaussian beam.

    Parameters
    ----------
    shape : tuple[int, int, int]
        Grid shape (nx, ny, nz)
    dx : float, default=1.0
        Grid spacing in x direction (m)
    dy : float, default=1.0
        Grid spacing in y direction (m)
    dz : float, default=1.0
        Grid spacing in z direction (m)
    wavelength : float, default=1.0
        Central wavelength (m)
    sigma0 : float, optional
        RMS transverse beam size at the waist in meters (round beam).
        Related to waist size w0 via: w0 = 2·σ₀
        Related to Rayleigh length via: z_R = 4π·σ₀² / λ = π·w0² / λ
    z0 : float, default=0.0
        Propagation distance from the waist (m), the same quantity `drift`
        advances. Sets the transverse size, curvature and Gouy phase. It does
        not shift the pulse along the mesh `z` axis (set `zmid`), nor set the
        beamline position (set `s_position`).
    x0 : float, default=0.0
        Beam center position in x (m)
    y0 : float, default=0.0
        Beam center position in y (m)
    sigma_z : float or None, default=None
        Longitudinal Gaussian width (m). If None, use constant longitudinal profile.
    energy : float, default=1.0
        Total beam energy in Joules
    phase : float, default=0.0
        Global phase offset (radians)
    polarization : str, default='x'
        Polarization direction ('x' or 'y')

    Returns
    -------
    Wavefront
        New Wavefront instance with Gaussian beam

    Raises
    ------
    ValueError
        If sigma0 is not specified
        If polarization is not 'x' or 'y'
        If sigma_z is negative

    Notes
    -----
    The Gaussian beam uses the complex beam parameter formalism in the
    exp(-i·ω·t + i·k·z) convention shared by `drift` and openPMD EXT_Wavefront:

    q(z) = z - i·z_R  where z_R = 4π·σ₀² / λ = π·w0² / λ

    u_xy(x,y,z) = (1/q) · exp(i·k·(x² + y²)/(2q))

    For sigma_z > 0:
        u_z(z) = exp(-z²/(4σ_z²)) / √(√(2π)·σ_z)
    For sigma_z = None:
        u_z(z) = 1  (constant profile, normalized over grid)

    The field is normalized such that:
    ∫∫∫ |u|² dx dy dz = energy / (2·Z₀·c)

    The relationship between sigma0 and common Gaussian beam parameters:
    - Waist size (1/e² intensity radius): w0 = 2·σ₀
    - Rayleigh length: z_R = π·w0²/λ = 4π·σ₀²/λ

    References
    ----------
    - Wikipedia: Gaussian beam, Complex beam parameter
    - Siegman "Lasers" 1986, Chapter 16.3
    """
    if sigma0 is None:
        raise ValueError("sigma0 must be specified")

    if polarization not in ("x", "y"):
        raise ValueError(f"polarization must be 'x' or 'y', got '{polarization}'")

    if sigma_z is not None and sigma_z < 0:
        raise ValueError(f"sigma_z must be non-negative or None, got {sigma_z}")

    # Convert sigma0 to Rayleigh length: z_R = π·w0²/λ = 4π·σ₀²/λ where w0 = 2·σ₀
    zR = 4 * pi * sigma0**2 / wavelength

    # Create coordinate vectors
    nx, ny, nz = shape
    xvec = np.linspace(-((nx - 1) * dx) / 2, ((nx - 1) * dx) / 2, nx)
    yvec = np.linspace(-((ny - 1) * dy) / 2, ((ny - 1) * dy) / 2, ny)
    zvec = np.linspace(-((nz - 1) * dz) / 2, ((nz - 1) * dz) / 2, nz)

    X, Y, Z = np.meshgrid(xvec, yvec, zvec, indexing="ij")

    k = 2 * pi / wavelength

    # Complex beam parameter at position z0 from waist (exp(-iwt + ikz) convention)
    q = z0 - 1j * zR

    # Transverse Gaussian beam profile
    uxy = (1 / q) * np.exp(0.5j * k * ((X - x0) ** 2 + (Y - y0) ** 2) / q)

    # Analytic integral of |u_xy|² over transverse plane
    integral_uxy_squared = pi / (k * zR)

    # Longitudinal profile
    if sigma_z is None or sigma_z == 0:
        # Constant longitudinal profile
        uz = np.ones_like(Z)
    else:
        # Gaussian longitudinal profile (square root of Gaussian in intensity)
        uz = np.sqrt(1 / np.sqrt(2 * pi) / sigma_z) * np.exp(
            -(Z**2) / (4 * sigma_z**2)
        )

    # Combined field, normalized by transverse integral
    u = uxy * uz / np.sqrt(integral_uxy_squared)

    # Exact normalization over the discrete grid
    integral2 = np.sum(np.abs(u) ** 2) * dx * dy * dz
    u = u / np.sqrt(integral2)

    # Scale for desired energy using c_light for consistency
    u = u * np.sqrt(energy * 2 * Z0 * c_light) * np.exp(1j * phase)

    # Create wavefront with appropriate polarization
    if polarization == "x":
        return cls(Ex=u, dx=dx, dy=dy, dz=dz, wavelength=wavelength)
    else:
        return cls(Ey=u, dx=dx, dy=dy, dz=dz, wavelength=wavelength)

from_genesis4 `classmethod`

from_genesis4(file: Union[Path, str, Group])

Create a Wavefront instance from a Genesis4 field file.

This method reads a Genesis4 field file (either from disk or an open HDF5 group) and initializes a Wavefront instance with the extracted field data.

Parameters:

file (Union[Path, str, Group]) –

The path to the Genesis4 HDF5 field file or an open h5py.Group.

Returns:

Wavefront –

A new instance of Wavefront initialized with field data from the Genesis4 file.

Raises:

ValueError –

If the provided file is not a valid path, string, or h5py.Group.

Notes

The field data is extracted and converted into an electric field representation.
The grid spacing (dx and dz) and wavelength are obtained from the file metadata.
Genesis4's refposition is the coordinate of the first slice, so it sets zmin. The transverse grid is centered on the axis.

Source code in beamphysics/wavefront/wavefront.py

@classmethod
def from_genesis4(
    cls,
    file: Union[pathlib.Path, str, h5py.Group],
):
    """
    Create a Wavefront instance from a Genesis4 field file.

    This method reads a Genesis4 field file (either from disk or an open HDF5 group)
    and initializes a `Wavefront` instance with the extracted field data.

    Parameters
    ----------
    file : Union[pathlib.Path, str, h5py.Group]
        The path to the Genesis4 HDF5 field file or an open `h5py.Group`.

    Returns
    -------
    Wavefront
        A new instance of `Wavefront` initialized with field data from the Genesis4 file.

    Raises
    ------
    ValueError
        If the provided file is not a valid path, string, or `h5py.Group`.

    Notes
    -----
    - The field data is extracted and converted into an electric field representation.
    - The grid spacing (`dx` and `dz`) and wavelength are obtained from the file metadata.
    - Genesis4's `refposition` is the coordinate of the first slice, so it sets
      `zmin`. The transverse grid is centered on the axis.
    """

    if isinstance(file, (str, pathlib.Path)):
        with h5py.File(file, "r") as h5:
            dfl, param = load_genesis4_fields(h5)

    elif isinstance(file, h5py.Group):
        dfl, param = load_genesis4_fields(file)
    else:
        raise ValueError(f"{file=} must be a str, pathlib.Path, or h5py.Group")

    dx = param["gridsize"]
    dz = param["slicespacing"]
    wavelength = param["wavelength"]

    Ex = dfl * np.sqrt(2 * Z0) / dx

    return cls(
        Ex=Ex,
        dx=float(dx),
        dy=float(dx),
        dz=float(dz),
        wavelength=float(wavelength),
        zmid=float(param["refposition"]) + (Ex.shape[2] - 1) * float(dz) / 2,
    )

from_openpmd `classmethod`

from_openpmd(file: Union[Path, str, Group], iteration: int | None = None)

Create a Wavefront from an openPMD file using the EXT_Wavefront extension.

Parameters:

file (Union[Path, str, Group]) –

Path to the openPMD HDF5 file, or an open h5py.Group to use as the series root.
iteration (int, default: None ) –

openPMD iteration to read. If None, the file must hold exactly one.

Returns:

Wavefront –

A new instance holding the field data, with the extension's zCoordinate in s_position and any other record attributes carried in .attrs.

Raises:

ValueError –

If the file is not a valid path, string, or h5py.Group; if a required attribute is missing; or if the file holds a frequency-domain or k-space field, neither of which this class can represent.

Notes

Any axisLabels permutation of (x, y, z) is honored on read, though the writer always stores (z, y, x).

Source code in beamphysics/wavefront/wavefront.py

@classmethod
def from_openpmd(
    cls,
    file: Union[pathlib.Path, str, h5py.Group],
    iteration: int | None = None,
):
    """
    Create a Wavefront from an openPMD file using the EXT_Wavefront extension.

    Parameters
    ----------
    file : Union[pathlib.Path, str, h5py.Group]
        Path to the openPMD HDF5 file, or an open `h5py.Group` to use as the
        series root.
    iteration : int, optional
        openPMD iteration to read. If None, the file must hold exactly one.

    Returns
    -------
    Wavefront
        A new instance holding the field data, with the extension's
        `zCoordinate` in `s_position` and any other record attributes carried
        in `.attrs`.

    Raises
    ------
    ValueError
        If the file is not a valid path, string, or `h5py.Group`; if a required
        attribute is missing; or if the file holds a frequency-domain or k-space
        field, neither of which this class can represent.

    Notes
    -----
    Any `axisLabels` permutation of (x, y, z) is honored on read, though the
    writer always stores (z, y, x).
    """
    if isinstance(file, (str, pathlib.Path)):
        with h5py.File(file, "r") as h5:
            kwargs = load_wavefront_openpmd(h5, iteration=iteration)
    elif isinstance(file, h5py.Group):
        kwargs = load_wavefront_openpmd(file, iteration=iteration)
    else:
        raise ValueError(
            f"file must be a str, pathlib.Path, or h5py.Group, got {type(file)}"
        )

    return cls(**kwargs)

pad

pad(nx=(0, 0), ny=(0, 0), nz=(0, 0))

zero-pad the field arrays.

If only a single number is given, the array will be symmetrically padded by this

Source code in beamphysics/wavefront/wavefront.py

def pad(self, nx=(0, 0), ny=(0, 0), nz=(0, 0)):
    """
    zero-pad the field arrays.

    If only a single number is given, the array will be symmetrically padded
    by this
    """

    nx = (nx, nx) if np.isscalar(nx) else nx
    ny = (ny, ny) if np.isscalar(ny) else ny
    nz = (nz, nz) if np.isscalar(nz) else nz

    Ex = np.pad(self.Ex, (nx, ny, nz)) if self.Ex is not None else None
    Ey = np.pad(self.Ey, (nx, ny, nz)) if self.Ey is not None else None

    # Padding adds elements, which is a removal of a negative count.
    mids = self._resized_mids((-nx[0], -nx[1]), (-ny[0], -ny[1]), (-nz[0], -nz[1]))
    return replace(self, Ex=Ex, Ey=Ey, attrs=self.attrs.copy(), **mids)

plot2

plot2(cmap='inferno', logscale=False)

Simple fluence plot

Notes

This is experimental.

Source code in beamphysics/wavefront/wavefront.py

def plot2(self, cmap="inferno", logscale=False):
    """
    Simple fluence plot

    Notes
    -----
    This is experimental.
    """

    # xlabel = r"$x$ (cm)"
    # ylabel = r"$y$ (cm)"
    xfactor = 100
    yfactor = 100
    zfactor = 1 / (100 * 100)  # 1/m^2 -> 1/cm^2
    # label = r"$F$ (J/cm$^2$)"
    F = self.fluence

    plot_2d_density_with_marginals(
        F * zfactor,
        dx=self.dx * xfactor,
        dy=self.dy * yfactor,
        xmin=self.xmin * xfactor,
        ymin=self.ymin * yfactor,
        x_name=r"$x$",
        x_units="cm",
        y_name=r"$y$",
        y_units="cm",
        z_name=r"$F$",
        z_units="J/cm$^2$",
        cmap=cmap,
        log_scale_marginals=logscale,
        log_scale_z=logscale,
    )

plot_fluence

plot_fluence(cmap='inferno', logscale=False)

Simple fluence plot

Source code in beamphysics/wavefront/wavefront.py

def plot_fluence(self, cmap="inferno", logscale=False):
    """
    Simple fluence plot

    """

    xlabel = r"$x$ (cm)"
    ylabel = r"$y$ (cm)"
    xfactor = 100
    yfactor = 100
    zfactor = 1 / (100 * 100)  # 1/m^2 -> 1/cm^2
    label = r"$F$ (J/cm$^2$)"
    extent = (
        self.xmin * xfactor,
        self.xmax * xfactor,
        self.ymin * yfactor,
        self.ymax * yfactor,
    )

    F = self.fluence * zfactor
    Fmax = np.max(F)

    fig, ax = plt.subplots()
    im = ax.imshow(
        F.T,
        cmap=cmap,
        extent=extent,
        origin="lower",
    )  # Note data.T and origin='lower' are required
    if logscale:
        norm = LogNorm(vmin=Fmax / 1e6, vmax=Fmax)
        im.set_norm(norm)

    fig.colorbar(im, ax=ax, label=label)

    ax.set_xlabel(xlabel)
    ax.set_ylabel(ylabel)

to_kspace

to_kspace(backend=numpy, *, inplace: bool = False) -> WavefrontK

Transform to k-space according to the Fourier transform convention:

Ẽ(kx,ky,kz) = 1/(2π)^(2/3) ∫∫∫ E(x,y,z) exp(-i kx x) exp(-i ky y) exp(-i kz z) dx dy dz

E(x,y,z) = 1/(2π)^(2/3) ∫∫∫ Ẽ(kx, ky, kz) exp(i kx x) exp(i ky y) exp(i kz z) dkx dky dkz

This ensures that that Plancherel's theorm applies without any addition factors ∫∫∫ |E(x, y, z)|^2 dx dy dz = ∫∫∫ |Ẽ(kx, ky, kz)|^2 dkx dky dkz

And therefore ϵ0/2 |Ẽ|^2 can be directly interpreted as the spectral energy density.

This is realized with factors so that: Ẽ = [dx sqrt(nx/(2π))][beamphysics.wavefront.wavefront] [...] fftn(E, norm = 'ortho') in units of V * m^2

E = [dkx sqrt(nx/(2π))][beamphysics.wavefront.wavefront] [...] fftn(Ẽ, norm = 'ortho') in units of V/m

Here [...] is similar for y, z, etc. We use norm='ortho' to simplify the symmetry in the code.

Source code in beamphysics/wavefront/wavefront.py

def to_kspace(self, backend=np, *, inplace: bool = False) -> WavefrontK:
    """
    Transform to k-space according to the Fourier transform convention:

    Ẽ(kx,ky,kz) = 1/(2π)^(2/3) ∫∫∫ E(x,y,z) exp(-i kx x) exp(-i ky y) exp(-i kz z) dx dy dz

    E(x,y,z)   =  1/(2π)^(2/3) ∫∫∫ Ẽ(kx, ky, kz) exp(i kx x) exp(i ky y) exp(i kz z) dkx dky dkz

    This ensures that that Plancherel's theorm applies without any addition factors
        ∫∫∫ |E(x, y, z)|^2 dx dy dz =  ∫∫∫ |Ẽ(kx, ky, kz)|^2 dkx dky dkz

    And therefore ϵ0/2 |Ẽ|^2 can be directly interpreted as the spectral energy density.

    This is realized with factors so that:
    Ẽ = [dx sqrt(nx/(2π))] [...] [...]  fftn(E, norm = 'ortho') in units of V * m^2

    E = [dkx sqrt(nx/(2π))] [...] [...]  fftn(Ẽ, norm = 'ortho') in units of V/m

    Here  [...] is similar for y, z, etc. We use `norm='ortho'` to simplify the symmetry in the code.

    """

    if inplace:
        raise NotImplementedError("inplace not yet implemented")

    fftn = backend.fft.fftn
    fftshift = backend.fft.fftshift
    ifftshift = backend.fft.ifftshift

    # Normalized for the Plancherel theorem (see energy def)
    norm = (
        self.dx
        * self.dy
        * self.dz
        * np.sqrt(self.nx * self.ny * self.nz / (2 * pi) ** 3)
    )
    Ex = (
        fftshift(fftn(ifftshift(self.Ex), norm="ortho")) * norm
        if self.Ex is not None
        else None
    )
    Ey = (
        fftshift(fftn(ifftshift(self.Ey), norm="ortho")) * norm
        if self.Ey is not None
        else None
    )

    return WavefrontK(
        Ex=Ex,
        Ey=Ey,
        dkx=self.dkx,
        dky=self.dky,
        dkz=self.dkz,
        wavelength=self.wavelength,
        attrs=self.attrs.copy(),
        xmid=self.xmid,
        ymid=self.ymid,
        zmid=self.zmid,
        s_position=self.s_position,
    )

write_genesis4

write_genesis4(file: Path | str | Group)

Write the Wavefront field data to a Genesis4-style HDF5 file.

This method saves the wavefront field data into a Genesis4-compatible HDF5 file, either by creating a new file on disk or writing into an existing HDF5 group.

Parameters:

file (Union[Path, str, Group]) –

The path to the output HDF5 file or an open h5py.Group for writing.

Raises:

ValueError –

If the provided file is not a valid path, string, or h5py.Group.

Notes

If file is a path or string, a new HDF5 file is created and written to.
If file is an h5py.Group, the data is written directly into the provided group.
zmin is written as Genesis4's refposition. s_position and attrs have nowhere to go in this format and are dropped.

Source code in beamphysics/wavefront/wavefront.py

def write_genesis4(
    self,
    file: pathlib.Path | str | h5py.Group,
):
    """
    Write the Wavefront field data to a Genesis4-style HDF5 file.

    This method saves the wavefront field data into a Genesis4-compatible HDF5 file,
    either by creating a new file on disk or writing into an existing HDF5 group.

    Parameters
    ----------
    file : Union[pathlib.Path, str, h5py.Group]
        The path to the output HDF5 file or an open `h5py.Group` for writing.

    Raises
    ------
    ValueError
        If the provided file is not a valid path, string, or `h5py.Group`.

    Notes
    -----
    - If `file` is a path or string, a new HDF5 file is created and written to.
    - If `file` is an `h5py.Group`, the data is written directly into the provided group.
    - `zmin` is written as Genesis4's `refposition`. `s_position` and `attrs`
      have nowhere to go in this format and are dropped.
    """
    if isinstance(file, (str, pathlib.Path)):
        with h5py.File(file, "w") as h5:
            wavefront_write_genesis4(self, h5)
        return
    elif isinstance(file, h5py.Group):
        wavefront_write_genesis4(self, file)
        return
    else:
        raise ValueError(
            f"file must be a str, pathlib.Path, or h5py.Group, got {type(file)}"
        )

write_openpmd

write_openpmd(file: Path | str | Group, iteration: int = 1, **extension_attrs)

Write the Wavefront as an openPMD file using the EXT_Wavefront extension.

Parameters:

file (Union[Path, str, Group]) –

Path for a new HDF5 file, or an open h5py.Group to use as the series root.
iteration (int, default: 1 ) –

openPMD iteration index to write.
**extension_attrs –

WavefrontAttrs field names (beamline, radius_of_curvature_x, radius_of_curvature_y, delta_radius_of_curvature_x, delta_radius_of_curvature_y), taking precedence over .attrs.

Raises:

ValueError –

If the file is not a valid path, string, or h5py.Group.
TypeError –

If an attribute name is not a WavefrontAttrs field.

Notes

Both transverse polarizations are written into one file as components x and y of the electricField record; the z component is never written.

The required zCoordinate is written from s_position. There is no write-time override: it is a coordinate the propagators maintain, so a keyword here could put a value in the file that disagrees with the wavefront it came from. Set s_position on the wavefront instead.

Source code in beamphysics/wavefront/wavefront.py

def write_openpmd(
    self,
    file: pathlib.Path | str | h5py.Group,
    iteration: int = 1,
    **extension_attrs,
):
    """
    Write the Wavefront as an openPMD file using the EXT_Wavefront extension.

    Parameters
    ----------
    file : Union[pathlib.Path, str, h5py.Group]
        Path for a new HDF5 file, or an open `h5py.Group` to use as the series
        root.
    iteration : int, default=1
        openPMD iteration index to write.
    **extension_attrs
        `WavefrontAttrs` field names (`beamline`, `radius_of_curvature_x`,
        `radius_of_curvature_y`, `delta_radius_of_curvature_x`,
        `delta_radius_of_curvature_y`), taking precedence over `.attrs`.

    Raises
    ------
    ValueError
        If the file is not a valid path, string, or `h5py.Group`.
    TypeError
        If an attribute name is not a `WavefrontAttrs` field.

    Notes
    -----
    Both transverse polarizations are written into one file as components `x`
    and `y` of the `electricField` record; the `z` component is never written.

    The required `zCoordinate` is written from `s_position`. There is no
    write-time override: it is a coordinate the propagators maintain, so a
    keyword here could put a value in the file that disagrees with the
    wavefront it came from. Set `s_position` on the wavefront instead.
    """
    if isinstance(file, (str, pathlib.Path)):
        with h5py.File(file, "w") as h5:
            write_wavefront_openpmd(
                self,
                h5,
                iteration=iteration,
                **extension_attrs,
            )
        return
    elif isinstance(file, h5py.Group):
        write_wavefront_openpmd(
            self,
            file,
            iteration=iteration,
            **extension_attrs,
        )
        return
    else:
        raise ValueError(
            f"file must be a str, pathlib.Path, or h5py.Group, got {type(file)}"
        )

WavefrontK `dataclass`

WavefrontK(Ex: ndarray | None = None, Ey: ndarray | None = None, dkx: float = 1.0, dky: float = 1.0, dkz: float = 1.0, wavelength: float = 1.0, attrs: WavefrontAttrs = WavefrontAttrs(), xmid: float = 0.0, ymid: float = 0.0, zmid: float = 0.0, s_position: float = 0.0)

Bases: WavefrontBase


              flowchart TD
              beamphysics.WavefrontK[WavefrontK]
              beamphysics.wavefront.wavefront.WavefrontBase[WavefrontBase]

                              beamphysics.wavefront.wavefront.WavefrontBase --> beamphysics.WavefrontK
                


              click beamphysics.WavefrontK href "" "beamphysics.WavefrontK"
              click beamphysics.wavefront.wavefront.WavefrontBase href "" "beamphysics.wavefront.wavefront.WavefrontBase"

K-space (Fourier) representation of electromagnetic wavefront fields.

This class represents electromagnetic fields in k-space (spatial frequency domain) for spectral analysis and propagation. It stores complex electric field Fourier components Ẽx and Ẽy on a 3D reciprocal-space grid with uniform spacing.

The Fourier transform convention follows: Ẽ(kx,ky,kz) = (2π)^(-3/2) ∫∫∫ E(x,y,z) exp(-i·k·r) dx dy dz

Parameters:

Ex (ndarray, default: None ) –

x-polarized electric field Fourier component in V·m², shape (nx, ny, nz)
Ey (ndarray, default: None ) –

y-polarized electric field Fourier component in V·m², shape (nx, ny, nz)
dkx (float, default: 1.0 ) –

Grid spacing in kx direction (rad/m)
dky (float, default: 1.0 ) –

Grid spacing in ky direction (rad/m)
dkz (float, default: 1.0 ) –

Grid spacing in kz direction (rad/m)
wavelength (float, default: 1.0 ) –

Central wavelength (m)
attrs (WavefrontAttrs or mapping, default: WavefrontAttrs() ) –

openPMD EXT_Wavefront attributes carried through file I/O, such as beamline. Not interpreted by this class. A mapping keyed by either openPMD or field names is coerced to WavefrontAttrs.
xmid (float, default: 0.0 ) –

Midpoint of the real-space grid in m, carried through from the real-space wavefront. Not interpreted here: a real-space shift is a linear phase in k-space, not a shift of the k grid.
ymid (float, default: 0.0 ) –

Midpoint of the real-space grid in m, carried through from the real-space wavefront. Not interpreted here: a real-space shift is a linear phase in k-space, not a shift of the k grid.
zmid (float, default: 0.0 ) –

Midpoint of the real-space grid in m, carried through from the real-space wavefront. Not interpreted here: a real-space shift is a linear phase in k-space, not a shift of the k grid.
s_position (float, default: 0.0 ) –

Beamline position of this plane in m, carried through from the real-space wavefront.

Attributes:

spatial_domain (SpatialDomain) –

Always returns SpatialDomain.K (k-space)
intensity (ndarray) –

Spectral intensity |Ẽ|² (arbitrary units)
spectral_energy_density (ndarray) –

Spectral energy density ε₀/2 |Ẽ|² in J·m³
energy (float) –

Total energy ε₀/2 ∫∫∫ |Ẽ|² dkx dky dkz in J
spectral_fluence (ndarray) –

K-space fluence F(kx,ky) = ε₀/2 ∫ |Ẽ|² dkz in J·m²
photon_energy_vec (ndarray) –

Photon energy axis from kz in eV
photon_energy_spectrum (ndarray) –

Photon spectral energy density dU/dE in J/eV
k0 (float) –

Central wavenumber 2π/λ in rad/m
photon_energy (float) –

Central photon energy ℏck₀ in eV

Methods:

to_rspace –

Transform to real-space representation (Wavefront)
plot_spectral_intensity –

Plot angular spectrum intensity F(θx, θy)
plot_photon_energy_spectrum –

Plot photon energy spectrum dU/dE
pad –

Zero-pad the field arrays
crop –

Crop the field arrays by removing elements from edges
auto_crop –

Determine and optionally apply symmetric crop based on intensity thresholds

Notes

At least one of Ex or Ey must be provided
All spacing parameters (dkx, dky, dkz, wavelength) must be positive
Field arrays must be complex dtype
The spectral energy density ε₀/2|Ẽ|² preserves total energy with real space
Derived real-space grid spacing: dx = 2π/(nx·dkx)
Angular coordinates: θx = kx/k₀, θy = ky/k₀
Statistical moments (mean, sigma) are intensity-weighted in k-space

Returned by:

API Wavefront
- Wavefront
  - crop
  - to_kspace
- WavefrontK crop

dthetax `property`

dthetax

Angular spacing in thetax direction in rad.

dthetax = wavelength / (nx * dx) = dkx / k0

dthetay `property`

dthetay

Angular spacing in thetay direction in rad.

dthetay = wavelength / (ny * dy) = dky / k0

dx `property`

dx: float

grid spacing in m dx = 2π / (nx * dkx)

dy `property`

dy

grid spacing in m dy = 2π / (nx * dky)

dz `property`

dz

grid spacing in m dz = 2π / (nx * dkz)

energy `property`

energy

Total (time-averaged) energy in J

energy = ϵ0/2 ∫∫∫ |E|^2 dx dy dz

   = ϵ0/2  ∫∫∫ |Ẽ|^2 dkx dky dkz

in_kspace `property`

in_kspace

True if the wavefront is in k-space (Fourier domain).

in_rspace `property`

in_rspace

True if the wavefront is in real space (position domain).

intensity `property`

intensity: ndarray | float

Total spectral intensity in V²·m⁴.

|Ẽx|² + |Ẽy|²

intensity_x `property`

intensity_x

x-polarization spectral intensity in V²·m⁴.

Spectral intensity |Ẽx|²

intensity_y `property`

intensity_y

y-polarization spectral intensity in V²·m⁴.

Spectral intensity |Ẽy|²

k0 `property`

k0

Central wavenumber in rad/m.

k0 = 2π / wavelength

kxmax `property`

kxmax

Maximum kx value in rad/m.

kxmin `property`

kxmin

Minimum kx value in rad/m.

kxvec `property`

kxvec

Array of kx (transverse wavenumber) values in rad/m.

kymax `property`

kymax

Maximum ky value in rad/m.

kymin `property`

kymin

Minimum ky value in rad/m.

kyvec `property`

kyvec

Array of ky (transverse wavenumber) values in rad/m.

kzmax `property`

kzmax

Maximum kz value in rad/m.

kzmin `property`

kzmin

Minimum kz value in rad/m.

kzvec `property`

kzvec

Array of kz (longitudinal wavenumber) values in rad/m.

mean_kx `property`

mean_kx

Intensity-weighted mean kx in rad/m.

= ∫ kx |Ẽ|² dkx dky dkz / ∫ |Ẽ|² dkx dky dkz

mean_ky `property`

mean_ky

Intensity-weighted mean ky in rad/m.

= ∫ ky |Ẽ|² dkx dky dkz / ∫ |Ẽ|² dkx dky dkz

mean_kz `property`

mean_kz

Intensity-weighted mean kz in rad/m.

= ∫ kz |Ẽ|² dkx dky dkz / ∫ |Ẽ|² dkx dky dkz

mean_thetax `property`

mean_thetax

Intensity-weighted mean thetax in rad.

<θx> = / k0

mean_thetay `property`

mean_thetay

Intensity-weighted mean thetay in rad.

<θy> = / k0

nx `property`

nx

Number of grid points in the x direction.

ny `property`

ny

Number of grid points in the y direction.

nz `property`

nz

Number of grid points in the z direction.

photon_energy `property`

photon_energy: float

Central photon energy in eV

photon_energy_spectrum `property`

photon_energy_spectrum: ndarray

Photon energy spectrum dU/dE in J/eV

dU/dE = (ε0/2) ∫∫ |Ẽ|² dkx dky / (ħc/e)

See .photon_energy_vec for the corresponding photon energies in eV.

photon_energy_vec `property`

photon_energy_vec: ndarray

Photon energy axis from kz in eV

shape `property`

shape

Shape of the field arrays as (nx, ny, nz).

sigma_kx `property`

sigma_kx

RMS width in kx in rad/m.

σ_kx = sqrt( - ²)

sigma_ky `property`

sigma_ky

RMS width in ky in rad/m.

σ_ky = sqrt( - ²)

sigma_kz `property`

sigma_kz

RMS width in kz in rad/m.

σ_kz = sqrt( - ²)

sigma_thetax `property`

sigma_thetax

RMS angular width in thetax in rad.

σ_θx = sqrt(<θx²> - <θx>²)

sigma_thetay `property`

sigma_thetay

RMS angular width in thetay in rad.

σ_θy = sqrt(<θy²> - <θy>²)

spectral_energy_density `property`

spectral_energy_density

Spectral energy density ϵ0/2 |Ẽ|^2 in J * m^3

3D real array

spectral_fluence `property`

spectral_fluence

K-space Fluence in the z direction:

F(kx,ky) = ϵ0/2 ∫∫∫ |Ẽ(kx,ky,kz)|^2 dkz in J * m^2

thetaxmax `property`

thetaxmax

Maximum thetax value in rad.

thetaxmin `property`

thetaxmin

Minimum thetax value in rad.

thetaxvec `property`

thetaxvec

Array of thetax (angular deviation) values in rad.

thetax = kx / k0

thetaymax `property`

thetaymax

Maximum thetay value in rad.

thetaymin `property`

thetaymin

Minimum thetay value in rad.

thetayvec `property`

thetayvec

Array of thetay (angular deviation) values in rad.

thetay = ky / k0

xmax `property`

xmax

Maximum x coordinate in m.

xmax = xmid + (nx-1) * dx / 2

xmin `property`

xmin

Minimum x coordinate in m.

xmin = xmid - (nx-1) * dx / 2

xvec `property`

xvec

Array of x coordinates in m.

ymax `property`

ymax

Maximum y coordinate in m.

ymax = ymid + (ny-1) * dy / 2

ymin `property`

ymin

Minimum y coordinate in m.

ymin = ymid - (ny-1) * dy / 2

yvec `property`

yvec

Array of y coordinates in m.

zmax `property`

zmax

Maximum z coordinate in m.

zmax = zmid + (nz-1) * dz / 2

zmin `property`

zmin

Minimum z coordinate in m.

zmin = zmid - (nz-1) * dz / 2

zvec `property`

zvec

Array of z coordinates in m.

__post_init__

__post_init__()

Validate inputs after dataclass initialization

Source code in beamphysics/wavefront/wavefront.py

def __post_init__(self):
    """
    Validate inputs after dataclass initialization
    """
    # Validate field shapes match
    if self.Ex is not None and self.Ey is not None:
        if self.Ex.shape != self.Ey.shape:
            raise ValueError(
                f"Ex shape {self.Ex.shape} != Ey shape {self.Ey.shape}"
            )

    # Validate at least one field exists
    if self.Ex is None and self.Ey is None:
        raise ValueError("At least one of Ex or Ey must be provided")

    # Validate positive spacing
    for attr in ["dx", "dy", "dz", "wavelength"]:
        val = getattr(self, attr)
        if val <= 0:
            raise ValueError(f"{attr} must be positive, got {val}")

setattr

__setattr__(name, value)

Keep attrs a WavefrontAttrs however it is set.

A mapping keyed by either openPMD or field names is accepted for convenience, but the attribute itself is always the typed object. Coercing only in __post_init__ would leave a plain assignment such as w.attrs = {...} holding a raw dict, so that w.attrs.beamline raises AttributeError and the mapping's contents go unvalidated until write time.

Source code in beamphysics/wavefront/wavefront.py

def __setattr__(self, name, value):
    """
    Keep `attrs` a `WavefrontAttrs` however it is set.

    A mapping keyed by either openPMD or field names is accepted for
    convenience, but the attribute itself is always the typed object. Coercing
    only in `__post_init__` would leave a plain assignment such as
    ``w.attrs = {...}`` holding a raw dict, so that `w.attrs.beamline`
    raises `AttributeError` and the mapping's contents go unvalidated until
    write time.
    """
    if name == "attrs" and not isinstance(value, WavefrontAttrs):
        value = WavefrontAttrs.from_pmd(value)
    object.__setattr__(self, name, value)

auto_crop

auto_crop(threshold: float = 1e-06, apply: bool = True)

Determine symmetric crop values based on intensity profile thresholds.

Analyzes the 1D intensity projections along each axis to find regions where the signal falls below a relative threshold, then computes symmetric crop amounts for each axis.

Parameters:

threshold (float, default: 1e-6 ) –

Relative threshold (fraction of maximum) below which data is considered negligible and can be cropped.
apply (bool, default: True ) –

If True, apply the crop and return the cropped wavefront. If False, return only the crop amounts as a dictionary.

Returns:

Wavefront, WavefrontK, or dict –

If apply=True: New wavefront with cropped field arrays. If apply=False: Dictionary with keys 'nx', 'ny', 'nz' containing the symmetric crop amounts for each axis.

Examples:

>>> from beamphysics.wavefront import Wavefront
>>> w = Wavefront.from_gaussian(
...     shape=(64, 64, 16), dx=1e-6, dy=1e-6, dz=1e-6,
...     wavelength=1e-9, sigma0=3e-6, sigma_z=2e-6,
... )

Get the crop amounts without applying them:

>>> crop_info = w.auto_crop(threshold=1e-4, apply=False)
>>> crop_info
{'nx': 19, 'ny': 19, 'nz': 0}

Apply the crop directly, removing that many elements from each end:

>>> w.auto_crop(threshold=1e-4).shape
(26, 26, 16)

The two steps are equivalent, so the amounts can be inspected or adjusted first:

>>> w.crop(**crop_info).shape
(26, 26, 16)

Source code in beamphysics/wavefront/wavefront.py

def auto_crop(self, threshold: float = 1e-6, apply: bool = True):
    """
    Determine symmetric crop values based on intensity profile thresholds.

    Analyzes the 1D intensity projections along each axis to find
    regions where the signal falls below a relative threshold, then
    computes symmetric crop amounts for each axis.

    Parameters
    ----------
    threshold : float, default=1e-6
        Relative threshold (fraction of maximum) below which data
        is considered negligible and can be cropped.
    apply : bool, default=True
        If True, apply the crop and return the cropped wavefront.
        If False, return only the crop amounts as a dictionary.

    Returns
    -------
    Wavefront, WavefrontK, or dict
        If apply=True: New wavefront with cropped field arrays.
        If apply=False: Dictionary with keys 'nx', 'ny', 'nz' containing
        the symmetric crop amounts for each axis.

    Examples
    --------
    >>> from beamphysics.wavefront import Wavefront
    >>> w = Wavefront.from_gaussian(
    ...     shape=(64, 64, 16), dx=1e-6, dy=1e-6, dz=1e-6,
    ...     wavelength=1e-9, sigma0=3e-6, sigma_z=2e-6,
    ... )

    Get the crop amounts without applying them:

    >>> crop_info = w.auto_crop(threshold=1e-4, apply=False)
    >>> crop_info
    {'nx': 19, 'ny': 19, 'nz': 0}

    Apply the crop directly, removing that many elements from each end:

    >>> w.auto_crop(threshold=1e-4).shape
    (26, 26, 16)

    The two steps are equivalent, so the amounts can be inspected or adjusted
    first:

    >>> w.crop(**crop_info).shape
    (26, 26, 16)
    """

    def find_symmetric_crop(profile: np.ndarray, threshold: float) -> int:
        """Find symmetric crop amount for a 1D profile."""
        if profile.size == 0:
            return 0

        max_val = np.max(profile)
        if max_val == 0:
            return 0

        # Normalize and find where signal is above threshold
        normalized = profile / max_val
        above_threshold = normalized > threshold

        if not np.any(above_threshold):
            return 0

        # Find first and last indices above threshold
        indices = np.where(above_threshold)[0]
        first_idx = indices[0]
        last_idx = indices[-1]

        # Calculate symmetric crop (use minimum of both sides)
        crop_start = first_idx
        crop_end = len(profile) - 1 - last_idx

        # A plain int, so that the returned dict reprs and serializes cleanly.
        return int(min(crop_start, crop_end))

    # Get 1D intensity projections along each axis
    intensity = self.intensity
    profile_x = np.sum(intensity, axis=_axis_for_sum[self.axis_labels[0]])
    profile_y = np.sum(intensity, axis=_axis_for_sum[self.axis_labels[1]])
    profile_z = np.sum(intensity, axis=_axis_for_sum[self.axis_labels[2]])

    # Analyze each axis
    nx = find_symmetric_crop(profile_x, threshold)
    ny = find_symmetric_crop(profile_y, threshold)
    nz = find_symmetric_crop(profile_z, threshold)

    crop_amounts = {"nx": nx, "ny": ny, "nz": nz}

    if apply:
        return self.crop(**crop_amounts)
    else:
        return crop_amounts

axis_index

axis_index(key)

Returns axis index for a named axis label key.

Examples:

.axis_labels == ('x', 'y', 'z')
.axis_index('z') returns 2

Source code in beamphysics/wavefront/wavefront.py

def axis_index(self, key):
    """
    Returns axis index for a named axis label key.

    Examples:

    - `.axis_labels == ('x', 'y', 'z')`
    - `.axis_index('z')` returns `2`
    """
    return self.axis_labels.index(key)

copy

copy()

Returns a deep copy

Source code in beamphysics/wavefront/wavefront.py

def copy(self):
    """Returns a deep copy"""
    return deepcopy(self)

crop

crop(nx=(0, 0), ny=(0, 0), nz=(0, 0))

Crop the field arrays by removing elements from edges.

Parameters:

nx (int or tuple of (int, int), default: (0, 0) ) –

Number of elements to remove from (start, end) of first axis. If a single int, removes symmetrically from both sides.
ny (int or tuple of (int, int), default: (0, 0) ) –

Number of elements to remove from (start, end) of second axis. If a single int, removes symmetrically from both sides.
nz (int or tuple of (int, int), default: (0, 0) ) –

Number of elements to remove from (start, end) of third axis. If a single int, removes symmetrically from both sides.

Returns:

Wavefront or WavefrontK –

New wavefront with cropped field arrays.

Raises:

ValueError –

If crop amounts exceed array dimensions.

Notes

The surviving samples keep their physical coordinates: an asymmetric crop moves the grid origin rather than re-centering the grid.

Source code in beamphysics/wavefront/wavefront.py

def crop(self, nx=(0, 0), ny=(0, 0), nz=(0, 0)):
    """
    Crop the field arrays by removing elements from edges.

    Parameters
    ----------
    nx : int or tuple of (int, int), default=(0, 0)
        Number of elements to remove from (start, end) of first axis.
        If a single int, removes symmetrically from both sides.
    ny : int or tuple of (int, int), default=(0, 0)
        Number of elements to remove from (start, end) of second axis.
        If a single int, removes symmetrically from both sides.
    nz : int or tuple of (int, int), default=(0, 0)
        Number of elements to remove from (start, end) of third axis.
        If a single int, removes symmetrically from both sides.

    Returns
    -------
    Wavefront or WavefrontK
        New wavefront with cropped field arrays.

    Raises
    ------
    ValueError
        If crop amounts exceed array dimensions.

    Notes
    -----
    The surviving samples keep their physical coordinates: an asymmetric crop
    moves the grid origin rather than re-centering the grid.
    """
    nx = (nx, nx) if np.isscalar(nx) else nx
    ny = (ny, ny) if np.isscalar(ny) else ny
    nz = (nz, nz) if np.isscalar(nz) else nz

    # Calculate slice indices
    x_start = nx[0]
    x_end = self.shape[0] - nx[1] if nx[1] > 0 else None
    y_start = ny[0]
    y_end = self.shape[1] - ny[1] if ny[1] > 0 else None
    z_start = nz[0]
    z_end = self.shape[2] - nz[1] if nz[1] > 0 else None

    # Validate crop amounts
    if nx[0] + nx[1] >= self.shape[0]:
        raise ValueError(
            f"Crop in x ({nx}) exceeds array dimension ({self.shape[0]})"
        )
    if ny[0] + ny[1] >= self.shape[1]:
        raise ValueError(
            f"Crop in y ({ny}) exceeds array dimension ({self.shape[1]})"
        )
    if nz[0] + nz[1] >= self.shape[2]:
        raise ValueError(
            f"Crop in z ({nz}) exceeds array dimension ({self.shape[2]})"
        )

    Ex = (
        self.Ex[x_start:x_end, y_start:y_end, z_start:z_end]
        if self.Ex is not None
        else None
    )
    Ey = (
        self.Ey[x_start:x_end, y_start:y_end, z_start:z_end]
        if self.Ey is not None
        else None
    )

    return replace(
        self,
        Ex=Ex,
        Ey=Ey,
        attrs=self.attrs.copy(),
        **self._resized_mids(nx, ny, nz),
    )

pad

pad(nx=(0, 0), ny=(0, 0), nz=(0, 0))

zero-pad the field arrays.

If only a single number is given, the array will be symmetrically padded by this

Source code in beamphysics/wavefront/wavefront.py

def pad(self, nx=(0, 0), ny=(0, 0), nz=(0, 0)):
    """
    zero-pad the field arrays.

    If only a single number is given, the array will be symmetrically padded
    by this
    """

    nx = (nx, nx) if np.isscalar(nx) else nx
    ny = (ny, ny) if np.isscalar(ny) else ny
    nz = (nz, nz) if np.isscalar(nz) else nz

    Ex = np.pad(self.Ex, (nx, ny, nz)) if self.Ex is not None else None
    Ey = np.pad(self.Ey, (nx, ny, nz)) if self.Ey is not None else None

    # Padding adds elements, which is a removal of a negative count.
    mids = self._resized_mids((-nx[0], -nx[1]), (-ny[0], -ny[1]), (-nz[0], -nz[1]))
    return replace(self, Ex=Ex, Ey=Ey, attrs=self.attrs.copy(), **mids)

plot_spectral_intensity

plot_spectral_intensity(cmap='inferno', logscale=False)

Simple projected intensity plot

Source code in beamphysics/wavefront/wavefront.py

def plot_spectral_intensity(self, cmap="inferno", logscale=False):
    """
    Simple projected intensity plot

    """

    xlabel = r"$\theta_x$ (µrad)"
    ylabel = r"$\theta_y$ (µrad)"
    xfactor = 1e6
    yfactor = 1e6
    zfactor = self.k0**2 / (1e6 * 1e6)
    label = r"Spectral $F$ (J/µrad$^2$)"

    extent = (
        self.thetaxmin * xfactor,
        self.thetaxmax * xfactor,
        self.thetaymin * yfactor,
        self.thetaymax * yfactor,
    )

    # Alternatively:
    # extent = (self.kxmin, self.kxmax, self.kymin, self.kymax)
    # xlabel = r'$k_x$ (rad/m)'
    # ylabel = r'$k_y$ (rad/m)'
    # zfactor = 1
    # label = r"Spectral $F$ (J$\cdot$m$^2$)"

    F = zfactor * self.spectral_fluence
    Fmax = np.max(F)

    fig, ax = plt.subplots()
    im = ax.imshow(
        F.T,
        cmap=cmap,
        extent=extent,
        origin="lower",
    )  # Note data.T and origin='lower' are required
    if logscale:
        norm = LogNorm(vmin=Fmax / 1e6, vmax=Fmax)
        im.set_norm(norm)

    fig.colorbar(im, ax=ax, label=label)

    ax.set_xlabel(xlabel)
    ax.set_ylabel(ylabel)

to_rspace

to_rspace(*, inplace: bool = False) -> Wavefront

See Wavefront.to_kspace()

Source code in beamphysics/wavefront/wavefront.py

def to_rspace(self, *, inplace: bool = False) -> Wavefront:
    """
    See Wavefront.to_kspace()
    """
    if inplace:
        raise NotImplementedError("inplace not yet implemented")

    # Normalized for the Plancherel theorem (see def energy)
    norm = (
        self.dx
        * self.dy
        * self.dz
        * np.sqrt(self.nx * self.ny * self.nz / (2 * pi) ** 3)
    )
    Ex = (
        fftshift(ifftn(ifftshift(self.Ex), norm="ortho")) / norm
        if self.Ex is not None
        else None
    )
    Ey = (
        fftshift(ifftn(ifftshift(self.Ey), norm="ortho")) / norm
        if self.Ey is not None
        else None
    )

    return Wavefront(
        Ex=Ex,
        Ey=Ey,
        dx=self.dx,
        dy=self.dy,
        dz=self.dz,
        wavelength=self.wavelength,
        attrs=self.attrs.copy(),
        xmid=self.xmid,
        ymid=self.ymid,
        zmid=self.zmid,
        s_position=self.s_position,
    )

openpmd

openPMD EXT_Wavefront reading and writing for the Wavefront class.

The layout below follows the Wavefront extension of the openPMD standard (branch upcoming-2.0.0), which is authoritative::

/                             openPMD "2.0.0", openPMDextension "Wavefront",
                              basePath "/data/%T/", meshesPath "meshes/",
                              iterationEncoding "groupBased"
/data/<iteration>/            one iteration per file. Slices of a single pulse are
                              simultaneous, so the slice axis is a *mesh* axis and
                              never the openPMD iteration.
.../meshes/electricField      the mesh record. All required attributes live on the
                              record: the base standard's geometry, axisLabels,
                              gridSpacing, gridGlobalOffset, gridUnitSI,
                              gridUnitDimension, unitDimension and timeOffset, plus
                              the extension's photonEnergy [J], temporalDomain,
                              spatialDomain and zCoordinate.
.../electricField/x, y        complex compound ``{r, i}`` datasets in V/m, which
                              h5py maps to complex dtypes natively. The ``z``
                              component is never written: a paraxial field has none.

Datasets are stored in (z, y, x) order -- declared by axisLabels -- so that each transverse slice is one contiguous block in the file. The class's (nx, ny, nz) convention is recovered by a transpose that is applied one slice at a time, keeping peak memory at one field plus one transverse slice rather than two full copies.

Only the real-space, time-domain case is implemented; spatialDomain='k' and temporalDomain='frequency' are refused by name.

Classes:

WavefrontAttrs –

openPMD EXT_Wavefront attributes carried alongside a Wavefront.

Functions:

load_wavefront_openpmd –

Read an openPMD EXT_Wavefront series into Wavefront constructor arguments.
write_wavefront_openpmd –

Write a Wavefront into an open HDF5 group as an openPMD EXT_Wavefront series.

WavefrontAttrs `dataclass`

WavefrontAttrs(beamline: str | None = None, radius_of_curvature_x: float | None = None, radius_of_curvature_y: float | None = None, delta_radius_of_curvature_x: float | None = None, delta_radius_of_curvature_y: float | None = None, other: dict = dict())

openPMD EXT_Wavefront attributes carried alongside a Wavefront.

Field names are Python-style; the corresponding openPMD attribute names, used in the file, are camelCase and recorded in each field's metadata. Attributes the class derives from its own state -- photonEnergy, temporalDomain, spatialDomain, zCoordinate and the grid attributes -- are deliberately absent. zCoordinate in particular is Wavefront.s_position, which is a coordinate the propagators advance rather than provenance carried through I/O.

Parameters:

beamline (str, default: None ) –

Name of the beamline this wavefront belongs to.
radius_of_curvature_x (float, default: None ) –

Radius of curvature in x, in m.
radius_of_curvature_y (float, default: None ) –

Radius of curvature in y, in m.
delta_radius_of_curvature_x (float, default: None ) –

Uncertainty in radius_of_curvature_x, in m.
delta_radius_of_curvature_y (float, default: None ) –

Uncertainty in radius_of_curvature_y, in m.
other (dict, default: dict() ) –

Record attributes that are not part of the extension as this module knows it, keyed by their openPMD name. Populated on read and written back verbatim, so that a file using a newer revision of the extension survives a round trip. Names this module computes are refused. Named after the other dict that readers.load_field_attrs uses for the same purpose.

Raises:

ValueError –

If other holds a name that has a field of its own, or one that the writer computes.

Examples:

>>> attrs = WavefrontAttrs(beamline="SXR", radius_of_curvature_x=12.5)
>>> attrs.to_pmd()["radiusOfCurvatureX"]
12.5

A misspelled attribute is a TypeError at construction rather than a silent omission at write time:

>>> WavefrontAttrs(radius_of_curvature_z=1.0)
Traceback (most recent call last):
    ...
TypeError: ...

Returned by:

API Wavefront openpmd WavefrontAttrs
- copy
- from_pmd

Used by:

API Wavefront
- Wavefront
- WavefrontK

Methods:

copy –

Return an independent copy, including other.
from_pmd –

Build from a mapping keyed by openPMD attribute names.
pmd_keys –

Map field names to their openPMD attribute names.
to_pmd –

Render as a mapping keyed by openPMD attribute names.

copy

copy()

Return an independent copy, including other.

Returns:

WavefrontAttrs –

Source code in beamphysics/wavefront/openpmd.py

def copy(self):
    """
    Return an independent copy, including `other`.

    Returns
    -------
    WavefrontAttrs
    """
    return type(self)(
        **{name: getattr(self, name) for name in self.pmd_keys()},
        other=dict(self.other),
    )

from_pmd `classmethod`

from_pmd(attrs)

Build from a mapping keyed by openPMD attribute names.

Parameters:

attrs (mapping) –

Attribute names to values. Field names are also accepted, so that a dict written in Python style round-trips. A nested other mapping is merged into other rather than nested inside it. Anything unrecognized is kept in other.

Returns:

WavefrontAttrs –

Raises:

ValueError –

If an attribute is given under both its openPMD name and its field name. Silently preferring one would write the other to the file as a nonstandard attribute holding a conflicting value.

Source code in beamphysics/wavefront/openpmd.py

@classmethod
def from_pmd(cls, attrs):
    """
    Build from a mapping keyed by openPMD attribute names.

    Parameters
    ----------
    attrs : mapping
        Attribute names to values. Field names are also accepted, so that a
        dict written in Python style round-trips. A nested `other` mapping is
        merged into `other` rather than nested inside it. Anything unrecognized
        is kept in `other`.

    Returns
    -------
    WavefrontAttrs

    Raises
    ------
    ValueError
        If an attribute is given under both its openPMD name and its field
        name. Silently preferring one would write the other to the file as a
        nonstandard attribute holding a conflicting value.
    """
    if isinstance(attrs, cls):
        return attrs.copy()

    remaining = dict(attrs)
    kwargs = {}
    for name, key in cls.pmd_keys().items():
        has_key = key in remaining
        # Some fields spell the two the same way, which is not a collision.
        has_name = name != key and name in remaining
        if has_key and has_name:
            raise ValueError(
                f"{key!r} and {name!r} are two spellings of the same attribute "
                f"and both were given, with values {remaining[key]!r} and "
                f"{remaining[name]!r}"
            )
        if has_key:
            kwargs[name] = remaining.pop(key)
        elif has_name:
            kwargs[name] = remaining.pop(name)

    # `other` is a field of this class, not a record attribute, so a mapping
    # carrying one means Python-style input rather than something read from a
    # file. Merge it instead of nesting it, which would produce an `other`
    # entry whose value is a dict and fail at write time.
    other = dict(remaining.pop("other", {}))
    other.update(remaining)

    return cls(**kwargs, other=other)

pmd_keys `classmethod`

pmd_keys()

Map field names to their openPMD attribute names.

Returns:

dict –

{field_name: pmd_key}, excluding other.

Source code in beamphysics/wavefront/openpmd.py

@classmethod
def pmd_keys(cls):
    """
    Map field names to their openPMD attribute names.

    Returns
    -------
    dict
        ``{field_name: pmd_key}``, excluding `other`.
    """
    return {
        fld.name: fld.metadata["pmd_key"]
        for fld in fields(cls)
        if "pmd_key" in fld.metadata
    }

to_pmd

to_pmd()

Render as a mapping keyed by openPMD attribute names.

Returns:

dict –

Set fields plus other. Fields left as None are omitted, since the extension treats them as absent rather than zero.

Source code in beamphysics/wavefront/openpmd.py

def to_pmd(self):
    """
    Render as a mapping keyed by openPMD attribute names.

    Returns
    -------
    dict
        Set fields plus `other`. Fields left as None are omitted, since the
        extension treats them as absent rather than zero.
    """
    out = {}
    for name, key in self.pmd_keys().items():
        value = getattr(self, name)
        if value is not None:
            out[key] = value
    out.update(self.other)
    return out

load_wavefront_openpmd

load_wavefront_openpmd(h5, iteration=None)

Read an openPMD EXT_Wavefront series into Wavefront constructor arguments.

Parameters:

h5 (Group) –

Series root.
iteration (int, default: None ) –

Iteration to read. If None, the file must hold exactly one.

Returns:

dict –

Keyword arguments for Wavefront: Ex, Ey, dx, dy, dz, wavelength, attrs (a WavefrontAttrs) and, when the file declares a gridGlobalOffset, xmid, ymid and zmid.

Raises:

ValueError –

If a required attribute is missing, or if the file holds something this class cannot represent: a frequency-domain field, a k-space field, axisLabels that are not a permutation of x, y and z, or x and y components with different position values, or a nonzero z component.

Notes

Any axisLabels permutation of (x, y, z) is honored, with gridSpacing permuted alongside it. Scalar attributes stored as length-1 arrays and strings stored as bytes are both tolerated. gridUnitSI, gridGlobalOffset, the per-component position and the per-component unitSI are applied, so a file written in non-SI units or with cell-centered samples reads correctly; gridUnitSI is accepted either as the openPMD 2.0 per-axis array or as the 1.x scalar. A constant component may take its shape from a sibling.

Source code in beamphysics/wavefront/openpmd.py

def load_wavefront_openpmd(h5, iteration=None):
    """
    Read an openPMD EXT_Wavefront series into `Wavefront` constructor arguments.

    Parameters
    ----------
    h5 : h5py.Group
        Series root.
    iteration : int, optional
        Iteration to read. If None, the file must hold exactly one.

    Returns
    -------
    dict
        Keyword arguments for `Wavefront`: `Ex`, `Ey`, `dx`, `dy`, `dz`,
        `wavelength`, `attrs` (a `WavefrontAttrs`) and, when the file declares a
        `gridGlobalOffset`, `xmid`, `ymid` and `zmid`.

    Raises
    ------
    ValueError
        If a required attribute is missing, or if the file holds something this
        class cannot represent: a frequency-domain field, a k-space field,
        `axisLabels` that are not a permutation of x, y and z, or x and y
        components with different `position` values, or a nonzero z component.

    Notes
    -----
    Any `axisLabels` permutation of ``(x, y, z)`` is honored, with `gridSpacing`
    permuted alongside it. Scalar attributes stored as length-1 arrays and strings
    stored as bytes are both tolerated. `gridUnitSI`, `gridGlobalOffset`, the
    per-component `position` and the per-component `unitSI` are applied, so a file
    written in non-SI units or with cell-centered samples reads correctly;
    `gridUnitSI` is accepted either as the openPMD 2.0 per-axis array or as the 1.x
    scalar. A constant component may take its `shape` from a sibling.
    """
    mesh = _mesh_record(h5, iteration)

    if "z" in mesh:
        longitudinal = mesh["z"]
        if isinstance(longitudinal, h5py.Group):
            if "value" not in longitudinal.attrs:
                raise ValueError(
                    f"component {longitudinal.name} is a group without a 'value'"
                )
            nonzero = np.any(np.asarray(longitudinal.attrs["value"]) != 0)
        else:
            nonzero = any(
                np.any(longitudinal[index] != 0)
                for index in range(longitudinal.shape[0])
            )
        if nonzero:
            raise ValueError(
                "nonzero electricField/z is not supported; Wavefront represents "
                "only transverse fields"
            )

    temporal_domain = _required_attr(mesh, "temporalDomain")
    if temporal_domain != "time":
        raise ValueError(
            f"temporalDomain {temporal_domain!r}: only 'time' (a field in V/m) "
            "is implemented"
        )

    spatial_domain = _required_attr(mesh, "spatialDomain")
    if spatial_domain != "r":
        raise ValueError(
            f"spatialDomain {spatial_domain!r}: only 'r' (cartesian space) "
            "is implemented"
        )

    photon_energy = float(_required_attr(mesh, "photonEnergy", fallback=h5))

    raw_labels = _required_attr(mesh, "axisLabels")
    labels = tuple(decode_attr(label) for label in np.atleast_1d(raw_labels))
    if sorted(labels) != sorted(_CLASS_AXIS_LABELS):
        raise ValueError(
            f"axisLabels {labels}: only a permutation of "
            f"{_CLASS_AXIS_LABELS} is implemented"
        )

    raw_spacing = np.atleast_1d(_required_attr(mesh, "gridSpacing"))
    if len(raw_spacing) != len(labels):
        raise ValueError(
            f"gridSpacing has {len(raw_spacing)} values but axisLabels has "
            f"{len(labels)}"
        )

    # gridSpacing is in the file's own units; gridUnitSI converts it to meters.
    # openPMD 2.0 makes gridUnitSI one value per axis, but implementations written
    # against 1.x emit a single scalar, so accept either.
    raw_grid_unit = np.atleast_1d(decode_attr(mesh.attrs.get("gridUnitSI", 1.0)))
    if raw_grid_unit.size == 1:
        grid_unit = np.full(len(labels), float(raw_grid_unit[0]))
    elif raw_grid_unit.size == len(labels):
        grid_unit = raw_grid_unit.astype(float)
    else:
        raise ValueError(
            f"gridUnitSI has {raw_grid_unit.size} values but axisLabels has "
            f"{len(labels)}"
        )

    spacing = {
        label: float(value) * float(unit)
        for label, value, unit in zip(labels, raw_spacing, grid_unit)
    }

    # gridGlobalOffset is the position of the first cell, in the same units as
    # gridSpacing. A file that omits it is read as a centered grid, which is this
    # class's own default rather than an assertion about the file.
    global_offset = None
    if "gridGlobalOffset" in mesh.attrs:
        raw_offset = np.atleast_1d(decode_attr(mesh.attrs["gridGlobalOffset"]))
        if len(raw_offset) != len(labels):
            raise ValueError(
                f"gridGlobalOffset has {len(raw_offset)} values but axisLabels has "
                f"{len(labels)}"
            )
        global_offset = {
            label: float(value) * float(unit)
            for label, value, unit in zip(labels, raw_offset, grid_unit)
        }

    # Everything on the record that this module does not write from the wavefront
    # itself, including names it does not know: those land in `other` so that a file
    # written against a newer revision of the extension survives a round trip.
    pmd_attrs = {
        name: decode_attr(value)
        for name, value in mesh.attrs.items()
        if name not in _COMPUTED_RECORD_ATTRS
    }
    if "zCoordinate" not in mesh.attrs:
        raise ValueError("required EXT_Wavefront attribute 'zCoordinate' is missing")
    attrs = WavefrontAttrs.from_pmd(pmd_attrs)

    kwargs = {
        "dx": spacing["x"],
        "dy": spacing["y"],
        "dz": spacing["z"],
        "wavelength": wavelength_from_photon_energy(photon_energy),
        "s_position": float(decode_attr(mesh.attrs["zCoordinate"])),
        "Ex": None,
        "Ey": None,
        "attrs": attrs,
    }

    # Permutation of the class's (x, y, z) axes into the file's stored order, so the
    # transposed view of the output array lines up with the dataset slice by slice.
    class_to_stored = _axis_permutation(_CLASS_AXIS_LABELS, labels)
    positions = []
    for name, key in (("x", "Ex"), ("y", "Ey")):
        if name not in mesh:
            continue
        component = mesh[name]

        # Sample i sits at gridGlobalOffset + (i + position) * gridSpacing.
        raw_position = np.atleast_1d(
            component.attrs.get("position", np.zeros(len(labels)))
        ).astype(float)
        if len(raw_position) != len(labels):
            raise ValueError(
                f"component {name!r} position has {len(raw_position)} values but "
                f"axisLabels has {len(labels)}"
            )
        positions.append(dict(zip(labels, raw_position)))

        # Component values are in the file's own units; unitSI converts to V/m.
        unit_si = float(
            np.atleast_1d(decode_attr(component.attrs.get("unitSI", 1.0)))[0]
        )

        # openPMD allows a uniform component to be stored as a group carrying
        # `value` instead of a dataset. It has no axes to permute.
        if isinstance(component, h5py.Group):
            if "value" not in component.attrs:
                raise ValueError(
                    f"component {component.name} is a group without a 'value'"
                )
            stored_shape = _constant_shape(mesh, component)
            value = decode_attr(component.attrs["value"])
            out = np.full(
                tuple(
                    stored_shape[labels.index(label)] for label in _CLASS_AXIS_LABELS
                ),
                value,
                dtype=np.result_type(value, np.float64),
            )
        else:
            shape = tuple(
                component.shape[labels.index(label)] for label in _CLASS_AXIS_LABELS
            )
            out = np.empty(shape, dtype=component.dtype)

            # A view of `out` in the file's axis order. Filling it slice by slice
            # avoids allocating a second full-size array for the transpose.
            stored_view = out.transpose(class_to_stored)
            for islice in range(component.shape[0]):
                stored_view[islice] = component[islice]

        # Multiply in place so that a complex64 field is not silently widened.
        if unit_si != 1.0:
            out *= unit_si

        kwargs[key] = out

    if len(positions) == 2 and positions[0] != positions[1]:
        raise ValueError(
            f"components x and y have different positions {positions[0]} and "
            f"{positions[1]}; a staggered field cannot share one grid"
        )

    # The class's grid fields are the axis *midpoint*, while gridGlobalOffset is
    # the first cell, so converting between them needs the sample count.
    sample = kwargs["Ex"] if kwargs["Ex"] is not None else kwargs["Ey"]
    if global_offset is not None and sample is not None:
        counts = dict(zip(_CLASS_AXIS_LABELS, sample.shape))
        for label in _CLASS_AXIS_LABELS:
            kwargs[f"{label}mid"] = global_offset[label] + spacing[label] * (
                positions[0][label] + (counts[label] - 1) / 2
            )

    return kwargs

write_wavefront_openpmd

write_wavefront_openpmd(wavefront, h5, iteration=1, **extension_attrs)

Write a Wavefront into an open HDF5 group as an openPMD EXT_Wavefront series.

Parameters:

wavefront (Wavefront) –

Real-space, time-domain wavefront. Ex and Ey are in V/m with shape (nx, ny, nz); an absent polarization is not written. Its s_position is written as the extension's required zCoordinate.
h5 (Group) –

Group to use as the openPMD series root.
iteration (int, default: 1 ) –

openPMD iteration index to write under basePath.
**extension_attrs –

WavefrontAttrs field names (beamline, radius_of_curvature_x, ...), taking precedence over the same fields in wavefront.attrs.

Raises:

TypeError –

If an attribute name is not a WavefrontAttrs field.

Notes

Datasets are stored in (z, y, x) order and written one transverse slice at a time, so no full transposed copy of the field is ever materialized.

Source code in beamphysics/wavefront/openpmd.py

def write_wavefront_openpmd(
    wavefront,
    h5,
    iteration=1,
    **extension_attrs,
):
    """
    Write a `Wavefront` into an open HDF5 group as an openPMD EXT_Wavefront series.

    Parameters
    ----------
    wavefront : Wavefront
        Real-space, time-domain wavefront. `Ex` and `Ey` are in V/m with shape
        ``(nx, ny, nz)``; an absent polarization is not written. Its `s_position`
        is written as the extension's required `zCoordinate`.
    h5 : h5py.Group
        Group to use as the openPMD series root.
    iteration : int, default=1
        openPMD iteration index to write under `basePath`.
    **extension_attrs
        `WavefrontAttrs` field names (`beamline`, `radius_of_curvature_x`, ...),
        taking precedence over the same fields in ``wavefront.attrs``.

    Raises
    ------
    TypeError
        If an attribute name is not a `WavefrontAttrs` field.

    Notes
    -----
    Datasets are stored in ``(z, y, x)`` order and written one transverse slice at a
    time, so no full transposed copy of the field is ever materialized.
    """
    attrs = WavefrontAttrs.from_pmd(wavefront.attrs)
    for name, value in extension_attrs.items():
        # Assigning an unknown name would silently stick a new instance attribute on
        # the dataclass, so check first and let the constructor raise TypeError.
        if name not in attrs.pmd_keys():
            WavefrontAttrs(**{name: value})
        setattr(attrs, name, value)

    carried = attrs.to_pmd()

    # Series (root) attributes.
    _write_attr(h5, "openPMD", "2.0.0")
    _write_attr(h5, "openPMDextension", "Wavefront")
    _write_attr(h5, "basePath", _DEFAULT_BASE_PATH)
    _write_attr(h5, "meshesPath", _DEFAULT_MESHES_PATH)
    _write_attr(h5, "iterationEncoding", "groupBased")
    _write_attr(h5, "iterationFormat", _DEFAULT_BASE_PATH)

    base = _DEFAULT_BASE_PATH.replace("%T", str(int(iteration))).strip("/")
    iteration_group = h5.require_group(base)
    iteration_group.attrs["time"] = 0.0
    iteration_group.attrs["dt"] = 0.0
    iteration_group.attrs["timeUnitSI"] = 1.0

    mesh = h5.require_group(f"{base}/{_DEFAULT_MESHES_PATH}{_RECORD_NAME}")

    # Grid quantities, ordered like the stored axes. `gridGlobalOffset` is the
    # position of the beginning of the first cell; the components declare
    # `position = 0`, so that is exactly the first sample of each axis.
    spacing = {"x": wavefront.dx, "y": wavefront.dy, "z": wavefront.dz}
    offset = {"x": wavefront.xmin, "y": wavefront.ymin, "z": wavefront.zmin}

    # Base standard, required on the mesh record.
    _write_attr(mesh, "geometry", "cartesian")
    mesh.attrs["axisLabels"] = encode_attr(_STORED_AXIS_LABELS)
    mesh.attrs["gridSpacing"] = np.array(
        [spacing[label] for label in _STORED_AXIS_LABELS], dtype=float
    )
    mesh.attrs["gridGlobalOffset"] = np.array(
        [offset[label] for label in _STORED_AXIS_LABELS], dtype=float
    )
    # openPMD 2.0 makes gridUnitSI one value per axis.
    mesh.attrs["gridUnitSI"] = np.ones(len(_STORED_AXIS_LABELS), dtype=float)
    mesh.attrs["gridUnitDimension"] = np.array(
        _UNIT_DIMENSION_LENGTH * len(_STORED_AXIS_LABELS), dtype=float
    )
    mesh.attrs["unitDimension"] = np.array(_UNIT_DIMENSION_E_FIELD, dtype=float)
    mesh.attrs["timeOffset"] = 0.0

    # Extension, required on the mesh record.
    mesh.attrs["photonEnergy"] = photon_energy_joules(wavefront.wavelength)
    _write_attr(mesh, "temporalDomain", "time")
    _write_attr(mesh, "spatialDomain", "r")
    mesh.attrs["zCoordinate"] = float(wavefront.s_position)

    # Extension, optional.
    for name, value in carried.items():
        _write_attr(mesh, name, value)

    # Components. `z` is never written.
    to_stored = _axis_permutation(_CLASS_AXIS_LABELS, _STORED_AXIS_LABELS)
    for name, field_array in (("x", wavefront.Ex), ("y", wavefront.Ey)):
        if field_array is None:
            continue

        # A view, not a copy: the transpose is realized one slice at a time below.
        stored_view = np.asarray(field_array).transpose(to_stored)

        # EXT_Wavefront components are complex, so a real input is promoted.
        dtype = np.result_type(field_array.dtype, np.complex64)
        dataset = mesh.create_dataset(name, shape=stored_view.shape, dtype=dtype)
        for islice in range(stored_view.shape[0]):
            dataset[islice] = stored_view[islice].astype(dtype, copy=False)

        dataset.attrs["unitSI"] = 1.0
        dataset.attrs["position"] = np.zeros(len(_STORED_AXIS_LABELS), dtype=float)

Wavefront

Wavefront dataclass

dkx property

dky property

dkz property

dthetax property

dthetay property

energy property

energy_density property

fluence property

fluence_profile_x property

fluence_profile_y property

in_kspace property

in_rspace property

intensity property

intensity_x property

intensity_y property

k0 property

kxmax property

kxmin property

kxvec property

kymax property

kymin property

kyvec property

kzmax property

kzmin property

kzvec property

mean_x property

mean_y property

mean_z property

nx property

ny property

nz property

photon_energy property

power property

shape property

sigma_x property

sigma_y property

sigma_z property

thetaxmax property

thetaxmin property

thetaxvec property

thetaymax property

thetaymin property

thetayvec property

xmax property

xmin property

xvec property

ymax property

ymin property

yvec property

zmax property

zmin property

zvec property

__post_init__

__setattr__

auto_crop

axis_index

copy

crop

estimate_curvature

from_gaussian classmethod

from_genesis4 classmethod

from_openpmd classmethod

pad

plot2

plot_fluence

to_kspace

write_genesis4

write_openpmd

WavefrontK dataclass

dthetax property

dthetay property

dx property

dy property

dz property

energy property

in_kspace property

in_rspace property

intensity property

Wavefront `dataclass`

dkx `property`

dky `property`

dkz `property`

dthetax `property`

dthetay `property`

energy `property`

energy_density `property`

fluence `property`

fluence_profile_x `property`

fluence_profile_y `property`

in_kspace `property`

in_rspace `property`

intensity `property`

intensity_x `property`

intensity_y `property`

k0 `property`

kxmax `property`

kxmin `property`

kxvec `property`

kymax `property`

kymin `property`

kyvec `property`

kzmax `property`

kzmin `property`

kzvec `property`

mean_x `property`

mean_y `property`

mean_z `property`

nx `property`

ny `property`

nz `property`

photon_energy `property`

power `property`

shape `property`

sigma_x `property`

sigma_y `property`

sigma_z `property`

thetaxmax `property`

thetaxmin `property`

thetaxvec `property`

thetaymax `property`

thetaymin `property`

thetayvec `property`

xmax `property`

xmin `property`

xvec `property`

ymax `property`

ymin `property`

yvec `property`

zmax `property`

zmin `property`

zvec `property`

setattr

from_gaussian `classmethod`

from_genesis4 `classmethod`

from_openpmd `classmethod`

WavefrontK `dataclass`

dthetax `property`

dthetay `property`

dx `property`

dy `property`

dz `property`

energy `property`

in_kspace `property`

in_rspace `property`

intensity `property`

intensity_x `property`

intensity_y `property`

k0 `property`

kxmax `property`

kxmin `property`

kxvec `property`

kymax `property`

kymin `property`

kyvec `property`

kzmax `property`

kzmin `property`

kzvec `property`

mean_kx `property`