<jupyter_start><jupyter_text>GFA Zero Calibration GFA calibrations should normally be updated in the following sequence: zeros, flats, darks.This notebook should be run using a DESI kernel, e.g. `DESI master`.<jupyter_code>%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
import os
import sys
import json
import collections
from pathlib import Path
import scipy.interpolate
import scipy.stats
import fitsio<jupyter_output><empty_output><jupyter_text>Install / upgrade the `desietcimg` package:<jupyter_code>try:
    import desietcimg
    print('desietcimg already installed')
except ImportError:
    print('Installing desietcimg...')
    !{sys.executable} -m pip install --user git+https://github.com/dkirkby/desietcimg
upgrade = False
if upgrade:
    print('Upgrading desietcimg...')
    !{sys.executable} -m pip install --upgrade --user git+https://github.com/dkirkby/desietcimg
import desietcimg.util
import desietcimg.plot
import desietcimg.gfa<jupyter_output><empty_output><jupyter_text>NERSC configuration:<jupyter_code>assert os.getenv('NERSC_HOST', False)
ROOT = Path('/project/projectdirs/desi/spectro/data/')
assert ROOT.exists()<jupyter_output><empty_output><jupyter_text>Initial GFA calibration:<jupyter_code>CALIB = Path('/global/cscratch1/sd/dkirkby/GFA_calib.fits')
assert CALIB.exists()<jupyter_output><empty_output><jupyter_text>Directory for saving plots:<jupyter_code>plotdir = Path('zerocal')
plotdir.mkdir(exist_ok=True)<jupyter_output><empty_output><jupyter_text>Process Zero Sequences Use a sequence of 200 zeros from [20191027](http://desi-www.kpno.noao.edu:8090/nightsum/nightsum-2019-10-27/nightsum.html).**Since this data has not yet been staged to its final location, we fetch it from the `lost+found` directory** (by overriding the definition of `ROOT` above):<jupyter_code>ROOT = Path('/global/project/projectdirs/desi/spectro/staging/lost+found/')
files = desietcimg.util.find_files(ROOT / '20191027' / '{N}/gfa-{N}.fits.fz', min=21968, max=22167)<jupyter_output><empty_output><jupyter_text>Build master zero images:<jupyter_code>def build_master_zero():
    master_zero = {}
    GFA = desietcimg.gfa.GFACamera(calib_name=str(CALIB))
    for k, gfa in enumerate(GFA.gfa_names):
        raw, meta = desietcimg.util.load_raw(files, 'EXPTIME', hdu=gfa)
        assert np.all(np.array(meta['EXPTIME']) == 0)
        GFA.setraw(raw, name=gfa, subtract_master_zero=False, apply_gain=False)
        master_zero[gfa] = np.median(GFA.data, axis=0)
    return master_zero
%time master_zero = build_master_zero()<jupyter_output>WARNING:root:Ignoring 503 bad overscan pixels for GUIDE8-E.<jupyter_text>Estimate the readnoise in ADU for each amplifier, using the new master zero:<jupyter_code>desietcimg.gfa.GFACamera.master_zero = master_zero
def get_readnoise(hrange=70, hbins=141, nsig=6, save=None):
    GFA = desietcimg.gfa.GFACamera(calib_name=str(CALIB))
    fig, axes = plt.subplots(5, 2, sharex=True, figsize=(18, 11))
    bins = np.linspace(-hrange, +hrange, hbins)
    noise = {}
    for k, gfa in enumerate(GFA.gfa_names):
        GFA.name = gfa
        ax = axes[k // 2, k % 2]
        raw, meta = desietcimg.util.load_raw(files, 'EXPTIME', hdu=gfa)
        assert np.all(np.array(meta['EXPTIME']) == 0)
        GFA.setraw(raw, name=gfa, subtract_master_zero=True, apply_gain=False)
        noise[gfa] = {}
        for j, amp in enumerate(GFA.amp_names):
            # Extract data for this quadrant.
            qdata = GFA.data[GFA.quad[amp]]
            X = qdata.reshape(-1)
            # Clip for std dev calculation.
            Xclipped, lo, hi = scipy.stats.sigmaclip(X, low=nsig, high=nsig)
            noise[gfa][amp] = np.std(Xclipped)
            label = f'{amp} {noise[gfa][amp]:.2f}'
            c = plt.rcParams['axes.prop_cycle'].by_key()['color'][j]
            ax.hist(X, bins=bins, label=label, color=c, histtype='step')
            for x in lo, hi:
                ax.axvline(x, ls='-', c=c, alpha=0.5)
        ax.set_yscale('log')
        ax.set_yticks([])
        if k in (8, 9):
            ax.set_xlabel('Zero Residual [ADU]')
            ax.set_xlim(bins[0], bins[-1])
        ax.legend(ncol=2, title=f'{gfa}', loc='upper left')
    plt.subplots_adjust(left=0.03, right=0.99, bottom=0.04, top=0.99, wspace=0.07, hspace=0.04)
    if save:
        plt.savefig(save)
    return noise
%time readnoise = get_readnoise(save=str(plotdir / 'GFA_readnoise.png'))
repr(readnoise)<jupyter_output><empty_output><jupyter_text>Save Updated Calibrations<jupyter_code>desietcimg.gfa.save_calib_data('GFA_calib_zero.fits', master_zero=master_zero, readnoise=readnoise)<jupyter_output>Using default master_dark
Using default pixel_mask
Using default gain
Using default tempfit
Saved GFA calib data to GFA_calib_zero.fits.<jupyter_text>Use this for subsequent flat and dark calibrations:<jupyter_code>!cp GFA_calib_zero.fits {CALIB}<jupyter_output><empty_output><jupyter_text>Comparisons Compare with the read noise values from the lab studies and <NAME>'s [independent analysis](https://desi.lbl.gov/trac/wiki/Commissioning/Planning/gfachar/bias_readnoise_20191027):<jupyter_code>ameisner_rdnoise = {
    'GUIDE0': { 'E': 5.56, 'F': 5.46, 'G': 5.12, 'H': 5.24},
    'FOCUS1': { 'E': 5.21, 'F': 5.11, 'G': 4.88, 'H': 4.90},
    'GUIDE2': { 'E': 7.11, 'F': 6.23, 'G': 5.04, 'H': 5.29},
    'GUIDE3': { 'E': 5.28, 'F': 5.16, 'G': 4.89, 'H': 5.00},
    'FOCUS4': { 'E': 5.23, 'F': 5.12, 'G': 5.01, 'H': 5.11},
    'GUIDE5': { 'E': 5.11, 'F': 5.00, 'G': 4.80, 'H': 4.86},
    'FOCUS6': { 'E': 5.12, 'F': 5.09, 'G': 4.85, 'H': 5.07},
    'GUIDE7': { 'E': 5.00, 'F': 4.96, 'G': 4.63, 'H': 4.79},
    'GUIDE8': { 'E': 6.51, 'F': 5.58, 'G': 5.12, 'H': 5.47},
    'FOCUS9': { 'E': 6.85, 'F': 5.53, 'G': 5.07, 'H': 5.57},
}
def compare_rdnoise(label='20191027', save=None):
    # Use the new calibrations written above.
    desietcimg.gfa.GFACamera.calib_data = None
    GFA = desietcimg.gfa.GFACamera(calib_name='GFA_calib_zero.fits')
    markers = '+xo.'
    fig, ax = plt.subplots(1, 2, figsize=(12, 5))
    for k, gfa in enumerate(GFA.gfa_names):
        color = plt.rcParams['axes.prop_cycle'].by_key()['color'][k]
        ax[1].scatter([], [], marker='o', c=color, label=gfa)
        for j, amp in enumerate(desietcimg.gfa.GFACamera.amp_names):
            marker = markers[j]
            measured = GFA.calib_data[gfa][amp]['RDNOISE']
            # Lab results are given in elec so use lab gains to convert back to ADU
            lab = GFA.lab_data[gfa][amp]['RDNOISE'] / GFA.lab_data[gfa][amp]['GAIN']
            ax[0].scatter(lab, measured, marker=marker, c=color)
            ax[1].scatter(ameisner_rdnoise[gfa][amp], measured, marker=marker, c=color)
    for j, amp in enumerate(GFA.amp_names):
        ax[1].scatter([], [], marker=markers[j], c='k', label=amp)
    xylim = (4.3, 5.5)
    for axis in ax:
        axis.plot(xylim, xylim, 'k-', zorder=-10, alpha=0.25)
        axis.set_ylabel(f'{label} Read Noise [ADU]')
        axis.set_xlim(*xylim)
        axis.set_ylim(*xylim)
    ax[1].legend(ncol=3)
    ax[0].set_xlabel('Lab Data Read Noise [ADU]')
    ax[1].set_xlabel('ameisner Read Noise [ADU]')
    plt.tight_layout()
    if save:
        plt.savefig(save)

compare_rdnoise(save=str(plotdir / 'rdnoise_compare.png'))<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.15812491458247915, 'fraction_numerical': 0.035260352603526036, 'mean_word_length': 3.9445945945945944, 'pattern_counts': {'":': 0, '<': 44, '<?xml version=': 0, '>': 44, 'https://': 3, 'lorem ipsum': 0, 'www.': 1, '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': '11', 'dir': 'jupyter-structured-clean-dedup', 'id': '4884', 'license': 'BSD-3-Clause', 'n_tokens_mistral': 3047, 'n_tokens_neox': 2830, 'n_words': 568, 'path': 'analysis/gfa/GFA-Zero-Calibration.ipynb', 'repo_name': 'tanveerkarim/desicmx', 'size': '7324'}