From 98445271b80bf4042823067996de4d44d21f45d5 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Beno=C3=AEt=20Sierro?= Date: Wed, 22 Sep 2021 11:32:18 +0200 Subject: [PATCH] many changes --- setup.cfg | 4 +- src/scgenerator/cli/cli.py | 2 +- src/scgenerator/data/materials.toml | 252 +++++++++++++++++++++++++++ src/scgenerator/data/silica.toml | 31 +++- src/scgenerator/physics/__init__.py | 63 ++++--- src/scgenerator/physics/materials.py | 15 +- src/scgenerator/physics/simulate.py | 11 +- src/scgenerator/physics/units.py | 1 + src/scgenerator/plotting.py | 68 +++++++- src/scgenerator/scripts/__init__.py | 2 + src/scgenerator/spectra.py | 44 ++--- src/scgenerator/utils/__init__.py | 24 ++- src/scgenerator/utils/parameter.py | 57 ++++-- 13 files changed, 488 insertions(+), 86 deletions(-) create mode 100644 src/scgenerator/data/materials.toml diff --git a/setup.cfg b/setup.cfg index fa612cb..86b849d 100644 --- a/setup.cfg +++ b/setup.cfg @@ -30,9 +30,7 @@ install_requires = [options.package_data] scgenerator = data/hr_t.npz - data/default_params.json - data/gas.json - data/silica.json + data/materials.toml [options.packages.find] where = src diff --git a/src/scgenerator/cli/cli.py b/src/scgenerator/cli/cli.py index 23c45dd..2532fdf 100644 --- a/src/scgenerator/cli/cli.py +++ b/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) diff --git a/src/scgenerator/data/materials.toml b/src/scgenerator/data/materials.toml new file mode 100644 index 0000000..1eab050 --- /dev/null +++ b/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" diff --git a/src/scgenerator/data/silica.toml b/src/scgenerator/data/silica.toml index 70a9871..c847f3a 100644 --- a/src/scgenerator/data/silica.toml +++ b/src/scgenerator/data/silica.toml @@ -1,4 +1,29 @@ -[sellmeier] -B = [ 0.6961663, 0.4079426, 0.8974794,] -C = [ 4.67914826e-15, 1.35120631e-14, 9.79340025e-11,] +[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)", +] + +[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)"] diff --git a/src/scgenerator/physics/__init__.py b/src/scgenerator/physics/__init__.py index 0ae8a77..82bc01e 100644 --- a/src/scgenerator/physics/__init__.py +++ b/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)) diff --git a/src/scgenerator/physics/materials.py b/src/scgenerator/physics/materials.py index 9820649..812e60e 100644 --- a/src/scgenerator/physics/materials.py +++ b/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] - kind = 1 + ind = lambda_ < WL_THRESHOLD + 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.") diff --git a/src/scgenerator/physics/simulate.py b/src/scgenerator/physics/simulate.py index 6878261..4d34743 100644 --- a/src/scgenerator/physics/simulate.py +++ b/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) diff --git a/src/scgenerator/physics/units.py b/src/scgenerator/physics/units.py index fe21838..80089fe 100644 --- a/src/scgenerator/physics/units.py +++ b/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") diff --git a/src/scgenerator/plotting.py b/src/scgenerator/plotting.py index 2aa6e55..fbaa903 100644 --- a/src/scgenerator/plotting.py +++ b/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") - return + raise ValueError("Can only plot 1D values") 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 diff --git a/src/scgenerator/scripts/__init__.py b/src/scgenerator/scripts/__init__.py index 88bd0d1..ede8211 100644 --- a/src/scgenerator/scripts/__init__.py +++ b/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] diff --git a/src/scgenerator/spectra.py b/src/scgenerator/spectra.py index 2401b16..1ec018c 100644 --- a/src/scgenerator/spectra.py +++ b/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: - return self.cache[i] - spec = np.load(self.path / SPECN_FN.format(i)) + spec = load_spectrum(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) - return propagation_plot(vals, plt_range, self.params, ax, **kwargs) + plot_range = PlotRange(left, right, unit) + 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) - return single_position_plot(vals, plt_range, self.params, ax, **kwargs) + plot_range = PlotRange(left, right, unit) + 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)) - return mean_values_plot(vals, plt_range, self.params, ax, **kwargs) + plot_range = PlotRange(left, right, unit) + 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): - plt_range = PlotRange(left, right, unit) - if plt_range.unit.type == "TIME": + def retrieve_plot_values( + self, plot_range: PlotRange, z_pos: Optional[Union[int, float]], sim_ind: Optional[int] + ): + + 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 plt_range, vals + return vals[sim_ind] def rin_propagation( self, left: float, right: float, unit: str diff --git a/src/scgenerator/utils/__init__.py b/src/scgenerator/utils/__init__.py index 42f7e9e..b3296ac 100644 --- a/src/scgenerator/utils/__init__.py +++ b/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", - "gas.toml", + "materials.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("silica")) - else: - return toml.loads(Paths.gets("gas"))[name] + return toml.loads(Paths.gets("materials"))[name] def update_appended_params(source: Path, destination: Path, z: Sequence): diff --git a/src/scgenerator/utils/parameter.py b/src/scgenerator/utils/parameter.py index b017f9c..53bb8b5 100644 --- a/src/scgenerator/utils/parameter.py +++ b/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) - # instance.__dict__[self.name] = defaut - instance.__dict__[self.name] = None + defaut = None if self.default is None else copy(self.default) + instance.__dict__[self.name] = defaut + # 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]: - return {float(k): v for k, v in d.items()} +def fiber_map_converter(d: dict[str, str]) -> list[tuple[float, str]]: + 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]) - ), format_variable_list(reduce_all_variable(l)), unique_vary + reduce_all_variable(l[:index]), add_iden=True + ), 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(