many changes

setup.cfg

@@ -30,9 +30,7 @@ install_requires =
 [options.package_data]
 scgenerator =
     data/hr_t.npz
-    data/default_params.json
+    data/materials.toml
-    data/gas.json
-    data/silica.json
 
 [options.packages.find]
 where = src

src/scgenerator/cli/cli.py

@@ -186,7 +186,7 @@ def prep_ray():
 def plot_all(args):
     opts = {}
     if args.options is not None:
-        opts |= dict([o.split("=")[:2] for o in re.split("[, ]", args.options)])
+        opts |= dict([o.split("=")[:2] for o in re.split("[, ]+", args.options)])
     root = Path(args.sim_dir).resolve()
     scripts.plot_all(root, args.spectrum_limits, show=args.show, **opts)
 

src/scgenerator/data/materials.toml

@@ -0,0 +1,252 @@
+[silica.sellmeier]
+B = [0.6961663, 0.4079426, 0.8974794]
+C = [4.67914826e-15, 1.35120631e-14, 9.79340025e-11]
+kind = 1
+reference = [
+    "I. H. Malitson. Interspecimen comparison of the refractive index of fused silica, J. Opt. Soc. Am. 55, 1205-1208 (1965)",
+    "C. Z. Tan. Determination of refractive index of silica glass for infrared wavelengths by IR spectroscopy, J. Non-Cryst. Solids 223, 158-163 (1998)",
+]
+
+[nbk7.sellmeier]
+B = [1.03961212, 0.231792344, 1.01046945]
+C = [6.00069867e-15, 2.00179144e-14, 1.03560653e-10]
+kind = 1
+reference = [
+    "SCHOTT Zemax catalog 2017-01-20b (obtained from http://www.schott.com)",
+]
+
+[e7.sellmeier]
+B = [0.0, 7.2e-15, 3e-28]
+C = []
+kind = 3
+reference = [
+    "J. Li, C. H. Wen, S. Gauza, R. Lu and S. T. Wu. Refractive indices of liquid crystals for display applications, J. Disp. Technol. 1, 51-61, 2005",
+]
+
+[d-zlaf52la.sellmeier]
+B = [
+    -8.48750283e9,
+    3.1753525,
+    3.71301651e-14,
+    -1.09609062e-27,
+    2.4418582e-40,
+    -8.94583294e-54,
+]
+C = []
+kind = 3
+reference = ["CDGM Zemax catalog 2017-09 (obtained from http://www.cdgmgd.com)"]
+
+# -----
+# GASES
+# -----
+
+[air]
+a = 0.1358
+b = 3.64e-5
+
+[air.sellmeier]
+B = [57921050000.0, 1679170000.0]
+C = [238018500000000.0, 57362000000000.0]
+P0 = 101325
+T0 = 288.15
+kind = 2
+
+[air.kerr]
+P0 = 101325
+T0 = 293.15
+n2 = 4.0e-23
+source = "Pigeon, J. J., Tochitsky, S. Y., Welch, E. C., & Joshi, C. (2016). Measurements of the nonlinear refractive index of air, N 2, and O 2 at 10 μm using four-wave mixing. Optics letters, 41(17), 3924-3927."
+
+[nitrogen]
+a = 0.137
+b = 1.709e-5
+
+[nitrogen.sellmeier]
+B = [32431570000.0]
+C = [144000000000000.0]
+P0 = 101325
+T0 = 273.15
+const = 6.8552e-5
+kind = 2
+
+[nitrogen.kerr]
+P0 = 30400.0
+T0 = 273.15
+n2 = 2.2e-23
+source = "Wahlstrand, J. K., Cheng, Y. H., & Milchberg, H. M. (2012). High field optical nonlinearity and the Kramers-Kronig relations. Physical review letters, 109(11), 113904."
+
+[helium]
+a = 0.00346
+b = 2.38e-5
+
+[helium.sellmeier]
+B = [2.16463842e-5, 2.10561127e-7, 4.7509272e-5]
+C = [-6.80769781e-16, 5.13251289e-15, 3.18621354e-15]
+P0 = 101325
+T0 = 273.15
+kind = 1
+source = "A. Ermolov, K. F. Mak, M. H. Frosz, J. C. Travers, P. St. J. Russell, Supercontinuum generation in the vacuum ultraviolet through dispersive-wave and soliton-plasma interaction in a noble-gas-filled hollow-core photonic crystal fiber, Phys. Rev. A 92, 033821 (2015)"
+
+[helium.kerr]
+P0 = 30400.0
+T0 = 273.15
+n2 = 3.1e-25
+source = "Wahlstrand, J. K., Cheng, Y. H., & Milchberg, H. M. (2012). High field optical nonlinearity and the Kramers-Kronig relations. Physical review letters, 109(11), 113904."
+
+[helium_alt]
+a = 0.00346
+b = 2.38e-5
+
+[helium_alt.sellmeier]
+B = [14755297000.0]
+C = [426297400000000.0]
+P0 = 101325
+T0 = 273.15
+kind = 2
+source = " C. Cuthbertson and M. Cuthbertson. The refraction and dispersion of neon and helium. Proc. R. Soc. London A 135, 40-47 (1936)"
+
+[helium_alt.kerr]
+P0 = 30400.0
+T0 = 273.15
+n2 = 3.1e-25
+
+[hydrogen]
+a = 0.02453
+b = 2.651e-5
+
+[hydrogen.sellmeier]
+B = [0.0148956, 0.0049037]
+C = [1.807e-10, 9.2e-11]
+P0 = 101325
+T0 = 273.15
+kind = 2
+source = "E. R. Peck and S. Hung. Refractivity and dispersion of hydrogen in the visible and near infrared, J. Opt. Soc. Am. 67, 1550-1554 (1977)"
+
+[hydrogen.kerr]
+P0 = 30400.0
+T0 = 273.15
+n2 = 6.36e-24
+source = "Shelton, D. P., & Rice, J. E. (1994). Measurements and calculations of the hyperpolarizabilities of atoms and small molecules in the gas phase. Chemical Reviews, 94(1), 3-29"
+
+[neon]
+a = 0.02135
+b = 1.709e-5
+
+[neon.sellmeier]
+B = [1281450000.0, 22048600000.0]
+C = [184661000000000.0, 376840000000000.0]
+P0 = 101325
+T0 = 273.15
+kind = 2
+
+[neon.kerr]
+P0 = 30400.0
+T0 = 273.15
+n2 = 8.7e-25
+source = "Wahlstrand, J. K., Cheng, Y. H., & Milchberg, H. M. (2012). High field optical nonlinearity and the Kramers-Kronig relations. Physical review letters, 109(11), 113904."
+
+[argon]
+a = 0.1355
+b = 3.201e-5
+
+[argon.sellmeier]
+B = [2501410000.0, 500283000.0, 52234300000.0]
+C = [91012000000000.0, 87892000000000.0, 214020000000000.0]
+P0 = 101325
+T0 = 273.15
+kind = 2
+source = "A. Bideau-Mehu, Y. Guern, R. Abjean, A. Johannin-Gilles. Measurement of refractive indices of neon, argon, krypton and xenon in the 253.7-140.4 nm wavelength range. Dispersion relations and estimated oscillator strengths of the resonance lines. J. Quant. Spectrosc. Rad. Transfer 25, 395-402 (1981)"
+
+[argon.kerr]
+P0 = 30400.0
+T0 = 273.15
+n2 = 9.7e-24
+source = "Wahlstrand, J. K., Cheng, Y. H., & Milchberg, H. M. (2012). High field optical nonlinearity and the Kramers-Kronig relations. Physical review letters, 109(11), 113904."
+
+[argon_alt]
+a = 0.1355
+b = 3.201e-5
+
+[argon_alt.sellmeier]
+B = [0.0002033229, 0.0003445831]
+C = [2.0612e-16, 8.066e-15]
+P0 = 101325
+T0 = 273.15
+kind = 1
+source = "A. Börzsönyi, Z. Heiner, M. P. Kalashnikov, A. P. Kovács, and K. Osvay, Dispersion measurement of inert gases and gas mixtures at 800 nm, Appl. Opt. 47, 4856-4863 (2008)"
+
+[argon_alt.kerr]
+P0 = 30400.0
+T0 = 273.15
+n2 = 9.7e-24
+source = "Wahlstrand, J. K., Cheng, Y. H., & Milchberg, H. M. (2012). High field optical nonlinearity and the Kramers-Kronig relations. Physical review letters, 109(11), 113904."
+
+[argon_alt2]
+a = 0.1355
+b = 3.201e-5
+
+[argon_alt2.sellmeier]
+B = [0.030182943]
+C = [144000000000000.0]
+P0 = 101325
+T0 = 273.15
+const = 6.7867e-5
+kind = 2
+source = "E. R. Peck and D. J. Fisher. Dispersion of argon, J. Opt. Soc. Am. 54, 1362-1364 (1964)"
+
+[argon_alt2.kerr]
+P0 = 30400.0
+T0 = 273.15
+n2 = 9.7e-24
+source = "Wahlstrand, J. K., Cheng, Y. H., & Milchberg, H. M. (2012). High field optical nonlinearity and the Kramers-Kronig relations. Physical review letters, 109(11), 113904."
+
+[krypton]
+a = 0.2349
+b = 3.978e-5
+
+[krypton.sellmeier]
+B = [2536370000.0, 2736490000.0, 62080200000.0]
+C = [65474200000000.0, 73698000000000.0, 181080000000000.0]
+P0 = 101325
+T0 = 273.15
+kind = 2
+
+[krypton.kerr]
+P0 = 30400.0
+T0 = 273.15
+n2 = 2.2e-23
+source = "Wahlstrand, J. K., Cheng, Y. H., & Milchberg, H. M. (2012). High field optical nonlinearity and the Kramers-Kronig relations. Physical review letters, 109(11), 113904."
+
+[xenon]
+a = 0.425
+b = 5.105e-5
+
+[xenon.sellmeier]
+B = [3228690000.0, 3553930000.0, 60676400000.0]
+C = [46301000000000.0, 59578000000000.0, 112740000000000.0]
+P0 = 101325
+T0 = 273.15
+kind = 2
+
+[xenon.kerr]
+P0 = 30400.0
+T0 = 273.15
+n2 = 5.8e-23
+source = "Wahlstrand, J. K., Cheng, Y. H., & Milchberg, H. M. (2012). High field optical nonlinearity and the Kramers-Kronig relations. Physical review letters, 109(11), 113904."
+
+[vacuum]
+a = 0
+b = 0
+
+[vacuum.sellmeier]
+B = []
+C = []
+P0 = 101325
+T0 = 273.15
+kind = 1
+
+[vacuum.kerr]
+P0 = 30400.0
+T0 = 273.15
+n2 = 0
+source = "none"

src/scgenerator/data/silica.toml

@@ -1,4 +1,29 @@
-[sellmeier]
+[silica.sellmeier]
-B = [ 0.6961663, 0.4079426, 0.8974794,]
+B = [0.6961663, 0.4079426, 0.8974794]
-C = [ 4.67914826e-15, 1.35120631e-14, 9.79340025e-11,]
+C = [4.67914826e-15, 1.35120631e-14, 9.79340025e-11]
 kind = 1
+reference = [
+    "I. H. Malitson. Interspecimen comparison of the refractive index of fused silica, J. Opt. Soc. Am. 55, 1205-1208 (1965)",
+    "C. Z. Tan. Determination of refractive index of silica glass for infrared wavelengths by IR spectroscopy, J. Non-Cryst. Solids 223, 158-163 (1998)",
+]
+
+[nbk7.sellmeier]
+B = [1.03961212, 0.231792344, 1.01046945]
+C = [6.00069867e-15, 2.00179144e-14, 1.03560653e-10]
+kind = 1
+reference = [
+    "SCHOTT Zemax catalog 2017-01-20b (obtained from http://www.schott.com)",
+]
+
+[d-zlaf52la.sellmeier]
+B = [
+    -8.48750283e9,
+    3.1753525,
+    3.71301651e-14,
+    -1.09609062e-27,
+    2.4418582e-40,
+    -8.94583294e-54,
+]
+C = []
+kind = 3
+reference = ["CDGM Zemax catalog 2017-09 (obtained from http://www.cdgmgd.com)"]

src/scgenerator/physics/__init__.py

@@ -11,6 +11,7 @@ from scipy.optimize import minimize_scalar
 from .. import math
 from . import fiber, materials, units, pulse
 from .. import utils
+from ..utils import cache
 
 T = TypeVar("T")
 
@@ -21,8 +22,9 @@ def group_delay_to_gdd(wavelength: np.ndarray, group_delay: np.ndarray) -> np.nd
     return gdd
 
 
+@cache.np_cache
 def material_dispersion(
-    wavelengths,
+    wavelengths: np.ndarray,
     material: str,
     pressure=None,
     temperature=None,
@@ -52,17 +54,6 @@ def material_dispersion(
 
     w = units.m(wavelengths)
 
-    order = np.argsort(w)
-
-    material_dico = utils.load_material_dico(material)
-    if ideal:
-        n_gas_2 = materials.sellmeier(wavelengths, material_dico, pressure, temperature) + 1
-    else:
-        N_1 = materials.number_density_van_der_waals(
-            pressure=pressure, temperature=temperature, material_dico=material_dico
-        )
-        N_0 = materials.number_density_van_der_waals(material_dico=material_dico)
-        n_gas_2 = materials.sellmeier(wavelengths, material_dico) * N_1 / N_0 + 1
     if safe:
         disp = np.zeros(len(w))
         ind = w > 0
@@ -71,7 +62,18 @@ def material_dispersion(
         )
         return disp
     else:
-        return fiber.beta2(w[order], np.sqrt(n_gas_2[order]))[order]
+        material_dico = utils.load_material_dico(material)
+        if ideal:
+            n_gas_2 = materials.sellmeier(wavelengths, material_dico, pressure, temperature) + 1
+        else:
+            N_1 = materials.number_density_van_der_waals(
+                pressure=pressure, temperature=temperature, material_dico=material_dico
+            )
+            N_0 = materials.number_density_van_der_waals(material_dico=material_dico)
+            n_gas_2 = materials.sellmeier(wavelengths, material_dico) * N_1 / N_0 + 1
+        order = np.argsort(w)
+        unorder = np.argsort(order)
+        return fiber.beta2(w[order], np.sqrt(n_gas_2[order]))[unorder]
 
 
 def find_optimal_depth(
@@ -108,15 +110,30 @@ def find_optimal_depth(
     def score(z):
         return -np.nansum(integrate(z) ** 6)
 
-    # import matplotlib.pyplot as plt
-
-    # to_test = np.linspace(0, max_z, 200)
-    # scores = [score(z) for z in to_test]
-    # fig, ax = plt.subplots()
-    # ax.plot(to_test, scores / np.min(scores))
-    # plt.show()
-    # plt.close(fig)
-    # ama = np.argmin(scores)
-
     opti = minimize_scalar(score, method="bounded", bounds=(0, max_z))
     return propagate(opti.x), opti
+
+
+def propagate_field(t: np.ndarray, field: np.ndarray, z: float, material: str) -> np.ndarray:
+    """propagates a field through bulk material
+
+    Parameters
+    ----------
+    t : np.ndarray, shape (n,)
+        time grid
+    field : np.ndarray, shape (n,)
+        corresponding complex field
+    z : float
+        distance to propagate in m
+    material : str
+        material name
+
+    Returns
+    -------
+    np.ndarray, shape (n,)
+        propagated field
+    """
+    w_c = math.wspace(t)
+    l = units.m(w_c + units.nm(1540))
+    disp = material_dispersion(l, material)
+    return np.fft.ifft(np.fft.fft(field) * np.exp(0.5j * disp * w_c ** 2 * z))

src/scgenerator/physics/materials.py

@@ -113,8 +113,8 @@ def sellmeier(lambda_, material_dico, pressure=None, temperature=None):
     logger = get_logger(__name__)
 
     WL_THRESHOLD = 8.285e-6
-    temp_l = lambda_[lambda_ < WL_THRESHOLD]
+    ind = lambda_ < WL_THRESHOLD
-    kind = 1
+    temp_l = lambda_[ind]
 
     B = material_dico["sellmeier"]["B"]
     C = material_dico["sellmeier"]["C"]
@@ -125,13 +125,20 @@ def sellmeier(lambda_, material_dico, pressure=None, temperature=None):
 
     chi = np.zeros_like(lambda_)  # = n^2 - 1
     if kind == 1:
+        logger.debug("materials : using Sellmeier 1st kind equation")
         for b, c in zip(B, C):
-            chi[lambda_ < WL_THRESHOLD] += temp_l ** 2 * b / (temp_l ** 2 - c)
+            chi[ind] += temp_l ** 2 * b / (temp_l ** 2 - c)
     elif kind == 2:  # gives n-1
+        logger.debug("materials : using Sellmeier 2nd kind equation")
         for b, c in zip(B, C):
-            chi[lambda_ < WL_THRESHOLD] += b / (c - 1 / temp_l ** 2)
+            chi[ind] += b / (c - 1 / temp_l ** 2)
         chi += const
         chi = (chi + 1) ** 2 - 1
+    elif kind == 3:  # Schott formula
+        logger.debug("materials : using Schott equation")
+        for i, b in reversed(list(enumerate(B))):
+            chi[ind] += b * temp_l ** (-2 * (i - 1))
+        chi[ind] = chi[ind] - 1
     else:
         raise ValueError(f"kind {kind} is not recognized.")
 

src/scgenerator/physics/simulate.py

@@ -44,6 +44,12 @@ class RK4IP:
         """
         self.set(params, save_data, job_identifier, task_id)
 
+    def __iter__(self) -> Generator[tuple[int, int, np.ndarray], None, None]:
+        yield from self.irun()
+
+    def __len__(self) -> int:
+        return self.len
+
     def set(
         self,
         params: Parameters,
@@ -73,6 +79,7 @@ class RK4IP:
         self.alpha = params.alpha_arr
         self.spec_0 = np.sqrt(params.input_transmission) * params.spec_0
         self.z_targets = params.z_targets
+        self.len = len(params.z_targets)
         self.z_final = params.length
         self.beta2_coefficients = (
             params.beta_func if params.beta_func is not None else params.beta2_coefficients
@@ -189,7 +196,7 @@ class RK4IP:
         """
         utils.save_data(data, self.data_dir, name)
 
-    def run(self):
+    def run(self) -> list[np.ndarray]:
         time_start = datetime.today()
 
         for step, num, _ in self.irun():
@@ -480,7 +487,7 @@ class Simulations:
 
     def _run_available(self):
         for variable, params in self.configuration:
-            v_list_str = format_variable_list(variable)
+            v_list_str = format_variable_list(variable, add_iden=True)
             utils.save_parameters(params.prepare_for_dump(), Path(params.output_path))
 
             self.new_sim(v_list_str, params)

src/scgenerator/physics/units.py

@@ -48,6 +48,7 @@ def unit(tpe: str, label: str, inv: Callable = None):
         setattr(To, name, inv.__call__)
         func.type = tpe
         func.label = label
+        func.name = name
         func.inv = inv
         if name in units_map:
             raise NameError(f"Two unit functions with the same name {name!r} were defined")

src/scgenerator/plotting.py

@@ -543,7 +543,7 @@ def transform_mean_values(
         values = abs2(values)
     new_axis, ind, ext = sort_axis(x_axis, plt_range)
     values = values[:, ind]
-    if plt_range.unit.type == "WL":
+    if plt_range.unit.type == "WL" and plt_range.conserved_quantity:
         values = np.apply_along_axis(units.to_WL, -1, values, new_axis)
 
     if isinstance(spacing, (float, np.floating)):
@@ -735,14 +735,13 @@ def transform_1D_values(
         x axis and values
     """
     if len(values.shape) != 1:
-        print(f"Shape was {values.shape}. Can only plot 1D arrays")
+        raise ValueError("Can only plot 1D values")
-        return
     is_complex, x_axis, plt_range = prep_plot_axis(values, plt_range, params)
     if is_complex:
         values = abs2(values)
     new_axis, ind, ext = sort_axis(x_axis, plt_range)
     values = values[ind]
-    if plt_range.unit.type == "WL":
+    if plt_range.unit.type == "WL" and plt_range.conserved_quantity:
         values = units.to_WL(values, new_axis)
 
     if isinstance(spacing, (float, np.floating)):
@@ -1000,3 +999,64 @@ def arrowstyle(direction=1, color="white"):
         color=color,
         backgroundcolor=(0.5, 0.5, 0.5, 0.5),
     )
+
+
+def measure_and_annotate_fwhm(
+    ax: plt.Axes,
+    t: np.ndarray,
+    field: np.ndarray,
+    side: Literal["left", "right"] = "right",
+    unit="fs",
+    arrow_length_pts: float = 20.0,
+    arrow_props: dict[str, Any] = None,
+) -> float:
+    """measured the FWHM of a pulse and plots it
+
+    Parameters
+    ----------
+    ax : plt.Axes
+        ax on which to plot
+    t : np.ndarray, shape (n,)
+        time in s
+    field : np.ndarray, shape (n,)
+        complex field
+    side : Literal["left", "right"]
+        whether to write the text on the right or left side
+    unit : str, optional
+        units of the plot, by default "fs"
+    arrow_length_pts : float, optional
+        length of the arrows in pts, by default 20.0
+    arrow_props : dict[str, Any], optional
+        style of the arrow to be passed to plt.annotate, by default None
+
+    Returns
+    -------
+    float
+        FWHM in units
+    """
+    unit = units.get_unit(unit)
+    if np.iscomplexobj(field):
+        field = abs2(field)
+    _, (left, right), *_ = pulse.find_lobe_limits(unit.inv(t), field)
+    arrow_label = f"{right - left:.1f} {unit.name}"
+    arrow_dict = dict(arrowstyle="->")
+    if arrow_props is not None:
+        arrow_dict |= arrow_props
+    ax.annotate(
+        "" if side == "right" else arrow_label,
+        (left, field.max() / 2),
+        xytext=(-arrow_length_pts, 0),
+        ha="right",
+        va="center",
+        textcoords="offset points",
+        arrowprops=arrow_dict,
+    )
+    ax.annotate(
+        "" if side == "left" else arrow_label,
+        (right, field.max() / 2),
+        xytext=(arrow_length_pts, 0),
+        textcoords="offset points",
+        arrowprops=arrow_dict,
+        va="center",
+    )
+    return right - left

src/scgenerator/scripts/__init__.py

@@ -31,6 +31,8 @@ def plot_all(sim_dir: Path, limits: list[str], show=False, **opts):
     for k, v in opts.items():
         if k in ["skip"]:
             opts[k] = int(v)
+        if k in {"log", "renormalize"}:
+            opts[k] = True if v == "True" else False
     dir_list = list(p for p in sim_dir.glob("*") if p.is_dir())
     if len(dir_list) == 0:
         dir_list = [sim_dir]

src/scgenerator/spectra.py

@@ -1,7 +1,7 @@
 import os
 from collections.abc import Sequence
 from pathlib import Path
-from typing import Callable, Dict, Iterable, Union
+from typing import Callable, Dict, Iterable, Optional, Union
 
 import matplotlib.pyplot as plt
 import numpy as np
@@ -17,6 +17,7 @@ from .plotting import (
     transform_2D_propagation,
 )
 from .utils.parameter import Parameters, PlotRange
+from .utils import load_spectrum
 
 
 class Spectrum(np.ndarray):
@@ -152,7 +153,7 @@ class Pulse(Sequence):
         self.params = Parameters.load(self.path / "params.toml")
         self.params.compute(["name", "t", "l", "w_c", "w0", "z_targets"])
         if self.params.fiber_map is None:
-            self.params.fiber_map = {0.0: self.params.name}
+            self.params.fiber_map = [(0.0, self.params.name)]
 
         try:
             self.z = np.load(os.path.join(path, "z.npy"))
@@ -161,7 +162,6 @@ class Pulse(Sequence):
                 self.z = self.params.z_targets
             else:
                 raise
-        self.cache: Dict[int, Spectrum] = {}
         self.nmax = len(list(self.path.glob("spectra_*.npy")))
         if self.nmax <= 0:
             raise FileNotFoundError(f"No appropriate file in specified folder {self.path}")
@@ -306,12 +306,9 @@ class Pulse(Sequence):
     def _load1(self, i: int):
         if i < 0:
             i = self.nmax + i
-        if i in self.cache:
+        spec = load_spectrum(self.path / SPECN_FN.format(i))
-            return self.cache[i]
-        spec = np.load(self.path / SPECN_FN.format(i))
         spec = np.atleast_2d(spec)
         spec = Spectrum(spec, self.params)
-        self.cache[i] = spec
         return spec
 
     def plot_2D(
@@ -324,8 +321,9 @@ class Pulse(Sequence):
         sim_ind: int = 0,
         **kwargs,
     ):
-        plt_range, vals = self.retrieve_plot_values(left, right, unit, z_pos, sim_ind)
+        plot_range = PlotRange(left, right, unit)
-        return propagation_plot(vals, plt_range, self.params, ax, **kwargs)
+        vals = self.retrieve_plot_values(plot_range, z_pos, sim_ind)
+        return propagation_plot(vals, plot_range, self.params, ax, **kwargs)
 
     def plot_1D(
         self,
@@ -337,8 +335,9 @@ class Pulse(Sequence):
         sim_ind: int = 0,
         **kwargs,
     ):
-        plt_range, vals = self.retrieve_plot_values(left, right, unit, z_pos, sim_ind)
+        plot_range = PlotRange(left, right, unit)
-        return single_position_plot(vals, plt_range, self.params, ax, **kwargs)
+        vals = self.retrieve_plot_values(plot_range, z_pos, sim_ind)
+        return single_position_plot(vals, plot_range, self.params, ax, **kwargs)
 
     def plot_mean(
         self,
@@ -349,20 +348,25 @@ class Pulse(Sequence):
         z_pos: int,
         **kwargs,
     ):
-        plt_range, vals = self.retrieve_plot_values(left, right, unit, z_pos, slice(None))
+        plot_range = PlotRange(left, right, unit)
-        return mean_values_plot(vals, plt_range, self.params, ax, **kwargs)
+        vals = self.retrieve_plot_values(plot_range, z_pos, slice(None))
+        return mean_values_plot(vals, plot_range, self.params, ax, **kwargs)
 
-    def retrieve_plot_values(self, left, right, unit, z_pos, sim_ind):
+    def retrieve_plot_values(
-        plt_range = PlotRange(left, right, unit)
+        self, plot_range: PlotRange, z_pos: Optional[Union[int, float]], sim_ind: Optional[int]
-        if plt_range.unit.type == "TIME":
+    ):
+
+        if plot_range.unit.type == "TIME":
             vals = self.all_fields(ind=z_pos)
         else:
             vals = self.all_spectra(ind=z_pos)
-        if vals.ndim == 3:
+
-            vals = vals[:, sim_ind]
+        if sim_ind is None:
+            return vals
+        elif z_pos is None:
+            return vals[:, sim_ind]
         else:
-            vals = vals[sim_ind]
+            return vals[sim_ind]
-        return plt_range, vals
 
     def rin_propagation(
         self, left: float, right: float, unit: str

src/scgenerator/utils/__init__.py

@@ -17,11 +17,10 @@ from collections import abc
 from io import StringIO
 from pathlib import Path
 from string import printable as str_printable
+from functools import cache
 from typing import Any, Callable, Generator, Iterable, MutableMapping, Sequence, TypeVar, Union
 
 import numpy as np
-from numpy.lib.arraysetops import isin
-from numpy.lib.function_base import insert
 import pkg_resources as pkg
 import toml
 from tqdm import tqdm
@@ -37,8 +36,7 @@ PathTree = list[tuple[Path, ...]]
 
 class Paths:
     _data_files = [
-        "silica.toml",
+        "materials.toml",
-        "gas.toml",
         "hr_t.npz",
         "submit_job_template.txt",
         "start_worker.sh",
@@ -87,7 +85,12 @@ class Paths:
 def load_previous_spectrum(prev_data_dir: str) -> np.ndarray:
     prev_data_dir = Path(prev_data_dir)
     num = find_last_spectrum_num(prev_data_dir)
-    return np.load(prev_data_dir / SPEC1_FN.format(num))
+    return load_spectrum(prev_data_dir / SPEC1_FN.format(num))
+
+
+@cache
+def load_spectrum(folder: os.PathLike) -> np.ndarray:
+    return np.load(folder)
 
 
 def conform_toml_path(path: os.PathLike) -> str:
@@ -97,6 +100,12 @@ def conform_toml_path(path: os.PathLike) -> str:
     return path
 
 
+def open_single_config(path: os.PathLike) -> dict[str, Any]:
+    d = open_config(path)
+    f = d.pop("Fiber")[0]
+    return d | f
+
+
 def open_config(path: os.PathLike):
     """returns a dictionary parsed from the specified toml file
     This also handle having a 'INCLUDE' argument that will fill
@@ -215,10 +224,7 @@ def load_material_dico(name: str) -> dict[str, Any]:
     ----------
         material_dico : dict
     """
-    if name == "silica":
+    return toml.loads(Paths.gets("materials"))[name]
-        return toml.loads(Paths.gets("silica"))
-    else:
-        return toml.loads(Paths.gets("gas"))[name]
 
 
 def update_appended_params(source: Path, destination: Path, z: Sequence):

src/scgenerator/utils/parameter.py

@@ -10,8 +10,9 @@ from copy import copy, deepcopy
 from dataclasses import asdict, dataclass, fields
 from functools import cache, lru_cache
 from pathlib import Path
-from typing import Any, Callable, Generator, Iterable, Literal, Optional, TypeVar, Union
+from typing import Any, Callable, Generator, Iterable, Literal, Optional, Sequence, TypeVar, Union
 import numpy as np
+from numpy.lib import isin
 
 from .. import math, utils
 from ..const import PARAM_FN, PARAM_SEPARATOR, __version__
@@ -287,14 +288,14 @@ class Parameter:
 
     def __set__(self, instance, value):
         if isinstance(value, Parameter):
-            # defaut = None if self.default is None else copy(self.default)
+            defaut = None if self.default is None else copy(self.default)
-            # instance.__dict__[self.name] = defaut
+            instance.__dict__[self.name] = defaut
-            instance.__dict__[self.name] = None
+            # instance.__dict__[self.name] = None
         else:
             if value is not None:
-                self.validator(self.name, value)
                 if self.converter is not None:
                     value = self.converter(value)
+                self.validator(self.name, value)
             instance.__dict__[self.name] = value
 
     def display(self, num: float):
@@ -308,8 +309,11 @@ class Parameter:
             return f"{num_str} {unit}"
 
 
-def fiber_map_converter(d: dict[str, str]) -> dict[float, str]:
+def fiber_map_converter(d: dict[str, str]) -> list[tuple[float, str]]:
-    return {float(k): v for k, v in d.items()}
+    if isinstance(d, dict):
+        return [(float(k), v) for k, v in d.items()]
+    else:
+        return [(float(k), v) for k, v in d]
 
 
 @dataclass
@@ -340,7 +344,7 @@ class Parameters:
     pitch_ratio: float = Parameter(in_range_excl(0, 1))
     core_radius: float = Parameter(in_range_excl(0, 1e-3))
     he_mode: tuple[int, int] = Parameter(int_pair, default=(1, 1))
-    fit_parameters: tuple[int, int] = Parameter(int_pair, default=(0.08, 200e-9))
+    fit_parameters: tuple[int, int] = Parameter(float_pair, default=(0.08, 200e-9))
     beta2_coefficients: Iterable[float] = Parameter(num_list)
     dispersion_file: str = Parameter(string)
     model: str = Parameter(
@@ -425,13 +429,15 @@ class Parameters:
     const_qty: np.ndarray = Parameter(type_checker(np.ndarray))
     beta_func: Callable[[float], list[float]] = Parameter(func_validator)
     gamma_func: Callable[[float], float] = Parameter(func_validator)
-    fiber_map: dict[float, str] = Parameter(type_checker(dict), converter=fiber_map_converter)
+    fiber_map: list[tuple[float, str]] = Parameter(
+        validator_list(type_checker(tuple)), converter=fiber_map_converter
+    )
     datetime: datetime_module.datetime = Parameter(type_checker(datetime_module.datetime))
     version: str = Parameter(string)
 
     def prepare_for_dump(self) -> dict[str, Any]:
         param = asdict(self)
-        param["fiber_map"] = {str(z): n for z, n in param.get("fiber_map", {}).items()}
+        param["fiber_map"] = [(str(z), n) for z, n in param.get("fiber_map", [])]
         param = Parameters.strip_params_dict(param)
         param["datetime"] = datetime_module.datetime.now()
         param["version"] = __version__
@@ -896,7 +902,7 @@ class Configuration:
                         length = 0.0
                     fiber_map.append((length, this_conf["name"]))
                 params.output_path = str(sim_config.output_path)
-                params.fiber_map = dict(fiber_map)
+                params.fiber_map = fiber_map
                 yield sim_config.vary_list, params
 
     def __iter_1_sim(
@@ -1020,7 +1026,7 @@ class DataPather:
     def __init__(self, dl: list[dict[str, Any]]):
         self.dict_list = dl
 
-    def dico_iterator(
+    def vary_list_iterator(
         self, index: int
     ) -> Generator[tuple[tuple[tuple[int, ...]], list[list[tuple[str, Any]]]], None, None]:
         """iterates through every possible combination of a list of dict of lists
@@ -1050,6 +1056,8 @@ class DataPather:
             [[(a, 57), (b, "!")], [(c, 1)]],
         ]
         """
+        if index < 0:
+            index = len(self.dict_list) - index
         d_tem_list = [el for d in self.dict_list[: index + 1] for el in d.items()]
         dict_pos = np.cumsum([0] + [len(d) for d in self.dict_list[: index + 1]])
         ranges = [range(len(l)) for _, l in d_tem_list]
@@ -1062,7 +1070,7 @@ class DataPather:
             yield pos, out
 
     def all_vary_list(self, index):
-        for sim_index, l in self.dico_iterator(index):
+        for sim_index, l in self.vary_list_iterator(index):
             unique_vary: list[tuple[str, Any]] = []
             for ll in l[: index + 1]:
                 for pname, pval in ll:
@@ -1072,14 +1080,14 @@ class DataPather:
                             break
                     unique_vary.append((pname, pval))
             yield sim_index, format_variable_list(
-                reduce_all_variable(l[:index])
+                reduce_all_variable(l[:index]), add_iden=True
-            ), format_variable_list(reduce_all_variable(l)), unique_vary
+            ), format_variable_list(reduce_all_variable(l), add_iden=True), unique_vary
 
     def __repr__(self):
         return f"DataPather([{', '.join(repr(d) for d in self.dict_list)}])"
 
 
-@dataclass
+@dataclass(frozen=True)
 class PlotRange:
     left: float = Parameter(type_checker(int, float))
     right: float = Parameter(type_checker(int, float))
@@ -1194,7 +1202,7 @@ def _mock_function(num_args: int, num_returns: int) -> Callable:
     return out_func
 
 
-def format_variable_list(l: list[tuple[str, Any]]) -> str:
+def format_variable_list(l: list[tuple[str, Any]], add_iden=False) -> str:
     """formats a variable list into a str such that each simulation has a unique
     directory name. A u_XXX unique identifier and b_XXX (ignoring repeat simulations)
     branch identifier are added at the beginning.
@@ -1203,6 +1211,8 @@ def format_variable_list(l: list[tuple[str, Any]]) -> str:
     ----------
     l : list[tuple[str, Any]]
         list of variable parameters
+    add_iden : bool
+        add unique simulation and parameter-set identifiers
 
     Returns
     -------
@@ -1215,6 +1225,8 @@ def format_variable_list(l: list[tuple[str, Any]]) -> str:
         vs = format_value(p_name, p_value).replace("/", "").replace(PARAM_SEPARATOR, "")
         str_list.append(ps + PARAM_SEPARATOR + vs)
     tmp_name = PARAM_SEPARATOR.join(str_list)
+    if not add_iden:
+        return tmp_name
     unique_id = "u_" + utils.to_62(hash(str(l)))
     branch_id = "b_" + utils.to_62(hash(str([el for el in l if el[0] != "num"])))
     return unique_id + PARAM_SEPARATOR + branch_id + PARAM_SEPARATOR + tmp_name
@@ -1324,6 +1336,17 @@ def reduce_all_variable(all_variable: list[list[tuple[str, Any]]]) -> list[tuple
     return out
 
 
+def strip_vary_list(all_variable: T) -> T:
+    if len(all_variable) == 0:
+        return all_variable
+    elif isinstance(all_variable[0], Sequence) and (
+        len(all_variable[0]) == 0 or not isinstance(all_variable[0][0], str)
+    ):
+        return [strip_vary_list(el) for el in all_variable]
+    else:
+        return [el for el in all_variable if el[0] != "num"]
+
+
 default_rules: list[Rule] = [
     # Grid
     *Rule.deduce(