diff --git a/amical/_cli/main.py b/amical/_cli/main.py index 20f34d6..5286f3e 100644 --- a/amical/_cli/main.py +++ b/amical/_cli/main.py @@ -59,13 +59,13 @@ def main(argv: Optional[List[str]] = None) -> int: action="store_true", help="Perform apodisation using a super-gaussian function " + "(known as windowing)." - + " The gaussian FWHM is set by the parameter `window`.", + + " The gaussian HWHM is set by the parameter `window`.", ) clean_parser.add_argument( "--window", default=65, type=int, - help="FWHM used for windowing (used with --apod)(default: %(default)s)", + help="HWHM used for windowing (used with --apod)(default: %(default)s)", ) clean_parser.add_argument( "--sky", diff --git a/amical/data_processing.py b/amical/data_processing.py index 1cbf2f0..c05e672 100644 --- a/amical/data_processing.py +++ b/amical/data_processing.py @@ -18,7 +18,7 @@ from rich import print as rprint from rich.progress import track -from amical.tools import apply_windowing, crop_max, find_max +from amical.tools import apply_windowing, crop_max, find_max, super_gaussian def _apply_patch_ghost(cube, xc, yc, radius=20, dx=0, dy=-200, method="bg"): @@ -369,11 +369,12 @@ def show_clean_params( f_kernel=3, offx=0, offy=0, - apod=False, + apod=True, window=None, *, ifu=False, mask=None, + window_contours=False, ): """Display the input parameters for the cleaning. @@ -390,6 +391,7 @@ def show_clean_params( `remove_bad` {bool}: If True, the bad pixels are removed using a gaussian interpolation,\n `nframe` {int}: Frame number to be shown (default: 0),\n `ihdu` {int}: Hdu number of the fits file. Normally 1 for NIRISS and 0 for SPHERE (default: 0). + `window_contours` {bool}: shown contours of the super-Gaussian windowing (default: False). """ import matplotlib.pyplot as plt from astropy.io import fits @@ -462,10 +464,6 @@ def show_clean_params( bg_x = bg_coords[0] bg_y = bg_coords[1] sky_method = "mask" - if window is not None: - r3 = window - x3 = r3 * np.cos(theta) + x0 - y3 = r3 * np.sin(theta) + y0 xs1, ys1 = x0 + isz // 2, y0 + isz // 2 xs2, ys2 = x0 - isz // 2, y0 + isz // 2 @@ -492,9 +490,44 @@ def show_clean_params( s=20, label="Pixels used for sky subtraction", ) - if apod: - if window is not None: - plt.plot(x3, y3, "--", label="Super-gaussian windowing") + + if window is not None: + # The window parameter gives the HWHM of the super-Gaussian. + # The value is used as the radius for the circle in the plot. + x3 = window * np.cos(theta) + x0 + y3 = window * np.sin(theta) + y0 + plt.plot(x3, y3, "--", label="Super-Gaussian windowing (HWHM)") + + if window_contours: + # Create distance grid for windowing, + # relative to the new image center (x0, y0) + y_coord = np.arange(img1.shape[0]) - y0 + x_coord = np.arange(img1.shape[1]) - x0 + xx_grid, yy_grid = np.meshgrid(x_coord, y_coord) + distance = np.hypot(xx_grid, yy_grid) + + # Create the super-Gaussian window function + # Mutiply the window value with 2 to change from HWHM to FWHM + super_gauss = super_gaussian(distance, sigma=window * 2) + + # Plot contours of the window function + # Create a new meshgrid because the coordinate system + # in the plot is relative to the bottom left corner + y_coord = np.arange(img1.shape[0]) + x_coord = np.arange(img1.shape[1]) + xx_grid, yy_grid = np.meshgrid(x_coord, y_coord) + levels = [0.1, 0.25, 0.5, 0.75, 0.9] + contours = plt.contour( + xx_grid, + yy_grid, + super_gauss, + levels=levels, + linestyles=":", + linewidths=0.8, + colors="white", + ) + plt.clabel(contours, contours.levels, inline=True, fontsize=7.0) + plt.plot(x0, y0, "+", color="c", ms=10, label="Centering position") plt.plot( [xs1, xs2, xs3, xs4, xs1], @@ -670,7 +703,7 @@ def select_clean_data( offx=0, offy=0, clip_fact=0.5, - apod=True, + apod=False, sky=True, window=None, darkfile=None, @@ -696,18 +729,22 @@ def select_clean_data( `edge` {int}: Patch the edges of the image (VLT/SPHERE artifact, default: {0}),\n `clip` {bool}: If True, sigma-clipping is used to reject frames with low integrated flux,\n `clip_fact` {float}: Relative sigma if rejecting frames by sigma-clipping,\n - `apod` {bool}: If True, apodisation is performed in the image plan using a super-gaussian - function (known as windowing). The gaussian FWHM is set by the parameter `window`,\n - `window` {float}: FWHM of the super-gaussian to apodise the image (smoothly go to zero - on the edges),\n + `apod` {bool}: If True, apodisation is performed in the image plan using a super-Gaussian + function (known as windowing). The Gaussian HWHM is set by the parameter `window`. This + parameter is deprecated and will be removed in a future release. Instead, the apodization + is applied when providing a value to the `window` parameter. Setting the argument of + `window` to None will not apply the super-Gaussian windowing,\n + `window` {float}: Half width at half maximum (HWHM) of the super-Gaussian to apodise + the image (smoothly go to zero on the edges). The windowing is not applied when the + argument is set to None,\n `sky` {bool}: If True, the sky is remove using the annulus technique (computed between `r1` - and `r1` + `dr`), + and `r1` + `dr`),\n `darkfile` {str}: If specified (default: None), the input dark (master_dark averaged if multiple integrations) is substracted from the raw image,\n image,\n `f_kernel` {float}: kernel size used in the applied median filter (to find the center). `remove_bad` {bool}: If True, the bad pixels are removed in the cleaning parameter - plots using a gaussian interpolation (default: {True}),\n + plots using a Gaussian interpolation (default: {True}),\n `nframe` {int}: Frame number used to show cleaning parameters (default: {0}),\n Returns: diff --git a/amical/tools.py b/amical/tools.py index 421fb5f..c2a173a 100644 --- a/amical/tools.py +++ b/amical/tools.py @@ -277,22 +277,60 @@ def cov2cor(cov): return cor, sigma -def super_gaussian(x, sigma, m, amp=1, x0=0): - sigma = float(sigma) - m = float(m) +def super_gaussian( + x: np.ndarray, sigma: float, m: float = 3.0, amp: float = 1.0, x0: float = 0.0 +) -> np.ndarray: + """ + Function for creating a super-Gaussian window. + + Parameters + ---------- + x : np.ndarray + 2D array with the distances of each pixel to the image center. + sigma : float + Full width at half maximum (FWHM) of the super-Gaussian + window. It is therefore not the standard deviation, as the + parameter name would suggest. + m : float + Exponent used for the super-Gaussian function (default: 3.0). + amp : float + Amplitude of the Gaussian function (default: 1.0). + x0 : float + Offset applied to the distances (default: 0.0) + + Returns + ------- + np.ndarray + 2D array with the super-Gaussian window function. + """ + return amp * ( - ( - np.exp( - -(2 ** (2 * m - 1)) - * np.log(2) - * (((x - x0) ** 2) / ((sigma) ** 2)) ** (m) - ) - ) + (np.exp(-(2 ** (2 * m - 1)) * np.log(2) * (((x - x0) ** 2) / (sigma**2)) ** m)) ** 2 ) -def apply_windowing(img, window=80, m=3): +def apply_windowing( + img: np.ndarray, window: float = 80.0, m: float = 3.0 +) -> np.ndarray: + """ + Function for applying a super-Gaussian window to an image. + + Parameters + ---------- + img : np.ndarray + 2D array with the input image. + window : float + Half width at half maximum (HWHM) of the window function + (default: 80.0). + m : float + Exponent used for the super-Gaussian function (default: 3.0). + + Returns + ------- + np.ndarray + 2D array with the windowed input image. + """ isz = len(img) xx, yy = np.arange(isz), np.arange(isz) xx2 = xx - isz // 2 @@ -301,11 +339,11 @@ def apply_windowing(img, window=80, m=3): distance = np.sqrt(xx2**2 + yy2[:, np.newaxis] ** 2) # Super-gaussian windowing - window = super_gaussian(distance, sigma=window * 2, m=m) + # Mutiply the window value with 2 to change from HWHM to FWHM + super_gauss = super_gaussian(distance, sigma=window * 2, m=m) # Apply the windowing - img_apod = img * window - return img_apod + return img * super_gauss def sanitize_array(dic): # pragma: no cover