# ---
# jupyter:
#   jupytext:
#     text_representation:
#       extension: .py
#       format_name: light
#       format_version: '1.5'
#       jupytext_version: 1.14.4
#   kernelspec:
#     display_name: Python 3
#     language: python
#     name: python3
# ---

# # Prepares dataset of HSC galaxies, PSFs and masks

import h5py
from astropy.table import Table
import astropy.units as u
from unagi import hsc
from unagi import task

# Define the HSC archive
archive = hsc.Hsc(dr='pdr2', rerun='pdr2_wide')

# +
query_mask = ''
-- Merge forced photometry and spectroscopic sample from HSC PDR 2 wide
SELECT object_id, ra, dec, tract, patch,
	-- Absorption
	a_g, a_r, a_i, a_z, a_y,
	-- Extendedness
	g_extendedness_value, r_extendedness_value, i_extendedness_value, z_extendedness_value, y_extendedness_value,
  -- Background Information
  g_localbackground_flux, r_localbackground_flux, i_localbackground_flux, z_localbackground_flux, y_localbackground_flux,
	-- Magnitudes
	g_cmodel_mag, g_cmodel_magsigma, g_cmodel_exp_mag, g_cmodel_exp_magsigma, g_cmodel_dev_mag, g_cmodel_dev_magsigma,
	r_cmodel_mag, r_cmodel_magsigma, r_cmodel_exp_mag, r_cmodel_exp_magsigma, r_cmodel_dev_mag, r_cmodel_dev_magsigma,
	i_cmodel_mag, i_cmodel_magsigma, i_cmodel_exp_mag, i_cmodel_exp_magsigma, i_cmodel_dev_mag, i_cmodel_dev_magsigma,
	z_cmodel_mag, z_cmodel_magsigma, z_cmodel_exp_mag, z_cmodel_exp_magsigma, z_cmodel_dev_mag, z_cmodel_dev_magsigma,
	y_cmodel_mag, y_cmodel_magsigma, y_cmodel_exp_mag, y_cmodel_exp_magsigma, y_cmodel_dev_mag, y_cmodel_dev_magsigma

FROM pdr2_wide.forced forced
  LEFT JOIN pdr2_wide.forced2 USING (object_id)
  LEFT JOIN pdr2_wide.forced3 USING (object_id)

-- Applying some data quality cuts
WHERE forced.isprimary
AND forced.i_cmodel_mag < 23.5
AND forced.i_cmodel_mag > 21
-- Simple Full Depth Full Colour cuts: At least 3 exposures in each band
AND forced.g_inputcount_value >= 3
AND forced.r_inputcount_value >= 3
AND forced.i_inputcount_value >= 3
AND forced.z_inputcount_value >= 3
AND forced.y_inputcount_value >= 3
-- Remove objects affected by bright stars
AND NOT forced.g_pixelflags_bright_objectcenter
AND NOT forced.r_pixelflags_bright_objectcenter
AND NOT forced.i_pixelflags_bright_objectcenter
AND NOT forced.z_pixelflags_bright_objectcenter
AND NOT forced.y_pixelflags_bright_objectcenter
AND NOT forced.g_pixelflags_bright_object
AND NOT forced.r_pixelflags_bright_object
AND NOT forced.i_pixelflags_bright_object
AND NOT forced.z_pixelflags_bright_object
AND NOT forced.y_pixelflags_bright_object
-- Remove objects intersecting edges
AND NOT forced.g_pixelflags_edge
AND NOT forced.r_pixelflags_edge
AND NOT forced.i_pixelflags_edge
AND NOT forced.z_pixelflags_edge
AND NOT forced.y_pixelflags_edge
-- Remove objects with saturated pixels
AND NOT forced.g_pixelflags_saturatedcenter
AND NOT forced.r_pixelflags_saturatedcenter
AND NOT forced.i_pixelflags_saturatedcenter
AND NOT forced.z_pixelflags_saturatedcenter
AND NOT forced.y_pixelflags_saturatedcenter
-- But force objects with interpolated pixels
AND forced.i_pixelflags_interpolatedcenter
-- Remove objects with generic cmodel fit failures
AND NOT forced.g_cmodel_flag
AND NOT forced.r_cmodel_flag
AND NOT forced.i_cmodel_flag
AND NOT forced.z_cmodel_flag
AND NOT forced.y_cmodel_flag
-- Sort by tract and patch for faster cutout query
ORDER BY object_id
LIMIT 10000
;
'''
# -

catalog = archive.sql_query(query_mask, from_file=False, verbose=True)

catalog.write('catalog_masked_obj.fits')

# +
# Defining size of cutouts
img_len = 42 / 2 # Size of cutouts in pixels
cutout_size = 0.168*(img_len) # Size of cutouts in Arcsecs

# Which filter we care about
filters = ['HSC-I']

tmp_dir='tmp_dir'
out_dir='./'
# !mkdir -p tmp_dir
# -

# Extract the cutouts
cutouts_filename = task.hsc_bulk_cutout(catalog, 
                                        cutout_size=cutout_size * u.Unit('arcsec'), 
                                        filters=filters, 
                                        archive=archive,  
                                        nproc=2, # Download using 2 parallel jobs
                                        tmp_dir=tmp_dir, 
                                        mask=True, variance=True,
                                        output_dir=out_dir)

# Extract the PSFs for all these objects
psfs_filename = task.hsc_bulk_psf(catalog, filters=filters, 
                                  archive=archive,  
                                  nproc=2, tmp_dir=tmp_dir,
                                  output_dir=out_dir)

cutouts_filename, psfs_filename



StarcoderdataJupyterScriptsDedupFiltered