{ "cells": [ { "cell_type": "raw", "id": "bba652ce", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Level-1 Dataset\n", "===============\n", "The ESIS Level-1 data product includes applying the following steps to the Level-0 data:\n", "\n", "* Bias subtraction\n", "* Non-active (NAP) pixel removal\n", "* Conversion to electrons\n", "* Dark frame subtraction\n", "* Cosmic ray removal\n", "\n", "These steps consist of all the effects that need to be accounted for before we transform\n", "into skyplane coordinates.\n", "All of these steps are packaged into the :func:`esis.flights.f1.data.level_1` function,\n", "but this notebook starts from the Level-0 dataset and examines each indivdual step." ] }, { "cell_type": "code", "id": "6098fc9709ad1b15", "metadata": {}, "source": [ "import IPython.display\n", "import numpy as np\n", "import scipy.stats\n", "import matplotlib.pyplot as plt\n", "import astropy.units as u\n", "import astropy.visualization\n", "import named_arrays as na\n", "import esis" ], "outputs": [], "execution_count": null }, { "cell_type": "code", "id": "2fa28361", "metadata": {}, "source": [ "plt.rcParams[\"animation.embed_limit\"] = 50" ], "outputs": [], "execution_count": null }, { "cell_type": "code", "id": "abb9b243", "metadata": {}, "source": [ "astropy.visualization.quantity_support();" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "9bee6912", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Bias Subtraction\n", "---------------------------\n", "The bias (or pedestal) is an offset charge used to keep all the outputs positive.\n", "The ESIS senors have 50 blank columns and 2 overscan columns on each tap which are not sensitive to\n", "light and can be used to estimate the bias.\n", "In :doc:`msfc_ccd:reports/bias`,\n", "we found that using only the overscan columns resulted in the\n", "best bias subtraction.\n", "To check this result, we'll start by loading the Level-0 dataset." ] }, { "cell_type": "code", "id": "87473ea1f4ed5335", "metadata": {}, "source": [ "level_0 = esis.flights.f1.data.level_0()" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "019233c7", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Estimate and remove the bias from all the images using the :attr:`~esis.data.Level_0.unbiased` property." ] }, { "cell_type": "code", "id": "10818edb83ed5c84", "metadata": {}, "source": [ "unbiased = level_0.unbiased" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "818c23c9", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "For brevity, we'll also save the names of the relevant logical axes to local variables." ] }, { "cell_type": "code", "id": "18aee5c9", "metadata": {}, "source": [ "axis_time = level_0.axis_time\n", "axis_channel = level_0.axis_channel\n", "axis_x = level_0.axis_x\n", "axis_y = level_0.axis_y\n", "axis_xy = (axis_x, axis_y)\n", "axis_txy = (axis_time, axis_x, axis_y)" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "0e21f819", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Plot the average bias-subtracted signal as a function of time\n", "and the intervals containing the light and dark frames." ] }, { "cell_type": "code", "id": "2c2bcba6335f55a2", "metadata": {}, "source": [ "time = level_0.inputs.time\n", "time = time.replace(ndarray=time.ndarray.datetime)" ], "outputs": [], "execution_count": null }, { "cell_type": "code", "id": "e4cb1e6595e0757c", "metadata": {}, "source": [ "fig, ax = plt.subplots(\n", " figsize=(6, 2.5),\n", " constrained_layout=True,\n", ")\n", "na.plt.plot(\n", " time,\n", " unbiased.outputs.mean(axis_xy),\n", " axis=axis_time,\n", " label=unbiased.channel,\n", " drawstyle='steps-mid',\n", ")\n", "ax.axvspan(\n", " xmin=unbiased.lights.inputs.time_start.ndarray.min().datetime,\n", " xmax=unbiased.lights.inputs.time_end.ndarray.max().datetime,\n", " alpha=0.3,\n", " label=\"lights\",\n", ")\n", "ax.axvspan(\n", " xmin=unbiased.darks_up.inputs.time_start.ndarray.min().datetime,\n", " xmax=unbiased.darks_up.inputs.time_end.ndarray.max().datetime,\n", " alpha=0.3,\n", " label=\"darks\",\n", " color=\"gray\",\n", ")\n", "ax.axvspan(\n", " xmin=unbiased.darks_down.inputs.time_start.ndarray.min().datetime,\n", " xmax=unbiased.darks_down.inputs.time_end.ndarray.max().datetime,\n", " alpha=0.3,\n", " color=\"gray\",\n", ")\n", "ax.set_xlabel(\"time (UTC)\")\n", "ax.set_ylabel(\"mean intensity (DN)\")\n", "ax.legend();" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "197ff3ef", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Isolate the dark frames using the :attr:`~esis.data.Level_0.darks` property and display them as an animation." ] }, { "cell_type": "code", "id": "1e585047fe0b8172", "metadata": {}, "source": [ "unbiased.darks.to_jshtml()" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "fb313843", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Split the dark frames into separate images for each tap using the :attr:`~esis.data.Level_0.taps` property.\n", "This will allows us to compare the different taps more easily." ] }, { "cell_type": "code", "id": "160f06c2b273def4", "metadata": {}, "source": [ "taps = unbiased.darks.taps" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "fa9ea9ee", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Plot the time/row-averaged signal as a function of column for every tap of every dark frame." ] }, { "cell_type": "code", "id": "209f75cca9946af1", "metadata": {}, "source": [ "axes_tap_xy = (taps.axis_tap_x, taps.axis_tap_y)\n", "axis_tap_xy = \"tap_xy\"\n", "tap_labels = taps.label.combine_axes(axes_tap_xy, axis_tap_xy)\n", "taps = taps.combine_axes(axes_tap_xy, axis_tap_xy)\n", "\n", "fig, ax = na.plt.subplots(\n", " axis_rows=axis_channel,\n", " axis_cols=axis_tap_xy,\n", " nrows=taps.shape[axis_channel],\n", " ncols=taps.shape[axis_tap_xy],\n", " sharex=True,\n", " sharey=True,\n", " constrained_layout=True,\n", " origin=\"upper\",\n", ")\n", "na.plt.plot(\n", " taps.outputs.mean_trimmed(.01, (axis_time, axis_y)),\n", " axis=axis_x,\n", " ax=ax,\n", ")\n", "na.plt.set_ylim(\n", " bottom=-1,\n", " top=1,\n", " ax=ax,\n", ")\n", "na.plt.axvspan(\n", " xmin=0,\n", " xmax=taps.camera.sensor.num_blank,\n", " color=\"green\",\n", " alpha=0.2,\n", " ax=ax,\n", " label=\"blank\",\n", ")\n", "na.plt.axvspan(\n", " xmin=taps.num_x - taps.camera.sensor.num_overscan,\n", " xmax=taps.num_x,\n", " color=\"red\",\n", " alpha=0.2,\n", " ax=ax,\n", " label=\"overscan\",\n", ")\n", "na.plt.set_xlabel(\"columns\", ax=ax[{level_0.axis_channel: ~0}])\n", "na.plt.text(\n", " x=0.5,\n", " y=1.02,\n", " s=tap_labels[{axis_channel: 0}],\n", " ax=ax[{axis_channel: 0}],\n", " transform=na.plt.transAxes(ax[{level_0.axis_channel: 0}]),\n", " ha=\"center\",\n", " va=\"bottom\",\n", ")\n", "na.plt.text(\n", " x=1.02,\n", " y=0.5,\n", " s=level_0.channel,\n", " ax=ax[{axis_tap_xy: ~0}],\n", " transform=na.plt.transAxes(ax[{axis_tap_xy: ~0}]),\n", " ha=\"left\",\n", " va=\"center\",\n", ");" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "9a5897cb", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Cropping Non-Active Pixels\n", "--------------------------\n", "The blank and overscan columns are useful for estimating the bias,\n", "but they are not sensitive to light and need to be removed to\n", "make a continuous image.\n", "We can remove these pixels easily using the :attr:`~esis.data.Level_0.active`\n", "property." ] }, { "cell_type": "code", "id": "30c267ae4450663f", "metadata": {}, "source": [ "active = unbiased.active" ], "outputs": [], "execution_count": null }, { "cell_type": "code", "id": "0aad43b0", "metadata": {}, "source": [ "unbiased.shape" ], "outputs": [], "execution_count": null }, { "cell_type": "code", "id": "e8afd63b", "metadata": {}, "source": [ "active.shape" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "9ca2aee5", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "From the printed shapes, notice how the `detector_x` axis is smaller by 104 pixels." ] }, { "cell_type": "raw", "id": "bdc90879", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Conversion to Electrons\n", "-----------------------\n", "To interpret the images in terms of electrons,\n", "multiply the values of each image by a tap-dependent gain value \n", "using the :attr:`~esis.data.Level_0.electrons` property." ] }, { "cell_type": "code", "id": "2ba62cc5", "metadata": {}, "source": [ "electrons = active.electrons" ], "outputs": [], "execution_count": null }, { "cell_type": "code", "id": "c66f78b3", "metadata": {}, "source": [ "active.outputs.mean(axis_txy).ndarray" ], "outputs": [], "execution_count": null }, { "cell_type": "code", "id": "06a9df49", "metadata": {}, "source": [ "electrons.outputs.mean(axis_txy).ndarray" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "f7058cdf", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "By printing the means of the `active` and `electrons` variables,\n", "we can see that the units are now in terms of electrons." ] }, { "cell_type": "raw", "id": "dfe39c45", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Dark Subtraction\n", "----------------\n", "To calculate the master dark frame,\n", "we'll remove the cosmic rays from the dark frames using the\n", ":attr:`~esis.data.Level_0.despiked` property,\n", "and then take the mean along the time axis." ] }, { "metadata": {}, "cell_type": "code", "source": "darks = electrons.darks.despiked", "id": "7b482a566ba8ae28", "outputs": [], "execution_count": null }, { "metadata": {}, "cell_type": "code", "source": "dark = darks.mean(darks.axis_time)", "id": "1667176b2b7fa5e2", "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "6854be50", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Plot the dark frame for each channel to inspect" ] }, { "cell_type": "code", "id": "e02eff69", "metadata": { "scrolled": false }, "source": [ "fig, ax = na.plt.subplots(\n", " axis_rows=axis_channel,\n", " nrows=dark.shape[axis_channel],\n", " sharex=True,\n", " sharey=True,\n", " constrained_layout=True,\n", " figsize=(5, 8),\n", " origin=\"upper\",\n", ")\n", "na.plt.set_xlabel(\"detector $x$ (pix)\", ax=ax[{axis_channel: ~0}])\n", "na.plt.set_ylabel(\"detector $y$ (pix)\", ax=ax)\n", "\n", "norm = plt.Normalize(-5, 5)\n", "\n", "colorizer = plt.Colorizer(norm=norm)\n", "\n", "i = {darks.axis_time: 0}\n", "\n", "na.plt.pcolormesh(\n", " dark[i].inputs.pixel.x,\n", " dark[i].inputs.pixel.y,\n", " C=dark[i].outputs,\n", " ax=ax,\n", " colorizer=colorizer,\n", ")\n", "na.plt.text(\n", " x=0.5,\n", " y=1.01,\n", " s=dark[i].channel,\n", " transform=na.plt.transAxes(ax),\n", " ax=ax,\n", " ha=\"center\",\n", " va=\"bottom\",\n", ")\n", "na.plt.set_aspect(\"equal\", ax=ax)\n", "\n", "plt.colorbar(\n", " mappable=plt.cm.ScalarMappable(colorizer=colorizer),\n", " ax=ax.ndarray,\n", " label=f\"signal ({darks.outputs.unit:latex_inline})\",\n", ");" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "28a63ad9", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Subtract the master dark from each frame." ] }, { "cell_type": "code", "id": "bbade136", "metadata": {}, "source": [ "dark_subtracted = electrons - dark.outputs" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "585a78d7", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Also subtract the master dark from each despiked dark frame to compute the residual" ] }, { "cell_type": "code", "id": "08464613", "metadata": {}, "source": [ "residual = darks - dark.outputs" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "0b1010fe", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Create a histogram of the residual dark frames to inpect their distribution." ] }, { "cell_type": "code", "id": "94c2f8c1", "metadata": { "scrolled": false }, "source": [ "hist = na.histogram(\n", " residual.taps.outputs, \n", " axis=axis_txy,\n", " bins=dict(values=51),\n", " min=-20 * u.electron,\n", " max=+20 * u.electron,\n", " density=True,\n", ")" ], "outputs": [], "execution_count": null }, { "cell_type": "code", "id": "e39b0f54", "metadata": {}, "source": [ "fig, ax = na.plt.subplots(\n", " axis_rows=axis_channel,\n", " nrows=4,\n", " sharex=True,\n", " figsize=(7, 6),\n", " constrained_layout=True,\n", " origin=\"upper\"\n", ")\n", "na.plt.stairs(\n", " hist.inputs,\n", " hist.outputs,\n", " ax=ax,\n", " axis=\"values\",\n", " label=residual.taps.label,\n", ")\n", "na.plt.text(\n", " x=0.01,\n", " y=.98,\n", " s=residual.channel,\n", " ax=ax,\n", " transform=na.plt.transAxes(ax),\n", " va=\"top\",\n", ")\n", "ax_bottom = ax.ndarray[~0]\n", "na.plt.axvline(0, ax=ax, color=\"black\", linestyle=\"--\")\n", "na.plt.set_ylabel(\"probability density\", ax=ax)\n", "na.plt.set_xlabel(f\"residual ({hist.inputs.unit:latex_inline})\", ax=ax_bottom)\n", "ax.ndarray[0].legend();" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "cd6f65fc", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Compute the standard deviation of each channel.\n", "This is equivalent to the read noise for each channel." ] }, { "cell_type": "code", "id": "25b91dc6", "metadata": {}, "source": [ "residual.outputs.std(axis_txy).ndarray" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "00f7a57a", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Compute the standard deviation of each tap of Channel 3,\n", "since that channel appears to have one tap that's worse than the others." ] }, { "cell_type": "code", "id": "99b18ed4", "metadata": {}, "source": [ "residual[{axis_channel: 2}].taps.outputs.std(axis_txy).ndarray" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "a046aeab", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Isolate only the frames with light from the Sun since we no longer need the dark frames." ] }, { "cell_type": "code", "id": "48fc4ecf", "metadata": {}, "source": [ "lights = dark_subtracted.lights" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "40e3221f", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Cosmic Ray Removal\n", "------------------\n", "All of the images have cosmic ray spikes that must be removed before proceeding.\n", "We'll remove spikes from all the solar images using :func:`astroscrappy.detect_cosmics()`." ] }, { "cell_type": "code", "id": "8374badbba8acb80", "metadata": {}, "source": [ "despiked = lights.despiked" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "4a2c93ae", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Take the diffence of the original lights and the despiked lights to isolate the cosmic rays." ] }, { "cell_type": "code", "id": "f8d50a4d", "metadata": {}, "source": [ "spikes = lights - despiked" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "ac8493ef", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Blink the max of the original frames against the despiked frames and the identified spikes." ] }, { "cell_type": "code", "id": "fb4d5988", "metadata": {}, "source": [ "blink = na.stack(\n", " arrays=[\n", " lights.max(axis_time),\n", " despiked.max(axis_time),\n", " spikes.max(axis_time),\n", " ],\n", " axis=\"blink\",\n", ")[{axis_time: 0}]\n", "blink.axis_time = \"blink\"\n", "blink.to_jshtml(fps=1)" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "9d4a66db", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Compute the `Pearson correlation coefficient `_\n", "between the spike mask and the despiked images.\n", "This gives us an understanding of how many false positives may have been identified." ] }, { "cell_type": "code", "id": "67550b13", "metadata": {}, "source": [ "where_spike = spikes.outputs >= .5 * u.electron" ], "outputs": [], "execution_count": null }, { "cell_type": "code", "id": "51302e0b", "metadata": {}, "source": [ "scipy.stats.pearsonr(\n", " x=where_spike.ndarray,\n", " y=despiked.outputs.ndarray,\n", " axis=None,\n", ").statistic" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "1cda4f3f", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "As we can see, the spikes are slightly correlated with the image.\n", "This means that there are some small fraction of pixels wrongly identified as spikes." ] }, { "cell_type": "raw", "id": "5e090187", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "The last thing to do is to convert the output to an instance of :class:`esis.data.Level_1`.\n", "This class will have it's own specified methods for dealing with this type of data." ] }, { "cell_type": "code", "id": "baa1562d", "metadata": {}, "source": [ "level_1 = esis.data.Level_1(\n", " inputs=despiked.inputs,\n", " outputs=despiked.outputs,\n", " axis_time=axis_time,\n", " axis_channel=axis_channel,\n", " axis_x=axis_x,\n", " axis_y=axis_y,\n", ")" ], "outputs": [], "execution_count": null }, { "cell_type": "raw", "id": "d5537381", "metadata": { "raw_mimetype": "text/restructuredtext" }, "source": [ "Plot the final Level-1 images for inspection." ] }, { "cell_type": "code", "id": "bfa8e0ec", "metadata": {}, "source": [ "level_1.to_jshtml()" ], "outputs": [], "execution_count": null } ], "metadata": { "celltoolbar": "Raw Cell Format", "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.11.1" } }, "nbformat": 4, "nbformat_minor": 5 }