<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Color-color plots for LRG targets
The goal of this notebook is to compare the colors of LRG targets against those of the LRG templates.
Step2:  Read the targets catalog
Step3:  Read the templates and compute colors on a redshift grid.
Step5:  Generate some plots
<ASSISTANT_TASK:>
Python Code:
import os
from time import time
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.ticker as ticker
import fitsio
import seaborn as sns
from speclite import filters
from desitarget import desi_mask
from desisim.io import read_basis_templates
%pylab inline
sns.set(style='white', font_scale=1.5, font='sans-serif', palette='Set2')
setcolors = sns.color_palette()
Explanation: Color-color plots for LRG targets
The goal of this notebook is to compare the colors of LRG targets against those of the LRG templates.
End of explanation
def flux2colors(cat):
    Convert DECam/WISE fluxes to magnitudes and colors.
    colors = dict()
    with warnings.catch_warnings(): # ignore missing fluxes (e.g., for QSOs)
        warnings.simplefilter('ignore')
        for ii, band in zip((1, 2, 4), ('g', 'r', 'z')):
            colors[band] = 22.5 - 2.5 * np.log10(cat['DECAM_FLUX'][..., ii].data)
        for ii, band in zip((0, 1), ('W1', 'W2')):
            colors[band] = 22.5 - 2.5 * np.log10(cat['WISE_FLUX'][..., ii].data)
    colors['gr'] = colors['g'] - colors['r']
    colors['rz'] = colors['r'] - colors['z']
    colors['rW1'] = colors['r'] - colors['W1']
    colors['W1W2'] = colors['W1'] - colors['W2']
    return colors
lrgfile = os.path.join( os.getenv('DESI_ROOT'), 'data', 'targets-dr3.1-EisDawLRG.fits' )
# Select just LRG targets.
print('Reading {}'.format(lrgfile))
cat = fitsio.read(lrgfile, ext=1, upper=True, columns=['DESI_TARGET'])
these = np.where( (cat['DESI_TARGET'] & desi_mask.LRG) != 0 )[0]
print('Number of LRG targets = {}'.format(len(these)))
cat = fitsio.read(lrgfile, ext=1, upper=True, rows=these)
data = flux2colors(cat)
Explanation: Read the targets catalog
End of explanation
filt = filters.load_filters('decam2014-g', 'decam2014-r', 'decam2014-z', 'wise2010-W1')
flux, wave, meta = read_basis_templates(objtype='LRG')
nt = len(meta)
print('Number of templates = {}'.format(nt))
zmin, zmax, dz = 0.0, 2.0, 0.1
nz = np.round( (zmax - zmin) / dz ).astype('i2')
print('Number of redshift points = {}'.format(nz))
cc = dict(
    redshift = np.linspace(0.0, 2.0, nz),
    gr = np.zeros( (nt, nz) ),
    rz = np.zeros( (nt, nz) ),
    rW1 = np.zeros( (nt, nz), )
)
t0 = time()
for iz, red in enumerate(cc['redshift']):
    zwave = wave.astype('float') * (1 + red)
    phot = filt.get_ab_maggies(flux, zwave, mask_invalid=False)
    cc['gr'][:, iz] = -2.5 * np.log10( phot['decam2014-g'] / phot['decam2014-r'] )
    cc['rz'][:, iz] = -2.5 * np.log10( phot['decam2014-r'] / phot['decam2014-z'] )
    cc['rW1'][:, iz] = -2.5 * np.log10( phot['decam2014-r'] / phot['wise2010-W1'] )
print('Total time = {:.2f} sec.'.format(time() - t0))
Explanation: Read the templates and compute colors on a redshift grid.
End of explanation
figsize = (8, 6)
grrange = (0.0, 3.0)
rzrange = (0.0, 2.5)
rW1range = (-1, 5)
mzrange = (17.5, 20.5)
ntspace = 5 # spacing between model curves
def rzz(pngfile=None):
    r-z vs apparent magnitude z
    fig, ax = plt.subplots(figsize=figsize)
    hb = ax.hexbin(data['z'], data['rz'], bins='log', cmap='Blues_r', 
                   mincnt=100, extent=mzrange+rzrange)
    ax.set_xlabel('z')
    ax.set_ylabel('r - z')
    ax.set_xlim(mzrange)
    ax.set_ylim(rzrange)
    cb = fig.colorbar(hb, ax=ax)
    cb.set_label(r'log$_{10}$ (Number of Galaxies per Bin)')
    if pngfile:
        fig.savefig(pngfile)
def grz(models=False, pngfile=None):
    fig, ax = plt.subplots(figsize=figsize)
    hb = ax.hexbin(data['rz'], data['gr'], bins='log', cmap='Blues_r', 
                   mincnt=100, extent=rzrange+grrange)
    ax.set_xlabel('r - z')
    ax.set_ylabel('g - r')
    ax.set_xlim(rzrange)
    ax.set_ylim(grrange)
    cb = fig.colorbar(hb, ax=ax)
    cb.set_label(r'log$_{10}$ (Number of Galaxies per Bin)')
    if models:
        for tt in np.arange(0, nt, ntspace):
            ax.scatter(cc['rz'][tt, 0], cc['gr'][tt, 0], marker='o', 
                       facecolors='none', s=80, edgecolors='k',
                       linewidth=1)
            ax.plot(cc['rz'][tt, :], cc['gr'][tt, :], marker='s', 
                    markersize=5, ls='-', alpha=0.5)
        ax.text(0.1, 0.05, 'z=0', ha='left', va='bottom',
                transform=ax.transAxes, fontsize=14)
    if pngfile:
        fig.savefig(pngfile)
def rzW1(models=False, pngfile=None):
    fig, ax = plt.subplots(figsize=figsize)
    hb = ax.hexbin(data['rz'], data['rW1'], bins='log', cmap='Blues_r', 
                   mincnt=100, extent=rzrange+grrange)
    ax.set_xlabel('r - z')
    ax.set_ylabel('r - W1')
    ax.set_xlim(rzrange)
    ax.set_ylim(rW1range)
    cb = fig.colorbar(hb, ax=ax)
    cb.set_label(r'log$_{10}$ (Number of Galaxies per Bin)')
    if models:
        for tt in np.arange(0, nt, ntspace):
            ax.scatter(cc['rz'][tt, 0], cc['rW1'][tt, 0], marker='o', 
                       facecolors='none', s=80, edgecolors='k',
                       linewidth=1)
            ax.plot(cc['rz'][tt, :], cc['rW1'][tt, :], marker='s', 
                    markersize=5, ls='-', alpha=0.5)
        ax.text(0.1, 0.05, 'z=0', ha='left', va='bottom',
                transform=ax.transAxes, fontsize=14)
    if pngfile:
        fig.savefig(pngfile)
grz(models=True)
rzW1(models=True)
rzz()
Explanation: Generate some plots
End of explanation
<END_TASK>