From 9b86b13b364747d90284455e708ef285c76099d3 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Beno=C3=AEt=20Sierro?= Date: Fri, 8 Oct 2021 13:02:09 +0200 Subject: [PATCH] better sort_axis, extinction distance --- src/scgenerator/physics/fiber.py | 20 +++++------ src/scgenerator/physics/pulse.py | 2 +- src/scgenerator/utils/parameter.py | 57 ++++++++++++++---------------- 3 files changed, 37 insertions(+), 42 deletions(-) diff --git a/src/scgenerator/physics/fiber.py b/src/scgenerator/physics/fiber.py index ccddbb1..07653c8 100644 --- a/src/scgenerator/physics/fiber.py +++ b/src/scgenerator/physics/fiber.py @@ -2,6 +2,7 @@ from typing import Any, Iterable, Literal, TypeVar import numpy as np from numpy.fft import fft, ifft +from numpy import e from numpy.polynomial.chebyshev import Chebyshev, cheb2poly from scipy.interpolate import interp1d @@ -792,7 +793,7 @@ def compute_capillary_loss( interpolation_range: tuple[float, float], he_mode: tuple[int, int], ) -> np.ndarray: - mask = l < interpolation_range[1] + mask = (l < interpolation_range[1]) & (l > 0) alpha = capillary_loss(l[mask], he_mode, core_radius) out = np.zeros_like(l) out[mask] = alpha @@ -1099,14 +1100,12 @@ def effective_radius_HCARF(core_radius, t, f1, f2, wl_for_disp): return f1 * core_radius * (1 - f2 * wl_for_disp ** 2 / (core_radius * t)) -def capillary_loss( - wl_for_disp: np.ndarray, he_mode: tuple[int, int], core_radius: float -) -> np.ndarray: +def capillary_loss(wl: np.ndarray, he_mode: tuple[int, int], core_radius: float) -> np.ndarray: """computes the wavelenth dependent capillary loss according to Marcatili Parameters ---------- - wl_for_disp : np.ndarray, shape (n, ) + wl : np.ndarray, shape (n, ) wavelength array he_mode : tuple[int, int] tuple of mode (n, m) @@ -1118,9 +1117,10 @@ def capillary_loss( np.ndarray loss in 1/m """ - alpha = np.zeros_like(wl_for_disp) - mask = wl_for_disp > 0 - chi_silica = mat.sellmeier(wl_for_disp[mask], utils.load_material_dico("silica")) + chi_silica = mat.sellmeier(wl, utils.load_material_dico("silica")) nu_n = 0.5 * (chi_silica + 2) / np.sqrt(chi_silica) - alpha[mask] = nu_n * (u_nm(*he_mode) * wl_for_disp[mask] / pipi) ** 2 * core_radius ** -3 - return alpha + return nu_n * (u_nm(*he_mode) * wl / pipi) ** 2 * core_radius ** -3 + + +def extinction_distance(loss: T, ratio=1 / e) -> T: + return np.log(ratio) / -loss diff --git a/src/scgenerator/physics/pulse.py b/src/scgenerator/physics/pulse.py index f44191e..ae570f1 100644 --- a/src/scgenerator/physics/pulse.py +++ b/src/scgenerator/physics/pulse.py @@ -275,7 +275,7 @@ def conform_pulse_params( return width, t0, peak_power, energy else: if soliton_num is None: - soliton_num = np.sqrt(peak_power * gamma * t0 ** 2 / abs(beta2)) + soliton_num = np.sqrt(peak_power * gamma * L_D(t0, beta2)) return width, t0, peak_power, energy, soliton_num diff --git a/src/scgenerator/utils/parameter.py b/src/scgenerator/utils/parameter.py index b23450d..ff161ea 100644 --- a/src/scgenerator/utils/parameter.py +++ b/src/scgenerator/utils/parameter.py @@ -1095,6 +1095,12 @@ class PlotRange: unit: Callable[[float], float] = Parameter(units.is_unit, converter=units.get_unit) conserved_quantity: bool = Parameter(boolean, default=True) + def __post_init__(self): + if self.left >= self.right: + raise ValueError( + f"left value {self.left!r} must be strictly smaller than right value {self.right!r}" + ) + def __str__(self): return f"{self.left:.1f}-{self.right:.1f} {self.unit.__name__}" @@ -1110,51 +1116,40 @@ class PlotRange: def sort_axis( axis: np.ndarray, plt_range: PlotRange ) -> tuple[np.ndarray, np.ndarray, tuple[float, float]]: - """ - given an axis, returns this axis cropped according to the given range, converted and sorted + """convert an array according to the given range Parameters ---------- - axis : 1D array containing the original axis (usual the w or t array) - plt_range : tupple (min, max, conversion_function) used to crop the axis + axis : np.ndarray, shape (n,) + array + plt_range : PlotRange + range to crop in Returns ------- - cropped : the axis cropped, converted and sorted - indices : indices to use to slice and sort other array in the same fashion - extent : tupple with min and max of cropped + np.ndarray + new array converted to the desired unit and cropped in the given range + np.ndarray + indices of the concerved values + tuple[float, float] + actual minimum and maximum of the new axis Example ------- - w = np.append(np.linspace(0, -10, 20), np.linspace(0, 10, 20)) - t = np.linspace(-10, 10, 400) - W, T = np.meshgrid(w, t) - y = np.exp(-W**2 - T**2) - - # Define ranges - rw = (-4, 4, s) - rt = (-2, 6, s) - - w, cw = sort_axis(w, rw) - t, ct = sort_axis(t, rt) - - # slice y according to the given ranges - y = y[ct][:, cw] + >> sort_axis([18.0, 19.0, 20.0, 13.0, 15.2], PlotRange(1400, 1900, "cm")) + ([1520.0, 1800.0, 1900.0], [4, 0, 1], (1520.0, 1900.0)) """ if isinstance(plt_range, tuple): plt_range = PlotRange(*plt_range) - r = np.array((plt_range.left, plt_range.right), dtype="float") - indices = np.arange(len(axis))[ - (axis <= np.max(plt_range.unit(r))) & (axis >= np.min(plt_range.unit(r))) - ] - cropped = axis[indices] - order = np.argsort(plt_range.unit.inv(cropped)) - indices = indices[order] - cropped = cropped[order] - out_ax = plt_range.unit.inv(cropped) + masked = np.ma.array(axis, mask=~np.isfinite(axis)) + converted = plt_range.unit.inv(masked) + converted[(converted < plt_range.left) | (converted > plt_range.right)] = np.ma.masked + indices = np.arange(len(axis))[~converted.mask] + cropped = converted.compressed() + order = cropped.argsort() - return out_ax, indices, (out_ax[0], out_ax[-1]) + return cropped[order], indices[order], (cropped.min(), cropped.max()) def get_arg_names(func: Callable) -> list[str]: