<jupyter_start><jupyter_text>Check the input imagesI here want to investigate the input files to check which patches they overlap and what other relevant files are reuired.The first run on exposures was tried without performing dark or flat correction. We hope that these might allow the removal of the large negative pixel values.<jupyter_code>from astropy.io import fits
from astropy.wcs import WCS
import glob
import numpy as np

%matplotlib inline
#%config InlineBackend.figure_format = 'svg'

import matplotlib.pyplot as plt
inputs = []
minimal_stacks = open('minimal_stacks.lis', 'r')
for line in minimal_stacks.readlines():
    #print(line.strip())
    inputs.append(line.strip())
    
    
inputs_exposures = []
minimal_exposures = open('minimal_exposures.lis', 'r')
for line in minimal_exposures.readlines():
    #print(line.strip())
    inputs_exposures.append(line.strip())
DATA_LOC = '/Volumes/Raph500/lsst-ir-fusion'
inputs[0]
inputs[0][5:]
f = fits.open(DATA_LOC + inputs[0][5:])
f[0].header
f[1].header
f = fits.open(DATA_LOC + inputs_exposures[0][5:])
f[0].header
f[1].header
flats = []
darks = []
conf = []
lincor = []
for f in inputs:
    header = fits.open(DATA_LOC + f[5:])[1].header
    flats.append(header['FLATCOR'].split('[')[0])
    darks.append(header['DARKCOR'].split('[')[0])
    conf.append(header['CIR_CPM'].split('[')[0])
    lincor.append(header['LINCOR'].split('[')[0])
set(flats) , set(darks), set(lincor)
flats_exposures = []
darks_exposures = []
conf_exposures = []
lincor_exposures = []
for f in inputs_exposures:
    header = fits.open(DATA_LOC + f[5:])[1].header
    flats_exposures.append(header['FLATCOR'].split('[')[0])
    darks_exposures.append(header['DARKCOR'].split('[')[0])
    conf_exposures.append(header['CIR_CPM'].split('[')[0])
    lincor_exposures.append(header['LINCOR'].split('[')[0])
set(flats_exposures) ,set(darks_exposures),set(conf_exposures),set(lincor_exposures)<jupyter_output><empty_output><jupyter_text>Patches<jupyter_code>#We are using a rings skymap
from lsst.geom import SpherePoint 
from lsst.geom import degrees
from lsst.skymap.ringsSkyMap import RingsSkyMap, RingsSkyMapConfig

#Following taken from dmu1/2_Survey_comparisons.ipynb test choice
coord = SpherePoint(35.429025*degrees,-4.90853*degrees)
config = RingsSkyMapConfig()
#These config options are chose to be the same as HSC:
#https://github.com/lsst/obs_subaru/blob/master/config/hsc/makeSkyMap.py
#and copied for obs_vista
#https://github.com/lsst-uk/obs_vista/blob/master/config/makeSkyMap.py
config.numRings = 120
config.projection = "TAN"
config.tractOverlap = 1.0/60 # Overlap between tracts (degrees)
config.pixelScale = 0.168
sm = RingsSkyMap(config)
sm.findTract(coord)
sxds_tracts = [8282,8283,8284,8523,8524,8525,8765,8766,8767] 
fig, ax = plt.subplots(figsize=(4, 4), dpi=140)
for tract in sxds_tracts:
    vertices = sm.generateTract(tract).getVertexList()
    ra = [vertices[n][0].asDegrees() for n in np.mod(np.arange(5),4)]
    dec = [vertices[n][1].asDegrees() for n in np.mod(np.arange(5),4)]
    ax.plot(ra, dec, c = 'k')
    for p in np.nditer(np.meshgrid(np.arange(0,9), np.arange(0,9))):
        #print(p[0], p[1])
        t = sm.generateTract(tract)
        vertices = t.getPatchInfo([p[0], p[1]]).getInnerSkyPolygon(t.getWcs()).getVertices()
        ra = [np.arctan(vertices[n][1]/vertices[n][0])* 180/np.pi for n in np.mod(np.arange(5),4)]
        #print(tract, ra)
        dec = [(vertices[n][2])* 180/np.pi for n in np.mod(np.arange(5),4)]
        ax.plot(ra, dec, c = 'k', alpha=0.5, linewidth=0.5)

#Plot the 9 test patches
for p in np.nditer(np.meshgrid(np.arange(3,6), np.arange(3,6))):
    #print(p[0], p[1])
    t = sm.generateTract(8524)
    vertices = t.getPatchInfo([p[0], p[1]]).getInnerSkyPolygon(t.getWcs()).getVertices()
    ra = [np.arctan(vertices[n][1]/vertices[n][0])* 180/np.pi for n in np.mod(np.arange(5),4)]
    #print(tract, ra)
    dec = [(vertices[n][2])* 180/np.pi for n in np.mod(np.arange(5),4)]
    ax.plot(ra, dec, c = 'r')
    
ax.axis('scaled')
ax.set_xlim([39, 33])
ax.set_xlabel('Right Ascension [deg]')
ax.set_ylabel('Declination[deg]')
fig.savefig('./figs/sxds_tracts_patches_tiles.pdf', overwrite=True)
fig.savefig('./figs/sxds_tracts_patches_tiles.png', overwrite=True)
tiles = glob.glob(DATA_LOC + '/dmu0/dmu0_VISTA/dmu0_VIDEO/data/*/*_tl.fit')
fits.open(tiles[0])[1].header
tiles
hdu = fits.open(tiles[0])
wcs = WCS(hdu[1].header)
ra, dec = wcs.all_pix2world([0,0,hdu[1].header['NAXIS1'],hdu[1].header['NAXIS1']], 
                            [0,hdu[1].header['NAXIS2'],0,hdu[1].header['NAXIS2']], 1) 
ra
fig, ax = plt.subplots(figsize=(4, 4), dpi=140)
for tract in sxds_tracts:
    vertices = sm.generateTract(tract).getVertexList()
    ra = [vertices[n][0].asDegrees() for n in np.mod(np.arange(5),4)]
    dec = [vertices[n][1].asDegrees() for n in np.mod(np.arange(5),4)]
    ax.plot(ra, dec, c = 'k', linewidth=0.7)
    for p in np.nditer(np.meshgrid(np.arange(0,9), np.arange(0,9))):
        #print(p[0], p[1])
        t = sm.generateTract(tract)
        vertices = t.getPatchInfo([p[0], p[1]]).getInnerSkyPolygon(t.getWcs()).getVertices()
        ra = [np.arctan(vertices[n][1]/vertices[n][0])* 180/np.pi for n in np.mod(np.arange(5),4)]
        #print(tract, ra)
        dec = [(vertices[n][2])* 180/np.pi for n in np.mod(np.arange(5),4)]
        ax.plot(ra, dec, c = 'k', alpha=0.5, linewidth=0.3)

#Plot the 9 test patches
for p in np.nditer(np.meshgrid(np.arange(3,6), np.arange(3,6))):
    #print(p[0], p[1])
    t = sm.generateTract(8524)
    vertices = t.getPatchInfo([p[0], p[1]]).getInnerSkyPolygon(t.getWcs()).getVertices()
    ra = [np.arctan(vertices[n][1]/vertices[n][0])* 180/np.pi for n in np.mod(np.arange(5),4)]
    #print(tract, ra)
    dec = [(vertices[n][2])* 180/np.pi for n in np.mod(np.arange(5),4)]
    ax.plot(ra, dec, c = 'r')
    
    
for t in tiles:
    hdu = fits.open(t)
    wcs = WCS(hdu[1].header)
    ra, dec = wcs.all_pix2world([0,0,hdu[1].header['NAXIS1'],hdu[1].header['NAXIS1']], 
                            [0,hdu[1].header['NAXIS2'],0,hdu[1].header['NAXIS2']], 1) 
    c = [0,1,3,2,0]
    ax.plot([ra[n] for n in c],[dec[n] for n in c], c='b')

ax.axis('scaled')
ax.set_xlim([39, 33])
ax.set_xlabel('Right Ascension [deg]')
ax.set_ylabel('Declination [deg]')

plt.grid(linestyle="--", color='black', linewidth=0.5, alpha=0.5)

fig.savefig('./figs/sxds_tracts_patches_tiles.pdf', overwrite=True)
fig.savefig('./figs/sxds_tracts_patches_tiles.png', overwrite=True)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.18831268165825302, 'fraction_numerical': 0.04099739941869359, 'mean_word_length': 4.554800339847069, 'pattern_counts': {'":': 0, '<': 9, '<?xml version=': 0, '>': 9, 'https://': 2, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '2', 'dir': 'jupyter-structured-clean-dedup', 'id': '324978', 'license': 'Apache-2.0', 'n_tokens_mistral': 2830, 'n_tokens_neox': 2704, 'n_words': 467, 'path': 'dmu4/dmu4_Example/4_Check_inputs.ipynb', 'repo_name': 'raphaelshirley/lsst-ir-fusion', 'size': '6537'}