aboutsummaryrefslogtreecommitdiffstats
path: root/scripts/ciao_bkg_spectra.sh
diff options
context:
space:
mode:
Diffstat (limited to 'scripts/ciao_bkg_spectra.sh')
-rwxr-xr-xscripts/ciao_bkg_spectra.sh156
1 files changed, 109 insertions, 47 deletions
diff --git a/scripts/ciao_bkg_spectra.sh b/scripts/ciao_bkg_spectra.sh
index 4cb8475..96560e6 100755
--- a/scripts/ciao_bkg_spectra.sh
+++ b/scripts/ciao_bkg_spectra.sh
@@ -1,8 +1,6 @@
-#!/bin/sh -
+#!/bin/sh
#
-trap date INT
unalias -a
-export GREP_OPTIONS=""
export LC_COLLATE=C
###########################################################
## extract background spectra from src and blanksky ##
@@ -15,14 +13,22 @@ export LC_COLLATE=C
## Ref: CIAO v4.4 region bugs ##
## http://cxc.harvard.edu/ciao/bugs/regions.html#bug-12187
## ##
-## LIweitiaNux, July 24, 2012 ##
+## Weitian LI ##
+## 2012/07/24 ##
###########################################################
-###########################################################
+VERSION="v5.0"
+UPDATED="2015/06/02"
# ChangeLogs:
-# v3, 2012/08/09
-# fix `scientific notation' for `bc'
-# change `spec group' method to `min 15'
+# v5.0, 2015/06/02, Aaron LI
+# * Removed 'GREP_OPTIONS' and replace 'grep' with '\grep'
+# * Removed 'trap INT date'
+# * Copy needed pfiles to current working directory, and
+# set environment variable $PFILES to use these first.
+# * replaced 'grppha' with 'dmgroup' to group spectra
+# (dmgroup will add history to fits file, while grppha NOT)
+# v4.1, 2014/07/29, Weitian LI
+# fix 'pbkfile' parameters for CIAO-4.6
# v4, 2012/08/13
# add `clobber=yes'
# improve error code
@@ -31,22 +37,22 @@ export LC_COLLATE=C
# (through cmdline which similar to CIAO,
# and default filename match patterns)
# add simple `logging' function
-# v4.1, 2014/07/29, Weitian LI
-# fix 'pbkfile' parameters for CIAO-4.6
-###########################################################
+# v3, 2012/08/09
+# fix `scientific notation' for `bc'
+# change `spec group' method to `min 15'
+#
-## about, used in `usage' {{{
-VERSION="v4"
-UPDATE="2012-08-14"
-## about }}}
## usage, help {{{
case "$1" in
-[hH]*|--[hH]*)
printf "usage:\n"
- printf " `basename $0` evt=<evt2_clean> reg=<reglist> blank=<blanksky_evt> basedir=<base_dir> nh=<nH> z=<redshift> [ grpcmd=<grppha_cmd> log=<log_file> ]\n"
+ printf " `basename $0` evt=<evt2_clean> reg=<reglist> blank=<blanksky_evt> basedir=<base_dir> nh=<nH> z=<redshift> [ grouptype=<NUM_CTS|BIN> grouptypeval=<number> binspec=<binspec> log=<log_file> ]\n"
+ printf "\nNotes:\n"
+ printf " If grouptype=NUM_CTS, then grouptypeval required.\n"
+ printf " If grouptype=BIN, then binspec required.\n"
printf "\nversion:\n"
- printf "${VERSION}, ${UPDATE}\n"
+ printf " ${VERSION}, ${UPDATED}\n"
exit ${ERR_USG}
;;
esac
@@ -55,16 +61,18 @@ esac
## default parameters {{{
# default `event file' which used to match `blanksky' files
#DFT_EVT="_NOT_EXIST_"
-DFT_EVT="`ls evt2*_clean.fits`"
+DFT_EVT="`\ls evt2*_clean.fits`"
# default `blanksky file'
#DFT_BLANK="_NOT_EXIST_"
-DFT_BLANK="`ls blanksky*.fits`"
+DFT_BLANK="`\ls blanksky*.fits`"
# default dir which contains `asols, asol.lis, ...' files
#DFT_BASEDIR="_NOT_EXIST_"
DFT_BASEDIR=".."
-# default `group command' for `grppha'
-#DFT_GRP_CMD="group 1 128 2 129 256 4 257 512 8 513 1024 16"
-DFT_GRP_CMD="group min 20"
+# default parameters for 'dmgroup'
+DFT_GROUPTYPE="NUM_CTS"
+DFT_GROUPTYPEVAL="20"
+#DFT_GROUPTYPE="BIN"
+DFT_BINSPEC="1:128:2,129:256:4,257:512:8,513:1024:16"
# default `log file'
DFT_LOGFILE="bkg_spectra_`date '+%Y%m%d'`.log"
@@ -115,7 +123,7 @@ CH_LOW=651
CH_HI=822
pb_flux() {
punlearn dmstat
- COUNTS=`dmstat "$1[channel=${CH_LOW}:${CH_HI}][cols COUNTS]" | grep -i 'sum:' | awk '{ print $2 }'`
+ COUNTS=`dmstat "$1[channel=${CH_LOW}:${CH_HI}][cols COUNTS]" | \grep -i 'sum:' | awk '{ print $2 }'`
punlearn dmkeypar
EXPTIME=`dmkeypar $1 EXPOSURE echo=yes`
BACK=`dmkeypar $1 BACKSCAL echo=yes`
@@ -222,13 +230,28 @@ fi
# remove the trailing '/'
BASEDIR=`echo ${BASEDIR} | sed 's/\/*$//'`
printf "## use basedir: \`${BASEDIR}'\n" | ${TOLOG}
-# check given `grpcmd'
-if [ ! -z "${grpcmd}" ]; then
- GRP_CMD="${grpcmd}"
+# check given dmgroup parameters: grouptype, grouptypeval, binspec
+if [ -z "${grouptype}" ]; then
+ GROUPTYPE="${DFT_GROUPTYPE}"
+elif [ "x${grouptype}" = "xNUM_CTS" ] || [ "x${grouptype}" = "xBIN" ]; then
+ GROUPTYPE="${grouptype}"
+else
+ printf "ERROR: given grouptype \`${grouptype}' invalid.\n"
+ exit ${ERR_GRPTYPE}
+fi
+printf "## use grouptype: \`${GROUPTYPE}'\n" | ${TOLOG}
+if [ ! -z "${grouptypeval}" ]; then
+ GROUPTYPEVAL="${grouptypeval}"
else
- GRP_CMD="${DFT_GRP_CMD}"
+ GROUPTYPEVAL="${DFT_GROUPTYPEVAL}"
fi
-printf "## use grppha cmd: \`${GRP_CMD}'\n" | ${TOLOG}
+printf "## use grouptypeval: \`${GROUPTYPEVAL}'\n" | ${TOLOG}
+if [ ! -z "${binspec}" ]; then
+ BINSPEC="${binspec}"
+else
+ BINSPEC="${DFT_BINSPEC}"
+fi
+printf "## use binspec: \`${BINSPEC}'\n" | ${TOLOG}
## parameters }}}
## check needed files {{{
@@ -242,7 +265,7 @@ for reg_f in ${REGLIST}; do
done
# check the validity of *pie* regions
printf "check pie reg validity ...\n"
-INVALID=`cat ${REGLIST} | grep -i 'pie' | awk -F, '{ print $6 }' | tr -d ')' | awk '$1 > 360'`
+INVALID=`cat ${REGLIST} | \grep -i 'pie' | awk -F, '{ print $6 }' | tr -d ')' | awk '$1 > 360'`
if [ "x${INVALID}" != "x" ]; then
printf "WARNING: some pie region's END_ANGLE > 360\n" | ${TOLOG}
printf " CIAO v4.4 tools may run into trouble\n"
@@ -251,28 +274,28 @@ fi
# check files in `basedir'
printf "check needed files in basedir \`${BASEDIR}' ...\n"
# check asolis files
-ASOLIS=`ls -1 ${BASEDIR}/${DFT_ASOLIS_PAT} | head -n 1`
+ASOLIS=`\ls -1 ${BASEDIR}/${DFT_ASOLIS_PAT} | head -n 1`
if [ -z "${ASOLIS}" ]; then
printf "ERROR: cannot find \"${DFT_ASOLIS_PAT}\" in dir \`${BASEDIR}'\n"
exit ${ERR_ASOL}
fi
printf "## use asolis: \`${ASOLIS}'\n" | ${TOLOG}
# check badpixel file
-BPIX=`ls -1 ${BASEDIR}/${DFT_BPIX_PAT} | head -n 1`
+BPIX=`\ls -1 ${BASEDIR}/${DFT_BPIX_PAT} | head -n 1`
if [ -z "${BPIX}" ]; then
printf "ERROR: cannot find \"${DFT_BPIX_PAT}\" in dir \`${BASEDIR}'\n"
exit ${ERR_BPIX}
fi
printf "## use badpixel: \`${BPIX}'\n" | ${TOLOG}
# check pbk file
-PBK=`ls -1 ${BASEDIR}/${DFT_PBK_PAT} | head -n 1`
+PBK=`\ls -1 ${BASEDIR}/${DFT_PBK_PAT} | head -n 1`
if [ -z "${PBK}" ]; then
printf "ERROR: cannot find \"${DFT_PBK_PAT}\" in dir \`${BASEDIR}'\n"
exit ${ERR_PBK}
fi
printf "## use pbk: \`${PBK}'\n" | ${TOLOG}
# check msk file
-MSK=`ls -1 ${BASEDIR}/${DFT_MSK_PAT} | head -n 1`
+MSK=`\ls -1 ${BASEDIR}/${DFT_MSK_PAT} | head -n 1`
if [ -z "${MSK}" ]; then
printf "ERROR: cannot find \"${DFT_MSK_PAT}\" in dir \`${BASEDIR}'\n"
exit ${ERR_MSK}
@@ -280,6 +303,19 @@ fi
printf "## use msk: \`${MSK}'\n" | ${TOLOG}
## check files }}}
+## prepare parameter files (pfiles) {{{
+CIAO_TOOLS="dmstat dmkeypar dmhedit specextract dmextract dmgroup"
+
+# Copy necessary pfiles for localized usage
+for tool in ${CIAO_TOOLS}; do
+ pfile=`paccess ${tool}`
+ [ -n "${pfile}" ] && punlearn ${tool} && cp -Lvf ${pfile} .
+done
+
+# Modify environment variable 'PFILES' to use local pfiles first
+export PFILES="./:${PFILES}"
+## pfiles }}}
+
## main part {{{
## use 'for' loop to process every region file
for reg_i in ${REGLIST}; do
@@ -290,7 +326,7 @@ for reg_i in ${REGLIST}; do
[ -f "${REG_TMP}" ] && rm -fv ${REG_TMP} # remove tmp files
cp -fv ${reg_i} ${REG_TMP}
# check the validity of *pie* regions {{{
- INVALID=`grep -i 'pie' ${REG_TMP} | awk -F, '{ print $6 }' | tr -d ')' | awk '$1 > 360'`
+ INVALID=`\grep -i 'pie' ${REG_TMP} | awk -F, '{ print $6 }' | tr -d ')' | awk '$1 > 360'`
if [ "x${INVALID}" != "x" ]; then
printf "WARNING: fix for *pie* region in file \`${reg_i}'\n"
cat ${REG_TMP}
@@ -303,30 +339,52 @@ for reg_i in ${REGLIST}; do
## check pie region }}}
LBKG_PI="${reg_i%.reg}.pi"
+
+ printf "use \`specextract' to generate spectra and response ...\n"
## use `specextract' to extract local bkg spectrum {{{
# NOTE: set `binarfwmap=2' to save the time for generating `ARF'
# I have tested that this bin factor has little impact on the results.
# NO background response files
# NO background spectrum (generate by self)
- # NO spectrum grouping (group by self using `grppha')
+ # NO spectrum grouping (group by self using `dmgroup')
+ # Determine parameters for different versions of specextract {{{
# 'pbkfile' parameter deprecated in CIAO-4.6
if `pget specextract pbkfile >/dev/null 2>&1`; then
P_PBKFILE="pbkfile=${PBK}"
else
P_PBKFILE=""
fi
- #
- printf "use \`specextract' to generate spectra and response ...\n"
+ # specextract: revision 2013-06:
+ # 'correct' parameter renamed to 'correctpsf'
+ if `pget specextract correct >/dev/null 2>&1`; then
+ P_CORRECT="correct=no"
+ else
+ P_CORRECT="correctpsf=no"
+ fi
+ # specextract: revision 2013-12:
+ # 'weight' parameter controls whether ONLY ARFs are weighted.
+ # 'weight_rmf' added to control whether RMFs are weighted.
+ # NOTE:
+ # (1) only when 'weight=yes' will the 'weight_rmf' parameter be used.
+ # (2) no longer distingush between unweighted and weighted reponses
+ # in file extension; only .arf & .rmf are now used.
+ if `pget specextract weight_rmf >/dev/null 2>&1`; then
+ P_WEIGHTRMF="weight_rmf=yes"
+ else
+ P_WEIGHTRMF=""
+ fi
+ # }}}
punlearn specextract
specextract infile="${EVT}[sky=region(${REG_TMP})]" \
outroot=${LBKG_PI%.pi} bkgfile="" asp="@${ASOLIS}" \
- ${P_PBKFILE} mskfile="${MSK}" badpixfile="${BPIX}" \
- weight=yes correct=no bkgresp=no \
+ mskfile="${MSK}" badpixfile="${BPIX}" \
+ ${P_PBKFILE} ${P_CORRECT} \
+ weight=yes ${P_WEIGHTRMF} \
energy="0.3:11.0:0.01" channel="1:1024:1" \
- combine=no binarfwmap=2 \
+ bkgresp=no combine=no binarfwmap=2 \
grouptype=NONE binspec=NONE \
- verbose=2 clobber=yes
- # `specextract' }}}
+ clobber=yes verbose=2
+ # specextract }}}
## generate the blanksky bkg spectrum {{{
printf "generate the blanksky bkg spectrum ...\n"
@@ -342,11 +400,16 @@ for reg_i in ${REGLIST}; do
## bkg renorm }}}
## group spectrum {{{
- printf "group spectrum using \`grppha'\n"
+ # use 'dmgroup' instead of 'grppha', because 'dmgroup' will add
+ # command history to FITS header (maybe useful for later reference).
+ printf "group spectrum using \`dmgroup'\n"
LBKG_GRP_PI="${LBKG_PI%.pi}_grp.pi"
- grppha infile="${LBKG_PI}" outfile="${LBKG_GRP_PI}" \
- comm="${GRP_CMD} & exit" clobber=yes > /dev/null
- ## group spectra }}}
+ punlearn dmgroup
+ dmgroup infile="${LBKG_PI}" outfile="${LBKG_GRP_PI}" \
+ grouptype="${GROUPTYPE}" grouptypeval=${GROUPTYPEVAL} \
+ binspec="${BINSPEC}" xcolumn="CHANNEL" ycolumn="COUNTS" \
+ clobber=yes
+ ## group }}}
## generate a script for XSPEC {{{
XSPEC_XCM="xspec_${LBKG_PI%.pi}_model.xcm"
@@ -420,6 +483,5 @@ printf "clean ...\n"
rm -f ${REG_TMP}
printf "DONE\n"
-###########################################################
# vim: set ts=8 sw=4 tw=0 fenc=utf-8 ft=sh: #