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:
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
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)
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]:
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]:
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)
In [43]:
print('Writing {}'.format(figfile))
sample.write(figfile, overwrite=True)
In [44]:
sample
Out[44]:
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))
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)
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)
In [ ]: