diff --git a/src/vip_hci/fm/fakecomp.py b/src/vip_hci/fm/fakecomp.py index 7f8d7cec..2b8e18c0 100644 --- a/src/vip_hci/fm/fakecomp.py +++ b/src/vip_hci/fm/fakecomp.py @@ -824,7 +824,7 @@ def psf_norm_2d(psf, fwhm, threshold, mask_core, full_output, verbose): def cube_planet_free(planet_parameter, cube, angs, psfn, imlib='vip-fft', interpolation='lanczos4', transmission=None, - radial_gradient=False): + radial_gradient=False, weights=None): """Return a cube in which we have injected negative fake companion at the\ position/flux given by planet_parameter. @@ -863,6 +863,10 @@ def cube_planet_free(planet_parameter, cube, angs, psfn, imlib='vip-fft', at the very edge of a physical mask (e.g. ALC) or the effect on the light distribution of a marginally extended source near the IWA of the coronagraph. + weights : 1d array, optional + If provided, the negative fake companion fluxes will be scaled according + to these weights before injection in the cube. Can reflect changes in + the observing conditions throughout the sequence. Returns ------- @@ -906,8 +910,12 @@ def cube_planet_free(planet_parameter, cube, angs, psfn, imlib='vip-fft', transmission=transmission, radial_gradient=radial_gradient) else: + if weights is None: + flevel = -planet_parameter[i, 2] + else: + flevel = -planet_parameter[i, 2]*weights cpf = cube_inject_companions(cube_temp, psfn, angs, n_branches=1, - flevel=-planet_parameter[i, 2], + flevel=flevel, rad_dists=[planet_parameter[i, 0]], theta=planet_parameter[i, 1], imlib=imlib, verbose=False, diff --git a/src/vip_hci/metrics/snr_source.py b/src/vip_hci/metrics/snr_source.py index 1ff67088..996d4c6c 100644 --- a/src/vip_hci/metrics/snr_source.py +++ b/src/vip_hci/metrics/snr_source.py @@ -66,6 +66,9 @@ def snrmap(array, fwhm, approximated=False, plot=False, known_sources=None, use2alone: bool, optional Whether to use array2 alone to estimate the noise (might be useful to estimate the snr of extended disk features). + exclude_negative_lobes : bool, opt + Whether to include the adjacent aperture lobes to the tested location + or not. Can be set to True if the image shows significant neg lobes. verbose: bool, optional Whether to print timing or not. **kwargs : dictionary, optional diff --git a/src/vip_hci/preproc/badframes.py b/src/vip_hci/preproc/badframes.py index 87d3344a..09429219 100644 --- a/src/vip_hci/preproc/badframes.py +++ b/src/vip_hci/preproc/badframes.py @@ -80,7 +80,8 @@ def cube_detect_badfr_pxstats(array, mode='annulus', in_radius=10, width=10, n = array.shape[0] res = cube_basic_stats(array, mode, radius=in_radius, - inner_radius=in_radius, size=width, full_output=True) + inner_radius=in_radius, size=width, + full_output=True) if method == 'mean': mean_values = res[0] else: @@ -97,12 +98,7 @@ def cube_detect_badfr_pxstats(array, mode='annulus', in_radius=10, width=10, top_boundary = np.empty([n]) bot_boundary = np.empty([n]) for i in range(n): - if mode == 'annulus': - i_mean_value = get_annulus_segments(array[i], width=width, - inner_radius=in_radius, - mode="val")[0].mean() - elif mode == 'circle': - i_mean_value = mean_values[i] + i_mean_value = mean_values[i] top_boundary[i] = mean_smooth[i] + top_sigma*sigma bot_boundary[i] = mean_smooth[i] - low_sigma*sigma if i_mean_value > top_boundary[i] or i_mean_value < bot_boundary[i]: diff --git a/src/vip_hci/preproc/recentering.py b/src/vip_hci/preproc/recentering.py index 5b7fff73..6719a62f 100644 --- a/src/vip_hci/preproc/recentering.py +++ b/src/vip_hci/preproc/recentering.py @@ -1251,6 +1251,8 @@ def cube_recenter_dft_upsampling(array, upsample_factor=100, subi_size=None, start_time = time_ini() check_array(array, dim=3) + if array.shape[0] < 2: + raise TypeError("Input array cube should have more than 1 frame.") if mask is not None: if mask.shape != array.shape[-2:]: msg = "If provided, mask should have same shape as frames" @@ -1296,7 +1298,7 @@ def cube_recenter_dft_upsampling(array, upsample_factor=100, subi_size=None, if nproc is None: nproc = cpu_count() // 2 # Hyper-threading doubles the # of cores - if nproc == 1: + if nproc == 1 or n_frames < nproc: for i in Progressbar(range(1, n_frames), desc="frames", verbose=verbose): y[i], x[i], array_rec[i] = _shift_dft(array_rec, array_rec, i, diff --git a/src/vip_hci/preproc/rescaling.py b/src/vip_hci/preproc/rescaling.py index d338f94f..5af4b51c 100644 --- a/src/vip_hci/preproc/rescaling.py +++ b/src/vip_hci/preproc/rescaling.py @@ -807,6 +807,7 @@ def find_scal_vector( interpolation="lanczos4", hpf=False, fwhm_max=5, + check_vec=True, **kwargs ): """ @@ -845,6 +846,10 @@ def find_scal_vector( fwhm_max: float, optional Maximum FWHM of the PSF across all wavelengths, in pixels. Only used if hpf is set to True. + check_vec: True, optional + Whether to check (and force) the output vector to be >=1. This will set + the smallest scaling factor to 1, and scale accordingly the other + factors. **kwargs: optional Optional arguments to the scipy.optimize.minimize function @@ -908,7 +913,8 @@ def find_scal_vector( scal_vec[z] = scal_fac flux_vec[z] = flux_fac - scal_vec = check_scal_vector(scal_vec) + if check_vec: + scal_vec = check_scal_vector(scal_vec) return scal_vec, flux_vec diff --git a/src/vip_hci/psfsub/pca_fullfr.py b/src/vip_hci/psfsub/pca_fullfr.py index f1f16175..9f04b556 100644 --- a/src/vip_hci/psfsub/pca_fullfr.py +++ b/src/vip_hci/psfsub/pca_fullfr.py @@ -1032,6 +1032,22 @@ def _adi_rdi_pca( interpolation=interpolation, **rot_options, ) + if smooth is not None: + nout = len(gridre) + pca_res = cube_filter_lowpass(gridre[0], mode='gauss', + fwhm_size=smooth, verbose=False) + if full_output==True and source_xy is not None: + frame = frame_filter_lowpass(gridre[1], mode='gauss', + fwhm_size=smooth) + else: + frame = gridre[1] + gridre_new = [pca_res, frame] + if nout>2: + diff = nout-2 + for d in range(diff): + gridre_new.append(gridre[2+d]) + gridre = tuple(gridre_new) + return gridre diff --git a/src/vip_hci/var/filters.py b/src/vip_hci/var/filters.py index 6db42364..16ac1717 100644 --- a/src/vip_hci/var/filters.py +++ b/src/vip_hci/var/filters.py @@ -52,6 +52,8 @@ no_opencv = True import numpy as np from scipy.ndimage import median_filter +from skimage.filters import median +from skimage.morphology import disk from skimage.restoration import richardson_lucy from astropy.convolution import (convolve_fft, convolve, Gaussian2DKernel) from astropy.convolution import interpolate_replace_nans as interp_nan @@ -62,7 +64,7 @@ def cube_filter_iuwt(cube, coeff=5, rel_coeff=1, full_output=False): """ - Isotropic Undecimated Wavelet Transform filtering, as implemented in + Isotropic Undecimated Wavelet Transform filtering, as implemented in\ [KEN15]_ and detailed in [DAB15]_. Parameters @@ -203,6 +205,8 @@ def frame_filter_highpass(array, mode, median_size=5, kernel_size=5, ``kernel_size``) and using the ``convolve_fft`` Astropy function. ``median-subt`` subtracts a median low-pass filtered version of the image. + ``cmedian-subt`` + subtracts a circular median low-pass filtered version of the image. ``gauss-subt`` subtracts a Gaussian low-pass filtered version of the image. ``fourier-butter`` @@ -210,8 +214,9 @@ def frame_filter_highpass(array, mode, median_size=5, kernel_size=5, ``hann`` uses a Hann window. - median_size : int, optional - Size of the median box for the ``median-subt`` filter. + median_size : int or float, optional + Size of the median box for the ``median-subt`` and ``cmedian-subt`` + filters. In the former case, input median_size should be int. kernel_size : int, optional Size of the Laplacian kernel used in ``laplacian`` mode. It must be an positive odd integer value. @@ -357,6 +362,12 @@ def round_away(x): median_size=median_size) filtered = array - medianed + elif mode == 'cmedian-subt': + # Subtracting the low_pass filtered (median) image from the image itself + medianed = frame_filter_lowpass(array, 'cmedian', + median_size=median_size) + filtered = array - medianed + elif mode == 'gauss-subt': # Subtracting the low_pass filtered (median) image from the image itself gaussed = frame_filter_lowpass(array, 'gauss', fwhm_size=fwhm_size, @@ -412,10 +423,11 @@ def frame_filter_lowpass(array, mode='gauss', median_size=5, fwhm_size=5, ---------- array : numpy ndarray Input array, 2d frame. - mode : {'median', 'gauss', 'psf'}, str optional - Type of low-pass filtering. - median_size : int, optional - Size of the median box for filtering the low-pass median filter. + mode : {'median', 'cmedian', 'gauss', 'psf'}, str optional + Type of low-pass filtering. 'cmedian' stands for circular median. + median_size : int or float, optional + Size of the median box for the low-pass median filter ('median' and + 'cmedian' modes). Should be an int for the 'median' mode. fwhm_size : float or tuple of 2 floats, optional Size of the Gaussian kernel for the low-pass Gaussian filter. If a tuple is provided, it should correspond to y and x kernel sizes, @@ -457,8 +469,6 @@ def frame_filter_lowpass(array, mode='gauss', median_size=5, fwhm_size=5, """ if array.ndim != 2: raise TypeError('Input array is not a frame or 2d array.') - if not isinstance(median_size, int): - raise ValueError('`Median_size` must be integer') if mask is not None: if mode == 'median': @@ -468,8 +478,14 @@ def frame_filter_lowpass(array, mode='gauss', median_size=5, fwhm_size=5, raise TypeError(msg) if mode == 'median': + if not isinstance(median_size, int): + raise ValueError('`Median_size` must be integer') # creating the low_pass filtered (median) image filtered = median_filter(array, median_size, mode='nearest') + elif mode == 'cmedian': + circle = disk(median_size) + # apply median filter with given footprint = structuring element + filtered = median(array, footprint=circle) elif mode == 'gauss': # 2d Gaussian filter kernel_sz_y = kernel_sz @@ -615,7 +631,7 @@ def cube_filter_lowpass(array, mode='gauss', median_size=5, fwhm_size=5, def frame_deconvolution(array, psf, n_it=30): """ - Iterative image deconvolution following the scikit-image implementation + Deconvolve image iteratively following the scikit-image implementation\ of the Richardson-Lucy algorithm, described in [RIC72]_ and [LUC74]_. Considering an image that has been convolved by the point spread function