diff --git a/amical/data_processing.py b/amical/data_processing.py index a1b014c0..c9ddece3 100644 --- a/amical/data_processing.py +++ b/amical/data_processing.py @@ -117,13 +117,22 @@ def select_data(cube, clip_fact=0.5, clip=False, verbose=True, display=True): ind_clip = [] cube_cleaned_checked = np.array(good_fram) + cube_cleaned_checked = np.array(cube_cleaned_checked) + ind_clip2 = np.where(fluxes <= limit_flux)[0] if ((worst_fr in ind_clip2) and clip) or (worst_fr in flag_fram): ext = "(rejected)" else: ext = "" - diffmm = 100 * abs(np.max(fluxes) - np.min(fluxes)) / med_flux + fluxes_check = np.array([x.sum() for x in cube_cleaned_checked]) + + best_fr = np.argmax(fluxes_check) + worst_fr = np.argmin(fluxes_check) + + med_flux = np.median(fluxes_check) + std_flux = np.std(fluxes_check) + diffmm = 100 * abs(np.max(fluxes_check) - np.min(fluxes_check)) / med_flux if display: import matplotlib.pyplot as plt from matplotlib.colors import PowerNorm @@ -429,7 +438,7 @@ def show_clean_params( noBadPixel = False bad_pix_x, bad_pix_y = [], [] - if np.any(bmap0): + if np.any(bmap0) or np.any(ab0): if len(ab0) != 0: for j in range(len(ab0)): bmap0[ab0[j][1], ab0[j][0]] = 1 @@ -544,6 +553,39 @@ def _remove_dark(img1, darkfile=None, ihdu=0, verbose=False): return img1 +def _cosmic_bad_frames_finder(data, f_kernel): + """ + Summary + ------------- + Review the datacube and attempt centering. If centering fails, mark the frame as + a "bad frame" and exclude it from further cleaning steps. + + Parameters + ---------- + `data` : {numpy.array} + datacube containing the NRM data.\n + + Returns: + -------- + `bad_frames` : {list} + List of bad frames index.\n + """ + n_im = data.shape[0] + bad_frames = [] + for i in range(n_im): + filtmed = f_kernel is not None + try: + crop_max(data[i], 64, iframe=i, filtmed=filtmed, f=f_kernel)[0] + except ValueError: + bad_frames.append(i) + print( + "[AMICAL] %i unusable frames have been identified in the data cube," + % (len(bad_frames)) + + " primarily due to cosmic rays or persistent bad pixels." + ) + return bad_frames + + def clean_data( data, isz=None, @@ -585,6 +627,8 @@ def clean_data( bad_map, add_bad = _get_3d_bad_pixels(bad_map, add_bad, data) + bad_frames = _cosmic_bad_frames_finder(data, f_kernel) + for i in track(range(n_im), description="Cleaning"): img0 = data[i] img0 = _apply_edge_correction(img0, edge=edge) @@ -617,34 +661,42 @@ def clean_data( img_biased = img1.copy() img_biased[img_biased < 0] = 0 # Remove negative pixels - if isz is not None: - # Get expected center for sky correction - filtmed = f_kernel is not None - im_rec_max = crop_max( - img_biased, isz, offx=offx, offy=offy, filtmed=filtmed, f=f_kernel - )[0] - else: - im_rec_max = img_biased.copy() - - if ( - (im_rec_max.shape[0] != im_rec_max.shape[1]) - or (isz is not None and im_rec_max.shape[0] != isz) - or (isz is None and im_rec_max.shape[0] != img0.shape[0]) - ): - l_bad_frame.append(i) - else: - if apod and window is not None: - img = apply_windowing(im_rec_max, window=window) - elif apod: - warnings.warn( - "apod is set to True, but window is None. Skipping apodisation", - RuntimeWarning, - stacklevel=2, - ) - img = im_rec_max.copy() + if i not in bad_frames: + if isz is not None: + # Get expected center for sky correction + filtmed = f_kernel is not None + im_rec_max = crop_max( + img_biased, + isz, + iframe=i, + offx=offx, + offy=offy, + filtmed=filtmed, + f=f_kernel, + )[0] + else: + im_rec_max = img_biased.copy() + + if ( + (im_rec_max.shape[0] != im_rec_max.shape[1]) + or (isz is not None and im_rec_max.shape[0] != isz) + or (isz is None and im_rec_max.shape[0] != img0.shape[0]) + ): + l_bad_frame.append(i) else: - img = im_rec_max.copy() + if apod and window is not None: + img = apply_windowing(im_rec_max, window=window) + elif apod: + warnings.warn( + "apod is set to True, but window is None. Skipping apodisation", + RuntimeWarning, + stacklevel=2, + ) + img = im_rec_max.copy() + else: + img = im_rec_max.copy() cube_cleaned.append(img) + if verbose: print("Bad centering frame number:", l_bad_frame) cube_cleaned = np.array(cube_cleaned) diff --git a/amical/tools.py b/amical/tools.py index 90dd0895..f2c29bf0 100644 --- a/amical/tools.py +++ b/amical/tools.py @@ -79,7 +79,7 @@ def find_max(img, filtmed=True, f=3): return X, Y -def crop_max(img, dim, offx=0, offy=0, filtmed=True, f=3): +def crop_max(img, dim, iframe=0, offx=0, offy=0, filtmed=True, f=3): """ Summary ------------- @@ -112,7 +112,7 @@ def crop_max(img, dim, offx=0, offy=0, filtmed=True, f=3): if isz_max < dim: size_msg = ( f"The specified cropped image size, {dim}, is greater than the distance to" - " the PSF center in at least one dimension. The max size for this image is" + f" the PSF center in at least one dimension (frame {iframe}). The max size for this image is" f" {isz_max}" ) raise ValueError(size_msg)