Note

The following shows the code used to run this marx test. You can inspect it and adapt it to your needs, but you cannot copy and paste it directly because it depends on local $PATH and other environment variables. For example, we use a python function to manage the directory structure for all the images generated by all the tests instead of giving the file name directly to save images.

Simulated off-axis PSF¶

CIAO : Set up default marx.par file¶

cp /nfs/melkor/d1/guenther/marx/installed/dev//share/marx/pfiles/marx.par marx.par
pset marx.par GratingType=NONE
pset marx.par DetectorType=HRC-I
pset marx.par DitherModel=INTERNAL
pset marx.par ExposureTime=10000
pset marx.par SourceFlux=0.2
pset marx.par SourceRA=0.
pset marx.par SourceDEC=0.
pset marx.par RA_Nom=0.
pset marx.par Dec_Nom=0.
pset marx.par Roll_Nom=0.
pset marx.par MinEnergy=4
pset marx.par MaxEnergy=4

marx : run marx for different energies¶

marx SourceRA=0.0 OutputDir=marx0
marx SourceRA=0.008333333333333333 OutputDir=marx0.5
marx SourceRA=0.016666666666666666 OutputDir=marx1
marx SourceRA=0.03333333333333333 OutputDir=marx2
marx SourceRA=0.05 OutputDir=marx3
marx SourceRA=0.08333333333333333 OutputDir=marx5
marx SourceRA=0.11666666666666667 OutputDir=marx7
marx SourceRA=0.16666666666666666 OutputDir=marx10
marx SourceRA=0.25 OutputDir=marx15

marxasp : Generate asol file¶

Since all simulations use the same pointing and exposure time, it is enough to run marxasp once.

marxasp MarxDir=marx0

marx2fits : Fits files from marx runs¶

marx2fits --pixadj=NONE marx0 marx0.fits
marx2fits --pixadj=NONE marx0.5 marx0.5.fits
marx2fits --pixadj=NONE marx1 marx1.fits
marx2fits --pixadj=NONE marx2 marx2.fits
marx2fits --pixadj=NONE marx3 marx3.fits
marx2fits --pixadj=NONE marx5 marx5.fits
marx2fits --pixadj=NONE marx7 marx7.fits
marx2fits --pixadj=NONE marx10 marx10.fits
marx2fits --pixadj=NONE marx15 marx15.fits

SAOTrace : Now use SAOTrace ¶

trace-nest tag=sao0 srcpars="point{ position = { ra = 0.0, dec = 0., ra_aimpt=0., dec_aimpt=0. }, spectrum = { { 4., 0.2 } } } roll(0) dither_asol_marx{ file = 'marx_asol1.fits', ra = 0., dec = 0., roll = 0. }" tstart=362912400.01  limit=9999.983202617168 limit_type=sec
trace-nest tag=sao0.5 srcpars="point{ position = { ra = 0.008333333333333333, dec = 0., ra_aimpt=0., dec_aimpt=0. }, spectrum = { { 4., 0.2 } } } roll(0) dither_asol_marx{ file = 'marx_asol1.fits', ra = 0., dec = 0., roll = 0. }" tstart=362912400.01  limit=9999.983202617168 limit_type=sec
trace-nest tag=sao1 srcpars="point{ position = { ra = 0.016666666666666666, dec = 0., ra_aimpt=0., dec_aimpt=0. }, spectrum = { { 4., 0.2 } } } roll(0) dither_asol_marx{ file = 'marx_asol1.fits', ra = 0., dec = 0., roll = 0. }" tstart=362912400.01  limit=9999.983202617168 limit_type=sec
trace-nest tag=sao2 srcpars="point{ position = { ra = 0.03333333333333333, dec = 0., ra_aimpt=0., dec_aimpt=0. }, spectrum = { { 4., 0.2 } } } roll(0) dither_asol_marx{ file = 'marx_asol1.fits', ra = 0., dec = 0., roll = 0. }" tstart=362912400.01  limit=9999.983202617168 limit_type=sec
trace-nest tag=sao3 srcpars="point{ position = { ra = 0.05, dec = 0., ra_aimpt=0., dec_aimpt=0. }, spectrum = { { 4., 0.2 } } } roll(0) dither_asol_marx{ file = 'marx_asol1.fits', ra = 0., dec = 0., roll = 0. }" tstart=362912400.01  limit=9999.983202617168 limit_type=sec
trace-nest tag=sao5 srcpars="point{ position = { ra = 0.08333333333333333, dec = 0., ra_aimpt=0., dec_aimpt=0. }, spectrum = { { 4., 0.2 } } } roll(0) dither_asol_marx{ file = 'marx_asol1.fits', ra = 0., dec = 0., roll = 0. }" tstart=362912400.01  limit=9999.983202617168 limit_type=sec
trace-nest tag=sao7 srcpars="point{ position = { ra = 0.11666666666666667, dec = 0., ra_aimpt=0., dec_aimpt=0. }, spectrum = { { 4., 0.2 } } } roll(0) dither_asol_marx{ file = 'marx_asol1.fits', ra = 0., dec = 0., roll = 0. }" tstart=362912400.01  limit=9999.983202617168 limit_type=sec
trace-nest tag=sao10 srcpars="point{ position = { ra = 0.16666666666666666, dec = 0., ra_aimpt=0., dec_aimpt=0. }, spectrum = { { 4., 0.2 } } } roll(0) dither_asol_marx{ file = 'marx_asol1.fits', ra = 0., dec = 0., roll = 0. }" tstart=362912400.01  limit=9999.983202617168 limit_type=sec
trace-nest tag=sao15 srcpars="point{ position = { ra = 0.25, dec = 0., ra_aimpt=0., dec_aimpt=0. }, spectrum = { { 4., 0.2 } } } roll(0) dither_asol_marx{ file = 'marx_asol1.fits', ra = 0., dec = 0., roll = 0. }" tstart=362912400.01  limit=9999.983202617168 limit_type=sec

marx : run marx for all SAOTrace runs¶

marx SourceType=SAOSAC SAOSACFile=sao0.fits OutputDir=saomarx0 DitherModel=FILE DitherFile=marx_asol1.fits
marx SourceType=SAOSAC SAOSACFile=sao0.5.fits OutputDir=saomarx0.5 DitherModel=FILE DitherFile=marx_asol1.fits
marx SourceType=SAOSAC SAOSACFile=sao1.fits OutputDir=saomarx1 DitherModel=FILE DitherFile=marx_asol1.fits
marx SourceType=SAOSAC SAOSACFile=sao2.fits OutputDir=saomarx2 DitherModel=FILE DitherFile=marx_asol1.fits
marx SourceType=SAOSAC SAOSACFile=sao3.fits OutputDir=saomarx3 DitherModel=FILE DitherFile=marx_asol1.fits
marx SourceType=SAOSAC SAOSACFile=sao5.fits OutputDir=saomarx5 DitherModel=FILE DitherFile=marx_asol1.fits
marx SourceType=SAOSAC SAOSACFile=sao7.fits OutputDir=saomarx7 DitherModel=FILE DitherFile=marx_asol1.fits
marx SourceType=SAOSAC SAOSACFile=sao10.fits OutputDir=saomarx10 DitherModel=FILE DitherFile=marx_asol1.fits
marx SourceType=SAOSAC SAOSACFile=sao15.fits OutputDir=saomarx15 DitherModel=FILE DitherFile=marx_asol1.fits

marx2fits : Fits files from marx + SAOTrace runs¶

marx2fits --pixadj=NONE saomarx0 saomarx0.fits
marx2fits --pixadj=NONE saomarx0.5 saomarx0.5.fits
marx2fits --pixadj=NONE saomarx1 saomarx1.fits
marx2fits --pixadj=NONE saomarx2 saomarx2.fits
marx2fits --pixadj=NONE saomarx3 saomarx3.fits
marx2fits --pixadj=NONE saomarx5 saomarx5.fits
marx2fits --pixadj=NONE saomarx7 saomarx7.fits
marx2fits --pixadj=NONE saomarx10 saomarx10.fits
marx2fits --pixadj=NONE saomarx15 saomarx15.fits

Python : No image gallery needed¶

# Note that this code might not run if you directly copy and paste it:
# - Not all import statements are shown here
# - `self` is a reference to a test instance, which allows access to
#   parameters such as the directory where the test is run etc.

'''No image gallery needed'''
pass

Python : Plots of radial distribution¶

# Note that this code might not run if you directly copy and paste it:
# - Not all import statements are shown here
# - `self` is a reference to a test instance, which allows access to
#   parameters such as the directory where the test is run etc.

def radial_distribution(x, y):
    '''Calculate the radial distance from the mean position for events.

    Parameters
    ----------
    x, y : np.array
        x and y coordinates of events

    Returns
    -------
    r : np.array
        radial distance from mean position of the events.
    '''
    centx = x - astropy.stats.sigma_clipped_stats(x, sigma=1.5, maxiters=20)[0]
    centy = y - astropy.stats.sigma_clipped_stats(y, sigma=1.5, maxiters=20)[0]

    xy = np.vstack([centx, centy])
    return np.linalg.norm(xy, axis=0)

'''Plots of radial distribution'''
import os
import numpy as np
from matplotlib import pyplot as plt
from astropy.table import Table
from astropy.io import fits

ecf = [0.5, 0.7, 0.9]
fig = plt.figure()
axecf = fig.add_subplot(111)
out = np.zeros((len(self.parameter), 2, len(ecf)))

for i, e in enumerate(self.parameter):

    for prog in ['marx', 'saomarx']:
        f = os.path.join(self.basepath, '{0}{1}.fits'.format(prog, e))
        with fits.open(f) as hdus:
            header = hdus[1].header
        tab = Table.read(f, hdu=1)
        # The old question: Is an index the center of a pixel of the corner?
        # Differs by 0.5...
        d_ra = (tab['X'] - header['TCRPX9'] - 0.5) * header['TCDLT9'] * 3600.
        # Count from the left or right (RA is reversed on the sky)?
        d_dec = (tab['Y'] - header['TCRPX10'] + 0.5) * header['TCDLT10'] * 3600.

        # The simulation is set up to make this simple: RA=DEC=0
        # so cos(DEC) = 1 and we can approximate with Euklidian distance
        r = radial_distribution(d_ra, d_dec)
        if prog == 'marx':
            j = 0
        else:
            j = 1
        out[i, j, :] = np.percentile(r, np.array(ecf) * 100.)

for i, e in enumerate(ecf):
    line, = axecf.plot(self.parameter, out[:, 0, i],
                       label='{0} %'.format(e*100))
    axecf.plot(self.parameter, out[:, 1, i], color=line.get_color(),
               ls=':')
axecf.legend(loc='lower right')
axecf.set_ylabel('radius [arcsec]')
axecf.set_xlabel('off-axis angle [arcmin]')
axecf.grid()
fig.savefig(self.figpath('ECF'), bbox_inches='tight')

Simulated off-axis PSF¶

CIAO : Set up default marx.par file¶

marx : run marx for different energies¶

marxasp : Generate asol file¶

marx2fits : Fits files from marx runs¶

SAOTrace : Now use SAOTrace ¶

marx : run marx for all SAOTrace runs¶

marx2fits : Fits files from marx + SAOTrace runs¶

Python : No image gallery needed¶

Python : Plots of radial distribution¶

Page Content

Previous topic

Next topic

Simulated off-axis PSF¶

CIAO : Set up default marx.par file¶

marx : run marx for different energies¶

marxasp : Generate asol file¶

marx2fits : Fits files from marx runs¶

SAOTrace : Now use SAOTrace¶

marx : run marx for all SAOTrace runs¶

marx2fits : Fits files from marx + SAOTrace runs¶

Python : No image gallery needed¶

Python : Plots of radial distribution¶

SAOTrace : Now use SAOTrace ¶