diff --git a/examples/ssi/config_ssi_depth.yml b/examples/ssi/config_ssi_depth.yml new file mode 100755 index 00000000..b09b45fe --- /dev/null +++ b/examples/ssi/config_ssi_depth.yml @@ -0,0 +1,21 @@ +# Values in this section are accessible to all the different stages. +# They can be overridden by individual stages though. +global: + # This is read by many stages that read complete + # catalog data, and tells them how many rows to read + # at once + chunk_rows: 100_000 + # These mapping options are also read by a range of stages + pixelization: healpix + nside: 256 + sparse: true # Generate sparse maps - faster if using small areas + + +TXIngestSSIDetectionDESBalrog: + name: TXIngestSSIDetectionDESBalrog + +TXMapPlots: + name: TXMapPlots + projection: cart + rot180: True + diff --git a/examples/ssi/pipeline_ssi_depth.yml b/examples/ssi/pipeline_ssi_depth.yml new file mode 100755 index 00000000..15ab0069 --- /dev/null +++ b/examples/ssi/pipeline_ssi_depth.yml @@ -0,0 +1,41 @@ +# Stages to run +stages: + - name: TXIngestSSIDetectionDESBalrog + - name: TXIngestSSIMatchedDESBalrog + - name: TXAuxiliarySSIMaps + - name: TXMapPlotsSSI + +modules: txpipe + +# Where to put outputs +output_dir: data/example/outputs_ssi_des_depth/ + +# How to run the pipeline: mini, parsl, or cwl +launcher: + name: mini + interval: 1.0 + +# Where to run the pipeline: cori-interactive, cori-batch, or local +site: + name: local + max_threads: 2 + +# configuration settings +config: examples/ssi/config_ssi_depth.yml + +inputs: + # See README for paths to download these files + + #NERSC locations of public DES balrog data + balrog_detection_catalog: /global/cfs/cdirs/lsst/groups/WL/users/jelvinpo/balrog_txpipe_example/balrog_detection_catalog_sof_run2_v1.4.fits + balrog_matched_catalog: /global/cfs/cdirs/lsst/groups/WL/users/jelvinpo/balrog_txpipe_example/balrog_matched_catalog_sof_run2_v1.4.fits + + fiducial_cosmology: data/fiducial_cosmology.yml + +# if supported by the launcher, restart the pipeline where it left off +# if interrupted +resume: false +# where to put output logs for individual stages +log_dir: data/example/outputs_ssi_des_mag/logs/ +# where to put an overall parsl pipeline log +pipeline_log: data/example/outputs_ssi_des_mag/log.txt diff --git a/txpipe/auxiliary_maps.py b/txpipe/auxiliary_maps.py index 4c6afa8f..17391e04 100644 --- a/txpipe/auxiliary_maps.py +++ b/txpipe/auxiliary_maps.py @@ -256,3 +256,100 @@ def run(self): print(len(pix)) print(len(depth[pix])) f.write_map("depth/depth", pix, depth[pix], metadata) + +class TXAuxiliarySSIMaps(TXBaseMaps): + """ + Generate auxiliary maps from SSI catalogs + + This class generates maps of: + - depth (measured magnitude) + - depth (true magnitude) + """ + name = "TXAuxiliarySSIMaps" + dask_parallel = True + inputs = [ + ("injection_catalog", HDFFile), # for injection locations + ("matched_ssi_photometry_catalog", HDFFile), + ] + outputs = [ + ("aux_ssi_maps", MapsFile), + ] + + ################### + ################## + + config_options = { + "block_size": 0, + "depth_band": "i", # Make depth maps for this band + "snr_threshold": 10.0, # The S/N value to generate maps for (e.g. 5 for 5-sigma depth) + "snr_delta": 1.0, # The range threshold +/- delta is used for finding objects at the boundary + } + + def run(self): + # Import dask and alias it as 'da' + _, da = import_dask() + + + # Retrieve configuration parameters + block_size = self.config["block_size"] + if block_size == 0: + block_size = "auto" + band = self.config["depth_band"] + + # Open the input photometry catalog file. + # We can't use a "with" statement because we need to keep the file open + # while we're using dask. + f = self.open_input("matched_ssi_photometry_catalog", wrapper=True) + + # Load photometry data into dask arrays. + # This is lazy in dask, so we're not actually loading the data here. + ra = da.from_array(f.file["photometry/ra"], block_size) + block_size = ra.chunksize + dec = da.from_array(f.file["photometry/dec"], block_size) + snr = da.from_array(f.file[f"photometry/snr_{band}"], block_size) + mag_meas = da.from_array(f.file[f"photometry/mag_{band}"], block_size) + mag_true = da.from_array(f.file[f"photometry/inj_mag_{band}"], block_size) + + # Choose the pixelization scheme based on the configuration. + # Might need to review this to make sure we use the same scheme everywhere + pixel_scheme = choose_pixelization(**self.config) + + # Initialize a dictionary to store the maps. + # To start with this is all lazy too, until we call compute + maps = {} + + # Create depth maps using dask and measured magnitudes + pix2, count_map, depth_map, depth_var = make_dask_depth_map( + ra, dec, mag_meas, snr, self.config["snr_threshold"], self.config["snr_delta"], pixel_scheme) + maps["depth/depth_meas"] = (pix2, depth_map[pix2]) + maps["depth/depth_meas_count"] = (pix2, count_map[pix2]) + maps["depth/depth_meas_var"] = (pix2, depth_var[pix2]) + + + # Create depth maps using daskand true magnitudes + pix2, count_map, depth_map, depth_var = make_dask_depth_map( + ra, dec, mag_true, snr, self.config["snr_threshold"], self.config["snr_delta"], pixel_scheme) + maps["depth/depth_true"] = (pix2, depth_map[pix2]) + maps["depth/depth_true_count"] = (pix2, count_map[pix2]) + maps["depth/depth_true_var"] = (pix2, depth_var[pix2]) + + maps, = da.compute(maps) + + # Prepare metadata for the maps. Copy the pixelization-related + # configuration options only here + metadata = { + key: self.config[key] + for key in map_config_options + if key in self.config + } + # Then add the other configuration options + metadata["depth_band"] = band + metadata["depth_snr_threshold"] = self.config["snr_threshold"] + metadata["depth_snr_delta"] = self.config["snr_delta"] + metadata.update(pixel_scheme.metadata) + + # Write the output maps to the output file + with self.open_output("aux_ssi_maps", wrapper=True) as out: + for map_name, (pix, m) in maps.items(): + out.write_map(map_name, pix, m, metadata) + out.file['maps'].attrs.update(metadata) diff --git a/txpipe/data_types/types.py b/txpipe/data_types/types.py index c5d6d6bb..c012c9e5 100755 --- a/txpipe/data_types/types.py +++ b/txpipe/data_types/types.py @@ -274,12 +274,20 @@ def write_map(self, map_name, pixel, value, metadata): subgroup.create_dataset("pixel", data=pixel) subgroup.create_dataset("value", data=value) - def plot_healpix(self, map_name, view="cart", **kwargs): + def plot_healpix(self, map_name, view="cart", rot180=False, **kwargs): import healpy import numpy as np m, pix, nside = self.read_healpix(map_name, return_all=True) lon, lat = healpy.pix2ang(nside, pix, lonlat=True) + if rot180: #(optional) rotate 180 degrees in the lon direction + lon += 180 + lon[lon > 360.] -= 360. + pix_rot = healpy.ang2pix(nside, lon, lat, lonlat=True) + m_rot = np.ones(healpy.nside2npix(nside))*healpy.UNSEEN + m_rot[pix_rot] = m[pix] + m = m_rot + pix = pix_rot npix = healpy.nside2npix(nside) if len(pix) == 0: print(f"Empty map {map_name}") @@ -291,6 +299,7 @@ def plot_healpix(self, map_name, view="cart", **kwargs): lon_range = [lon[w].min() - 0.1, lon[w].max() + 0.1] lat_range = [lat[w].min() - 0.1, lat[w].max() + 0.1] lat_range = np.clip(lat_range, -90, 90) + lon_range = np.clip(lon_range, 0, 360.) m[m == 0] = healpy.UNSEEN title = kwargs.pop("title", map_name) if view == "cart": diff --git a/txpipe/ingest/__init__.py b/txpipe/ingest/__init__.py index 64b85659..016c0c5a 100644 --- a/txpipe/ingest/__init__.py +++ b/txpipe/ingest/__init__.py @@ -2,4 +2,4 @@ from .mocks import TXCosmoDC2Mock, TXGaussianSimsMock from .gcr import TXMetacalGCRInput, TXIngestStars from .redmagic import TXIngestRedmagic -from .ssi import TXIngestSSIGCR, TXMatchSSI, TXIngestSSIMatched, TXIngestSSIMatchedDESBalrog \ No newline at end of file +from .ssi import TXIngestSSIGCR, TXMatchSSI, TXIngestSSIDESBalrog, TXIngestSSIMatchedDESBalrog, TXIngestSSIDetectionDESBalrog \ No newline at end of file diff --git a/txpipe/ingest/ssi.py b/txpipe/ingest/ssi.py index 546496b3..cfad6a1d 100644 --- a/txpipe/ingest/ssi.py +++ b/txpipe/ingest/ssi.py @@ -16,7 +16,126 @@ "coord_dec": "dec", } -class TXIngestSSIGCR(PipelineStage): +class TXIngestSSI(PipelineStage): + """ + Base-Class for ingesting SSI catalogs + + These can be GCR catalogs, or precursor surveys + """ + + name = "TXIngestSSI" + + def create_photometry(self, input_name, output_name, column_names, dummy_columns): + """ + Create and save photometry dataset from an input catalog. + + This method sets up the output HDF5 file structure and processes photometry + data from the input catalog. It maps the input columns to the desired output + column names, creates datasets for photometric magnitudes and their errors, + and iterates over chunks of data to save them in the output file. + + Parameters: + ---------- + input_name : str + name of (FITS) input + + output_name : str + name of (HDF5) output + + column_names : dict + A dictionary mapping the input column names to the desired output column names. + Keys are input column names, and values are the corresponding output names. + + dummy_columns : dict + A dictionary of columns to be added to the output with fixed values. + Keys are the names of the dummy columns, and values are the constant values + to fill those columns. + + """ + # get some basic onfo about the input file + f = self.open_input(input_name) + n = f[1].get_nrows() + dtypes = f[1].get_rec_dtype()[0] + f.close() + + chunk_rows = self.config["chunk_rows"] + + # set up the output file columns + output, g = self.setup_output(output_name, column_names, dtypes, n) + + # iterate over the input file and save to the output columns + self.add_columns(g, input_name, column_names, chunk_rows, n) + + # Set up any dummy columns with sentinal values + # that were not in the original files + for col_name in dummy_columns.keys(): + g.create_dataset(col_name, data=np.full(n, dummy_columns[col_name])) + + output.close() + + def setup_output(self, output_name, column_names, dtypes, n): + """ + Set up the output HDF5 file structure. + + Parameters + ---------- + output_name : str + The name of the output HDF5 file. + + column_names : dict + dict of column names to include in the output file. + + dtypes : dict + A dictionary mapping column names to their corresponding data types. + + n : int + The total number of rows in the dataset. + + Returns + ------- + tuple + A tuple containing the open HDF5 output file object and the "photometry" group. + """ + cols = list(column_names.keys()) + output = self.open_output(output_name) + g = output.create_group("photometry") + for col in cols: + dtype = dtypes[col] + g.create_dataset(column_names[col], (n,), dtype=dtype) + + return output, g + + def add_columns(self, g, input_name, column_names, chunk_rows, n): + """ + Add data to the HDF5 output file in chunks. + + This method reads chunks of data from the input file and writes them + to the corresponding datasets in the output file. + + Parameters + ---------- + g : h5py.Group + The HDF5 group where data will be written. + + input_name : str + The name of the input file (e.g., a FITS file). + + column_names : dict + Dict of column names to read from the input file. + + chunk_rows : int + Number of rows to process in each chunk. + + n : int + Total number of rows in the dataset. + """ + cols = list(column_names.keys()) + for s, e, data in self.iterate_fits(input_name, 1, cols, chunk_rows): + print(s, e, n) + for col in cols: + g[column_names[col]][s:e] = data[col] + +class TXIngestSSIGCR(TXIngestSSI): """ Class for ingesting SSI catalogs using GCR @@ -57,7 +176,7 @@ def run(self): sys.path.insert(0, self.config["GCRcatalog_path"]) import GCRCatalogs - #attempt to ingest three types of catalog + # attempt to ingest three types of catalog catalogs_labels = [ "injection_catalog", "ssi_photometry_catalog", @@ -94,11 +213,15 @@ def run(self): group[column_names[q]] = group[q] except KeyError: # skip quantities that are missing - warnings.warn(f"quantity {q} was missing from the GCRCatalog object") + warnings.warn( + f"quantity {q} was missing from the GCRCatalog object" + ) continue except TypeError: - warnings.warn(f"Quantity {q} coud not be saved as it has a data type not recognised by hdf5") + warnings.warn( + f"Quantity {q} coud not be saved as it has a data type not recognised by hdf5" + ) # convert fluxes to mags using txpipe/utils/conversion.py bands = "ugrizy" @@ -116,17 +239,13 @@ def run(self): output_file.close() - class TXMatchSSI(PipelineStage): """ - Class for ingesting SSI injection and photometry catalogs + Class for matching an SSI injection catalog and a photometry catalog Default inputs are in TXPipe photometry catalog format - Will perform its own matching between the catalogs to output a - matched SSI catalog for further use - - TO DO: make a separate stage to ingest the matched catalog directly + Outpus a matched SSI catalog for further use """ name = "TXMatchSSI" @@ -192,10 +311,14 @@ def match_cats(self): ) out_start = 0 - for ichunk, (in_start,in_end,data) in enumerate(self.iterate_hdf("ssi_photometry_catalog", "photometry", ["ra","dec"], batch_size)): + for ichunk, (in_start, in_end, data) in enumerate( + self.iterate_hdf( + "ssi_photometry_catalog", "photometry", ["ra", "dec"], batch_size + ) + ): phot_coord = SkyCoord( - ra=data["ra"]*u.degree, - dec=data["dec"]*u.degree, + ra=data["ra"] * u.degree, + dec=data["dec"] * u.degree, ) idx, d2d, d3d = phot_coord.match_to_catalog_sky(inj_coord) @@ -276,30 +399,62 @@ def finalize_output(self, outfile, group, ntot): outfile.close() return - - -class TXIngestSSIMatched(PipelineStage): +class TXIngestSSIDESBalrog(TXIngestSSI): + """ + Base-stage for ingesting a DES SSI catalog AKA "Balrog" """ - Base-stage for ingesting a matched SSI catalog - This stage will just read in a file in a given format and output to a - HDF file that TXPIPE can use + name = "TXIngestSSIDESBalrog" - """ + def setup_output(self, output_name, column_names, dtypes, n): + """ + For balrog, we need to include if statements to catch the 2D data entries + and possibly add an extendedness column + """ + cols = list(column_names.keys()) + output = self.open_output(output_name) + g = output.create_group("photometry") + for col in cols: + dtype = dtypes[col] - name = "TXIngestSSIMatched" + if dtype.subdtype is not None: #this is a multi-dimentional column + assert dtype.subdtype[1]==(4,) #We are assuming this entry is a 2D array with 4 columns (corresponding to griz) + dtype = dtype.subdtype[0] + for b in "griz": + g.create_dataset(column_names[col] + f"_{b}", (n,), dtype=dtype) + else: + g.create_dataset(column_names[col], (n,), dtype=dtype) - outputs = [ - ("matched_ssi_photometry_catalog", HDFFile), - ] + if col == "meas_EXTENDED_CLASS_SOF": + #also create an "extendedness" column + g.create_dataset("extendedness", (n,), dtype=dtype) + + return output, g - config_options = { - "chunk_rows": 100_000, - "magnification":0, # magnification label - } + def add_columns(self, g, input_name, column_names, chunk_rows, n): + """ + For balrog, we need to include if statements to catch the 2D data entries + and possibly add a extendedness column + """ + cols = list(column_names.keys()) + for s, e, data in self.iterate_fits(input_name, 1, cols, chunk_rows): + print(s, e, n) + for col in cols: + if len(data[col].shape) == 2: + assert data[col].shape[1] == 4 + for iband, b in enumerate("griz"): + g[column_names[col] + f"_{b}"][s:e] = data[col][:, iband] + else: + g[column_names[col]][s:e] = data[col] + + if col == "meas_EXTENDED_CLASS_SOF": + # "meas_EXTENDED_CLASS_SOF" is (0 or 1) for star, (2 or 3) for galaxy + # extendedness is 0 for star, 1 for galaxy + extendedness = np.where((data[col] == 2) | (data[col] == 3), 1, 0) + g[column_names[col]][s:e] = extendedness -class TXIngestSSIMatchedDESBalrog(TXIngestSSIMatched): +class TXIngestSSIMatchedDESBalrog(TXIngestSSIDESBalrog): """ Class for ingesting a matched "SSI" catalog from DES (AKA Balrog) """ @@ -309,76 +464,97 @@ class TXIngestSSIMatchedDESBalrog(TXIngestSSIMatched): inputs = [ ("balrog_matched_catalog", FitsFile), ] - + + outputs = [ + ("matched_ssi_photometry_catalog", HDFFile), + ] + def run(self): """ Run the analysis for this stage. """ - print('Ingesting DES Balrog matched catalog') + print("Ingesting DES Balrog matched catalog") - #get some basic onfo about the input file - f = self.open_input("balrog_matched_catalog") - n = f[1].get_nrows() - dtypes = f[1].get_rec_dtype()[0] - f.close() - - print(f'{n} objects in matched catalog') - - chunk_rows = self.config["chunk_rows"] - - #we will only load a subset of columns to save space + # we will only load a subset of columns to save space column_names = { - "bal_id": "bal_id", # Unique identifier for object (created during balrog process) - "true_bdf_mag_deredden": "inj_mag", # Magnitude of the original deep field object, dereddened - "true_id": "inj_id", # Original coadd_obj_id of deep field object - "meas_id": "id", # Coadd_object_id of injection - "meas_ra": "ra", # measured RA of the injection - "meas_dec": "dec", # measured DEC of the injection - "meas_cm_mag_deredden": "mag", # measured magnitude of the injection - "meas_cm_T": "cm_T", # measured size parameter T (x^2+y^2) - "meas_EXTENDED_CLASS_SOF": "EXTENDED_CLASS_SOF", # Star galaxy classifier (0,1=star, 2,3=Galaxy) - "meas_FLAGS_GOLD_SOF_ONLY": "FLAGS_GOLD", # Measured flags (short version) - } - cols = list(column_names.keys()) - - #set up the output file columns - output = self.open_output("matched_ssi_photometry_catalog") - g = output.create_group("photometry") - for col in cols: - dtype = dtypes[col] - - if "_mag" in col: - #per band - dtype = dtype.subdtype[0] - for b in "griz": - g.create_dataset(column_names[col]+f"_{b}", (n,), dtype=dtype) - - #also create an empty array for the mag errors - #TO DO: compute mag errors from flux errors - g.create_dataset(column_names[col]+f"_err_{b}", (n,), dtype=dtype) - else: - g.create_dataset(column_names[col], (n,), dtype=dtype) - - #iterate over the input file and save to the output columns - for (s, e, data) in self.iterate_fits("balrog_matched_catalog", 1, cols, chunk_rows): - print(s,e,n) - for col in cols: - if "_mag" in col: - for iband,b in enumerate("griz"): - g[column_names[col]+f"_{b}"][s:e] = data[col][:,iband] - else: - g[column_names[col]][s:e] = data[col] - - # Set up any dummy columns with sentinal values - # that were not in the original files + "bal_id": "bal_id", # Unique identifier for object (created during balrog process) + "true_bdf_mag_deredden": "inj_mag", # Magnitude of the original deep field object, dereddened + "true_id": "inj_id", # Original coadd_obj_id of deep field object + "meas_id": "id", # Coadd_object_id of injection + "meas_ra": "ra", # measured RA of the injection + "meas_dec": "dec", # measured DEC of the injection + "meas_cm_mag_deredden": "mag", # measured magnitude of the injection + "meas_cm_max_flux_s2n": "snr", # measured S2N of the injection + "meas_cm_T": "cm_T", # measured size parameter T (x^2+y^2) + "meas_EXTENDED_CLASS_SOF": "EXTENDED_CLASS_SOF", # Star galaxy classifier (0,1=star, 2,3=Galaxy) + "meas_FLAGS_GOLD_SOF_ONLY": "FLAGS_GOLD", # Measured flags (short version) + } dummy_columns = { - "redshift_true":10.0, - } - for col_name in dummy_columns.keys(): - g.create_dataset(col_name, data=np.full(n,dummy_columns[col_name])) + "redshift_true": 10.0, + "meas_err_g":-99., + "meas_err_r":-99., + "meas_err_i":-99., + "meas_err_z":-99., + } + + self.create_photometry( + "balrog_matched_catalog", + "matched_ssi_photometry_catalog", + column_names, + dummy_columns, + ) + +class TXIngestSSIDetectionDESBalrog(TXIngestSSIDESBalrog): + """ + Class for ingesting an "SSI" "detection" catalog from DES (AKA Balrog) + """ - output.close() + name = "TXIngestSSIDetectionDESBalrog" + inputs = [ + ("balrog_detection_catalog", FitsFile), + ] + outputs = [ + ("injection_catalog", HDFFile), + ("ssi_detection_catalog", HDFFile), + ] + def run(self): + """ + We will split the Balrog "detection" catalog into two catalogs + One containing the injected catalog + The other contains info on whether a given injection has been detected + """ + ## Extract the injection catalog + column_names_inj = { + "bal_id": "bal_id", # Unique identifier for object (created during balrog process) + "true_ra": "ra", # *injected* ra + "true_dec": "dec", # *injected* dec + "true_bdf_mag_deredden": "inj_mag", # true magnitude (de-reddened) + "meas_FLAGS_GOLD_SOF_ONLY": "flags", # measured data flags + "meas_EXTENDED_CLASS_SOF": "meas_EXTENDED_CLASS_SOF", # star galaxy separator + } + + self.create_photometry( + "balrog_detection_catalog", "injection_catalog", column_names_inj, {} + ) + + # Extract the "detection" file + # We will only load a subset of columns to save space + column_names_det = { + "bal_id": "bal_id", # Unique identifier for object (created during balrog process) + "detected": "detected", # 0 or 1. Is there a match between this object in injection_catalog and the ssi_photometry_catalog (search radius of 0.5'') + "match_flag_0.5_asec": "match_flag_0.5_asec", # 0,1 or 2. Is there a match between this object in injection_catalog and the ssi_uninjected_photometry_catalog (search radius 0.5''). 1 if detection is lower flux than injection, 2 if brighter + "match_flag_0.75_asec": "match_flag_0.75_asec", # as above but with search radius of 0.75'' + "match_flag_1.0_asec": "match_flag_1.0_asec", # etc + "match_flag_1.25_asec": "match_flag_1.25_asec", + "match_flag_1.5_asec": "match_flag_1.5_asec", + "match_flag_1.75_asec": "match_flag_1.75_asec", + "match_flag_2.0_asec": "match_flag_2.0_asec", + } + + self.create_photometry( + "balrog_detection_catalog", "ssi_detection_catalog", column_names_det, {} + ) diff --git a/txpipe/map_plots.py b/txpipe/map_plots.py index ee3ac433..0c30686d 100644 --- a/txpipe/map_plots.py +++ b/txpipe/map_plots.py @@ -42,6 +42,7 @@ class TXMapPlots(PipelineStage): config_options = { # can also set moll "projection": "cart", + "rot180": False, "debug": False, } @@ -96,7 +97,7 @@ def aux_source_plots(self): for i in range(flag_max): plt.subplot(1, flag_max, i + 1) f = 2**i - m.plot(f"flags/flag_{f}", view=self.config["projection"]) + m.plot(f"flags/flag_{f}", view=self.config["projection"], rot180=self.config["rot180"]) fig.close() # PSF plots - 2 x n, for g1 and g2 @@ -104,9 +105,9 @@ def aux_source_plots(self): _, axes = plt.subplots(2, nbin_source, squeeze=False, num=fig.file.number) for i in range(nbin_source): plt.sca(axes[0, i]) - m.plot(f"psf/g1_{i}", view=self.config["projection"]) + m.plot(f"psf/g1_{i}", view=self.config["projection"], rot180=self.config["rot180"]) plt.sca(axes[1, i]) - m.plot(f"psf/g2_{i}", view=self.config["projection"]) + m.plot(f"psf/g2_{i}", view=self.config["projection"], rot180=self.config["rot180"]) fig.close() def aux_lens_plots(self): @@ -122,11 +123,11 @@ def aux_lens_plots(self): # Depth plots with self.open_output("depth_map", wrapper=True, figsize=(5, 5)) as fig: - m.plot("depth/depth", view=self.config["projection"]) + m.plot("depth/depth", view=self.config["projection"], rot180=self.config["rot180"]) # Bright objects with self.open_output("bright_object_map", wrapper=True, figsize=(5, 5)) as fig: - m.plot("bright_objects/count", view=self.config["projection"]) + m.plot("bright_objects/count", view=self.config["projection"], rot180=self.config["rot180"]) def source_plots(self): import matplotlib.pyplot as plt @@ -148,11 +149,11 @@ def source_plots(self): for i in range(nbin_source): # g1 plt.sca(axes[0, i]) - m.plot(f"g1_{i}", view=self.config["projection"], min=-0.1, max=0.1) + m.plot(f"g1_{i}", view=self.config["projection"], rot180=self.config["rot180"], min=-0.1, max=0.1) # g2 plt.sca(axes[1, i]) - m.plot(f"g2_{i}", view=self.config["projection"], min=-0.1, max=0.1) + m.plot(f"g2_{i}", view=self.config["projection"], rot180=self.config["rot180"], min=-0.1, max=0.1) fig.close() def lens_plots(self): @@ -173,9 +174,9 @@ def lens_plots(self): for i in range(nbin_lens): plt.sca(axes[0, i]) - m.plot(f"ngal_{i}", view=self.config["projection"]) + m.plot(f"ngal_{i}", view=self.config["projection"], rot180=self.config["rot180"]) plt.sca(axes[1, i]) - rho.plot(f"delta_{i}", view=self.config["projection"]) + rho.plot(f"delta_{i}", view=self.config["projection"], rot180=self.config["rot180"]) fig.close() def mask_plots(self): @@ -184,5 +185,66 @@ def mask_plots(self): m = self.open_input("mask", wrapper=True) fig = self.open_output("mask_map", wrapper=True, figsize=(5, 5)) - m.plot("mask", view=self.config["projection"]) + m.plot("mask", view=self.config["projection"], rot180=self.config["rot180"]) fig.close() + +class TXMapPlotsSSI(TXMapPlots): + """ + Make plots of all the available maps that use SSI inputs + + This makes plots of: + - depth (meas mag) + - depth (true mag) + """ + name = "TXMapPlotsSSI" + + inputs = [ + ("aux_ssi_maps", MapsFile), + ] + + outputs = [ + ("depth_ssi_meas_map", PNGFile), + ("depth_ssi_true_map", PNGFile), + ] + + def run(self): + # PSF tests + import matplotlib + + matplotlib.use("agg") + import matplotlib.pyplot as plt + + # Plot from each file separately, just + # to organize this file a bit + methods = [ + self.aux_ssi_plots, + ] + + # We don't want this to fail if some maps are missing. + for m in methods: + try: + m() + except: + if self.config["debug"]: + raise + sys.stderr.write(f"Failed to make maps with method {m.__name__}") + + + def aux_ssi_plots(self): + import matplotlib.pyplot as plt + if self.get_input("aux_ssi_maps") == "none": + with self.open_output("depth_ssi_meas_map", wrapper=True) as f: + pass + with self.open_output("depth_ssi_true_map", wrapper=True) as f: + pass + return + + m = self.open_input("aux_ssi_maps", wrapper=True) + + # Depth plots (measured magnitude) + with self.open_output("depth_ssi_meas_map", wrapper=True, figsize=(5, 5)) as fig: + m.plot("depth/depth_meas", view=self.config["projection"], rot180=self.config["rot180"]) + + # Depth plots (true magnitude) + with self.open_output("depth_ssi_true_map", wrapper=True, figsize=(5, 5)) as fig: + m.plot("depth/depth_true", view=self.config["projection"], rot180=self.config["rot180"])