Load extragalactic, galactic and CMB maps

Sky emission

load_sky_emission[source]

load_sky_emission(config_sky_emission, site, channel)

Load foreground maps for a channel

Parameters
----------
config_sky_emission : dict
    CMB-S4 sky emission configuration,
    generally config['sky_emission'] for a configuration
    loaded from a TOML file
    dictionary with standard emission names and their weights
site : str
    'Pole' or 'Chile', case doesn't matter
channel : str
    Channel tag, e.g. 'MFHS1'

Returns
-------
output_map : numpy array
    Output map with all emissions combined, uses nan for missing pixels
def load_sky_emission(config_sky_emission, site, channel):
    """Load foreground maps for a channel

    Parameters
    ----------
    config_sky_emission : dict
        CMB-S4 sky emission configuration,
        generally config['sky_emission'] for a configuration
        loaded from a TOML file
        dictionary with standard emission names and their weights
    site : str
        'Pole' or 'Chile', case doesn't matter
    channel : str
        Channel tag, e.g. 'MFHS1'

    Returns
    -------
    output_map : numpy array
        Output map with all emissions combined, uses nan for missing pixels
    """
    log.info("Configuration %s", str(config_sky_emission))
    telescope = get_telescope(channel)
    nside = 512 if telescope == "SAT" else 4096
    npix = hp.nside2npix(nside)
    output_map = np.zeros((3, npix), dtype=np.float32)
    realization = 0  # foregrounds are deterministic
    channel_noP = channel.replace("P", "")
    for emission, weight in config_sky_emission.items():
        if weight == 0:
            log.info("Skip %s", emission)
            continue
        log.info("Processing %s", emission)

        emission_map = hp.read_map(
            Path(base_folder)
            / f"{realization:08d}"
            / f"{site.lower()}_{emission_naming[emission]}_{telescope}_{channel_noP}_{mapmaking_naming[telescope]}",
            (0, 1, 2),
            dtype=None,
            verbose=False,
        )
        emission_map[emission_map == 0] = np.nan
        emission_map[emission_map == hp.UNSEEN] = np.nan

        if emission == "CMB_tensor_to_scalar_ratio":
            # input maps are simulated with r=3e-3
            # tensor-to-scalar ratio is defined in the power spectrum, we need to weight the maps with the square root
            weight = np.sqrt(weight / 3e-3)
        output_map += emission_map * weight

    output_map /= 1e6  # uK -> K
    return output_map
import toml

config = toml.load("s4_design.toml")
s4 = read_instrument_model()
s4
QTable length=21
bandtelescopecenter_frequencyfwhmdetectors_per_tubeNETchile_temp_sensitivity_DSRchile_pol_sensitivity_DSRpole_pol_sensitivity_DSR
GHzarcmins(1/2) uKarcmin uKarcmin uKarcmin uK
str6str3float64float64int64float64float64float64float64
ULFPL1LAT20.011.454332.20.00.08.4
LFL1LAT25.757.496275.021.830.85.0
LFPL1LAT25.758.496278.621.830.85.0
LFL2LAT38.755.196215.012.417.64.5
LFPL2LAT38.755.896268.612.417.64.5
MFPL1LAT91.52.5864285.22.02.90.68
MFL1LAT91.52.2864251.02.02.90.68
MFL2LAT148.51.4864280.02.02.80.96
MFPL2LAT148.51.6864264.62.02.80.96
HFL1LAT227.01.0864625.06.99.85.7
HFPL1LAT227.01.1864534.96.99.85.7
HFL2LAT285.50.98641528.016.723.69.8
HFPL2LAT285.51.08641228.316.723.69.8
LFS1SAT25.7572.8288169.20.00.03.5
LFS2SAT38.7572.8288204.30.00.04.5
MFLS1SAT85.025.53528290.00.00.00.88
MFHS1SAT95.022.73528248.20.00.00.78
MFLS2SAT145.025.53528280.40.00.01.2
MFHS2SAT155.122.73528297.40.00.01.3
HFS1SAT227.013.08592609.10.00.03.5
HFS2SAT285.513.085921434.90.00.06.0
log.basicConfig(level=log.INFO)
ch = "LFS1"
site = "Pole"
output_map = load_sky_emission(config["sky_emission"], site, ch)
INFO:root:Configuration {'foreground_emission': 1, 'CMB_unlensed': 1, 'CMB_lensing_signal': 1, 'CMB_tensor_to_scalar_ratio': 0.005}
INFO:root:Processing foreground_emission
INFO:root:Processing CMB_unlensed
INFO:root:Processing CMB_lensing_signal
INFO:root:Processing CMB_tensor_to_scalar_ratio
assert len(output_map[output_map == hp.UNSEEN]) == 0, \
    "We should be using nan not UNSEEN to simplify summing"
np.nanmin(output_map), np.nanmax(output_map)
(-0.0015811161, 0.0018167086)
assert np.nanmax(output_map) > 1e-4 and np.nanmax(output_map) < 1e-2, "Rough unit check"
%matplotlib inline
hp.mollview(output_map[0], min=-1e-4, max=1e-4, unit="K", title="Foregrounds I")
/global/homes/z/zonca/condajupynersc/lib/python3.7/site-packages/healpy/projaxes.py:920: MatplotlibDeprecationWarning: You are modifying the state of a globally registered colormap. In future versions, you will not be able to modify a registered colormap in-place. To remove this warning, you can make a copy of the colormap first. cmap = copy.copy(mpl.cm.get_cmap("viridis"))
  newcm.set_over(newcm(1.0))
/global/homes/z/zonca/condajupynersc/lib/python3.7/site-packages/healpy/projaxes.py:921: MatplotlibDeprecationWarning: You are modifying the state of a globally registered colormap. In future versions, you will not be able to modify a registered colormap in-place. To remove this warning, you can make a copy of the colormap first. cmap = copy.copy(mpl.cm.get_cmap("viridis"))
  newcm.set_under(bgcolor)
/global/homes/z/zonca/condajupynersc/lib/python3.7/site-packages/healpy/projaxes.py:922: MatplotlibDeprecationWarning: You are modifying the state of a globally registered colormap. In future versions, you will not be able to modify a registered colormap in-place. To remove this warning, you can make a copy of the colormap first. cmap = copy.copy(mpl.cm.get_cmap("viridis"))
  newcm.set_bad(badcolor)
/global/homes/z/zonca/condajupynersc/lib/python3.7/site-packages/healpy/projaxes.py:211: MatplotlibDeprecationWarning: Passing parameters norm and vmin/vmax simultaneously is deprecated since 3.3 and will become an error two minor releases later. Please pass vmin/vmax directly to the norm when creating it.
  **kwds
hp.mollview(output_map[1], min=-1e-5, max=1e-5, unit="K", title="Foregrounds Q")
hp.mollview(output_map[2], min=-1e-5, max=1e-5, unit="K", title="Foregrounds U")