Some plots of Arcus ray-trace outputs¶

In [1]:

from nbtemplate import display_header, get_path
display_header('Plots.ipynb', status='Arcus 2020')

Revision status: Arcus 2020¶

Last revision in version control:

Author: Hans Moritz Guenther hgunther@mit.edu
commit cec12fa383d7ab490123c6ca02cd0a13286196af
Date: Thu Oct 22 15:35:06 2020 -0400

This document is git version controlled. The repository is available at https://github.com/hamogu/arcus. See git commit log for full revision history.

Code was last run with:

MARXS ray-trace code version 1.2.dev769 (commit hash: 45bf7363f2 from 2020-11-22)
ARCUS python code version 0.2.dev186 (commit hash: 8c6308a4ec from 2020-10-19)
ARCUS CALDB version hash: b'b1a9c0f' (commited on b'2020-10-22 ')

Purpose¶

The goal of this notebook is to show a few example outputs from ray-trace simulations to give a feel of how Arcus data will look. The goal of these plots is not to make the best possible prediction including all sorts of effects like sky background, instrumental background, bad pixels, read-out streaks, and charge transfer inefficiences - that sort of work is better done with the SIXTE simulation tool. Instead, the goal is to present a few results that show off the layout of the CCDs and the shape of the spectral traces.

Emission line source¶

For this example, I chose the spectrum of EQ Peg A (an active star) with the parameters from Liefke et al. (2008). In practice, there is a close companion EQ Peg B and the two spectra will overlap in ARCUS, but for the purposes of this example, EQ Peg A alone will be good enough to represent an emission line spectrum.

The code below reads the input spectrum from table of energy and flux density that I generated with Sherpa (could have used XSPEC or ISIS, too).

In [2]:

import os
import numpy as np
from matplotlib import pyplot as plt
from astropy.table import Table
from matplotlib import ticker

%matplotlib inline

In [3]:

obs = Table.read(os.path.join(get_path('rays'), 'EQPegA100ks_evt.fits'))

In [4]:

from arcus.arcus import defaultconf, DetCamera
det = DetCamera(defaultconf)

In [5]:

hist, xedges, yedges = np.histogram2d(obs['proj_x'], obs['proj_y'], weights=obs['probability'],
                            bins=[np.arange(-700, 700, 2), np.arange(-15, 15.01, .2)])

from astropy.visualization import (MinMaxInterval, AsymmetricPercentileInterval, SqrtStretch,
                                   ImageNormalize)
import matplotlib.patches as mpatches
from matplotlib.collections import PatchCollection

fig, (ax) = plt.subplots(figsize=(12, 4))

# Create an ImageNormalize object
norm = ImageNormalize(hist, interval=AsymmetricPercentileInterval(0, 99.9),
                      stretch=SqrtStretch())


cax = ax.imshow(hist.T, extent=[xedges.min(), xedges.max(), yedges.min(), yedges.max()],
                origin='lower', aspect='auto', 
                cmap=plt.get_cmap('inferno'),
                norm=norm)

for ccd in det.elements:
    # add a rectangle
    width = ccd.geometry['v_y'][0]
    height = ccd.geometry['v_z'][1]
    ax.add_patch(mpatches.Rectangle(ccd.geometry['center'][:2] - [width, height], 
                              width * 2, height * 2, fill=False, ec='w', fc="none"))

ax.set_xlabel("dispersion coordinate [mm]")
ax.set_ylabel('cross-dispersion [mm]')
ax.set_title('Detector image for emission line source\nNote vastly different scale in x and y axis!')


ax.annotate('zeroth order',
            xy=(302.5, 10), xycoords='data',
            xytext=(-70, 13), textcoords='data',
            bbox=dict(boxstyle="round", fc="0.8"),
            arrowprops=dict(arrowstyle="->",
                            connectionstyle="angle,angleA=0,angleB=90,rad=10",
                            color='w'),
           )

ax.annotate('zeroth order',
            xy=(-297.5, 5), xycoords='data',
            xytext=(-70, 13), textcoords='data',
            bbox=dict(boxstyle="round", fc="0.8"),
            arrowprops=dict(arrowstyle="->",
                            connectionstyle="angle,angleA=0,angleB=90,rad=10",
                            color='w'),
           )

ax.annotate('zeroth order',
            xy=(297.5, -5), xycoords='data',
            xytext=(-70, -14), textcoords='data',
            bbox=dict(boxstyle="round", fc="0.8"),
            arrowprops=dict(arrowstyle="->",
                            connectionstyle="angle,angleA=0,angleB=90,rad=10",
                            color='w'),
           )

ax.annotate('zeroth order',
            xy=(-302.5, -10), xycoords='data',
            xytext=(-70, -14), textcoords='data',
            bbox=dict(boxstyle="round", fc="0.8"),
            arrowprops=dict(arrowstyle="->",
                            connectionstyle="angle,angleA=0,angleB=90,rad=10",
                            color='w'),
           )

for chan in ['1', '2', '1m', '2m']:
    x, y = defaultconf['pos_opt_ax'][chan][:2]
    sig = -1 if 'm' in chan else 1
    ax.annotate("",
                xy=(x + sig * 500, y + 0), 
                xycoords='data',
                xytext=(x + sig * 10, y), textcoords='data',
                arrowprops=dict(arrowstyle="->",
                                connectionstyle="arc3",
                                ls=':', lw=3, ec='0.5',
                               mutation_scale=20),
               )

cbar = fig.colorbar(cax, ax=ax)
cbar.set_label('counts / bin')
#ax.set_xlim(-320, -280)
#fig.savefig(os.path.join(get_path("figures"), 'EQFullDet.png'), bbox_inches='tight', dpi=300)
#fig.savefig(os.path.join(get_path("figures"), 'EQFullDet.pdf'), bbox_inches='tight', dpi=300)

This plot shows an extremely distorted view of the CCD plane. The $x$-axis extends over 1.2 m, while the $y$-axis covers only 2.5 cm. Each CCD is surrounded by a white frame. Due to the dimensions of the image, CCDs look high and narrow, but in reality they are twice as long in $x$ direction as in $y$ direction. Arcus has four separate channels, each with its own zero order location and dispersed trace. The positions of the zero orders are marked in the figure and gray arrows indicate the dispersion direction.

Each spectral trace is slightly curved. This is the the result of astigmatism coupled with sub-aperturing. If each channel had a full 360$^\circ$ mirror, we would see an hour-glass shape that is narrowest at the position of the zeroth order and then widens as the spectral focus diverges from the imaging focus. In Arcus' design with a tilted torus, there is a second point where the spectral focus coincides with the imaging focus and thus the spectral trace would become narrower in cross-dispersion direction again. However, the SPOs for each channel cover only a very small fraction of the full circle, and thus the spectral trace fills only one "edge" of the hourglass shape, which results in the curved image. The dispersion direction is exactly parallel to the $x$ axis at all locations, but the $y$ extend of the extraction region can be chosen to follow the curve to reduce the number of background events in the spatial extraction region.

This particular simulation does not include any background. The photons that are seen in between the spectra are cross-dispersed photons. The CAT gratings disperse photons in one direction. However, the main CAT grating bars are held in place by a support structure (L1) that also acts as a grating and disperses some of the photons perpendicular to the main dispersion direction. This is visible in the plot for the zeroth orders and the strongest emission lines.

In [6]:

fig, axes = plt.subplots(4, 4, sharex=True, sharey=True, subplot_kw={'aspect': 'equal'}, figsize=(12, 7))

for i in range(16):
    ind = obs['CCD'] == i + 1
    hist, xedges, yedges = np.histogram2d(obs['detpix_x'][ind], obs['detpix_y'][ind], weights=obs['probability'][ind],
                            bins=[np.arange(0, 2048, 8), np.arange(0, 1024, 8)])

    norm = ImageNormalize(hist, interval=AsymmetricPercentileInterval(0, 99.9),
                          stretch=SqrtStretch())

    ax = axes[np.unravel_index(i, axes.shape)]
    cax = ax.imshow(hist.T, extent=[0, 2048, 0, 1024],
                origin='lower', aspect='equal', 
                cmap=plt.get_cmap('inferno'),
                norm=norm)
    ax.set_xlim(0, 2048)
    ax.set_ylim(0, 1024)
    ax.set_title(i+1)

#cbar = fig.colorbar(cax, ax=ax)
#cbar.set_label('counts / bin')
fig.savefig(os.path.join(get_path("figures"), 'EQAllDet.png'), bbox_inches='tight', dpi=300)
fig.savefig(os.path.join(get_path("figures"), 'EQAllDet.pdf'), bbox_inches='tight', dpi=300)

This figure shows a more faithful representation of the Arcus data. Each panel shows a binned CCD image. The zeroth orders are seen on CCD 7 and 10 and the dispersed photons from the zeroth orders on CCD 7 run towards higher CCD numbers, those starting on CCD 10 towards lower CCD numbers. The color scale is different for each CCD to bring out faint features in regions of low signal. The high count number in the zeroth orders means that the dispersed spectrum is displayed relatively weak in comparison.

In [7]:

ind = obs['CCD'] == 7
hist, xedges, yedges = np.histogram2d(obs['detpix_x'][ind], obs['detpix_y'][ind], weights=obs['probability'][ind],
                            bins=[np.arange(850, 950), np.arange(100, 270)])
fig = plt.figure()
ax = fig.add_subplot(111, aspect='equal')
# Create an ImageNormalize object
norm = ImageNormalize(hist, interval=AsymmetricPercentileInterval(0, 99.9),
                      stretch=SqrtStretch())


cax = ax.imshow(hist.T, extent=[xedges.min(), xedges.max(), 
                                yedges.min(), yedges.max()],
                origin='lower', aspect='equal', 
                interpolation='nearest',
                cmap=plt.get_cmap('inferno'),
                norm=norm)
cbar = fig.colorbar(cax, ax=ax)
cbar.set_label('counts / bin')
ax.set_title("Zero order image")
ax.set_ylabel("cross-dispersion direction (pixel)")
ax.set_xlabel("dispersion direction (pixel)")
ax.xaxis.set_major_locator(ticker.MaxNLocator(4))
fig.savefig(os.path.join(get_path("figures"), 'EQZeroOrder.png'), bbox_inches='tight', dpi=300)
fig.savefig(os.path.join(get_path("figures"), 'EQZeroOrder.pdf'), bbox_inches='tight', dpi=300)

This picture shows a zoom into one on the zeroth orders which is roughly bow-tie shaped with most photons concentrated in the middle, where the structure is narrowest in the dispersion direction (left to right). There are some photons detected above and below the main bow-tie. Those are the result of L1 cross-dispersion. The distribution becomes wider for $y$ coordinates further away from the center. Dispersed orders look similar in shape and this opens an interesting trade: The spectral resolution of the extracted spectra is higher, if only data from the center of the order is extracted, at the cost of the some effective area.

In [8]:

ind = obs['CCD'] == 12
hist, xedges, yedges = np.histogram2d(obs['detpix_x'][ind], obs['detpix_y'][ind], weights=obs['probability'][ind],
                                      bins=100) 
fig = plt.figure()
ax = fig.add_subplot(111, aspect='equal')
# Create an ImageNormalize object
norm = ImageNormalize(hist, interval=AsymmetricPercentileInterval(0, 99.9),
                      stretch=SqrtStretch())


cax = ax.imshow(hist.T, extent=[xedges.min(), xedges.max(), 
                                yedges.min(), yedges.max()],
                origin='lower', aspect='equal', 
                interpolation='nearest',
                cmap=plt.get_cmap('inferno'),
                norm=norm)
cbar = fig.colorbar(cax, ax=ax)
cbar.set_label('counts / bin')

In [9]:

ind = obs['CCD'] == 12
hist, xedges, yedges = np.histogram2d(obs['detpix_x'][ind], obs['detpix_y'][ind], weights=obs['probability'][ind],
                            bins=[np.arange(1350, 1600), np.arange(120, 270)])
fig = plt.figure()
ax = fig.add_subplot(111, aspect='equal')
# Create an ImageNormalize object
norm = ImageNormalize(hist, interval=AsymmetricPercentileInterval(0, 99.9),
                      stretch=SqrtStretch())


cax = ax.imshow(hist.T, extent=[xedges.min(), xedges.max(), 
                                yedges.min(), yedges.max()],
                origin='lower', aspect='equal', 
                interpolation='nearest',
                cmap=plt.get_cmap('inferno'),
                norm=norm)
cbar = fig.colorbar(cax, ax=ax)
cbar.set_label('counts / bin')
ax.set_title("Continuum and lines in a spectrum")
ax.set_ylabel("cross-dispersion direction (pixel)")
ax.set_xlabel("dispersion direction (pixel)")
fig.savefig(os.path.join(get_path("figures"), 'EQline.png'), bbox_inches='tight', dpi=300)
fig.savefig(os.path.join(get_path("figures"), 'EQline.pdf'), bbox_inches='tight', dpi=300)

This plot shows a zoom into a spectral region with a few emission lines and some weaker signal (continuum and weak lines) in between. The bow-tie shape of the dipersed lines is not easy to see given the limited number of photons in each line.