<jupyter_start><jupyter_text>Deep dive on Overscan    v1 -- Commissioning data    v2 -- minisv2<jupyter_code>%matplotlib inline
# imports
import os
import numpy as np
from matplotlib import pyplot as plt
from importlib import reload
import glob

from scipy import signal

from astropy.io import fits
from astropy.stats import sigma_clipped_stats

from desispec import preproc

from bokeh import plotting
from bokeh import models
#exp_path = '/home/xavier/DESI/Commissioning/exposures'
exp_path = '/home/xavier/DESI/DESI_SCRATCH/'
night = '20191024'<jupyter_output><empty_output><jupyter_text>Methods<jupyter_code>def add_gain(header, SP='sp2'):
    if SP=='sp2':
        header['GAINA'] = 1.21763880607
        header['GAINB'] = 1.26599257375
        header['GAINC'] = 1.37207665312
        header['GAIND'] = 1.28793077884
def get_oscan_col(rawimage, header, amp='A'):
    ov_col = preproc._parse_sec_keyword(header['BIASSEC'+amp])
    raw_overscan_col = rawimage[ov_col].copy()
    vals = []
    for row in range(raw_overscan_col.shape[0]):
        o, _ = preproc._overscan(raw_overscan_col[row,:])
        vals.append(o)
    overscan_col = np.array(vals)
    return overscan_col
def get_oscan_row(rawimage, header, amp='A'):
    ov_row_idx = preproc._parse_sec_keyword(header['ORSEC'+amp])
    raw_overscan_row = rawimage[ov_row_idx].copy()
    vals = []
    for col in range(raw_overscan_row.shape[1]):
        o, _ = preproc._overscan(raw_overscan_row[:,col])
        vals.append(o)
    overscan_row = np.array(vals)
    return overscan_row<jupyter_output><empty_output><jupyter_text>Find the right exposure<jupyter_code>files = glob.glob(os.path.join(exp_path, night, '*','desi-*.fits.fz'))
files
for ifile in files:
    head = fits.open(ifile)[1].header
    #import pdb; pdb.set_trace()
    print(ifile, head['OBSTYPE'], head['EXPTIME'])<jupyter_output>/home/xavier/DESI/DESI_SCRATCH/20191024/00020665/desi-00020665.fits.fz ARC 10.0
/home/xavier/DESI/DESI_SCRATCH/20191024/00020664/desi-00020664.fits.fz ARC 10.0
/home/xavier/DESI/DESI_SCRATCH/20191024/00020655/desi-00020655.fits.fz ZERO 0.0
/home/xavier/DESI/DESI_SCRATCH/20191024/00020675/desi-00020675.fits.fz FLAT 120.0
/home/xavier/DESI/DESI_SCRATCH/20191024/00020668/desi-00020668.fits.fz ARC 5.0
/home/xavier/DESI/DESI_SCRATCH/20191024/00020669/desi-00020669.fits.fz ARC 5.0
/home/xavier/DESI/DESI_SCRATCH/20191024/00020657/desi-00020657.fits.fz ZERO 0.0
/home/xavier/DESI/DESI_SCRATCH/20191024/00020652/desi-00020652.fits.fz ZERO 0.0
/home/xavier/DESI/DESI_SCRATCH/20191024/00020671/desi-00020671.fits.fz ARC 5.0
/home/xavier/DESI/DESI_SCRATCH/20191024/00020651/desi-00020651.fits.fz ZERO 0.0
/home/xavier/DESI/DESI_SCRATCH/20191024/00020676/desi-00020676.fits.fz FLAT 120.0
/home/xavier/DESI/DESI_SCRATCH/20191024/00020658/desi-00020658.fits.fz ZERO 0.0
/home/xavier/DESI/DESI_SCRATCH/20191024[...]<jupyter_text>120s exposure<jupyter_code>iid = 20675
ifile = os.path.join(exp_path, night, '000{}'.format(iid), 'desi-000{}.fits.fz'.format(iid))
ifile
hdulist = fits.open(ifile)
hdulist.info()<jupyter_output>Filename: /home/xavier/DESI/DESI_SCRATCH/20191024/00020675/desi-00020675.fits.fz
No.    Name      Ver    Type      Cards   Dimensions   Format
  0  PRIMARY       1 PrimaryHDU       6   ()      
  1  SPS           1 CompImageHDU     58   (10,)   int16   
  2  Z3            1 CompImageHDU    137   (4256, 4194)   int16   
  3  R3            1 CompImageHDU    137   (4256, 4194)   int16   
  4  B3            1 CompImageHDU    137   (4232, 4162)   int16   
  5  SPECTCONS     1 BinTableHDU    119   1R x 42C   [K, K, 1PA(26), 1PA(15), D, 1PA(12), 1PA(19), 1PA(5), 1PA(4), 1PA(3), 1PA(4), 1PA(3), 1PA(5), 1PA(4), 1PA(8), 1PA(2), 1PA(6), 1PA(8), 1PA(2), 1PA(8), D, D, D, D, D, D, D, D, D, D, 1PA(3), 1PA(3), 1PA(3), 1PA(3), 1PA(3), 1PA(3), 1PA(7), 1PA(7), 1PA(7), 1PA(3), 1PA(3), 1PA(3)]<jupyter_text>z3<jupyter_code>z3_hdu_120 = hdulist['z3']<jupyter_output><empty_output><jupyter_text>Examine overscan_column<jupyter_code>header = z3_hdu_120.header
rawimage = z3_hdu_120.data.astype(np.float64)
amp = 'A'
ov_col = preproc._parse_sec_keyword(header['BIASSEC'+amp])
raw_overscan_col = rawimage[ov_col].copy()
raw_overscan_col.shape
vals = []
for row in range(raw_overscan_col.shape[0]):
    o, _ = preproc._overscan(raw_overscan_col[row,:])
    vals.append(o)
overscan_col = np.array(vals)
overscan_col.shape
# Plot
plt.clf()
ax = plt.gca()
ax.plot(overscan_col)
plt.show()
plt.clf()
ax = plt.gca()
ax.hist(overscan_col)
plt.show()
mean, med_sn2, stddev = sigma_clipped_stats(overscan_col,sigma_lower=3.0,sigma_upper=3.0)
stddev<jupyter_output><empty_output><jupyter_text>Let's try a savgol<jupyter_code>oscol_sub = signal.savgol_filter(overscan_col, 65, 5)
# Plot
plt.clf()
ax = plt.gca()
ax.plot(overscan_col)
ax.plot(oscol_sub, 'r-')
ax.set_xlim(775, 1020)
ax.set_ylim(2006., 2007.25)
plt.show()
plt.clf()
ax = plt.gca()
ax.hist(overscan_col-oscol_sub)
plt.show()
mean, med_sn2, stddev_sub = sigma_clipped_stats(overscan_col-oscol_sub,sigma_lower=3.0,sigma_upper=3.0)
stddev_sub<jupyter_output><empty_output><jupyter_text>Overscan row<jupyter_code>ov_row_idx = preproc._parse_sec_keyword(header['ORSEC'+amp])
raw_overscan_row = rawimage[ov_row_idx].copy()
raw_overscan_row.shape
vals = []
for col in range(raw_overscan_row.shape[1]):
    o, _ = preproc._overscan(raw_overscan_row[:,col])
    vals.append(o)
overscan_row = np.array(vals)
overscan_row.shape
# Plot
plt.clf()
ax = plt.gca()
ax.plot(overscan_row)
plt.show()<jupyter_output><empty_output><jupyter_text>Sobolev<jupyter_code>osrow_sub = signal.savgol_filter(overscan_row, 65, 5)
# Plot
plt.clf()
ax = plt.gca()
ax.plot(overscan_row)
ax.plot(osrow_sub, 'r-')
plt.show()
plt.clf()
ax = plt.gca()
ax.hist(overscan_row-osrow_sub, bins=100)
plt.show()
# Stats
mean, med_sn2, stddev_sub_row = sigma_clipped_stats(overscan_row-osrow_sub,sigma_lower=3.0,sigma_upper=3.0)
stddev_sub_row<jupyter_output><empty_output><jupyter_text>Explore a zero image<jupyter_code>night = '20190620'
iid = 16394
ifile = os.path.join(exp_path, night, '000{}'.format(iid), 'desi-000{}.fits.fz'.format(iid))
ifile
hdulist = fits.open(ifile)
z2_zero = hdulist['z2']
header = z2_zero.header
rawimage = z2_zero.data.astype(np.float64)<jupyter_output><empty_output><jupyter_text>Overscan col<jupyter_code>oscan_col = get_oscan_col(rawimage, header)
plt.clf()
ax = plt.gca()
ax.plot(oscan_col)
plt.show()<jupyter_output><empty_output><jupyter_text>Overscan row<jupyter_code>oscan_row = get_oscan_row(rawimage, header)
plt.clf()
ax = plt.gca()
ax.plot(oscan_row)
plt.show()<jupyter_output><empty_output><jupyter_text>Let's process with savgol<jupyter_code>reload(preproc)
img_zero = preproc.preproc(rawimage, header, None, bias=False, pixflat=False, mask=False, 
                      nocrosstalk=True, ccd_calibration_filename=False, dark=False)
img_zero.pix.shape
#from pypeit import ginga
#ginga.show_image(img_zero.pix)
vmin, vmax = -9., 9.
plt.clf()
#plt.figure(figsize=(13,7))
plt.imshow(img_zero.pix[1800:2200,500:1000], vmin=vmin, vmax=vmax)
plt.show()<jupyter_output><empty_output><jupyter_text>----  Image 20675, r3<jupyter_code>#night = '20190627'
#iid = 16727
iid = 20675
ifile = os.path.join(exp_path, night, '000{}'.format(iid), 'desi-000{}.fits.fz'.format(iid))
ifile
hdulist = fits.open(ifile)
r3 = hdulist['r3']
reload(preproc)
img_r3 = preproc.preproc(r3.data, r3.header, None, bias=False, pixflat=False, mask=False, 
                      nocrosstalk=True, ccd_calibration_filename=False, dark=False)
from pypeit import ginga
ginga.show_image(img_r3.pix)
vmin, vmax = -9., 9.
plt.clf()
#plt.figure(figsize=(13,7))
plt.imshow(img_r3.pix, vmin=vmin, vmax=vmax, origin='lower')
plt.show()<jupyter_output><empty_output><jupyter_text>----  Bias<jupyter_code>iid = 20655
ifile = os.path.join(exp_path, night, '000{}'.format(iid), 'desi-000{}.fits.fz'.format(iid))
ifile
hdulist = fits.open(ifile)
#hdulist.info()
z3 = hdulist['z3']
reload(preproc)
img_z3 = preproc.preproc(z3.data, z3.header, None, bias=False, pixflat=False, mask=False, 
                      nocrosstalk=True, ccd_calibration_filename=False, dark=False)
ginga.show_image(img_z3.pix)
vmin, vmax = -9., 9.
plt.clf()
#plt.figure(figsize=(13,7))
plt.imshow(img_z3.pix, vmin=vmin, vmax=vmax, origin='lower')
plt.show()<jupyter_output><empty_output><jupyter_text>----  Exploring yet another bias frame : z3-00022217<jupyter_code>def set_fontsize(p, fsz):
    '''
    Parameters
    ----------
    ax : Bokeh plot class
    fsz : float
      Font size
    '''
    p.xaxis.axis_label_text_font_size = '{:d}pt'.format(fsz)
    p.xaxis.major_label_text_font_size = "{:d}pt".format(fsz)
    #
    p.yaxis.axis_label_text_font_size = '{:d}pt'.format(fsz)
    p.yaxis.major_label_text_font_size = "{:d}pt".format(fsz)<jupyter_output><empty_output><jupyter_text>Load em<jupyter_code>path = '/home/xavier/DESI/DESI_SCRATCH/Commissioning/Overscan/'
hdu_sav = fits.open(os.path.join(path, 'preproc-z3-00022217_savgol.fits'))
data_sav = hdu_sav[0].data
hdu_sav.info()
hdu_no = fits.open(os.path.join(path, 'preproc-z3-00022217_nosavgol.fits'))
data_no = hdu_no[0].data<jupyter_output><empty_output><jupyter_text>Compare em<jupyter_code>med_sav = np.median(data_sav[500:1000,:],axis=0)
med_no = np.median(data_no[500:1000,:],axis=0)
med_sav.size
# Plot
plotting.output_notebook()

p = plotting.figure(title='Median Overscan', x_axis_label='Column',
                   y_axis_label='Counts')

p.line(np.arange(4114), med_no, legend='No Sav Gol', line_width=2, color='red')
p.line(np.arange(4114), med_sav, legend='Sav Gol', line_width=2, color='blue')
#p.line(np.arange(4114), med_master, legend='Master', line_width=2, color='green')


p.legend.location = "bottom_left"

set_fontsize(p, 16)

plotting.show(p)<jupyter_output><empty_output><jupyter_text>----  MiniSV2<jupyter_code># ginga
from pypeit import ginga
os.environ['DESI_SPECTRO_CALIB'] = '/home/xavier/DESI/spectro/desi_spectro_calib'<jupyter_output><empty_output><jupyter_text>Load<jupyter_code>exp_path = '/scratch/DESI_SCRATCH/minisv2/data'
night = '20200304'
iid = 53018
ifile = os.path.join(exp_path, night, '000{}'.format(iid), 'desi-000{}.fits.fz'.format(iid))
ifile, os.path.isfile(ifile)
hdulist = fits.open(ifile)
hdulist.info()
camera = 'z3'
hdu = hdulist[camera.upper()]
#r4_calib_file = os.path.join(os.getenv('DESI_SPECTRO_CALIB'), 'spec', 'sm1', 'sm1-r.yaml')\n",
smnum = hdu.header['SPECGRPH']
smnum = hdu.header['SPECID']
#calib_file = os.path.join(os.getenv('DESI_SPECTRO_CALIB'), 'spec', 'sp{}'.format(smnum),
#                         '{}{}.yaml'.format(camera[0],smnum))
calib_file = os.path.join(os.getenv('DESI_SPECTRO_CALIB'), 'spec', 'sm{}'.format(smnum),
                         'sm{}-{}.yaml'.format(smnum, camera[0]))
hdu.header<jupyter_output><empty_output><jupyter_text>Default<jupyter_code>from desispec import calibfinder
reload(calibfinder)
reload(preproc)
img_def = preproc.preproc(hdu.data, hdu.header, hdulist[0].header, bias=True, pixflat=False, mask=False,
                        nocrosstalk=True, ccd_calibration_filename=calib_file, dark=False, flag_savgol=False)
ginga.show_image(img_def.pix, chname='default')<jupyter_output><empty_output><jupyter_text>Savgol recipe<jupyter_code>reload(preproc)
img_savgol = preproc.preproc(hdu.data, hdu.header, hdulist[0].header, bias=True, pixflat=False, mask=False,
                        nocrosstalk=True, ccd_calibration_filename=calib_file, dark=False, 
                          overscan_per_row=True, flag_savgol=True)
ginga.show_image(img_savgol.pix, chname='savgol')<jupyter_output><empty_output><jupyter_text>Compare<jupyter_code>slc = slice(3000,3500)
prof_def=np.median(img_def.pix[slc],axis=0)
#prof_both = np.median(img_full.pix[2070:4000], axis=0)
prof_gol = np.median(img_savgol.pix[slc], axis=0)

plotting.output_notebook()

p = plotting.figure(title='Overscan comparisons', x_axis_label='column',
                   y_axis_label='electrons', plot_width=800)#, x_range=(0,-90))

p.line(np.arange(prof_def.size), prof_def, legend='default')#, line_width=2)
#p.line(np.arange(prof_row.size), prof_both, legend='per_row, savgol', color='green')
p.line(np.arange(prof_gol.size), prof_gol, legend='savgol', color='red')
#p.line(z, v, legend='v', line_width=2, color='red')

p.legend.location = "top_right"
#set_fontsize(p, 16)\n",

plotting.show(p)
D = (slice(2100,4000), slice(2100,4000))
o_def, r_def = preproc._overscan(img_def.pix[D])
o_def, r_def
o_gol, r_gol = preproc._overscan(img_savgol.pix[D])
o_gol, r_gol
np.std(prof_def[:2000:]), np.std(prof_gol[:2000])
img_savgol.pix.shape<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.15678200966030564, 'fraction_numerical': 0.06888906485074035, 'mean_word_length': 5.36272040302267, 'pattern_counts': {'":': 0, '<': 91, '<?xml version=': 0, '>': 91, 'https://': 0, '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': '23', 'dir': 'jupyter-structured-clean-dedup', 'id': '553493', 'license': 'BSD-3-Clause', 'n_tokens_mistral': 6060, 'n_tokens_neox': 5375, 'n_words': 738, 'path': 'doc/nb/Overscan.ipynb', 'repo_name': 'echaussidon/desispec', 'size': '12629'}