import os
import gc
import subprocess
import numpy as np
from scipy.interpolate import RectBivariateSpline, splrep, splev
from astropy.io import fits
from astropy.table import Table
from astropy.modeling.models import Linear1D
from astropy.modeling.fitting import LinearLSQFitter
from astroquery.mast import Observations

import matplotlib
%matplotlib inline
from matplotlib import pyplot as plt

# Some plotting defaults:
matplotlib.rcParams['figure.figsize']      = [8, 8]  # inches
matplotlib.rcParams['image.origin']        = 'lower'
matplotlib.rcParams['image.cmap']          = 'inferno'  # https://matplotlib.org/stable/tutorials/colors/colormaps.html
matplotlib.rcParams['image.aspect']        = 'auto'
matplotlib.rcParams['image.interpolation'] = 'none'

# Dataset shown here:
# https://www.planetarylightshow.com/jupiter/prop_16193/v54-oef454vrq.html

obs = Observations.query_criteria(obs_id='OEF454VRQ')
pl = Observations.get_product_list(obs)

# Only download time-tag data product:
pl = pl[[x == 'TAG' for x in pl['productSubGroupDescription']]]

downloads = Observations.download_products(pl, mrp_only=False)  # 213 MB
downloads

INFO: Found cached file ./mastDownload/HST/oef454vrq/oef454vrq_tag.fits with expected size 223032960. [astroquery.query]

filename = downloads[0]['Local Path']
rootname = os.path.basename(filename).rsplit('_',1)[0]
filename

np.str_('./mastDownload/HST/oef454vrq/oef454vrq_tag.fits')

# STIS tag files contain these extensions:
fits.info(filename)

Filename: ./mastDownload/HST/oef454vrq/oef454vrq_tag.fits
No.    Name      Ver    Type      Cards   Dimensions   Format
  0  PRIMARY       1 PrimaryHDU     216   ()      
  1  EVENTS        1 BinTableHDU    127   22299384R x 4C   [1J, 1I, 1I, 1I]   
  2  GTI           1 BinTableHDU     22   1R x 2C   [1D, 1D]

t = Table.read(filename, hdu='EVENTS')
t

t['AXIS1'].min(), t['AXIS1'].max(), t['AXIS2'].min(), t['AXIS2'].max()

(np.int16(1), np.int16(2048), np.int16(1), np.int16(2048))

gti = Table.read(filename, hdu='GTI')
gti

max_time = int(np.ceil(gti[-1]['STOP']))
max_time

2351

fig, ax = plt.subplots(1, 1)
fig.set_size_inches(15, 5)

_ = ax.hist(t['TIME'], range=[0, max_time], bins=max_time, histtype='step')

ax.set_xlabel('Time (s)')
ax.set_ylabel('Total Flux (counts / s)')
ax.grid(True, alpha=0.2);

fig, axes = plt.subplots(1, 2)
fig.set_size_inches(15, 7)

for ax in axes:
    _ = ax.hist2d(t['AXIS1'], t['AXIS2'],
                  bins=[2048, 2048],
                  range=[[1, 2049], [1, 2049]],
                  vmin=0, vmax=75)  # counts/pixel used to scale the colormap

    ax.axis('equal')
    ax.set_xlabel('AXIS1 (high-res pix)')
    ax.set_ylabel('AXIS2 (high-res pix)')

axes[0].set_title('Native 2048x2048 high-res format')
axes[1].set_title('Zoomed in')

axes[0].set_xlim(0, 2048)
axes[0].set_ylim(0, 2048)
axes[1].set_xlim(1100, 1250)
axes[1].set_ylim(700, 850);

# Down-sample to low-res format:

fig, axes = plt.subplots(1, 2)
fig.set_size_inches(15, 7)

for ax in axes:
    _ = ax.hist2d(t['AXIS1'], t['AXIS2'],
                  bins=[1024, 1024],
                  range=[[1, 2049], [1, 2049]],
                  vmin=0, vmax=250)  # counts/pixel used to scale the colormap

    ax.axis('equal')
    ax.set_xlabel('AXIS1 (high-res pix)')
    ax.set_ylabel('AXIS2 (high-res pix)')

axes[0].set_title('Down-sampled to 1024x1024 low-res format')
axes[1].set_title('Zoomed in')

axes[0].set_xlim(0, 2048)
axes[0].set_ylim(0, 2048)
axes[1].set_xlim(1100, 1250)
axes[1].set_ylim(700, 850);

# Download a dataset useful for fitting the grid wire correction:
# https://www.planetarylightshow.com/europa/prop_14930/v02-odf702f3q.html

obs = Observations.query_criteria(obs_id='ODF702F3Q')
pl = Observations.get_product_list(obs)
# Only download 2D flat-fielded data product:
pl = pl[[x == 'FLT' for x in pl['productSubGroupDescription']]]
downloads = Observations.download_products(pl, mrp_only=False)
filename_gridwire = downloads[0]['Local Path']
print(filename_gridwire)

# An observation of Jupiter with relatively flat source flux near the grid wire:
gridwire = fits.getdata(filename_gridwire, ext=1)

# Manually-derived grid wire path (not shown here, but lots of clicking):
x_manual = np.array([
    11.973965, 14.053856, 54.265094, 84.076874, 86.850063, 119.43503,
    131.22108, 154.79319, 165.19265, 194.31113, 209.56367, 237.29556,
    262.25425, 293.45263, 321.87781, 335.74376, 375.26170, 405.76678,
    416.85953, 454.29758, 474.40320, 493.81552, 531.25357, 563.83854,
    572.85140, 616.52913, 658.82025, 698.33820, 751.02878, 796.09310,
    840.46412, 887.60833, 936.13914, 956.93805, 981.89675, 1003.3890,
    1020.0281, 1024.8812,])
y_manual = np.array([
    496.52994, 496.44585, 496.76175, 497.01404, 496.92994, 497.35043,
    497.85500, 498.52778, 499.03236, 499.87332, 499.95742, 500.79839,
    501.30296, 501.89164, 502.14393, 502.73261, 503.57357, 503.82586,
    504.07815, 505.67599, 506.09647, 506.68515, 507.94659, 508.70346,
    509.12395, 509.54443, 510.04901, 510.80588, 511.73094, 511.98323,
    512.31962, 513.07649, 514.08565, 515.34709, 516.35625, 518.29047,
    519.29963, 520.14060,])

fig, axes = plt.subplots(5, 1)
fig.set_size_inches(15, 5*len(axes))
fig.subplots_adjust(hspace=0.25)

axes[0].imshow(gridwire, vmin=0, vmax=100)
axes[0].set_title('Jupiter Observation ODF702F3Q, Grid Wire')

# Fit a smoothed cubic spline to the manually-derived coordinates along the grid wire path:
tck = splrep(x_manual, y_manual, k=3, s=1)
x_c = np.arange(0, 1024, 0.1, dtype=float)
axes[0].plot(x_c, splev(x_c, tck, der=0), color='c', alpha=0.5,
             label='smoothed spline of manually-fit positions')
axes[0].legend(loc='best')
#axes[0].plot(x_manual, y_manual, 'wx', ms=2)

# Remove the path of the grid wire to derive an average profile for all columns:
centered = splev(np.arange(1024, dtype=float), tck, der=0) - splev([512.], tck, der=0)[0]
rectified = np.zeros_like(gridwire, dtype=float)
for col in range(1024):
    new_coords = np.arange(1024) + centered[col]
    # Linear interpolate each column to rectify the gridwire:
    rectified[:, col] = np.interp(new_coords, np.arange(1024), gridwire[:,col])
# Normalize each column near the grid wire:
rectified /= rectified[475:540, :].mean(axis=0)

axes[1].imshow(rectified, vmin=0, vmax=1.5)
axes[1].set_title('Rectified and Normalized')

# Fit background of the profile, excluding left & right edges:
x_c = np.arange(1024)
g = ((x_c >= 475) & (x_c < 495)) | ((x_c >= 525) & (x_c < 540))
model = Linear1D(slope=0., intercept=1.)
fitter = LinearLSQFitter()
fit = fitter(model, np.arange(1024)[g], rectified[:, 20:970].mean(axis=1)[g])
# Extract an empirical profile, excluding the far left & right edges:
profile_raw = rectified[495:520, 15:-15].mean(axis=1)
profile = np.ones((1024), dtype=float)
# Renormalize the 1D profile using a linear fit to data on the edges:
# This helps avoid a discontinuity with the surrounding data.
profile[495:520] = profile_raw / fit(np.arange(495, 520))

axes[2].plot(np.arange(495, 520), profile_raw, label='profile')
axes[2].plot(profile, label='renormalized profile')
axes[2].set_title('Empirical Profile')
axes[2].set_xlim(485, 530)
axes[2].set_xlabel('Y')
axes[2].set_ylabel('Normalized Flux')
axes[2].legend(loc='best')
axes[2].grid(True, alpha=0.2)

# Apply the 1D profile back along the 2D grid wire path:
grid_flat = np.ones_like(gridwire, dtype=float)
for col in range(1024):
    new_coords = np.arange(1024) + centered[col]
    # Interpolate in the opposite direction to reintroduce the gridwire path:
    grid_flat[:,col] = np.interp(np.arange(1024), new_coords, profile)

# Derived grid wire flat:
axes[3].imshow(grid_flat, vmin=0.7, vmax=1.3, cmap='RdBu_r')
axes[3].set_title('2D Grid Wire Flat')

# Application of the grid wire flat to its source data:
axes[4].imshow(gridwire / grid_flat, vmin=0, vmax=100)
axes[4].set_title('Observation / Grid Wire Flat')

for ax in list(axes[:2]) + list(axes[3:]):
    ax.set_xlim(-0.5, 1024.5)
    ax.set_ylim(465, 545)
    ax.set_xlabel('X')
    ax.set_ylabel('Y')

INFO: Found cached file ./mastDownload/HST/odf702f3q/odf702f3q_flt.fits with expected size 10532160. [astroquery.query]
./mastDownload/HST/odf702f3q/odf702f3q_flt.fits

# Download STIS flat reference files, as specified in CRDS:
# https://hst-crds.stsci.edu
# (Or as specified in the ext=0 FITS headers, assuming they are up to date)

# pixel-to-pixel flat (PFLAT) for F25SRF2 aperture:
Observations.download_file(uri='mast:HST/reference/mbj1658ao_pfl.fits')
pflat = fits.getdata('mbj1658ao_pfl.fits', ext=1, ignore_missing_simple=True)

# low-order flat (LFLAT) for F25SRF2 aperture:
Observations.download_file(uri='mast:HST/reference/iaf1723io_lfl.fits')
lflat = fits.getdata('iaf1723io_lfl.fits', ext=1)

INFO: Found cached file /Users/.../mbj1658ao_pfl.fits with expected size 41976000. [astroquery.query]
INFO: Found cached file /Users/.../iaf1723io_lfl.fits with expected size 673920. [astroquery.query]

def rebin(arr, new_shape, method='sum'):
    '''Rebin 2D data down by integer steps.
    Derived from https://scipython.com/blog/binning-a-2d-array-in-numpy
    '''
    shape = (new_shape[0], arr.shape[0] // new_shape[0],
             new_shape[1], arr.shape[1] // new_shape[1])
    if method == 'sum':
        return arr.reshape(shape).sum(-1).sum(1)
    elif method == 'mean':
        return arr.reshape(shape).mean(-1).mean(1)
    else:
        raise IOError('"method" unrecognized:  {}'.format(method))

fig, axes = plt.subplots(2, 2)
fig.set_size_inches(15, 7.5*axes.shape[0])
axes = axes.ravel()

# Interpolate LFLAT from 256x256 to 2048x2048:
lflat_spline = RectBivariateSpline(
    np.arange(2048//256//2, 2048, 2048//256), np.arange(2048//256//2, 2048, 2048//256),
    lflat,
    bbox=[0, 2047, 0, 2047])
lflat_high_res = lflat_spline(np.arange(0, 2048), np.arange(0, 2048))

# Interpolate gridflat from 1024x1024 to 2048x2048:
grid_flat_spline = RectBivariateSpline(
    np.arange(2048//1024//2, 2048, 2048//1024), np.arange(2048//1024//2, 2048, 2048//1024),
    grid_flat,
    bbox=[0, 2047, 0, 2047])
grid_flat_high_res = grid_flat_spline(np.arange(0, 2048), np.arange(0, 2048))

# Multiply all flats together in 2048x2048 format:
flat = lflat_high_res * pflat * grid_flat_high_res
# Bin back down to 1024x1024 format:
flat_low_res = rebin(flat, (1024, 1024), method='mean')

axes[0].imshow(lflat_high_res,     vmin=0.5, vmax=1.5, cmap='RdBu_r')
axes[1].imshow(grid_flat_high_res, vmin=0.5, vmax=1.5, cmap='RdBu_r')
axes[2].imshow(flat,               vmin=0.5, vmax=1.5, cmap='RdBu_r')
axes[3].imshow(flat_low_res,       vmin=0.5, vmax=1.5, cmap='RdBu_r')

axes[0].set_title('Low-Order Flat')
axes[1].set_title('Gridwire Flat')
axes[2].set_title('Combined Flat (LFLAT * PFLAT * Gridwire Flat)')
axes[3].set_title('Combined Flat, Down-Sampled to 1024x1024');

# Note that AXIS1 & AXIS2 are given as one-indexed coordinates.

# Sample both the 2048x2048 and 1024x1024 flats for comparison:
t['flat_weight'] = flat[t['AXIS2']-1, t['AXIS1']-1]
t['flat_low_res_weight'] = flat_low_res[(t['AXIS2']-1)//2, (t['AXIS1']-1)//2]

fig, ax = plt.subplots(1, 1)
fig.set_size_inches(15, 5)

_ = ax.hist(t['flat_weight'],         bins=1000, range=[0, 3], histtype='step', alpha=0.8, label='High Res Flat')
_ = ax.hist(t['flat_low_res_weight'], bins=1000, range=[0, 3], histtype='step', alpha=0.8, label='Low Res Flat')

ax.legend(loc='best')
ax.set_xlabel(f"Weights of {rootname} Sampled from Combined Flats")
ax.set_ylabel(f"Photon Count (within {3/1000} cts bin)")
ax.grid(True, alpha=0.2)

t[:10]

fig, axes = plt.subplots(2, 3)
fig.set_size_inches(15, 10)

w = (t['TIME'] >= 0.) & (t['TIME'] < 100.)

# High-res uncorrected for flat-fielding:
_ = axes[0,0].hist2d(t['AXIS1'][w], t['AXIS2'][w],
                     bins=[2048, 2048],
                     range=[[1, 2049], [1, 2049]],
                     vmin=0, vmax=30/4,)

# Sampling the flat in high-res space:
_ = axes[0,1].hist2d(t['AXIS1'][w], t['AXIS2'][w],
                     bins=[2048, 2048],
                     range=[[1, 2049], [1, 2049]],
                     vmin=0, vmax=30/4,
                     weights=1./t['flat_weight'][w],)

# Sampling the flat in high-res space and then binning down the corrected product:
h, _, __ = np.histogram2d(t['AXIS2'][w], t['AXIS1'][w],
               bins=[2048, 2048],
               range=[[1, 2049], [1, 2049]],
               weights=1./t['flat_weight'][w],)
h_binned = rebin(h, (1024, 1024), method='sum')
_ = axes[0,2].imshow(h_binned, vmin=0, vmax=30)

# Low-res uncorrected for flat-fielding:
_ = axes[1,0].hist2d(t['AXIS1'][w], t['AXIS2'][w],
                     bins=[1024, 1024],
                     range=[[1, 2049], [1, 2049]],
                     vmin=0, vmax=30,)

# Sampling the flat in low-res space:
_ = axes[1,1].hist2d(t['AXIS1'][w], t['AXIS2'][w],
                     bins=[1024, 1024],
                     range=[[1, 2049], [1, 2049]],
                     vmin=0, vmax=30,
                     weights=1./t['flat_low_res_weight'][w],)

axes[1,2].set_visible(False)
axes[0,0].set_title('2048x2048')
axes[0,1].set_title('2048x2048, High-Res Flat-Fielded')
axes[0,2].set_title('High-Res Flat-Fielded, Rebinned to 1024x1024')
axes[1,0].set_title('1024x1024')
axes[1,1].set_title('1024x1024, Low-Res Flat-Fielded')

for ax in axes[:, :2].ravel():
    ax.set_xlim(900, 1100)  # high-res space
    ax.set_ylim(450,  650)  # high-res space
axes[0,2].set_xlim(900//2 - 1, 1100//2 - 1)   # low-res space
axes[0,2].set_ylim(450//2 - 1,  650//2 - 1);  # low-res space

# In full:
_ = plt.hist2d(t['AXIS1'][w], t['AXIS2'][w],
               bins=[1024, 1024],
               range=[[1, 2049], [1, 2049]],
               vmin=0, vmax=15,
               weights=1./t['flat_low_res_weight'][w],)

# Break a time-tag observation up into multiple frames:

step = 100  # s
starts = range(0, max_time, step)
stops = range(step, max_time + step, step)

fig, axes = plt.subplots(6, len(starts)//5)
fig.set_size_inches(15, 5*axes.shape[1])
axes = axes.ravel()  # Turn into 1D array

for i, (start, stop, ax) in enumerate(zip(starts, stops, axes)):
    w = (t['TIME'] >= start) & (t['TIME'] < stop)
    ax.hist2d(t['AXIS1'][w], t['AXIS2'][w],
              bins=[1024, 1024],
              range=[[1, 2049], [1, 2049]],
              vmin=0, vmax=15,
              weights=1./t['flat_low_res_weight'][w],)
    ax.set_title(f"{t['TIME'][w].min():.0f} – {t['TIME'][w].max():.0f} s")
    ax.axis('equal')
    ax.axis('off')

def determine_time_weights(timeline, starts, stops):
    if len(timeline) == 0:
        raise ValueError('GTI extension is empty!')

    starts = np.array(starts)
    stops = np.array(stops)

    tl_start, tl_stop = timeline['START'], timeline['STOP']
    time_delta = 100  # s

    step = 0.1
    on_time_x = np.arange(0, np.max(tl_stop) + step, step)
    on_time = np.zeros_like(on_time_x, dtype=int)
    for interval_start, interval_stop in zip(tl_start, tl_stop):
        on_time[(on_time_x >= interval_start) & (on_time_x < interval_stop)] = 1

    w = np.zeros_like(starts)
    for i, (time_start, time_stop) in enumerate(zip(starts, stops)):
        w[i] = np.sum(on_time[(on_time_x >= time_start) & (on_time_x < time_stop)]) * step

    return w / (stops - starts)


step = 100  # s
starts = range(0, max_time, step)
stops = range(step, max_time + step, step)
time_weights = determine_time_weights(gti, starts, stops)
print('Time weights for each frame:\n', time_weights)

fig, axes = plt.subplots(6, len(starts)//5)
fig.set_size_inches(15, 5*axes.shape[1])
axes = axes.ravel()  # Turn into 1D array

for i, (start, stop, time_weight, ax) in enumerate(zip(starts, stops, time_weights, axes)):
    w = (t['TIME'] >= start) & (t['TIME'] < stop)
    ax.hist2d(t['AXIS1'][w], t['AXIS2'][w],
              bins=[1024, 1024],
              range=[[1, 2049], [1, 2049]],
              vmin=0, vmax=15,
              weights=1./(time_weight * t['flat_low_res_weight'][w]),)
    ax.set_title(f"{t['TIME'][w].min():.0f} – {t['TIME'][w].max():.0f} s")
    ax.axis('equal')
    ax.axis('off')

Time weights for each frame:
 [0.99 1.   1.   1.   1.   1.   1.   1.   1.   1.   1.   1.   1.   1.
 1.   1.   1.   1.   1.   1.   1.   1.   1.   0.5 ]

# Zooming in:

fig, axes = plt.subplots(6, len(starts)//5)
fig.set_size_inches(15, 5*axes.shape[1])
axes = axes.ravel()  # Turn into 1D array

for i, (start, stop, time_weight, ax) in enumerate(zip(starts, stops, time_weights, axes)):
    w = (t['TIME'] >= start) & (t['TIME'] < stop)
    ax.hist2d(t['AXIS1'][w], t['AXIS2'][w],
              bins=[1024, 1024],
              range=[[1, 2049], [1, 2049]],
              vmin=0, vmax=15,
              weights=1./(time_weight * t['flat_low_res_weight'][w]),)
    ax.set_title(f"{t['TIME'][w].min():.0f} – {t['TIME'][w].max():.0f} s")
    ax.axis('equal')
    ax.axis('off')

    # Zoom in:
    ax.set_xlim(900, 1400)
    ax.set_ylim(350,  850)

step = 100  # s
starts = range(0, max_time, step)
stops = range(step, max_time + step, step)
time_weights = determine_time_weights(gti, starts, stops)

#ax.imshow(h[0], interpolation='none', vmin=0, vmax=vmax * bin**2, cmap='inferno', 
#    aspect='equal', origin='lower')

outnames = []
for i, (start, stop, time_weight) in enumerate(zip(starts, stops, time_weights)):
    # Create and save image:
    fig = plt.figure(frameon=False, figsize=(1024/96, 1024/96), dpi=96)
    ax = plt.Axes(fig, [0., 0., 1., 1.])
    ax.set_axis_off()
    fig.add_axes(ax)

    w = (t['TIME'] >= start) & (t['TIME'] < stop)
    ax.hist2d(t['AXIS1'][w], t['AXIS2'][w],
              bins=[1024, 1024],
              range=[[1, 2049], [1, 2049]],
              vmin=0, vmax=15,
              weights=1./(time_weight * t['flat_low_res_weight'][w]),)
    ax.axis('equal')
    ax.axis('off')

    outname = f'{rootname}_{i:03.0f}.png'
    fig.savefig(outname, dpi=96)
    outnames.append(outname)

    # Clean up memory:
    plt.close(fig)
    del fig
    if i % 5 == 0:
        gc.collect()
    
    print(f'{i + 1: 3.0f} / {len(starts)} frames written')

outnames

  1 / 24 frames written
  2 / 24 frames written
  3 / 24 frames written
  4 / 24 frames written
  5 / 24 frames written
  6 / 24 frames written
  7 / 24 frames written
  8 / 24 frames written
  9 / 24 frames written
 10 / 24 frames written
 11 / 24 frames written
 12 / 24 frames written
 13 / 24 frames written
 14 / 24 frames written
 15 / 24 frames written
 16 / 24 frames written
 17 / 24 frames written
 18 / 24 frames written
 19 / 24 frames written
 20 / 24 frames written
 21 / 24 frames written
 22 / 24 frames written
 23 / 24 frames written
 24 / 24 frames written

['oef454vrq_000.png',
 'oef454vrq_001.png',
 'oef454vrq_002.png',
 'oef454vrq_003.png',
 'oef454vrq_004.png',
 'oef454vrq_005.png',
 'oef454vrq_006.png',
 'oef454vrq_007.png',
 'oef454vrq_008.png',
 'oef454vrq_009.png',
 'oef454vrq_010.png',
 'oef454vrq_011.png',
 'oef454vrq_012.png',
 'oef454vrq_013.png',
 'oef454vrq_014.png',
 'oef454vrq_015.png',
 'oef454vrq_016.png',
 'oef454vrq_017.png',
 'oef454vrq_018.png',
 'oef454vrq_019.png',
 'oef454vrq_020.png',
 'oef454vrq_021.png',
 'oef454vrq_022.png',
 'oef454vrq_023.png']

def pngs_to_animated_gif(images, gif_outname, erase_pngs=True):
    '''Convert PNGs to an animated gif using ImageMagick.

    Parameters
    ----------
    images: iterable of str
        input PNG names

    gif_outname: str
        output name for animated gif

    erase_pngs: bool, default=True
        Delete the individual frames after generating the animated GIF?
    '''
    # Check inputs:
    try:
        subprocess.check_call(['which', 'convert'], stdout=subprocess.DEVNULL)
    except subprocess.CalledProcessError:
        raise ImportError('Please install ImageMagick:  https://imagemagick.org')

    if not gif_outname.lower().endswith('.gif'):
        raise ValueError('Output name must end in ".gif".')

    c = ['convert']
    c.extend(['-delay', '8'])
    c.extend(images)
    c.extend(['-loop', '0', gif_outname])

    # Call ImageMagick with appropriate parameters:
    subprocess.check_call(c, stdout=subprocess.DEVNULL)

    full_outname = os.path.abspath(gif_outname)
    if not os.access(full_outname, os.F_OK):
        raise IOError('Output file not written!')

    # Optionally remove the intermediate PNG files:
    if erase_pngs:
        for im in images:
            if os.access(im, os.F_OK):
                os.remove(im)

    return f"file://{full_outname}"


url = pngs_to_animated_gif(outnames, f"{rootname}.gif")
print(url)

file:///Users/.../oef454vrq.gif

Local Path	Status	Message	URL
str47	str8	object	object
./mastDownload/HST/oef454vrq/oef454vrq_tag.fits	COMPLETE	None	None

TIME	AXIS1	AXIS2	DETAXIS1
seconds	pixels	pixels	pixels
float64	int16	int16	int16
0.030625	929	553	929
0.03075	940	603	940
0.03075	1195	479	1195
0.031125	1658	595	1658
0.031125	1008	751	1008
0.031375	1426	713	1426
0.031375	1223	546	1223
0.0315	1643	203	1643
0.0315	836	304	836
...	...	...	...
2350.229625	155	24	155
2350.22975	1170	437	1170
2350.22975	1299	335	1299
2350.23	265	326	265
2350.23	1650	511	1650
2350.230125	1657	482	1657
2350.23025	1656	509	1656
2350.230375	1105	602	1105
2350.2305	1460	778	1460

START	STOP
seconds	seconds
float64	float64
0.030625	2350.2305

TIME	AXIS1	AXIS2	DETAXIS1	flat_weight	flat_low_res_weight
seconds	pixels	pixels	pixels
float64	int16	int16	int16	float64	float64
0.030625	929	553	929	0.8289870918519365	1.0290947868644418
0.03075	940	603	940	0.9701029209392998	1.0238008955158917
0.03075	1195	479	1195	1.0359367410774785	0.9940077875737259
0.031125	1658	595	1658	0.6521971189886022	0.9805573633767817
0.031125	1008	751	1008	1.0230089127952833	1.0316010004112917
0.031375	1426	713	1426	0.854371779665402	0.9811807111461933
0.031375	1223	546	1223	1.4374209428095486	1.009872424652115
0.0315	1643	203	1643	0.9967315362913256	0.9983080644330082
0.0315	836	304	836	0.684750936501488	0.9495118281281141
0.031625	1242	727	1242	0.8189160812808816	1.0657477948730696

Processing Hubble/STIS Aurora Data¶

Table of Contents

Download Data¶

Read in FITS Data¶

Time-Series¶

2D Image¶

Flat-Fielding¶

Grid Wire Shadow¶

Pipeline Flat Files¶

Sample the Flat-Field¶

Slice into Time Chunks¶

Effective Exposure Time Correction¶

Animated GIF¶

Steps Not Shown¶

Acknowledgements¶