diff --git a/frb/surveys/catalog_utils.py b/frb/surveys/catalog_utils.py index de99058d..9d993a04 100644 --- a/frb/surveys/catalog_utils.py +++ b/frb/surveys/catalog_utils.py @@ -160,9 +160,10 @@ def summarize_catalog(frbc, catalog, summary_radius, photom_column, magnitude): return summary_list -def xmatch_catalogs(cat1:Table, cat2:Table, skydist:units.Quantity = 5*units.arcsec, +def xmatch_catalogs(cat1:Table, cat2:Table, dist:units.Quantity = 5*units.arcsec, RACol1:str = "ra", DecCol1:str = "dec", RACol2:str = "ra", DecCol2:str = "dec", + distcol1:str = None, distcol2:str = None, return_match_idx:bool=False)->tuple: """ Cross matches two astronomical catalogs and returns @@ -171,15 +172,21 @@ def xmatch_catalogs(cat1:Table, cat2:Table, skydist:units.Quantity = 5*units.arc cat1, cat2: astropy Tables Two tables with sky coordinates to be matched. - skydist: astropy Quantity, optional + dist: astropy Quantity, optional Maximum separation for a valid match. - 5 arcsec by default. + 5 arcsec by default. Can be length units + if distcol1 and distcol2 are provided. RACol1, RACol2: str, optional Names of columns in cat1 and cat2 respectively that contain RA in degrees. DecCol1, DecCol2: str, optional Names of columns in cat1 and cat2 respectively that contain Dec in degrees. + distcol1, distcol2: str, optional + Names of columns in cat1 and cat2 + respectively that contain the radial distance + as floats in Mpc. If None, 2D cross-matches + are returned. return_match_idx: bool, optional Return the indices of the matched entries with with the distance instead? @@ -193,19 +200,31 @@ def xmatch_catalogs(cat1:Table, cat2:Table, skydist:units.Quantity = 5*units.arc assert isinstance(cat1, (Table, QTable))&isinstance(cat1, (Table, QTable)), "Catalogs must be astropy Table instances." assert (RACol1 in cat1.colnames)&(DecCol1 in cat1.colnames), " Could not find either {:s} or {:s} in cat1".format(RACol1, DecCol1) assert (RACol2 in cat2.colnames)&(DecCol2 in cat2.colnames), " Could not find either {:s} or {:s} in cat2".format(RACol2, DecCol2) + do_3d = (distcol1 is not None)&(distcol2 is not None) + if do_3d: + assert (distcol1 in cat1.colnames)&(distcol2 in cat2.colnames), "Could not find either {:s} or {:s} in cat1".format(distcol1, distcol2) + dist1 = cat1[distcol1]*units.Mpc + dist2 = cat2[distcol2]*units.Mpc + else: + dist1 = None + dist2 = None # Get corodinates - cat1_coord = SkyCoord(cat1[RACol1], cat1[DecCol1], unit = "deg") - cat2_coord = SkyCoord(cat2[RACol2], cat2[DecCol2], unit = "deg") + cat1_coord = SkyCoord(cat1[RACol1]*units.deg, cat1[DecCol1]*units.deg, distance=dist1) + cat2_coord = SkyCoord(cat2[RACol2]*units.deg, cat2[DecCol2]*units.deg, distance=dist2) # Match 2D - idx, d2d, _ = cat1_coord.match_to_catalog_sky(cat2_coord) + idx, d2d, d3d = cat1_coord.match_to_catalog_sky(cat2_coord) # Get matched tables - match1 = cat1[d2d < skydist] - match2 = cat2[idx[d2d < skydist]] + if do_3d: + separation = d3d + else: + separation = d2d + match1 = cat1[separation < dist] + match2 = cat2[idx[separation < dist]] if return_match_idx: - return idx, d2d + return idx, d2d, d3d else: return match1, match2 diff --git a/frb/surveys/cluster_search.py b/frb/surveys/cluster_search.py index 1b4738c4..2402ac1f 100644 --- a/frb/surveys/cluster_search.py +++ b/frb/surveys/cluster_search.py @@ -8,6 +8,7 @@ from frb.defs import frb_cosmo from astropy.coordinates import SkyCoord from astropy import units as u +from astropy.table import Table try: from astroquery.vizier import Vizier @@ -22,11 +23,11 @@ class VizierCatalogSearch(surveycoord.SurveyCoord): """ - def __init__(self, coord, radius = 90*u.deg, catalog = None, cosmo=None, **kwargs): + def __init__(self, coord, radius = 90*u.deg, survey=None, viziercatalog = None, cosmo=None, **kwargs): # Initialize a SurveyCoord object - surveycoord.SurveyCoord.__init__(self, coord, radius, **kwargs) - self.survey = None # Name - self.catalog = catalog # Name of the Vizier table to draw from. + super(VizierCatalogSearch, self).__init__(coord, radius, **kwargs) + self.survey = survey # Name + self.viziercatalog = viziercatalog # Name of the Vizier table to draw from. self.coord = coord # Location around which to perform the search self.radius = radius.to('deg').value # Radius of cone search if cosmo is None: # Use the same cosmology as elsewhere in this repository unless specified. @@ -34,13 +35,46 @@ def __init__(self, coord, radius = 90*u.deg, catalog = None, cosmo=None, **kwarg else: self.cosmo = cosmo - - def get_catalog(self, query_fields=None, transverse_distance_cut = 5*u.Mpc, **kwargs): + def clean_catalog(self, catalog): + """ + This will be survey specific. + """ pass - + def _transverse_distance_cut(self, catalog, transverse_distance_cut, distance_column='Dist'): + # Apply a transverse distance cut + angular_dist = self.coord.separation(SkyCoord(catalog['ra'], catalog['dec'], unit='deg')).to('rad').value + transverse_dist = catalog[distance_column]*np.sin(angular_dist) + catalog = catalog[transverse_dist 0: + result = self.clean_catalog(result) + + # Apply a transverse distance cut + if transverse_distance_cut=richness_cut] + self.catalog = result + + return self.catalog + +# Wen+2024 +class WenGroupCat(VizierCatalogSearch): + """ + A class to query sources within the Wen+2024 + group/cluster catalog. + """ + + + def __init__(self, coord, radius = 0.2*u.deg, cosmo=None, **kwargs): + # Initialize a SurveyCoord object + super(WenGroupCat, self).__init__(coord, radius, + survey="Wen+2024", + viziercatalog="J/ApJS/272/39/table2", + cosmo=cosmo, **kwargs) - # Convert distances from h^-1 Mpc to Mpc based on the cosmology being used. - result['Dist'] /=self.cosmo.h + def clean_catalog(self, catalog): + try: + catalog.rename_columns(['RAJ2000'],['ra']) + except KeyError: + assert 'ra' in catalog.keys() + try: + catalog.rename_columns(['DEJ2000'],['dec']) + except KeyError: + assert 'dec' in catalog.keys() + try: + catalog.rename_columns(['zCl'],['z']) + except KeyError: + assert 'z' in catalog.keys() + + + # Add a distance estimate in Mpc using the given cosmology + catalog['Dist'] = self.cosmo.angular_diameter_distance(catalog['z']).value + + return catalog + + def get_catalog(self, query_fields=None, + transverse_distance_cut = np.inf*u.Mpc, richness_cut = 5): + """ + Get the catalog of objects + Args: + z_lim (float): The maximum redshift of the objects to include in the catalog. + transverse_distance_cut (Quantity): The maximum impact parameter of the objects to include in the catalog. + richness_cut (int): The minimum number of members in any group/cluster returned. + query_fields (list): The fields to include in the catalog. If None, all fields are used. + Returns: + A table of objects within the given limits. + """ + result = super(WenGroupCat, self)._get_catalog(query_fields=query_fields) + if len(result) > 0: + result = self.clean_catalog(result) + + # Apply a transverse distance cut + if transverse_distance_cut=richness_cut] + self.catalog = result + return self.catalog + +# Bahk and Hwang 2024 (Updated Planck+2015) +class UPClusterSZCat(VizierCatalogSearch): + """ + A class to query sources within the Bahk and Hwang 2024 + group/cluster catalog. + """ + + + def __init__(self, coord, radius = 90*u.deg, cosmo=None, **kwargs): + # Initialize a SurveyCoord object + super(UPClusterSZCat, self).__init__(coord, radius, + survey="UPClusterSZ", + viziercatalog="J/ApJS/272/7/table2", + cosmo=cosmo, **kwargs) + + def clean_catalog(self, catalog): + if len(catalog) > 0: + try: + catalog.rename_columns(['RAJ2000', 'DEJ2000'], ['ra', 'dec']) # Rename the columns to match the SurveyCoord class + except KeyError: + print(catalog.keys()) + + assert 'ra' in catalog.keys() and 'dec' in catalog.keys() + + # Add a distance estimate in Mpc using the given cosmology + catalog['Dist'] = self.cosmo.angular_diameter_distance(catalog['z']).value + + return catalog + + def get_catalog(self, query_fields=None, + transverse_distance_cut = np.inf*u.Mpc): + """ + Get the catalog of objects + Args: + z_lim (float): The maximum redshift of the objects to include in the catalog. + transverse_distance_cut (Quantity): The maximum impact parameter of the objects to include in the catalog. + richness_cut (int): The minimum number of members in any group/cluster returned. + query_fields (list): The fields to include in the catalog. If None, all fields are used. + Returns: + A table of objects within the given limits. + """ + result = super(UPClusterSZCat, self)._get_catalog(query_fields=query_fields) + if len(result) > 0: + result = self.clean_catalog(result) + + # Apply a transverse distance cut + if transverse_distance_cut 0: + catalog.rename_columns(['RAJ2000', 'DEJ2000'], ['ra', 'dec']) # Rename the columns to match the SurveyCoord class + + # Add a distance estimate in Mpc using the given cosmology + catalog['Dist'] = self.cosmo.angular_diameter_distance(catalog['z']).value + + return catalog + + def get_catalog(self, query_fields=None, + transverse_distance_cut = np.inf*u.Mpc): + """ + Get the catalog of objects + Args: + z_lim (float): The maximum redshift of the objects to include in the catalog. + transverse_distance_cut (Quantity): The maximum impact parameter of the objects to include in the catalog. + richness_cut (int): The minimum number of members in any group/cluster returned. + query_fields (list): The fields to include in the catalog. If None, all fields are used. + Returns: + A table of objects within the given limits. + """ + result = super(ROSATXClusterCat, self)._get_catalog(query_fields=query_fields) + if len(result) > 0: + result = self.clean_catalog(result) + + # Apply a transverse distance cut + if transverse_distance_cut 0: + catalog.rename_columns(['RAJ2000', 'DEJ2000', 'zcmb'], ['ra', 'dec', 'z']) # Rename the columns to match the SurveyCoord class + + # Add a distance estimate in Mpc using the given cosmology + catalog['Dist'] = self.cosmo.angular_diameter_distance(catalog['z']).to('Mpc').value + + return catalog + + def get_catalog(self, query_fields=None, + transverse_distance_cut = np.inf*u.Mpc): + """ + Get the catalog of objects + Args: + z_lim (float): The maximum redshift of the objects to include in the catalog. + transverse_distance_cut (Quantity): The maximum impact parameter of the objects to include in the catalog. + richness_cut (int): The minimum number of members in any group/cluster returned. + query_fields (list): The fields to include in the catalog. If None, all fields are used. + Returns: + A table of objects within the given limits. + """ + result = super(TempelClusterCat, self)._get_catalog(query_fields=query_fields) + if len(result) > 0: + result = self.clean_catalog(result) + + # Apply a transverse distance cut + if transverse_distance_cut 0: + catalog.rename_columns(['RAJ2000', 'DEJ2000','zsp1'], ['ra', 'dec','z']) # Rename the columns to match the SurveyCoord class + + # Add a distance estimate in Mpc using the given cosmology + catalog['Dist'] = self.cosmo.angular_diameter_distance(catalog['z']).value + + return catalog + + def get_catalog(self, query_fields=None, + transverse_distance_cut = np.inf*u.Mpc): + """ + Get the catalog of objects + Args: + z_lim (float): The maximum redshift of the objects to include in the catalog. + transverse_distance_cut (Quantity): The maximum impact parameter of the objects to include in the catalog. + richness_cut (int): The minimum number of members in any group/cluster returned. + query_fields (list): The fields to include in the catalog. If None, all fields are used. + Returns: + A table of objects within the given limits. + """ + result = super(RASSClusterCat, self)._get_catalog(query_fields=query_fields) + if len(result) > 0: + result = self.clean_catalog(result) + + # Apply a transverse distance cut + if transverse_distance_cut 0: + catalog.rename_columns(['RAJ2000', 'DEJ2000','zspec'], ['ra', 'dec','z']) # Rename the columns to match the SurveyCoord class + + # Add a distance estimate in Mpc using the given cosmology + redshift = catalog['z'] + redshift[redshift<0] = catalog['zlambda'][redshift<0] + catalog['Dist'] = self.cosmo.angular_diameter_distance(redshift).value + + return catalog + + def get_catalog(self, query_fields=None, + transverse_distance_cut = np.inf*u.Mpc): + """ + Get the catalog of objects + Args: + z_lim (float): The maximum redshift of the objects to include in the catalog. + transverse_distance_cut (Quantity): The maximum impact parameter of the objects to include in the catalog. + richness_cut (int): The minimum number of members in any group/cluster returned. + query_fields (list): The fields to include in the catalog. If None, all fields are used. + Returns: + A table of objects within the given limits. + """ + result = super(RedMapperClusterCat, self)._get_catalog(query_fields=query_fields) + if len(result) > 0: + result = self.clean_catalog(result) # Apply a transverse distance cut - angular_dist = self.coord.separation(SkyCoord(result['ra'], result['dec'], unit='deg')).to('rad').value - transverse_dist = result['Dist']*np.sin(angular_dist) - result = result[transverse_dist=richness_cut] + if transverse_distance_cut 0: + catalog.rename_columns(['RAJ2000', 'DEJ2000', 'z1C'], ['ra', 'dec','z']) # Rename the columns to match the SurveyCoord class + + # Add a distance estimate in Mpc using the given cosmology + redshift = catalog['z'] + catalog['Dist'] = self.cosmo.angular_diameter_distance(redshift).value + + return catalog + + def get_catalog(self, query_fields=None, + transverse_distance_cut = np.inf*u.Mpc): + """ + Get the catalog of objects + Args: + z_lim (float): The maximum redshift of the objects to include in the catalog. + transverse_distance_cut (Quantity): The maximum impact parameter of the objects to include in the catalog. + richness_cut (int): The minimum number of members in any group/cluster returned. + query_fields (list): The fields to include in the catalog. If None, all fields are used. + Returns: + A table of objects within the given limits. + """ + result = super(ACTDR5ClusterCat, self)._get_catalog(query_fields=query_fields) + if len(result) > 0: + result = self.clean_catalog(result) + + # Apply a transverse distance cut + if transverse_distance_cut 0: + catalog.rename_columns(['RAJ2000', 'DEJ2000', 'Bestz'], ['ra', 'dec','z']) # Rename the columns to match the SurveyCoord class + + # Add a distance estimate in Mpc using the given cosmology + redshift = catalog['z'] + catalog['Dist'] = self.cosmo.angular_diameter_distance(redshift).value + + return catalog + + def get_catalog(self, query_fields=None, + transverse_distance_cut = np.inf*u.Mpc): + """ + Get the catalog of objects + Args: + z_lim (float): The maximum redshift of the objects to include in the catalog. + transverse_distance_cut (Quantity): The maximum impact parameter of the objects to include in the catalog. + richness_cut (int): The minimum number of members in any group/cluster returned. + query_fields (list): The fields to include in the catalog. If None, all fields are used. + Returns: + A table of objects within the given limits. + """ + result = super(ERASSClusterCat, self)._get_catalog(query_fields=query_fields) + if len(result) > 0: + result = self.clean_catalog(result) + + # Apply a transverse distance cut + if transverse_distance_cut