Gallery for the Overview Paper

The purpose of this notebook is to build a nice gallery of object images for the overview paper.

For future reference: The notebook must be run from https://jupyter-dev.nersc.gov with the following (approximate) activation script:

#!/bin/bash                                                                                                           
version=$1                                                                                                            
connection_file=$2                                                                                                    

desiconda_version=20170818-1.1.12-img                                                                                 
module use /global/common/${NERSC_HOST}/contrib/desi/desiconda/$desiconda_version/modulefiles                         
module load desiconda                                                                                                 

export LEGACYPIPE_DIR=$SCRATCH/repos/legacypipe                                                                       

export PATH=$LEGACYPIPE_DIR/bin:${PATH}                                                                               
export PATH=$SCRATCH//repos/build/bin:$PATH                                                                           
export PYTHONPATH=$LEGACYPIPE_DIR/py:${PYTHONPATH}                                                                    
export PYTHONPATH=$SCRATCH/repos/build/lib/python3.5/site-packages:$PYTHONPATH                                        

module use $LEGACYPIPE_DIR/bin/modulefiles/cori                                                                       
module load dust                                                                                                      

exec python -m ipykernel -f $connection_file

Some neat objects:

Imports and paths


In [1]:
import os, sys
import shutil, time, warnings
from contextlib import redirect_stdout
import numpy as np
import matplotlib.pyplot as plt

In [2]:
from astropy.table import Table, vstack
from PIL import Image, ImageDraw, ImageFont

In [3]:
import multiprocessing
nproc = multiprocessing.cpu_count() // 2

In [4]:
%matplotlib inline

Preliminaries

Define the data release and the various output directories.


In [17]:
PIXSCALE = 0.262

In [5]:
figdir = '/global/project/projectdirs/desi/users/ioannis/legacysurveys/overview-paper'
figfile = os.path.join(figdir, 'gallery.fits')

In [7]:
jpgdir = os.path.join(figdir, 'jpg')
if not os.path.isdir(jpgdir):
    os.mkdir(jpgdir)

In [8]:
pngdir = os.path.join(figdir, 'png')
if not os.path.isdir(pngdir):
    os.mkdir(pngdir)

Build a sample with the objects of interest.


In [35]:
cat = Table()
cat['name'] = (
    'NGC6742',
    'M92', 
    'Bow-Shock', 
    'NGC2782',
    'UGC10321', 
    'C4-2010'
)
cat['nicename'] = (
    'NGC 6742 Planetary Nebula',
    'Messier 92 Globular Cluster', 
    'Interstellar Bow Shock', 
    'NGC 2782',
    'UGC 10321 Galaxy Group', 
    'SDSS/C4 Galaxy Cluster 2010'
)
cat['viewer'] = (
    'http://legacysurvey.org/viewer/?layer=decals-dr6&ra=284.83291667&dec=48.46527778',
    'http://legacysurvey.org/viewer/?layer=decals-dr6&ra=259.28029167&dec=43.13652778&zoom=12',
    'http://legacysurvey.org/viewer?ra=325.6872&dec=1.0032&zoom=14&layer=decals-dr5',         
    'http://legacysurvey.org/viewer/?layer=decals-dr6&ra=138.52129167&dec=40.11369444&zoom=12',
    'http://legacysurvey.org/viewer?ra=244.5280&dec=21.5591&zoom=14&layer=decals-dr5',
    'http://legacysurvey.org/viewer?ra=29.0707&dec=1.0510&zoom=13&layer=decals-dr5'
)
cat['dr'] = (
    'dr6',
    'dr6', 
    'dr7',
    'dr6',
    'dr7', 
    'dr7'
)
cat['ra'] = (
    284.83291667,
    259.28029167, 
    325.6872, 
    138.52129167,
    244.5280, 
    29.070641492
)
cat['dec'] = (
    48.46527778,
    43.13652778, 
    1.0032, 
    40.11369444,
    21.5591, 
    1.050816667
)
cat['diam'] = np.array([
    1.5,
    20, 
    4, 
    7,
    3, 
    5 
]).astype('f4') # [arcmin]
cat


Out[35]:
<Table length=6>
namenicenameviewerdrradecdiam
str9str27str88str3float64float64float32
NGC6742NGC 6742 Planetary Nebulahttp://legacysurvey.org/viewer/?layer=decals-dr6&ra=284.83291667&dec=48.46527778dr6284.8329166748.465277781.5
M92Messier 92 Globular Clusterhttp://legacysurvey.org/viewer/?layer=decals-dr6&ra=259.28029167&dec=43.13652778&zoom=12dr6259.2802916743.1365277820.0
Bow-ShockInterstellar Bow Shockhttp://legacysurvey.org/viewer?ra=325.6872&dec=1.0032&zoom=14&layer=decals-dr5dr7325.68721.00324.0
NGC2782NGC 2782http://legacysurvey.org/viewer/?layer=decals-dr6&ra=138.52129167&dec=40.11369444&zoom=12dr6138.5212916740.113694447.0
UGC10321UGC 10321 Galaxy Grouphttp://legacysurvey.org/viewer?ra=244.5280&dec=21.5591&zoom=14&layer=decals-dr5dr7244.52821.55913.0
C4-2010SDSS/C4 Galaxy Cluster 2010http://legacysurvey.org/viewer?ra=29.0707&dec=1.0510&zoom=13&layer=decals-dr5dr729.0706414921.0508166675.0

Some rejected objects.


In [36]:
toss = Table()
toss['name'] = (
    'Abell383', 
    'NGC2874'
)
toss['nicename'] = (
    'Abell 383', 
    'NGC2874 Galaxy Group'
)
toss['viewer'] = (
    'http://legacysurvey.org/viewer?ra=42.0141&dec=-3.5291&zoom=15&layer=decals-dr5',
    'http://legacysurvey.org/viewer?ra=141.4373&dec=11.4284&zoom=13&layer=decals-dr5'
)
toss['dr'] = (
    'dr5', # Abell 383
    'dr5'  # C4 cluster
)
toss['ra'] = (
    42.0141, 
    141.44215000
)
toss['dec'] = (
    -3.5291, 
    11.43696000
)
toss['diam'] = np.array([
    6, 
    6
]).astype('f4') # [arcmin]
toss


Out[36]:
<Table length=2>
namenicenameviewerdrradecdiam
str8str20str79str3float64float64float32
Abell383Abell 383http://legacysurvey.org/viewer?ra=42.0141&dec=-3.5291&zoom=15&layer=decals-dr5dr542.0141-3.52916.0
NGC2874NGC2874 Galaxy Grouphttp://legacysurvey.org/viewer?ra=141.4373&dec=11.4284&zoom=13&layer=decals-dr5dr5141.4421511.436966.0

Ensure all objects are in the DR6+DR7 footprint before building coadds.


In [37]:
def init_survey(dr='dr7'):
    from legacypipe.survey import LegacySurveyData
    
    if dr == 'dr7':
        survey = LegacySurveyData(
            survey_dir='/global/project/projectdirs/cosmo/work/legacysurvey/dr7',
            output_dir=figdir)
    else:
        survey = LegacySurveyData(
            survey_dir='/global/project/projectdirs/cosmo/work/legacysurvey/dr6',
            output_dir=figdir)
    
    return survey

In [38]:
def simple_wcs(obj):
    """Build a simple WCS object for a single object."""
    from astrometry.util.util import Tan
    
    size = np.rint(obj['diam'] * 60 / PIXSCALE).astype('int') # [pixels]
    wcs = Tan(obj['ra'], obj['dec'], size/2+0.5, size/2+0.5,
                 -PIXSCALE/3600.0, 0.0, 0.0, PIXSCALE/3600.0, 
                 float(size), float(size))
    return wcs

In [39]:
def _build_sample_one(args):
    """Wrapper function for the multiprocessing."""
    return build_sample_one(*args)

In [40]:
def build_sample_one(obj, verbose=False):
    """Wrapper function to find overlapping grz CCDs for a given object.
    
    """
    survey = init_survey(dr=obj['dr'])
    
    print('Working on {}...'.format(obj['name']))
    wcs = simple_wcs(obj)
    try:
        ccds = survey.ccds_touching_wcs(wcs) # , ccdrad=2*diam/3600)
    except:
        return None
    
    if ccds:
        # Is there 3-band coverage?
        if 'g' in ccds.filter and 'r' in ccds.filter and 'z' in ccds.filter:
            if verbose:
                print('For {} found {} CCDs, RA = {:.5f}, Dec = {:.5f}, Diameter={:.4f} arcmin'.format(
                        obj['name'], len(ccds), obj['ra'], obj['dec'], obj['diam']))
            return obj
    return None

In [41]:
def build_sample(cat, use_nproc=nproc):
    """Build the full sample with grz coverage in DR6."""

    sampleargs = list()
    for cc in cat:
        sampleargs.append( (cc, True) ) # the False refers to verbose=False

    if use_nproc > 1:
        p = multiprocessing.Pool(nproc)
        result = p.map(_build_sample_one, sampleargs)
        p.close()
    else:
        result = list()
        for args in sampleargs:
            result.append(_build_sample_one(args))

    # Remove non-matching objects and write out the sample
    outcat = vstack(list(filter(None, result)))
    print('Found {}/{} objects in the DR6+DR7 footprint.'.format(len(outcat), len(cat)))
    
    return outcat

In [42]:
sample = build_sample(cat, use_nproc=1)


Working on NGC6742...
For NGC6742 found 7 CCDs, RA = 284.83292, Dec = 48.46528, Diameter=1.5000 arcmin
Working on M92...
For M92 found 33 CCDs, RA = 259.28029, Dec = 43.13653, Diameter=20.0000 arcmin
Working on Bow-Shock...
For Bow-Shock found 85 CCDs, RA = 325.68720, Dec = 1.00320, Diameter=4.0000 arcmin
Working on NGC2782...
For NGC2782 found 13 CCDs, RA = 138.52129, Dec = 40.11369, Diameter=7.0000 arcmin
Working on UGC10321...
For UGC10321 found 21 CCDs, RA = 244.52800, Dec = 21.55910, Diameter=3.0000 arcmin
Working on C4-2010...
For C4-2010 found 46 CCDs, RA = 29.07064, Dec = 1.05082, Diameter=5.0000 arcmin
Found 6/6 objects in the DR6+DR7 footprint.

In [43]:
print('Writing {}'.format(figfile))
sample.write(figfile, overwrite=True)


Writing /global/project/projectdirs/desi/users/ioannis/legacysurveys/overview-paper/gallery.fits

In [44]:
sample


Out[44]:
<Table length=6>
namenicenameviewerdrradecdiam
str9str27str88str3float64float64float32
NGC6742NGC 6742 Planetary Nebulahttp://legacysurvey.org/viewer/?layer=decals-dr6&ra=284.83291667&dec=48.46527778dr6284.8329166748.465277781.5
M92Messier 92 Globular Clusterhttp://legacysurvey.org/viewer/?layer=decals-dr6&ra=259.28029167&dec=43.13652778&zoom=12dr6259.2802916743.1365277820.0
Bow-ShockInterstellar Bow Shockhttp://legacysurvey.org/viewer?ra=325.6872&dec=1.0032&zoom=14&layer=decals-dr5dr7325.68721.00324.0
NGC2782NGC 2782http://legacysurvey.org/viewer/?layer=decals-dr6&ra=138.52129167&dec=40.11369444&zoom=12dr6138.5212916740.113694447.0
UGC10321UGC 10321 Galaxy Grouphttp://legacysurvey.org/viewer?ra=244.5280&dec=21.5591&zoom=14&layer=decals-dr5dr7244.52821.55913.0
C4-2010SDSS/C4 Galaxy Cluster 2010http://legacysurvey.org/viewer?ra=29.0707&dec=1.0510&zoom=13&layer=decals-dr5dr729.0706414921.0508166675.0

Generate the color mosaics for each object.


In [49]:
def custom_brickname(obj, prefix='custom-'): 
    brickname = 'custom-{:06d}{}{:05d}'.format(
        int(1000*obj['ra']), 'm' if obj['dec'] < 0 else 'p', 
        int(1000*np.abs(obj['dec'])))
    return brickname

In [50]:
def custom_coadds_one(obj, scale=PIXSCALE, clobber=False):
    from legacypipe.runbrick import run_brick
    #from astrometry.util.multiproc import multiproc
    #from legacypipe.runbrick import stage_tims, run_brick
    #from legacypipe.coadds import make_coadds

    name = obj['name']
    jpgfile = os.path.join(jpgdir, '{}.jpg'.format(name))
    if os.path.isfile(jpgfile) and not clobber:
        print('File {} exists...skipping.'.format(jpgfile))
    else:
        size = np.rint(obj['diam'] * 60 / scale).astype('int') # [pixels]
        print('Generating mosaic for {} with width={} pixels.'.format(name, size))
        
        bands = ('g', 'r', 'z')
        if 'Bow' in name:
            rgb_kwargs = dict({'Q': 200, 'm': 0.01})
        else:
            rgb_kwargs = dict({'Q': 20, 'm': 0.03})
            
        survey = init_survey(dr=obj['dr'])
        brickname = custom_brickname(obj, prefix='custom-')
        
        with warnings.catch_warnings():
            warnings.simplefilter("ignore")
            run_brick(None, survey, radec=(obj['ra'], obj['dec']), pixscale=scale, 
                      width=size, height=size, rgb_kwargs=rgb_kwargs, threads=nproc, 
                      stages=['image_coadds'], splinesky=True, early_coadds=True, pixPsf=True, 
                      hybridPsf=True, normalizePsf=True, write_pickles=False, depth_cut=False, 
                      apodize=True, do_calibs=False, ceres=False)

        sys.stdout.flush()    
        _jpgfile = os.path.join(survey.output_dir, 'coadd', 'cus', brickname, 
                               'legacysurvey-{}-image.jpg'.format(brickname))
        shutil.copy(_jpgfile, jpgfile)
        shutil.rmtree(os.path.join(survey.output_dir, 'coadd'))

In [51]:
#custom_coadds_one(sample[2], clobber=True)

In [52]:
def custom_coadds(sample, clobber=False):
    for obj in sample:
        custom_coadds_one(obj, clobber=clobber)

In [53]:
coaddslogfile = os.path.join(figdir, 'make-coadds.log')
print('Generating the coadds.')
print('Logging to {}'.format(coaddslogfile))
t0 = time.time()
with open(coaddslogfile, 'w') as log:
    with redirect_stdout(log):
        custom_coadds(sample, clobber=True)
print('Total time = {:.3f} minutes.'.format((time.time() - t0) / 60))


Generating the coadds.
Logging to /global/project/projectdirs/desi/users/ioannis/legacysurveys/overview-paper/make-coadds.log
Total time = 2.502 minutes.

Add labels and a scale bar.


In [54]:
barlen = np.round(60.0 / PIXSCALE).astype('int')
fonttype = os.path.join(figdir, 'Georgia.ttf')

In [55]:
def _add_labels_one(args):
    """Wrapper function for the multiprocessing."""
    return add_labels_one(*args)

In [56]:
def add_labels_one(obj, verbose=False):
    name = obj['name']
    nicename = obj['nicename']
    
    jpgfile = os.path.join(jpgdir, '{}.jpg'.format(name))
    pngfile = os.path.join(pngdir, '{}.png'.format(name))
    thumbfile = os.path.join(pngdir, 'thumb-{}.png'.format(name))
        
    im = Image.open(jpgfile)
    sz = im.size
    fntsize = np.round(sz[0]/28).astype('int')
    width = np.round(sz[0]/175).astype('int')
    font = ImageFont.truetype(fonttype, size=fntsize)
    draw = ImageDraw.Draw(im)
    # Label the object name--
    draw.text((0+fntsize*2, 0+fntsize*2), nicename, font=font)
    # Add a scale bar--
    x0, x1, yy = sz[1]-fntsize*2-barlen, sz[1]-fntsize*2, sz[0]-fntsize*2
    draw.line((x0, yy, x1, yy), fill='white', width=width)
    im.save(pngfile)    
        
    # Generate a thumbnail
    if False:
        cmd = '/usr/bin/convert -thumbnail 300x300 {} {}'.format(pngfile, thumbfile)
        os.system(cmd)

In [57]:
def add_labels(sample):
    labelargs = list()
    for obj in sample:
        labelargs.append((obj, False))

    if nproc > 1:
        p = multiprocessing.Pool(nproc)
        res = p.map(_add_labels_one, labelargs)
        p.close()
    else:
        for args in labelargs:
            res = _add_labels_one(args)

In [59]:
%time add_labels(sample)


CPU times: user 73.6 ms, sys: 940 ms, total: 1.01 s
Wall time: 4.86 s

Finally make a nice montage figure for the paper.


In [60]:
def make_montage(cat, clobber=False):
    montagefile = os.path.join(figdir, 'overview-gallery.png')

    ncol = 3
    nrow = np.ceil(len(sample) / ncol).astype('int')
    
    if not os.path.isfile(montagefile) or clobber:
        cmd = 'montage -bordercolor white -borderwidth 1 -tile {}x{} -geometry 512x512 '.format(ncol, nrow)
        cmd = cmd+' '.join([os.path.join(pngdir, '{}.png'.format(name)) for name in cat['name']])
        cmd = cmd+' {}'.format(montagefile)
        print(cmd)
        os.system(cmd)        
        print('Writing {}'.format(montagefile))

In [61]:
%time make_montage(cat, clobber=True)


montage -bordercolor white -borderwidth 1 -tile 3x2 -geometry 512x512 /global/project/projectdirs/desi/users/ioannis/legacysurveys/overview-paper/png/NGC6742.png /global/project/projectdirs/desi/users/ioannis/legacysurveys/overview-paper/png/M92.png /global/project/projectdirs/desi/users/ioannis/legacysurveys/overview-paper/png/Bow-Shock.png /global/project/projectdirs/desi/users/ioannis/legacysurveys/overview-paper/png/NGC2782.png /global/project/projectdirs/desi/users/ioannis/legacysurveys/overview-paper/png/UGC10321.png /global/project/projectdirs/desi/users/ioannis/legacysurveys/overview-paper/png/C4-2010.png /global/project/projectdirs/desi/users/ioannis/legacysurveys/overview-paper/overview-gallery.png
Writing /global/project/projectdirs/desi/users/ioannis/legacysurveys/overview-paper/overview-gallery.png
CPU times: user 36.2 ms, sys: 107 ms, total: 143 ms
Wall time: 3.92 s

In [ ]: