Source code for toolbox_scs.base.knife_edge

import numpy as np
from scipy import special
from scipy.optimize import curve_fit


__all__ = ['knife_edge', 'knife_edge_base']


[docs]def knife_edge(positions, intensities, axisRange=None, p0=None): """ Calculates the beam radius at 1/e^2 from a knife-edge scan by fitting with erfc function: f(a,b,u) = a*erfc(u) + b or where u = sqrt(2)*(x-x0)/w0 with w0 the beam radius at 1/e^2 and x0 the beam center. Parameters ---------- positions : np.ndarray Motor position values, typically 1D intensities : np.ndarray Intensity values, could be either 1D or 2D, with the number or rows equivalent with the position size axisRange : sequence of two floats or None Edges of the scanning axis between which to apply the fit. p0 : list of floats, numpy 1D array Initial parameters used for the fit: x0, w0, a, b. If None, a beam radius of 100 um is assumed. Returns ------- width : float The beam radius at 1/e^2 std : float The standard deviation of the width """ popt, pcov = knife_edge_base(positions, intensities, axisRange=axisRange, p0=p0) width, std = 0, 0 if popt is not None and pcov is not None: width, std = np.abs(popt[1]), pcov[1, 1]**0.5 return width, std
[docs]def knife_edge_base(positions, intensities, axisRange=None, p0=None): """ The base implementation of the knife-edge scan analysis. Calculates the beam radius at 1/e^2 from a knife-edge scan by fitting with erfc function: f(a,b,u) = a*erfc(u) + b or where u = sqrt(2)*(x-x0)/w0 with w0 the beam radius at 1/e^2 and x0 the beam center. Parameters ---------- positions : np.ndarray Motor position values, typically 1D intensities : np.ndarray Intensity values, could be either 1D or 2D, with the number or rows equivalent with the position size axisRange : sequence of two floats or None Edges of the scanning axis between which to apply the fit. p0 : list of floats, numpy 1D array Initial parameters used for the fit: x0, w0, a, b. If None, a beam radius of 100 um is assumed. Returns ------- popt : sequence of floats or None The parameters of the resulting fit. pcov : sequence of floats The covariance matrix of the resulting fit. """ # Prepare arrays positions, intensities = prepare_arrays(positions, intensities, xRange=axisRange) # Estimate initial fitting params if p0 is None: p0 = [np.mean(positions), 0.1, np.max(intensities) / 2, 0] # Fit popt, pcov = function_fit(erfc, positions, intensities, p0=p0) return popt, pcov
def function_fit(func, x, y, **kwargs): """A wrapper for scipy.optimize curve_fit() """ # Fit try: popt, pcov = curve_fit(func, x, y, **kwargs) except (TypeError, RuntimeError) as err: print("Fit did not converge:", err) popt, pcov = (None, None) return popt, pcov def prepare_arrays(arrX: np.ndarray, arrY: np.ndarray, xRange=None, yRange=None): """ Preprocessing of the input x and y arrays. This involves the following steps. 1. Converting the arrays to 1D of the same size 2. Select the ranges from the input x- and y-ranges 3. Retrieve finite values. """ # Convert both arrays to 1D of the same size arrX, arrY = arrays_to1d(arrX, arrY) # Select ranges if xRange is not None: low, high = xRange if low == high: raise ValueError('The supplied xRange is not a valid range.') mask_ = range_mask(arrX, low, high) arrX = arrX[mask_] arrY = arrY[mask_] if yRange is not None: low, high = yRange if low == high: raise ValueError('The supplied xRange is not a valid range.') mask_ = range_mask(arrY, low, high) arrX = arrX[mask_] arrY = arrY[mask_] # Clean both arrays by only getting finite values finite_idx = np.isfinite(arrX) & np.isfinite(arrY) arrX = arrX[finite_idx] arrY = arrY[finite_idx] return arrX, arrY def arrays_to1d(arrX: np.ndarray, arrY: np.ndarray): """Flatten two arrays and matches their sizes """ assert arrX.shape[0] == arrY.shape[0] arrX, arrY = arrX.flatten(), arrY.flatten() if len(arrX) > len(arrY): arrY = np.repeat(arrY, len(arrX) // len(arrY)) else: arrX = np.repeat(arrX, len(arrY) // len(arrX)) return arrX, arrY def range_mask(array, minimum=None, maximum=None): """Retrieve the resulting array from the given minimum and maximum """ default = np.ones(array.shape, dtype=bool) min_slice, max_slice = default, default if minimum is not None: if minimum > np.nanmax(array): raise ValueError('The range minimum is too large.') min_slice = array >= minimum if maximum is not None: if maximum < np.nanmin(array): raise ValueError('The range maximum is too small.') max_slice = array <= maximum return min_slice & max_slice def erfc(x, x0, w0, a, b): return a * special.erfc(np.sqrt(2) * (x - x0) / w0) + b