如何用Python（如healpy）将天文源目录与fits掩码交叉匹配？

核心思路 #

把星系的赤道坐标转成Healpy适配的球极坐标，匹配掩码对应的像素，筛选出位于有效掩码区域内的源。

步骤实现 #

1. 安装依赖库
   先装好所需工具：

   pip install healpy astropy pandas numpy
   

2. 加载数据

   • 读取星系目录DataFrame：

     import pandas as pd
     df = pd.read_csv("galaxy_catalog.csv")  # 替换成你的目录路径/格式
     

   • 读取FITS二进制掩码：

     import healpy as hp
     # verbose=False关闭加载日志
     mask = hp.read_map("mask.fits", verbose=False)
     # 把掩码里的NaN转成0（标记为无效区域）
     mask = np.nan_to_num(mask, nan=0)
     

3. 坐标转换
   Healpy默认用球极坐标（θ=余纬，φ=赤经，单位都是弧度），得把星系的RA/Dec转成对应格式：

   import numpy as np
   from astropy.coordinates import SkyCoord
   from astropy import units as u
   
   # 从DataFrame提取赤经(RA)、赤纬(Dec)
   ra = df["ra"].values
   dec = df["dec"].values
   
   # 转成天球坐标对象
   coords = SkyCoord(ra=ra*u.degree, dec=dec*u.degree)
   # 计算余纬θ（90°-Dec）和赤经φ，转成弧度
   theta = np.radians(90 - coords.dec.degree)
   phi = np.radians(coords.ra.degree)
   
   # 如果掩码用的是银河坐标，就转成银河系坐标
   # coords_gal = coords.galactic
   # theta = np.radians(90 - coords_gal.b.degree)
   # phi = np.radians(coords_gal.l.degree)
   

4. 匹配掩码并筛选

   • 获取掩码的像素分辨率参数nside：

     nside = hp.get_nside(mask)
     

   • 计算每个星系对应的掩码像素索引：

     pix_indices = hp.ang2pix(nside, theta, phi)
     

   • 提取对应像素的掩码值，筛选出在有效区域（掩码值为1）的星系：

     mask_values = mask[pix_indices]
     filtered_df = df[mask_values == 1]
     

注意事项 #

• 必须保证星系目录的坐标系统和掩码的坐标系统一致（赤道/银河），不然匹配会出错。
• 如果掩码里有非0/1的数值，得先归一化处理成二进制（1=有效，0=无效）。

内容的提问来源于stack exchange，提问作者NeStack