Retrieval across pulse shapes and interactions#

Compare SHG, PG and SD retrievals for Gaussian and sech pulses, then include a dispersive slab. This is a compact survey after your first retrieval.

Inputs and method#

The script generates its own traces on a matched Fourier grid. It uses a fixed noise seed and fixed retrieval seeds. Clean cases use COPRA; noisy cases use L-BFGS with its package defaults, so this is not a controlled solver ranking. The final case uses a 260 nm carrier and a 10 µm fused-silica slab.

Run from a source checkout:

uv run python examples/example_retrieval.py

Calculation and output#

The cell below loads and executes the runnable script itself. Expand it to inspect the source, or download example_retrieval.py. Generated files are isolated in a temporary directory for this documentation run.

Hide code cell source

import os
import sys
import tempfile
from pathlib import Path
from IPython.display import Image, display

# The notebook starts in docs/examples. Keep generated files outside the checkout.
repository = Path.cwd().parents[1]
workspace = tempfile.TemporaryDirectory(prefix="croak-gallery-")
(Path(workspace.name) / "examples").symlink_to(repository / "examples", target_is_directory=True)
sys.path.insert(0, str(repository / "examples"))
os.chdir(workspace.name)

Hide code cell source

"""End-to-end FROG retrieval examples for the croak package.

Demonstrates the uniform API across pulse shapes, interaction geometries, noise
levels and retrieval algorithms, including a dispersive (pseudo-propagation)
example. Run with::

    uv run python examples/example_retrieval.py

Prints the final FROG error and the recovered temporal FWHM for each case. This script
prints a numerical comparison and does not save a figure.
"""

from __future__ import annotations

import numpy as np

import croak
from croak.maths import fwhm as measure_fwhm
from croak.maths import wlfreq


def recovered_fwhm(res: croak.RetrievalResult) -> float:
    """Temporal intensity FWHM (fs) of a retrieved pulse, peak-centred."""
    intensity = res.intensity_t
    return measure_fwhm(res.t, intensity) * 1e15


def banner(title: str) -> None:
    print("\n" + "=" * 68)
    print(title)
    print("=" * 68)


def thin_examples() -> None:
    """Gaussian and sech pulses across all three interactions, clean + noisy."""
    banner("Thin-medium retrieval (Gaussian & sech, PG/SD/SHG)")
    # A 38.4 fs window: the ±18 fs scan stays inside the half-window, so for these
    # few-fs pulses the circular delay shift never wraps the gate onto the pulse.
    g = croak.Grid(128, dt=0.3e-15)
    delays = np.linspace(-18e-15, 18e-15, 87)
    rng = np.random.default_rng(0)

    pulses = {
        "Gaussian 2.5 fs, GDD": croak.gaussian_pulse(g, 2.5e-15, phases=[4e-30]),
        "sech 3 fs": croak.sech_pulse(g, 3.0e-15),
    }
    print(
        f"{'pulse':22s} {'inter':5s} {'algo':6s} {'noise':>6s} "
        f"{'R':>10s} {'FWHM/fs':>9s}"
    )
    for label, ew in pulses.items():
        true_fwhm = measure_fwhm(g.t, np.abs(g.ifft(ew)) ** 2) * 1e15
        for interaction in ("pg", "sd", "shg"):
            clean = croak.maketrace(g.omega, delays, ew, interaction)
            for sigma in (0.0, 0.01):
                trace = clean
                if sigma:
                    trace = np.clip(
                        clean + sigma * clean.max() * rng.standard_normal(clean.shape),
                        0.0,
                        None,
                    )
                algo = "copra" if sigma == 0 else "lbfgs"
                near = croak.gaussian_pulse(g, 1.8e-15)
                res = croak.retrieve(
                    trace,
                    g.omega,
                    delays,
                    interaction,
                    algorithm=algo,
                    guess=near,
                    maxiters=250,
                    rng=np.random.default_rng(0),
                    progress=False,
                )
                print(
                    f"{label:22s} {interaction:5s} {algo:6s} {sigma:6.0%} "
                    f"{res.error:10.2e} {recovered_fwhm(res):9.2f}"
                )
        print(f"{'  (true FWHM = ' + f'{true_fwhm:.2f} fs)':22s}")


def dispersive_example() -> None:
    """Few-fs UV pulse through a fused-silica slab: D-COPRA and dispersive L-BFGS."""
    banner("Dispersive retrieval (UV pulse, 10 um SiO2 slab)")
    g = croak.Grid(128, dt=0.5e-15)
    omega0 = float(wlfreq(260e-9))  # 260 nm carrier
    ew = croak.gaussian_pulse(g, 4.0e-15, phases=[6e-30])
    delays = np.linspace(-22e-15, 22e-15, 80)
    trace = croak.maketrace(
        g.omega,
        delays,
        ew,
        "pg",
        material="SiO2",
        thickness=10e-6,
        npoints=20,
        omega0=omega0,
    )
    true_fwhm = measure_fwhm(g.t, np.abs(g.ifft(ew)) ** 2) * 1e15
    near = croak.gaussian_pulse(g, 3.0e-15)

    print(f"true temporal FWHM = {true_fwhm:.2f} fs\n")
    print(f"{'algorithm':10s} {'R':>10s} {'FWHM/fs':>9s}")
    for algorithm in ("copra", "lbfgs"):
        res = croak.retrieve(
            trace,
            g.omega,
            delays,
            "pg",
            algorithm=algorithm,
            guess=near,
            maxiters=300,
            rng=np.random.default_rng(0),
            progress=False,
            material="SiO2",
            thickness=10e-6,
            npoints=20,
            omega0=omega0,
        )
        print(f"{algorithm:10s} {res.error:10.2e} {recovered_fwhm(res):9.2f}")


def main() -> None:
    thin_examples()
    dispersive_example()
    print("\nDone.")


if __name__ == "__main__":
    main()
====================================================================
Thin-medium retrieval (Gaussian & sech, PG/SD/SHG)
====================================================================
pulse                  inter algo    noise          R   FWHM/fs
Gaussian 2.5 fs, GDD   pg    copra      0%   8.12e-17      5.10
Gaussian 2.5 fs, GDD   pg    lbfgs      1%   6.31e-03      4.99
Gaussian 2.5 fs, GDD   sd    copra      0%   6.54e-17      5.10
Gaussian 2.5 fs, GDD   sd    lbfgs      1%   6.45e-03      4.84
Gaussian 2.5 fs, GDD   shg   copra      0%   1.13e-16      5.10
Gaussian 2.5 fs, GDD   shg   lbfgs      1%   3.15e-02      4.01
  (true FWHM = 5.10 fs)
sech 3 fs              pg    copra      0%   2.77e-10      3.00
sech 3 fs              pg    lbfgs      1%   6.44e-03      3.10
sech 3 fs              sd    copra      0%   2.26e-10      3.00
sech 3 fs              sd    lbfgs      1%   6.59e-03      2.81
sech 3 fs              shg   copra      0%   9.66e-13      3.00
sech 3 fs              shg   lbfgs      1%   6.19e-03      2.95
  (true FWHM = 3.00 fs)

====================================================================
Dispersive retrieval (UV pulse, 10 um SiO2 slab)
====================================================================
true temporal FWHM = 5.78 fs

algorithm           R   FWHM/fs
copra        3.53e-17      5.78
lbfgs        1.69e-05      5.75

Done.

Interpretation and limits#

Read the printed duration against the generating value for each pulse. In the seeded noisy SHG case, the retrieved duration differs appreciably from the generating Gaussian; retain that discrepancy when assessing the method. GDD stretches a pulse beyond its transform-limited input width. Noise can change the residual and field differently, and clipping negative noisy samples changes the noise distribution. This script is an API survey, not a coverage or noise-model validation. Follow solver selection and retrieval theory for controlled comparisons.