A Simple Single Observation, Spectroscopy and Imaging

By a simple observation, we mean an observation like a “blank-sky” observation, where the residual soft proton contamination can be fit from a spectrum extracted from the entire FOV. This section provides a recipe for extracting the spectrum and constructing an image of the entire FOV. If you have a big cluster in the middle of your FOV, you'll want to follow the recipe in this section to extract the raw image, the quiescent particle background image, and the spectrum of the entire FOV, but you will want to follow the recipe in the following section to extract the spectrum from the region least dominated by the cluster and, after fitting that, to construct the soft proton image.

One of the things that you will notice is that this script does a large amount of logging, so you can retain a record of what actually happened, and listing, so you know what files were produced by which routine. Given the number of files that SAS and ESAS produces, it's often a good thing to keep track of. This script also uses a number of environment variables to reduce the amount of typing (and errors).

######################################################################################################
###                     THIS IS THE XMM-ESAS SCRIPT KDK-RUN                                        ###
######################################################################################################
#  SAS V21
#
# This script is designed to run under MACOSX.  It assumes that SAS has been set up, but sets several
# of the SAS parameters for the specific observation. The required data files are in directories which
# are inputs to scripts and programs which use them. 
#
# All ODF data have been placed in the directory 
#     /path/ODF
# but nothing further has been done to them
#
# All processing takes place in the directory
#     /path/analysis
#
# Date: 11/09/22 but with constant modifications
# Purp: this script is for <insert your identifier> to be run in the "analysis" directory
# Note: this script is meant to be modified by the user for their own purpose
#       There are many things in here that may not be useful for you particular analysis,
#       Just delete what you don't need!
#       This script lists the files repeatedly. That is so you can compare the lists to
#       determine what was made in a particular step...in case you need to delete them.

# Modify the following lines before you do anything else to make sure that you are pointing to
# the correct directory

# Set a few SAS environment variables

setenv OBSID 0000000001
setenv SAS_CCFPATH /Users/kuntz/Soft/sas_ccf
setenv SAS_CCF /Users/kuntz/Data/xmm_proj/${OBSID}/analysis/ccf.cif
setenv SAS_ODF /Users/kuntz/Data/xmm_proj/${OBSID}/ODF

# Cleaning up the ODF directory - this only has to be done once
# IF
# the data downloaded from the HEASARC:

cd $SAS_ODF
gzip -d *.gz

# ENDIF ELSE downloaded from the SOC

tar -xvf *.tar # for an ODF extracted from XSA this may have to be run several times
rm *.SAS       # this is a "just in case" step which would be needed 
               # if you have moved the directory since the last odfingest

# ENDELSE
# End of cleaning up the ODF

cd ..
mkdir analysis
cd ../analysis

# Run the inital tasks to set up SAS for this ObsID - that is - set up the CCF
#   Here is where you want to make sure that there is no *.SAS file in the ODF directory

cifbuild withccfpath=no analysisdate=now category=XMMCCF calindexset=$SAS_CCF fullpath=yes verbosity=1

odfingest odfdir=$SAS_ODF outdir=$SAS_ODF verbosity=1

# At this point there should be only the ccf.cif in the analysis directory

# Run epchain and emchain to create calibrated photon event files using the latest calibrations.
#  emchain will process all sub-observations for both MOS detectors
#  epchain will process only the first sub-observation
#    and will need to be run for each pn sub-observation
#  epchain will need to be run a second time for each sub-observation
#    to produce the OOT event files

emchain verbosity=1 >& emchain.log

# While emchain processes all segments automatically, epchain only processes the first pn segment
# automatically. Further segements require specifying the segment number.
# To determine the number of pn segments to process:

epchain exposure=90

# This will end with an error since there is no exposure with index 90. However,
# the output will show you what exposures exist for the pn. 
# For each pn segment, where the argument to "exposure" is the index of the exposure

epchain exposure=1 verbosity=1 >& epchain_1.log
epchain exposure=1 withoutoftime=true verbosity=1 >& epchain_oot_1.log

# Clean up after the chains since the chains produce many files that you won't need

ls -l *.* > post_chain.list
mkdir chains
mkdir logs
mv *.FIT chains
cp chains/*IEVLI*.FIT .
cp chains/*OOEVLI*.FIT .
mv *.log logs
mv *.list logs

# At this point you should determine the filter and exposure time for each segment, that is, for each
# IEVLI file determine whether it should be retained or deleted. 
# You may find it useful to create a routine to extract information from all of the 
# *IEVLI* files the way I have:

#P0000000001M1S001MIEVLI0000.FIT 336.86800  61.02600  Thin1            PrimeFullWindow  35.51
#R0000000001M1U002MIEVLI0000.FIT 336.86800  61.02600 CalClosed         PrimeFullWindow   0.46
#P0000000001M1U003MIEVLI0000.FIT 336.86800  61.02600  Thin1            PrimeFullWindow  10.26
#P0000000001M1U004MIEVLI0000.FIT 336.86800  61.02600 CalClosed         PrimeFullWindow   0.64
#P0000000001M1U005MIEVLI0000.FIT 336.86800  61.02600  Thin1            PrimeFullWindow   1.26
#P0000000001M2S002MIEVLI0000.FIT 336.86800  61.02600  Thin1            PrimeFullWindow  47.32
#P0000000001M2U002MIEVLI0000.FIT 336.86800  61.02600 CalClosed         PrimeFullWindow   0.11
#P0000000001M2U005MIEVLI0000.FIT 336.86800  61.02600  Thin1.           PrimeFullWindow   1.21
#P0000000001PNS003PIEVLI0000.FIT 336.86800  61.02600  Thin1    PrimeFullWindowExtended  33.36

# The event files for the "CalClosed" segments are not useful and can be deleted.

rm P${OBSID}M1U002MIEVLI0000.FIT
rm P${OBSID}M1U004MIEVLI0000.FIT
rm P${OBSID}M2U002MIEVLI0000.FIT

# It is often a good idea to remove really short observations, usually <5 ks

rm P${OBSID}M1S002MIEVLI0000.FIT
rm P${OBSID}M2S003MIEVLI0000.FIT

# For this observations the useful segments are mos1S001, (possibly) mos1U003,
# mos2S002, and pnS003.

# If there are multiple segments for a given detector, you may want to merge them here, but it is
# probably better to proceed with the emanom and espfilt steps, and only merge if both segments
# should be retained

merge set1=P${OBSID}M1S001MIEVLI0000.FIT set2=P${OBSID}M1U003MIEVLI0000.FIT outset=P${OBSID}M1T000MIEVLI0000.FIT

# Rename the files that you need for ease and clarity
# These lines work if you used the "chains", but won't if you used the "procs"
# Please don't try to give segments new numbers...many of the ESAS programs
# look for the segment number in the header and use it to create output file names.

mv P${OBSID}M1T000MIEVLI0000.FIT mos1S001.fits
mv P${OBSID}M2S002MIEVLI0000.FIT mos2S002.fits
mv P${OBSID}PNS003OOEVLI0000.FIT pnS003-oot.fits
mv P${OBSID}PNS003PIEVLI0000.FIT pnS003.fits

# Set up some system variables so you don't have to go through this script,
# changing every "mos1S001" to "mos1U017", etc.!

setenv M1 "mos1S001"
setenv M2 "mos2S002"
setenv PN "pnS003"

ls -l *.* > post_chain_clean.list

# Take a look at the data to see what you've got

evselect table=${M1}.fits withimageset=yes imageset=${M1}-diag-det-unfilt.fits 
  expression="((PI in [300:1000])&&(PATTERN<=12)&&((FLAG & 0x766aa000)==0))" filtertype=expression
  ignorelegallimits=yes imagebinning=imageSize xcolumn=DETX ximagesize=780 ximagemax=19500
  ximagemin=-19499 ycolumn=DETY yimagesize=780 yimagemax=19500 yimagemin=-19499

evselect table=${M2}.fits withimageset=yes imageset=${M2}-diag-det-unfilt.fits
  expression="(PI in [300:1000])&&(PATTERN<=12)&&((FLAG & 0x766aa000)==0)" filtertype=expression
  ignorelegallimits=yes imagebinning=imageSize xcolumn=DETX ximagesize=780 ximagemax=19500
  ximagemin=-19499 ycolumn=DETY yimagesize=780 yimagemax=19500 yimagemin=-19499

evselect table=${PN}.fits withimageset=yes imageset=${PN}-diag-det-unfilt.fits
  expression="(PI in [300:1000])&&(PATTERN<=4)&&(#XMMEA_EP)" ignorelegallimits=yes imagebinning=imageSize
  xcolumn=DETX ximagesize=780 ximagemax=19500 ximagemin=-19499 ycolumn=DETY yimagesize=780
  yimagemax=19500 yimagemin=-19499

ds9 *diag-det-unfilt.fits &

# Some clean-up of the diagnostic files (of course you can also just delete them)

mkdir diags
mv *diag-det-unfilt.fits diags

# Get the anomalous chip information - you should look at the screen output or 
# the anom*.log files carefully so you will know which CCDs to deselect in some of the 
# future steps (mosspectra and mosback)

emanom eventfile=${M1}.fits keepcorner=FALSE

emanom eventfile=${M2}.fits keepcorner=FALSE

# Determine which chips should be on, and build environment variables for them
#   you'll be glad that you set these as environment variables later!
#   These are examples only! (These get rid of chips that were problems in my obsid.)

setenv M1ON "T T F F T F T" 
setenv M2ON "T T T T T T T"
setenv PNON "T T T T"            

ls -l *.* > post_emanom.list
mv *-anom.log logs
mv *.list logs

# Soft proton filtering - this will remove periods of obvious soft proton flares from 
# the event files. It will also produce QDP plot files showing the light curves and 
# the accepted time intervals. Note that while much of the SP contamination is removed,
# there is likely to be some residual contamination.
# Note that a different rangescale is used for the pn than for the MOS!
# The most important output is the cleaned event file: x-allevc.fits 
# For the pn, both a cleaned event file and a cleaned OOT event file are produced.
# Note that the energy range is chosen to maximize the soft proton/cosmic ratio.

espfilt eventfile=${M1}.fits method=histogram elow=2500 ehigh=8500 withsmoothing=yes smooth=51
  rangescale=6.0 allowsigma=3.0 keepinterfiles=true > espfilt_mos1.log

espfilt eventfile=${M2}.fits method=histogram elow=2500 ehigh=8500 withsmoothing=yes smooth=51
  rangescale=6.0 allowsigma=3.0 keepinterfiles=true > espfilt_mos2.log

espfilt eventfile=${PN}.fits method=histogram elow=2500 ehigh=8500 withsmoothing=yes smooth=51
  rangescale=15.0 allowsigma=3 withoot=Y ootfile=${PN}-oot.fits keepinterfiles=true > espfilt_pn.log

# Look at the QDP output from espfilt. The plots should be examined to get an idea 
# of the extent of any residual contamination.  Asymmetry of the peaks in the count 
# rate histograms or strong high count-rate tails often imply that there is going 
# to be residual soft proton contamination.

qdp ${M1}-hist.qdp
/cps
exit
mv pgplot.ps ${M1}-hist.ps
qdp ${M2}-hist.qdp
/cps
exit
mv pgplot.ps ${M2}-hist.ps
qdp ${PN}-hist.qdp
/cps
exit
mv pgplot.ps ${PN}-hist.ps

# At this point you might want to check the diagnostic images.
# See whether the edge of the FOV is still visible in the cleaned images. 
# It shouldn't be unless your extended source is brighter than the cosmic X-ray background.

ds9 *-allim.fits *-allimc.fits &

# IF
# the espfilt cleaning is not sufficient, you may edit the -gti.fits files by hand to produce a 
# -gti-edit.fits file that does the cleaning that you want.
# You will need to apply the new -gti-edit.fits files to the original event files as follows:

evselect table=${M1}.fits withfilteredset=yes expression="(GTI(${M1}-gti-edit.fits,TIME))&&
  ((FLAG & 0x766aa000)==0)&&(PATTERN<=12)" filtertype=expression keepfilteroutput=yes updateexposure=yes
  filterexposure=yes filteredset=${M1}-allevc.fits

evselect table=${M2}.fits withfilteredset=yes expression="(GTI(${M2}-gti-edit.fits,TIME))&&
  ((FLAG & 0x766aa000)==0)&&(PATTERN<=12)" filtertype=expression keepfilteroutput=yes updateexposure=yes
  filterexposure=yes filteredset=${M2}-allevc.fits

evselect table=${PN}-oot.fits withfilteredset=yes expression="(GTI(${PN}-gti-edit.fits,TIME))
  &&(#XMMEA_EP)&&(PATTERN<=4)" filtertype=expression keepfilteroutput=yes updateexposure=yes
  filterexposure=yes filteredset=${PN}-allevcoot.fits

evselect table=${PN}.fits withfilteredset=yes expression="(GTI(${PN}-gti-edit.fits,TIME))&&
  (#XMMEA_EP)&&(PATTERN<=4)" filtertype=expression keepfilteroutput=yes updateexposure=yes
  filterexposure=yes filteredset=${PN}-allevc.fits

# ENDIF 

# Note that the gti selection must be applied to both the pn and the pn-oot event files.
# Note that by default espfilt applies filters to both the flags and the patterns, which you must
# reproduce to create consistent files.

# Take a look at the data to see whether it's really cleaned.

evselect table=${M1}-allevc.fits withimageset=yes imageset=${M1}-diag-det-filt.fits 
  expression="((PI in [300:1000])&&(PATTERN<=12)&&((FLAG & 0x766aa000)==0))" filtertype=expression
  ignorelegallimits=yes imagebinning=imageSize xcolumn=DETX ximagesize=780 ximagemax=19500
  ximagemin=-19499 ycolumn=DETY yimagesize=780 yimagemax=19500 yimagemin=-19499

evselect table=${M2}-allevc.fits withimageset=yes imageset=${M2}-diag-det-filt.fits
  expression="(PI in [300:1000])&&(PATTERN<=12)&&((FLAG & 0x766aa000)==0)" filtertype=expression
  ignorelegallimits=yes imagebinning=imageSize xcolumn=DETX ximagesize=780 ximagemax=19500
  ximagemin=-19499 ycolumn=DETY yimagesize=780 yimagemax=19500 yimagemin=-19499

evselect table=${PN}-allevc.fits withimageset=yes imageset=${PN}-diag-det-filt.fits
  expression="(PI in [300:1000])&&(PATTERN<=4)&&(#XMMEA_EP)" ignorelegallimits=yes imagebinning=imageSize
  xcolumn=DETX ximagesize=780 ximagemax=19500 ximagemin=-19499 ycolumn=DETY yimagesize=780
  yimagemax=19500 yimagemin=-19499

ds9 *-diag-det-filt.fits
ls -l *.* > post_espfilt.list

# Clean up after espfilt

mv *-diag-* diags
mkdir espfilt
mv *.qdp espfilt
mv *-corev.fits espfilt
mv *-corimc.fits espfilt
mv *-corim.fits espfilt
mv *-corlc.fits espfilt
mv *-allim.fits espfilt
mv *-allimc.fits espfilt
mv *-fovlc.fits espfilt
mv *-hist.* espfilt
mv *-gti.fits espfilt
mv *-gti-edit.fits espfilt
mv *-gti.txt espfilt
mv *-corevc.fits espfilt
mv ${M1}.fits espfilt
mv ${M2}.fits espfilt
mv ${PN}.fits espfilt
mv ${PN}-oot.fits espfilt

ls -l *.* > post_espfilt_clean.list
mv *.log logs
mv *.list logs

# At this point all you need in your analysis directory are 
# ${M1}-allevc.fits, ${M2}-allevc.fits, ${PN}-allevc.fits, & ${PN}-allevcoot.fits

# Source detection - cheese is named such because it creates masks that look much like Swiss cheese.
# Here we are running it on two bands, but it can run for only one band as well.
# If run for two bands, the program produces three bands: soft (s), hard (h), and total (t = s+h)
# (though no total band cheese is produced, it can be produced as the product of the soft and hard)
# If run for a single band, it is considered the total (t) band.

cheese mos1file=${M1}-allevc.fits mos2file=${M2}-allevc.fits pnfile=${PN}-allevc.fits 
  pnootfile=${PN}-allevcoot.fits elowlist=350 ehighlist=1100 scale=0.5 mlmin=10 ratetotal=0.2 dist=50.
  keepinterfiles=no -V 7 >& cheese.log

# Inspect the cheese images for any anomalies: sources that didn't get removed, or spurious sources
# that did. See the cook book for advice on what to do to get rid of the right sources. One may need 
# to edit a copy of the emllist.fits file and rerun a few steps of cheese in order to produce better 
# source exclusion regions, x-bkregBdet.fits and x-bkregBsky.fits and the source mask x-cheeseB.fits

ds9 *fovim?.fits *cheese?.fits &

ls -l *.* > post-cheese.list

# Clean up after running cheese but before the ultimate background region files are produced

mkdir cheese
mv cheese.log cheese
mv eboxlist* cheese
cp mos?????-cheese?.fits cheese
cp pn????-cheese?.fits cheese

# IF
# You really need to reconstruct with a better source list then here is how to do it:
# If all "bad" sources are removed from emllist_fixd_filt.fits then you can set the selected 
# DET_ML arbitrarily low and reconstruct the region files (-bkgregtdet.fits and -bkgregtsky.fits) and
# cheese file. You will need the fovimN.fits and fovimNmask.fits files

# MOS1:

region eventset=${M1}-allevc.fits operationstyle=global srclisttab=emllist_fixd_filt.fits:SRCLIST
  expression="(ID_INST == 2)&&(DET_ML >= 10)" bkgregionset=${M1}-bkgregtdet.fits energyfraction=0.4
  radiusstyle=contour outunit=detxy verbosity=4
region eventset=${M1}-allevc.fits operationstyle=global srclisttab=emllist_fixd_filt.fits:SRCLIST
  expression="(ID_INST == 2)&&(DET_ML >= 10)" bkgregionset=${M1}-bkgregtsky.fits energyfraction=0.4
  radiusstyle=contour outunit=xy verbosity=4
makemask imagefile=${M1}-fovimt.fits maskfile=${M1}-fovimtmask.fits regionfile=${M1}-bkgregtsky.fits
  cheesefile=${M1}-cheeset.fits

# MOS2:

region eventset=${M2}-allevc.fits operationstyle=global srclisttab=emllist_fixd_filt.fits:SRCLIST
  expression="(ID_INST == 3)&&(DET_ML >= 10)" bkgregionset=${M2}-bkgregtdet.fits energyfraction=0.4
  radiusstyle=contour outunit=detxy verbosity=4
region eventset=${M2}-allevc.fits operationstyle=global srclisttab=emllist_fixd_filt.fits:SRCLIST
  expression="(ID_INST == 3)&&(DET_ML >= 10)" bkgregionset=${M2}-bkgregtsky.fits energyfraction=0.4
  radiusstyle=contour outunit=xy verbosity=4
makemask imagefile=${M2}-fovimt.fits maskfile=${M2}-fovimtmask.fits regionfile=${M2}-bkgregtsky.fits
  cheesefile=${M2}-cheeset.fits

# PN:

region eventset=${PN}-allevc.fits operationstyle=global srclisttab=emllist_fixd_filt.fits:SRCLIST
  expression="(ID_INST == 1)&&(DET_ML >= 10)" bkgregionset=${PN}-bkgregtdet.fits energyfraction=0.4
  radiusstyle=contour outunit=detxy verbosity=4
region eventset=${PN}-allevc.fits operationstyle=global srclisttab=emllist_fixd_filt.fits:SRCLIST
  expression="(ID_INST == 1)&&(DET_ML >= 10)" bkgregionset=${PN}-bkgregtsky.fits energyfraction=0.4
  radiusstyle=contour outunit=xy verbosity=4
makemask imagefile=${PN}-fovimt.fits maskfile=${PN}-fovimtmask.fits regionfile=${PN}-bkgregtsky.fits
  cheesefile=${PN}-cheeset.fits

ds9 *fovim?.fits *cheese?.fits &

# ENDIF

# Second post cheese cleanup

mv emllist*.* cheese
mv *-fovim?*.fits cheese
ls -l *.* > post_cheese_clean.list
mv *.log logs
mv *.list logs

# At this point in the main directory you should have only the files *-allevc, *-bkgregdet, 
# *-bkgregtsky, and *-cheeseB for each instrument, as well as the atthk.fit and ccf.cif files.
# At this point I find it useful to save the analysis directory so if something goes wrong with 
# mos/pn spectra and mos/pn back, I can quickly resurrect the directory with which I started!

cd ..
tar -cvf ffov_analysis.tar analysis
cd analysis

# Extracting and  constructing the background spectrum - this is the most time consuming step.
# For this step you will need to have constructed your region selection files (ffov_mos1.txt, etc.)
# I've called the region 'ffov' because it is the full FOV but without the point sources. 
# Use whatever descriptive name you like.
# This is the part of ESAS that produces the most files and can cause the most trouble if you are
# extracting more than one region. To reduce problems, I stick all of the region specific products
# into separate subdirectories. There is more to be done with these files, but I've set it up to be easy
# to follow rather than efficient.

# BREAK POINT 1 (see text in the following subsection of the appendix)

# MAKE SURE THAT THE REGION FILES ARE IN THE ANALYSIS DIRECTORY AT THIS POINT
# You need files containing expressions that describe the region that you want to extract. We assume
# that that region it is either the full FOV, or some large fraction of it. These files do not contain
# information about the sources to be excised (that information is already in the *-bkgregtdet.fits
# files).
# Model region files can be found [webpage to be determined].
# In the following steps the region files are assumed to be called 
#   regmos1-ffov.txt, regmos2-ffov.txt, and regpn-ffov.txt

# Set the elow and ehigh with environment variables to make doing multiple bands easier

setenv ELO 350
setenv EHI 1100

# MOS1 & MOS2

mosspectra eventfile=${M1}-allevc.fits keepinterfiles=yes withregion=yes regionfile=regmos1.txt
  pattern=12 withsrcrem=yes maskdet=${M1}-bkgregtdet.fits masksky=${M1}-bkgregtsky.fits elow=${ELO}
  ehigh=${EHI} ccds="${M1ON}" -V=7  >& mosspectra_1.log

mosspectra eventfile=${M2}-allevc.fits keepinterfiles=yes withregion=yes regionfile=regmos2.txt
  pattern=12 withsrcrem=yes maskdet=${M2}-bkgregtdet.fits masksky=${M2}-bkgregtsky.fits elow=${ELO}
  ehigh=${EHI} ccds="${M2ON}" -V=7 >& mosspectra_2.log

mosback inspecfile=${M1}-fovt.pi elow=${ELO} ehigh=${EHI} ccds="${M1ON}"  >& mosback_1.log

mosback inspecfile=${M2}-fovt.pi elow=${ELO} ehigh=${EHI} ccds="${M2ON}"  >& mosback_2.log

# Look at the diagnostics! There are many diagnostic files. I've gotten out of the habit.
# (Shame on me!) You can catch some problems early by looking at them.

qdp ${M1}-augindiv.qdp
/cps
exit
mv pgplot.ps ${M1}-augindiv.ps
qdp ${M1}-augspec.qdp
/cps
exit
mv pgplot.ps ${M1}-augspec.ps
qdp ${M1}-bkgaccum.qdp
/cps
exit
mv pgplot.ps ${M1}-bkgaccum.ps
qdp ${M1}-bkgindiv.qdp
/cps
exit
mv pgplot.ps ${M1}-bkgindiv.ps
qdp ${M1}-bkgspec.qdp
/cps
exit
mv pgplot.ps ${M1}-bkgspec.ps
qdp ${M1}-bridgefit.qdp
/cps
exit
mv pgplot.ps ${M1}-bridgefit.ps
qdp ${M1}-ratehard.qdp
/cps
exit
mv pgplot.ps ${M1}-ratehard.ps
qdp ${M1}-ratioindiv.qdp
/cps
exit
mv pgplot.ps ${M1}-ratioindiv.ps

qdp ${M2}-augindiv.qdp
/cps
exit
mv pgplot.ps ${M2}-augindiv.ps
qdp ${M2}-augspec.qdp
/cps
exit
mv pgplot.ps ${M2}-augspec.ps
qdp ${M2}-bkgaccum.qdp
/cps
exit
mv pgplot.ps ${M2}-bkgaccum.ps
qdp ${M2}-bkgindiv.qdp
/cps
exit
mv pgplot.ps ${M2}-bkgindiv.ps
qdp ${M2}-bkgspec.qdp
/cps
exit
mv pgplot.ps ${M2}-bkgspec.ps
qdp ${M2}-bridgefit.qdp
/cps
exit
mv pgplot.ps ${M2}-bridgefit.ps
qdp ${M2}-ratehard.qdp
/cps
exit
mv pgplot.ps ${M2}-ratehard.ps
qdp ${M2}-ratioindiv.qdp
/cps
exit
mv pgplot.ps ${M2}-ratioindiv.ps

# Now do some clean-up in order to keep your working directory small enough to find stuff!

ls -l *.* > post_mos.list

mkdir ffov_${M1}
mv mos1*-${ELO}-${EHI}* ffov_${M1}
mv mos1*.pi ffov_${M1}
mv mos1*.qdp ffov_${M1}
mv mos1*.psf fov_${M1}
mv mos1*imt* ffov_${M1}
mv mos1*.arf ffov_${M1}
mv mos1*.rmf ffov_${M1}
mv mos1*imspdet* ffov_${M1}

mkdirffov_${M2}
mv mos2*-${ELO}-${EHI}* ffov_${M2}
mv mos2*.pi ffov_${M2}
mv mos2*.qdp ffov_${M2}
mv mos2*.psf fov_${M2}
mv mos2*imt* ffov_${M2}
mv mos2*.arf ffov_${M2}
mv mos2*.rmf ffov_${M2}
mv mos2*imspdet* ffov_${M2}
ls -l *.* > post_mos_clean.list

# Now we want to repeat the QPB extraction for the pn. However we need to do it once for PATTERN==0 
# (for the soft emission) and once for PATTERN<=4 (for the hard emission).

# PN PATTERN==0

pnspectra eventfile=${PN}-allevc.fits ootevtfile=${PN}-allevcoot.fits keepinterfiles=yes withregion=yes
  regionfile=regpn.txt pattern=0 withsrcrem=yes maskdet=${PN}-bkgregtdet.fits 
  masksky=${PN}-bkgregtsky.fits elow=${ELO} ehigh=${EHI} quads="${PNON}" -V=7 >& pnspectra_0.log

pnback inspecfile=${PN}-fovt.pi inspecoot=${PN}-fovtoot.pi elow=${ELO} ehigh=${EHI} 
  quads="${PNON}" >& pnback_0.log

# Look at the diagnostics

qdp ${PN}-augratehrev.qdp
/cps
exit
mv pgplot.ps ${PN}-augratehrev.ps 
qdp ${PN}-augspec.qdp  
/cps
exit
mv pgplot.ps ${PN}-augspec.ps  
qdp ${PN}-bridgefit.qdp
/cps
exit
mv pgplot.ps ${PN}-bridgefit.ps    
qdp ${PN}-ratehard.qdp
/cps
exit
mv pgplot.ps ${PN}-ratehard.ps
qdp ${PN}-augraterev.qdp
/cps
exit
mv pgplot.ps ${PN}-augraterev.ps  
qdp ${PN}-bkgspec.qdp
/cps
exit
mv pgplot.ps ${PN}-bkgspec.ps    
qdp ${PN}-quadspec.qdp
/cps
exit
mv pgplot.ps ${PN}-quadspec.ps   
qdp ${PN}-spec.qdp
/cps
exit
mv pgplot.ps ${PN}-spec.ps

# and do some clean-up

ls -l *.* > post_pn0.list

mkdir ffov_${PN}_0
mv pn*-${ELO}-${EHI}* ffov_${PN}_0
mv pn*.pi ffov_${PN}_0
mv pn*.qdp ffov_${PN}_0
mv pn*.ps ffov_${PN}_0
mv pn*imt* ffov_${PN}_0
mv pn*.arf ffov_${PN}_0
mv pn*.rmf ffov_${PN}_0
mv pn*imspdet* ffov_${PN}_0

# PN PATTERN<=4

pnspectra eventfile=${PN}-allevc.fits ootevtfile=${PN}-allevcoot.fits keepinterfiles=yes withregion=yes
  regionfile="regpn.txt" pattern=4 withsrcrem=yes maskdet=${PN}-bkgregtdet.fits withimages=yes
  elow=${ELO} ehigh=${EHI} quads="${PNON}" -V=7 >& pnspectra_4.log

pnback inspecfile=${PN}-fovt.pi inspecoot=${PN}-fovtoot.pi elow=${ELO} ehigh=${EHI} 
  quads="${PNON}" >& pnback_4.log

# It is always useful to look at the diagnostics, copy and past the commands from PN_0

# and do the clean-up

ls -l *.* > post_pn4.list

mkdir ffov_${PN}_4
mv pn*${ELO}-${EHI}* ffov_${PN}_4
mv pn*.pi ffov_${PN}_4
mv pn*.qdp ffov_${PN}_4
mv pn*.psf fov_${PN}_4
mv pn*imt* ffov_${PN}_4
mv pn*.arf ffov_${PN}_4
mv pn*.rmf ffov_${PN}_4
mv pn*imspdet* ffov_${PN}_4
ls -l *.* > post_pn_clean.list

mv *.log logs
mv *.list logs

# We now need to get all of the spectral products in one place. 
# IF you are working with something like a blank sky observation, these are the products that
# you will need for the spectral fitting to determine the residual soft proton flare emission
# (and maybe swcx). 
# IF you have a large, bright, extended source in your FOV (such as a cluster of galaxies)
# you will need some of these products for the imaging reduction/analysis

# We create a directory just for the spectra. Then we visit each of the ffov_xxxx directories 
# to copy the right stuff to the spectrum directory. Grouping the spectra is necessary for 
# chi-squared fitting. Populating the header keywords is to save time and effort when we are 
# actually fitting. While we are at it, we can use protonscale to determine the BACKSCAL parameter
# (in square arcminutes) which we will need to use in the fitting.

mkdir ffov_spectrum

cd ffov_${M1}

# An easy way to get the BACKSCAL header keyword which you will need for the spectral fitting

protonscale mode=1 maskfile={$M1}-fovimspdet.fits specfile={$M1}-fovt.pi

# Grouping the spectrum and adding the ANCRFILE, BACKFILE, and RESPFILE header keywords (important)

grppha ${M1}-fovt.pi mos1-grp.pi 'chkey BACKFILE mos1-bkg.pi & chkey RESPFILE mos1.rmf & 
  chkey ANCRFILE mos1.arf & group min 50 & exit'
mv mos1-grp.pi ../ffov_spectrum/mos1-grp.pi
cp ${M1}-bkg.pi ../ffov_spectrum/mos1-bkg.pi
cp ${M1}.rmf ../ffov_spectrum/mos1.rmf
cp ${M1}.arf ../ffov_spectrum/mos1.arf
cd ..

cd ffov_${M2}

protonscale mode=1 maskfile={$M2}-fovimspdet.fits specfile={$M2}-fovt.pi

grppha ${M2}-fovt.pi mos2-grp.pi 'chkey BACKFILE mos2-bkg.pi & chkey RESPFILE mos2.rmf & 
  chkey ANCRFILE mos2.arf & group min 50 & exit'
mv mos2-grp.pi ../ffov_spectrum/mos2-grp.pi
cp ${M2}-bkg.pi ../ffov_spectrum/mos2-bkg.pi
cp ${M2}.rmf ../ffov_spectrum/mos2.rmf
cp ${M2}.arf ../ffov_spectrum/mos2.arf
cd ..

cd ffov_${PN}_0

protonscale mode=1 maskfile={$PN}-fovimspdet.fits specfile={$PN}-fovt.pi

grppha ${PN}-fovtootsub.pi pn0-grp.pi 'chkey BACKFILE pn0-bkg.pi & chkey RESPFILE pn0.rmf &
  chkey ANCRFILE pn0.arf & group min 50 & exit'
mv pn0-grp.pi ../ffov_spectrum/pn0-grp.pi
cp ${PN}-bkg.pi ../ffov_spectrum/pn0-bkg.pi
cp ${PN}.rmf ../ffov_spectrum/pn0.rmf
cp ${PN}.arf ../ffov_spectrum/pn0.arf
cd ..

cd ffov_${PN}_4

protonscale mode=1 maskfile={$PN}-fovimspdet.fits specfile={$PN}-fovt.pi

grppha ${PN}-fovtootsub.pi pn4-grp.pi 'chkey BACKFILE pn4-bkg.pi & chkey RESPFILE pn4.rmf &
  chkey ANCRFILE pn4.arf & group min 50 & exit'
mv pn4-grp.pi ../ffov_spectrum/pn4-grp.pi
cp ${PN}-bkg.pi ../ffov_spectrum/pn4-bkg.pi
cp ${PN}.rmf ../ffov_spectrum/pn4.rmf
cp ${PN}.arf ../ffov_spectrum/pn4.arf
cd ..

# BREAK POINT 2 (see text in the following subsection of the appendix)

# IF you have a large, bright, extended source in your FOV (such as a cluster of galaxies) 
# you now need to refer to the recipe in the following section. If you have a "blank-sky-like"
# observation, proceed to spectral fitting.

# Spectral fitting: this is up to you!

# We now assume that you have fit the soft proton flare spectrum and the swcx spectrum.
# We are now ready to create all of the components of our image from the individual detectors.
# As with the spectra, we create an imaging directory, and then visit all of the ffov_xxxx 
# directories to copy the right stuff to the image directory. There are also a number of other
# pieces to be created for each instrument.

# Building the background images:

# At this point we have the raw images, the *-fovimsky-elo-ehi.fits files, and the QPB images, the 
# *-bkgimdet-elo-ehi.fits files. However the QPB images must be rotated from detector coordinates
# into sky coordinates. In this section we will do that rotation, as well as create and rotate the 
# SPF and SWCX images, one detector at a time.
# As with the spectral components, we sequestion all the image components into one directory

mkdir ffov_image

# copy the cheese masks if you will be using them

cp *cheese?.fits ffov_image

# BREAK POINT 3 (see text in the following subsection of the appendix)

# MOS1

cd ffov_${M1}

# Rotate the QPB image from detector to sky coordinates

rotdet2sky intemplate=${M1}-fovimsky-${ELO}-${EHI}.fits inimage=${M1}-bkgimdet-${ELO}-${EHI}.fits
  outimage=${M1}-bkgimsky-${ELO}-${EHI}.fits withdetxy=false withskyxy=false

# Create and rotate the soft proton flare background image (if you need to). Note that the
# normalization (pnorm) is for the entire FOV, and will be different for each instrument

proton imagefile=${M1}-fovimdet-${ELO}-${EHI}.fits specfile=${M1}-fovt.pi ccds="${M1ON}" elow=${ELO}
  ehigh=${EHI} speccontrol=1 pindex=1.26191 pnorm=2.208e-2 > proton_mos1.log
rotdet2sky intemplate=${M1}-fovimsky-${ELO}-${EHI}.fits inimage=${M1}-protimdet-${ELO}-${EHI}.fits
  outimage=${M1}-protimsky-${ELO}-${EHI}.fits withdetxy=false withskyxy=false

# Create and rotate the swcx background image (if you need to). Note that the normalization here is
# per square arcminute, and is the same for all instruments

swcx imagefile=${M1}-fovimdet-${ELO}-${EHI}.fits specfile=${M1}-fovt.pi ccds="${M1ON}" elow=${ELO}
  ehigh=${EHI} rmffile=${M1}.rmf arffile=${M1}.arf lines="OVII OVIII" 
  gnorms="1.69938e-7 1.75262e-8" > swcx_mos1.log
rotdet2sky intemplate=${M1}-fovimsky-${ELO}-${EHI}.fits inimage=${M1}-swcximdet-${ELO}-${EHI}.fits
  outimage=${M1}-swcximsky-${ELO}-${EHI}.fits withdetxy=false withskyxy=false

# Copy all the useful stuff to the image directory

cp ${M1}-fovimsky-${ELO}-${EHI}.fits ../ffov_image
cp ${M1}-bkgimsky-${ELO}-${EHI}.fits ../ffov_image
cp ${M1}-maskimsky-${ELO}-${EHI}.fits ../ffov_image
cp ${M1}-expimsky-${ELO}-${EHI}.fits ../ffov_image
cp ${M1}-protimsky-${ELO}-${EHI}.fits ../ffov_image
cp ${M1}-swcximsky-${ELO}-${EHI}.fits ../ffov_image
mv *.log ../logs
cd ..

# MOS2

cd ffov_${M2}

# Rotate the QPB image made by mosback into sky coordinates

rotdet2sky intemplate=${M2}-fovimsky-${ELO}-${EHI}.fits inimage=${M2}-bkgimdet-${ELO}-${EHI}.fits
  outimage=${M2}-bkgimsky-${ELO}-${EHI}.fits withdetxy=false withskyxy=false

# Create and rotate the soft proton flare background image (if you need to).

proton imagefile=${M2}-fovimdet-${ELO}-${EHI}.fits specfile=${M2}-fovt.pi ccds="${M2ON}" elow=${ELO}
  ehigh=${EHI} speccontrol=1 pindex=0.1 pnorm=2.546e-2 > proton_mos2.log
rotdet2sky intemplate=${M2}-fovimsky-${ELO}-${EHI}.fits inimage=${M2}-protimdet-${ELO}-${EHI}.fits
  outimage=${M2}-protimsky-${ELO}-${EHI}.fits withdetxy=false withskyxy=false

# Create and rotate the swcx background image (if you need to).

swcx imagefile=${M2}-fovimdet-${ELO}-${EHI}.fits specfile=${M2}-fovt.pi ccds="${M2ON}" elow=${ELO}
  ehigh=${EHI} rmffile=${M2}.rmf arffile=${M2}.arf lines="OVII OVIII" 
  gnorms="1.69938e-7 1.75262e-8" > swcx_mos2.log
rotdet2sky intemplate=${M2}-fovimsky-${ELO}-${EHI}.fits inimage=${M2}-swcximdet-${ELO}-${EHI}.fits
  outimage=${M2}-swcximsky-${ELO}-${EHI}.fits withdetxy=false withskyxy=false

# Copy all the useful stuff to the image directory

cp ${M2}-fovimsky-${ELO}-${EHI}.fits ../ffov_image
cp ${M2}-bkgimsky-${ELO}-${EHI}.fits ../ffov_image
cp ${M2}-maskimsky-${ELO}-${EHI}.fits ../ffov_image
cp ${M2}-expimsky-${ELO}-${EHI}.fits ../ffov_image
cp ${M2}-protimsky-${ELO}-${EHI}.fits ../ffov_image
cp ${M2}-swcximsky-${ELO}-${EHI}.fits ../ffov_image
mv *.log ../logs
cd ..

# PN_0 (assuming that you want the soft image. Use PN_4 if you want the hard image.)

cd ffov_${PN}_0

# Rotate the QPB image made by mosback into sky coordinates

rotdet2sky intemplate=${PN}-fovimsky-${ELO}-${EHI}.fits inimage=${PN}-bkgimdet-${ELO}-${EHI}.fits
  outimage=${PN}-bkgimsky-${ELO}-${EHI}.fits withdetxy=false withskyxy=false

# Create and rotate the soft proton flare background image (if you need to).

proton imagefile=${PN}-fovimdet-${ELO}-${EHI}.fits specfile=${PN}-fovt.pi ccds="${PNON}" elow=${ELO}
  ehigh=${EHI} speccontrol=1 pindex=0.1 pnorm=5.479e-2 > proton_pn.log
rotdet2sky intemplate=${PN}-fovimsky-${ELO}-${EHI}.fits inimage=${PN}-protimdet-${ELO}-${EHI}.fits
  outimage=${PN}-protimsky-${ELO}-${EHI}.fits withdetxy=false withskyxy=false

# Create and rotate the swcx background image (if you need to).

swcx imagefile=${PN}-fovimdet-${ELO}-${EHI}.fits specfile=${PN}-fovt.pi ccds="${PNON}" elow=${ELO}
  ehigh=${EHI} rmffile=${PN}.rmf arffile=${PN}.arf lines="OVII OVIII" 
  gnorms="1.69938e-7 1.75262e-8" > swcx_mos2.log
rotdet2sky intemplate=${PN}-fovimsky-${ELO}-${EHI}.fits inimage=${PN}-swcximdet-${ELO}-${EHI}.fits
  outimage=${PN}-swcximsky-${ELO}-${EHI}.fits withdetxy=false withskyxy=false

# Copy all the useful stuff to the image directory

cp ${PN}-fovimsky-${ELO}-${EHI}.fits ../ffov_image
cp ${PN}-fovimootsky-${ELO}-${EHI}.fits ../ffov_image
cp ${PN}-bkgimsky-${ELO}-${EHI}.fits ../ffov_image
cp ${PN}-maskimsky-${ELO}-${EHI}.fits ../ffov_image
cp ${PN}-expimsky-${ELO}-${EHI}.fits ../ffov_image
cp ${PN}-protimsky-${ELO}-${EHI}.fits ../ffov_image
cp ${PN}-swcximsky-${ELO}-${EHI}.fits ../ffov_image
mv *.log ../logs
cd ..

# Now, combine all of the image components

cd ffov_image
ls -l *.* > pre_comb.list

combimage prefixlist='1S001 2S006 S005' withpartbkg=true withspbkg=true withswcxbkg=true 
  withcheese=true cheesetype=t elowlist=${ELO} ehighlist=${EHI} alpha=1.7 >& combimage.log

ls -l *.* > post_comb.list
mv *.list logs

#Adaptive smoothing or binning? You can do either or both

# Just adaptive smoothing...

binadapt prefix=comb elow=${ELO} ehigh=${EHI} withpartbkg=true withswcxbkg=false withspbkg=true
  withmask=false withbinning=false binfactor=1 withsmoothing=true smoothcounts=50

# ...just binning...

binadapt prefix=comb elow=${ELO} ehigh=${EHI} withpartbkg=true withswcxbkg=true withspbkg=true
  withmask=false withbinning=true binfactor=4 withsmoothing=false smoothcounts=50

# ...or both!

binadapt prefix=comb elow=${ELO} ehigh=${EHI} withpartbkg=true withswcxbkg=true withspbkg=true
  withmask=false withbinning=true binfactor=4 withsmoothing=true smoothcounts=50

ls -l *.* > post_adapt.list
mv *.list logs

# I will note that this leaves the main directory clean to extract another region if desired. 
# You just need to name it something other than "ffov".
# Note the ESAS does not require this type of directory structure, but this directory structure
# has developed with the new ESAS as a way of keeping things neat, and keeping things from
# becoming too confusing.