diff --git a/config/calibration/mask_v1.X.6.yaml b/config/calibration/mask_v1.X.6.yaml index c3c7b254..9dee1944 100644 --- a/config/calibration/mask_v1.X.6.yaml +++ b/config/calibration/mask_v1.X.6.yaml @@ -31,13 +31,13 @@ dat: # Number of epochs - col_name: N_EPOCH - label: r"$n_{\rm epoch}$" + label: $n_{\rm epoch}$ kind: greater_equal value: 2 # Magnitude range - col_name: mag - label: mag range + label: r kind: range value: [15, 30] @@ -86,7 +86,7 @@ dat_ext: # Rough pointing coverage - col_name: npoint3 - label: r"$n_{\rm pointing}$" + label: $n_{\rm pointing}$ kind: greater_equal value: 3 diff --git a/notebooks/cosmo_val/catalog_paper_plot/hist_mag.py b/notebooks/cosmo_val/catalog_paper_plot/hist_mag.py new file mode 100644 index 00000000..4ced3036 --- /dev/null +++ b/notebooks/cosmo_val/catalog_paper_plot/hist_mag.py @@ -0,0 +1,549 @@ +# %% +# hist_mag.py +# +# Plot magnitude histogram for various cuts and selection criteria + +# %% +import matplotlib +import matplotlib.pylab as plt + +# enable autoreload for interactive sessions +from IPython import get_ipython +ipython = get_ipython() +if ipython is not None: + ipython.run_line_magic("matplotlib", "inline") + ipython.run_line_magic("reload_ext", "autoreload") + ipython.run_line_magic("autoreload", "2") + ipython.run_line_magic("reload_ext", "log_cell_time") + +# %% +import sys +import os +import re +import numpy as np +from astropy.io import fits +from io import StringIO + +from sp_validation import run_joint_cat as sp_joint +from sp_validation import util +from sp_validation.basic import metacal +from sp_validation import calibration +import sp_validation.cat as cat + +# %% +# Initialize calibration class instance +obj = sp_joint.CalibrateCat() + +config = obj.read_config_set_params("config_mask.yaml") + +test_only = True + +# %% +# Funcitons +def get_data(obj, test_only=False): + """Get Data. + + Returns catalogue. + + Parameters + ---------- + obj : CalibrateCat instance + Instance of CalibrateCat class + test_only : bool, optional + If True, only load a subset of data for testing; + default is False. + + """ + # Get data. Set load_into_memory to False for very large files + dat, dat_ext = obj.read_cat(load_into_memory=False) + + if test_only: + n_max = 1_000_000 + print(f"MKDEBUG testing only first {n_max} objects") + dat = dat[:n_max] + dat_ext = dat_ext[:n_max] + + return dat, dat_ext + + +def read_hist_data(hist_data_path): + """ + Read Hist Data. + + Read histogram data from npz file. + + Parameters + ---------- + hist_data_path : str + Path to the npz file containing histogram data + + Returns + ------- + hist_data : dict + Dictionary with keys for each selection criterion containing: + - 'counts': histogram counts + - 'bins': bin edges + - 'label': label for the histogram + + """ + loaded = np.load(hist_data_path, allow_pickle=True) + hist_data = {} + + for key in loaded.files: + data = loaded[key] + hist_data[key] = { + 'counts': data[0], + 'bins': data[1], + 'label': str(data[2]) + } + + return hist_data + + +def get_mask(masks, col_name): + """Get Mask. + + Returns mask corresponding to col_name. + + Parameters + ---------- + masks : list + List of mask objects + col_name : str + Column name to identify the mask + Returns + ------- + list + Mask object + integer + Mask position in list + + """ + # Get mask fomr masks with col_name = col_name + for idx, mask in enumerate(masks): + if mask._col_name == col_name: + return mask, idx + + +def compute_hist(masks, col_name, mask_cumul, mag, bins): + """ + Compute histogram for given mask and magnitude data. + + Parameters + ---------- + masks : list + List of mask objects + col_name : str + Column name to identify the mask + mask_cumul : array or None + Cumulative mask array + mag : array + Magnitude data + bins : array + Bin edges for histogram + + Returns + ------- + counts : array + Histogram counts + bins : array + Bin edges + label : str + Label for the histogram + n_valid : int + Number of valid data points after masking + mask_cumul : array + Updated cumulative mask array + + """ + this_mask, _ = get_mask(masks, col_name) + + # First time: + if mask_cumul is None: + mask_cumul = this_mask._mask + else: + mask_cumul &= this_mask._mask + + # Data values + my_mag = mag[mask_cumul] + + # Data count + n_valid = np.sum(mask_cumul) + + # Label + string_buffer = StringIO() + this_mask.print_condition(string_buffer, latex=True) + label = string_buffer.getvalue().strip() + if label == "": + label = col_name + print("MKDEBUG", label) + + counts, bin_edges = np.histogram(my_mag, bins=bins) + + return counts, bin_edges, rf"{label}", n_valid, mask_cumul + + +def plot_hist(counts, bins, label, alpha=1, ax=None, color=None): + """ + Plot histogram from counts and bins using bar plot. + + Parameters + ---------- + counts : array + Histogram counts + bins : array + Bin edges + label : str + Label for the histogram + alpha : float + Transparency for plot + ax : matplotlib axis + Axis to plot on + color : str or None + Color for the histogram + + """ + bin_centers = 0.5 * (bins[:-1] + bins[1:]) + ax.bar( + bin_centers, + counts, + width=np.diff(bins), + alpha=alpha, + label=label, + align='center', + color=color + ) + + +def plot_all_hists( + hist_data, + col_names, + figsize=10, + alpha=0.5, + color_map=None, + fraction=False, + ax=None, + out_path=None, +): + + if ax is None: + plt.figure() + fig, (ax) = plt.subplots( + 1, + 1, + figsize=(figsize, figsize) + ) + + counts0 = None + for col_name in col_names: + if col_name in hist_data: + + data = hist_data[col_name] + if fraction: + if counts0 is None: + counts0 = data["counts"] + counts = data["counts"] / counts0 + else: + counts = data["counts"] + + plot_hist( + counts, + data['bins'], + data['label'], + alpha=alpha, + ax=ax, + color=color_map[col_name] + ) + #print(f"{col_name}: n_valid = {data['n_valid']}") + + ax.set_xlabel('$r$') + ylabel = "fraction" if fraction else "number" + ax.set_ylabel(ylabel) + ax.set_xlim(17.5, 26.5) + if not fraction: + ax.legend() + + if out_path: + plt.tight_layout() + plt.savefig( + out_path, + dpi=150, + bbox_inches='tight' + ) + + +# %% +# Main program +scenario = 1 +hist_data_path = f"magnitude_histograms_data_scenario-{scenario}.npz" + +if os.path.exists(hist_data_path): + print(f"Histogram data file {hist_data_path} found.") + print("Reading and plotting.") + + dat = dat_ext = None + hist_data = read_hist_data(hist_data_path) + +else: + print(f"Histogram data file {hist_data_path} not found.") + print("Reading UNIONS cat and computing.") + + dat, dat_ext = get_data(obj, test_only=test_only) + hist_data = None + + +# %% +# Masking +# Get all masks, with or without dat, dat_ext +masks, labels = sp_joint.get_masks_from_config( + config, + dat, + dat_ext, + verbose=obj._params["verbose"], +) + +# %% +# Combine mask according to scenario +# List of basic masks to apply to all cases + +masks_labels_basic = ["overlap", "mag", "64_r"] +col_names = ["basic masks"] + +if scenario == 0: + masks_labels_basic.extend([ + "FLAGS", + "IMAFLAGS_ISO", + "NGMIX_MOM_FAIL", + "NGMIX_ELL_PSFo_NOSHEAR_0", + "NGMIX_ELL_PSFo_NOSHEAR_1", + "4_Stars", + "8_Manual", + "1024_Maximask", + ]) + + col_names.extend(["N_EPOCH", "npoint3", "metacal"]) + +elif scenario == 1: + + col_names.extend([ + "IMAFLAGS_ISO", + "FLAGS", + "NGMIX_MOM_FAIL", + "NGMIX_ELL_PSFo_NOSHEAR_0", + "NGMIX_ELL_PSFo_NOSHEAR_1", + "4_Stars", + "8_Manual", + "1024_Maximask", + "N_EPOCH", + "npoint3", + "metacal", + ]) + + combine_cols = { + "ngmix failures": [ + "NGMIX_MOM_FAIL", + "NGMIX_ELL_PSFo_NOSHEAR_0", + "NGMIX_ELL_PSFo_NOSHEAR_1", + ] + } + +# %% +# Combine columns if specified. +# Remove old columns after combining. +if combine_cols is not None: + for new_col, old_cols in combine_cols.items(): + print(f"Combining columns {old_cols} into {new_col}") + # Create combined mask + old_masks = [] + idx_first = None + for col in old_cols: + mask, idx = get_mask(masks, col) + old_masks.append(mask) + if idx_first is None: + idx_first = idx + if dat is not None: + print(f"Creating combined mask for {new_col}") + masks_combined = sp_joint.Mask.from_list( + old_masks, + label=new_col, + verbose=obj._params["verbose"], + ) + else: + print(f"Creating dummy mask for {new_col} (for plot label)") + masks_combined = sp_joint.Mask( + new_col, + new_col, + kind="none", + ) + masks.insert(idx, masks_combined) + col_names.insert(idx, new_col) + + for old_mask, old_col in zip(old_masks, old_cols): + masks.remove(old_mask) + col_names.remove(old_col) + + print("After combining: masks =", [mask._col_name for mask in masks]) + +# %% +# Createe list of masks +masks_basic = [] +for mask in masks: + if mask._col_name in masks_labels_basic: + masks_basic.append(mask) + +if dat is not None: + print("Creating combined basic mask") + masks_basic_combined = sp_joint.Mask.from_list( + masks_basic, + label="basic masks", + verbose=obj._params["verbose"], + ) +else: + # Create dummy combined mask (for plot label) + print("Creating dummy combined basic mask") + masks_basic_combined = sp_joint.Mask( + "basic masks", + "basic masks", + kind="none", + ) + +masks.append(masks_basic_combined) + +# %% +# Metacal mask (cuts) +mask_tmp = sp_joint.Mask( + "metacal", + "metacal", + kind="none", +) + +# %% +# Call metacal if data is available +if dat is not None: + cm = config["metacal"] + gal_metacal = metacal( + dat, + masks_basic_combined._mask, + snr_min=cm["gal_snr_min"], + snr_max=cm["gal_snr_max"], + rel_size_min=cm["gal_rel_size_min"], + rel_size_max=cm["gal_rel_size_max"], + size_corr_ell=cm["gal_size_corr_ell"], + sigma_eps=cm["sigma_eps_prior"], + global_R_weight=cm["global_R_weight"], + col_2d=False, + verbose=True, + ) + + g_corr_mc, g_uncorr, w, mask_metacal, c, c_err = ( + calibration.get_calibrated_m_c(gal_metacal) + ) + + # Convert index array to boolean mask + mask_metacal_bool = np.zeros(len(dat), dtype=bool) + mask_metacal_bool[mask_metacal] = True + + mask_tmp._mask = mask_metacal_bool + +masks.append(mask_tmp) + + +# %% +# Define magnitude bins +mag_bins = np.arange(15, 30, 0.05) +mag_centers = 0.5 * (mag_bins[:-1] + mag_bins[1:]) + +# %% +# Create figure with multiple subplots +figsize = 10 +alpha = 0.5 + +# Define explicit colors for each histogram +colors = [f'C{i}' for i in range(len(col_names))] # Use matplotlib default color cycle +color_map = dict(zip(col_names, colors)) + + +# %% +# If hist_data not loaded, compute it +if hist_data is None: + hist_data = {} + +if dat is not None: + # Get magnitude column + mag = dat['mag'] + + mask_cumul = None + for col_name in col_names: + counts, bins, label, n_valid, mask_cumul = compute_hist( + masks=masks, + col_name=col_name, + mask_cumul=mask_cumul, + mag=mag, + bins=mag_bins + ) + hist_data[col_name] = { + 'counts': counts, + 'bins': bins, + 'label': label, + 'n_valid': n_valid, + } +# %% +# Create plots +fig, axes = plt.subplots( + 1, + 2, + figsize=(2 * figsize, figsize) +) +# Plot histogram data +plot_all_hists( + hist_data, + col_names, + alpha=alpha, + color_map=color_map, + ax=axes[0], +) +plot_all_hists( + hist_data, + col_names, + alpha=alpha, + color_map=color_map, + fraction=True, + ax=axes[1], +) +plt.tight_layout() +out_path = f"magnitude_histograms_scenario-{scenario}.png" +plt.savefig( + out_path, + dpi=150, + bbox_inches='tight' +) + +# %% +# Save histogram data to file (only if we computed it) +if dat is not None: + np.savez( + hist_data_path, + **{ + key: np.array( + [ + val['counts'], + val['bins'], + val['label'], + val["n_valid"], + ], + dtype=object + ) + for key, val in hist_data.items() + } + ) + print(f"Histogram data saved to {hist_data_path}") + +# %% +if dat is not None: + obj.close_hd5() + +# %% +for mask in masks: + mask.print_condition(sys.stdout, latex=True) + +# %% diff --git a/src/sp_validation/__init__.py b/src/sp_validation/__init__.py index b31411d6..62e2307b 100644 --- a/src/sp_validation/__init__.py +++ b/src/sp_validation/__init__.py @@ -16,14 +16,14 @@ ] # Explicit imports to avoid circular issues -from . import util -from . import io -from . import basic -from . import galaxy -from . import cosmology -from . import calibration -from . import cat -from . import plot_style -from . import plots -from . import run_joint_cat -from . import survey \ No newline at end of file +#from . import util +#from . import io +#from . import basic +#from . import galaxy +#from . import cosmology +#from . import calibration +#from . import cat +#from . import plot_style +#from . import plots +#from . import run_joint_cat +#from . import survey diff --git a/src/sp_validation/basic.py b/src/sp_validation/basic.py index 9c0972ee..22ba28fa 100644 --- a/src/sp_validation/basic.py +++ b/src/sp_validation/basic.py @@ -15,8 +15,6 @@ from scipy.spatial import cKDTree from scipy.special import gamma -import matplotlib.pyplot as plt - from tqdm import tqdm import operator as op import itertools as itools @@ -361,27 +359,23 @@ def _masking_gal(self): [self.ns, self.m1, self.p1, self.m2, self.p2], ['ns', 'm1', 'p1', 'm2', 'p2'] ): - Tr_tmp = data['T'] - if self._size_corr_ell: - Tr_tmp *= ( - (1 - (data['g1'] ** 2 + data['g2'] ** 2)) - / (1 + (data['g1'] ** 2 + data['g2'] ** 2)) - ) + # Add SNR if available from sextractor if hasattr(self, 'snr_sextractor'): - snr_flux = self.snr_sextractor - else: - snr_flux = data['flux'] / data['flux_err'] - - mask_tmp = ( - (data['flag'] == 0) - & (Tr_tmp / data['Tpsf'] > self._rel_size_min) - & (Tr_tmp / data['Tpsf'] < self._rel_size_max) - & (snr_flux > self._snr_min) - & (snr_flux < self._snr_max) + data['snr_flux'] = self.snr_sextractor + + # Get masks using standalone function + # Data dict already has the required keys: g1, g2, T, Tpsf, flag, flux, flux_err + masks = get_metacal_masks_separate( + data, + snr_min=self._snr_min, + snr_max=self._snr_max, + rel_size_min=self._rel_size_min, + rel_size_max=self._rel_size_max, + size_corr_ell=self._size_corr_ell, ) - # Take care of rotated version - ind_masked = np.where(mask_tmp == True)[0] + # Extract indices from combined mask + ind_masked = np.where(masks['combined'] == True)[0] self.mask_dict[name] = ind_masked @@ -391,6 +385,7 @@ def _masking_gal_mom(self): ... """ + print("Masking of galaxies with moments measurements is deprecated") self.mask_dict = {} for data, name in zip( [self.ns, self.m1, self.p1, self.m2, self.p2], @@ -564,6 +559,124 @@ def _return(): ) +def get_metacal_masks_separate( + data, + prefix='NGMIX', + snr_min=10, + snr_max=500, + rel_size_min=0.5, + rel_size_max=3.0, + size_corr_ell=True, + shear_variant='NOSHEAR', + col_2d=True, +): + """Get Metacal Masks Separate. + + Get separate masks for metacalibration quality cuts without running + full metacal computation. + + Parameters + ---------- + data : array or Table or dict + input galaxy catalogue, or dictionary with keys 'T', 'Tpsf', + 'flag', and either 'flux'/'flux_err' or 'snr_flux'. + If size_corr_ell=True, also requires 'g1' and 'g2'. + prefix : str, optional, default='NGMIX' + column name prefix in catalogue (ignored if data is dict) + snr_min : float, optional, default=10 + signal-to-noise minimum + snr_max : float, optional, default=500 + signal-to-noise maximum + rel_size_min : float, optional, default=0.5 + relative size minimum + rel_size_max : float, optional, default=3.0 + relative size maximum + size_corr_ell : bool, optional, default=True + if True, correct size for ellipticity + shear_variant : str, optional, default='NOSHEAR' + shear variant name, e.g. 'NOSHEAR', '1M', '1P', '2M', '2P' + (ignored if data is dict) + col_2d : bool, optional, default=True + if True, ellipticity in one 2D column; if False, in two columns + (ignored if data is dict) + + Returns + ------- + dict + Dictionary containing separate masks: + - 'flag': boolean mask for flag == 0 + - 'rel_size': boolean mask for relative size cuts + - 'snr': boolean mask for SNR cuts + - 'combined': boolean mask for all cuts combined + + """ + # Handle dictionary input (from metacal class) + if isinstance(data, dict): + T = data['T'] + Tpsf = data['Tpsf'] + flag = data['flag'] + if 'snr_flux' in data: + snr_flux = data['snr_flux'] + else: + snr_flux = data['flux'] / data['flux_err'] + + # Only extract g1, g2 if needed for size correction + if size_corr_ell: + g1 = data['g1'] + g2 = data['g2'] + else: + # Handle catalogue input (raw data) + T = data[f'{prefix}_T_{shear_variant}'] + Tpsf = data[f'{prefix}_Tpsf_{shear_variant}'] + flag = data[f'{prefix}_FLAGS_{shear_variant}'] + + # Get SNR + if f'SNR_WIN' in data.dtype.names or (hasattr(data, 'colnames') and 'SNR_WIN' in data.colnames): + snr_flux = data['SNR_WIN'] + else: + flux = data[f'{prefix}_FLUX_{shear_variant}'] + flux_err = data[f'{prefix}_FLUX_ERR_{shear_variant}'] + snr_flux = flux / flux_err + + # Only extract g1, g2 if needed for size correction + if size_corr_ell: + if col_2d: + g1 = data[f'{prefix}_ELL_{shear_variant}'][:, 0] + g2 = data[f'{prefix}_ELL_{shear_variant}'][:, 1] + else: + g1 = data[f'{prefix}_ELL_{shear_variant}_0'] + g2 = data[f'{prefix}_ELL_{shear_variant}_1'] + + # Apply size correction for ellipticity if requested + Tr_tmp = T.copy() if hasattr(T, 'copy') else np.array(T) + if size_corr_ell: + Tr_tmp *= ( + (1 - (g1 ** 2 + g2 ** 2)) + / (1 + (g1 ** 2 + g2 ** 2)) + ) + + # Create separate masks + mask_flag = (flag == 0) + mask_rel_size = ( + (Tr_tmp / Tpsf > rel_size_min) + & (Tr_tmp / Tpsf < rel_size_max) + ) + mask_snr = ( + (snr_flux > snr_min) + & (snr_flux < snr_max) + ) + + # Combined mask + mask_combined = mask_flag & mask_rel_size & mask_snr + + return { + 'flag': mask_flag, + 'rel_size': mask_rel_size, + 'snr': mask_snr, + 'combined': mask_combined, + } + + def jackknif_weighted_average2( data, weights, diff --git a/src/sp_validation/plots.py b/src/sp_validation/plots.py index d9c5d0e5..d565048c 100644 --- a/src/sp_validation/plots.py +++ b/src/sp_validation/plots.py @@ -669,3 +669,86 @@ def hsp_map_logical_or(maps, verbose=False): ) return map_comb + + +def plot_area_mask(ra, dec, zoom, mask=None): + """Plot Area Mask. + + Create sky plot of objects. + + Parameters + ---------- + ra : list + R.A. coordinates + dec : list + Dec. coordinates + zoom : TBD + mask: TBD, optional + + """ + if mask is None: + mask == np.ones_like(ra) + + fig, axes = plt.subplots(nrows=1, ncols=2, figsize=(30,15)) + axes[0].hexbin(ra[mask], dec[mask], gridsize=100) + axes[1].hexbin(ra[mask & zoom], dec[mask & zoom], gridsize=200) + for idx in (0, 1): + axes[idx].set_xlabel("R.A. [deg]") + axes[idx].set_ylabel("Dec [deg]") + + +def sky_plots(dat, masks, labels, zoom_ra, zoom_dec): + """Sky Plots. + + Plot sky regions with different masks. + + Parameters + ---------- + masks : list + masks to be applied + labels : dict + labels for masks + zoom_ra : list + min and max R.A. for zoom-in plot + zoom_dec : list + min and max Dec. for zoom-in plot + + """ + ra = dat["RA"][:] + dec = dat["Dec"][:] + + zoom_ra = (room_ra[0] < dat["RA"]) & (dat["RA"] < zoom_ra[1]) + zoom_dec = (zoom_dec[0] < dat["Dec"]) & (dat["Dec"] < zoom_dec[1]) + zoom = zoom_ra & zoom_dec + + # No mask + plot_area_mask(ra, dec, zoom) + + # SExtractor and SP flags + m_flags = masks[labels["FLAGS"]]._mask & masks[labels["IMAFLAGS_ISO"]]._mask + plot_area_mask(ra, dec, zoom, mask=m_flags) + + # Overlap regions + m_over = masks[labels["overlap"]]._mask & m_flags + plot_area_mask(ra, dec, zoom, mask=m_over) + + # Coverage mask + m_point = masks[labels["npoint3"]]._mask & m_over + plot_area_mask(ra, dec, zoom, mask=m_point) + + # Maximask + m_maxi = masks[labels["1024_Maximask"]]._mask & m_point + plot_area_mask(ra, dec, zoom, mask=m_maxi) + + m_comb = mask_combined._mask + plot_area_mask(ra, dec, zoom, mask=m_comb) + + m_man = m_maxi & masks[labels["8_Manual"]]._mask + plot_area_mask(ra, dec, zoom, mask=m_man) + + m_halos = ( + m_maxi + & masks[labels['1_Faint_star_halos']]._mask + & masks[labels['2_Bright_star_halos']]._mask + ) + plot_area_mask(ra, dec, zoom, mask=m_halos) diff --git a/src/sp_validation/run_joint_cat.py b/src/sp_validation/run_joint_cat.py index ef404a58..0b17d274 100644 --- a/src/sp_validation/run_joint_cat.py +++ b/src/sp_validation/run_joint_cat.py @@ -1183,89 +1183,6 @@ def run(self): """ -def sky_plots(dat, masks, labels, zoom_ra, zoom_dec): - """Sky Plots. - - Plot sky regions with different masks. - - Parameters - ---------- - masks : list - masks to be applied - labels : dict - labels for masks - zoom_ra : list - min and max R.A. for zoom-in plot - zoom_dec : list - min and max Dec. for zoom-in plot - - """ - ra = dat["RA"][:] - dec = dat["Dec"][:] - - zoom_ra = (room_ra[0] < dat["RA"]) & (dat["RA"] < zoom_ra[1]) - zoom_dec = (zoom_dec[0] < dat["Dec"]) & (dat["Dec"] < zoom_dec[1]) - zoom = zoom_ra & zoom_dec - - # No mask - plot_area_mask(ra, dec, zoom) - - # SExtractor and SP flags - m_flags = masks[labels["FLAGS"]]._mask & masks[labels["IMAFLAGS_ISO"]]._mask - plot_area_mask(ra, dec, zoom, mask=m_flags) - - # Overlap regions - m_over = masks[labels["overlap"]]._mask & m_flags - plot_area_mask(ra, dec, zoom, mask=m_over) - - # Coverage mask - m_point = masks[labels["npoint3"]]._mask & m_over - plot_area_mask(ra, dec, zoom, mask=m_point) - - # Maximask - m_maxi = masks[labels["1024_Maximask"]]._mask & m_point - plot_area_mask(ra, dec, zoom, mask=m_maxi) - - m_comb = mask_combined._mask - plot_area_mask(ra, dec, zoom, mask=m_comb) - - m_man = m_maxi & masks[labels["8_Manual"]]._mask - plot_area_mask(ra, dec, zoom, mask=m_man) - - m_halos = ( - m_maxi - & masks[labels['1_Faint_star_halos']]._mask - & masks[labels['2_Bright_star_halos']]._mask - ) - plot_area_mask(ra, dec, zoom, mask=m_halos) - - -def plot_area_mask(ra, dec, zoom, mask=None): - """Plot Area Mask. - - Create sky plot of objects. - - Parameters - ---------- - ra : list - R.A. coordinates - dec : list - Dec. coordinates - zoom : TBD - mask: TBD, optional - - """ - if mask is None: - mask == np.ones_like(ra) - - fig, axes = plt.subplots(nrows=1, ncols=2, figsize=(30,15)) - axes[0].hexbin(ra[mask], dec[mask], gridsize=100) - axes[1].hexbin(ra[mask & zoom], dec[mask & zoom], gridsize=200) - for idx in (0, 1): - axes[idx].set_xlabel("R.A. [deg]") - axes[idx].set_ylabel("Dec [deg]") - - def confusion_matrix(mask, confidence_level=0.9): n_key = len(mask) @@ -1363,7 +1280,7 @@ class Mask(): """ - def __init__(self, col_name, label, kind="equal", value=0, dat=None, verbose=False): + def __init__(self, col_name, label, kind=None, value=0, dat=None, verbose=False): self._col_name = col_name self._label = label @@ -1441,8 +1358,7 @@ def print_strings(cls, coln, lab, num, fnum, f_out=None): print(msg) if f_out: print(msg, file=f_out) - - + def print_stats(self, num_obj, f_out=None): if self._num_ok is None: self._num_ok = sum(self._mask) @@ -1451,29 +1367,37 @@ def print_stats(self, num_obj, f_out=None): sf = f"{self._num_ok/num_obj:10.2%}" self.print_strings(self._col_name, self._label, si, sf, f_out=f_out) - def get_sign(self): + def get_sign(self, latex=False): sign = None - if self._kind =="equal": - sign = "=" - elif self._kind =="not_equal": - sign = "!=" - elif self._kind =="greater_equal": - sign = ">=" - elif self._kind =="smaller_equal": - sign = "<=" + if self._kind == "equal": + sign = "$=$" if latex else "=" + elif self._kind == "not_equal": + sign = "$\ne$" if latex else "!=" + elif self._kind in ("greater_equal", "range"): + sign = "$\leq$" if latex else ">=" + elif self._kind == "smaller_equal": + sign = "$\geq$" if latex else "<=" return sign + + def print_condition(self, f_out, latex=False): + + if self._value is None: + return "" + + sign = self.get_sign(latex=latex) - def print_summary(self, f_out): - print(f"[{self._label}]\t\t\t", file=f_out, end="") - - sign = self.get_sign() + name = self._label if latex else self._col_name if sign is not None: - print(f"{self._col_name} {sign} {self._value}", file=f_out) - + print(f"{name} {sign} {self._value}", file=f_out) + if self._kind == "range": - print(f"{self._value[0]} <= {self._col_name} <= {self._value[1]}", file=f_out) + print(f"{self._value[0]} {sign} {name} {sign} {self._value[1]}", file=f_out) + + def print_summary(self, f_out): + print(f"[{self._label}]\t\t\t", file=f_out, end="") + self.print_condition(f_out) def create_descr(self): """Create Descr.