ciao_deproj_spectra.sh: Cleanup and simplification

Use newly developed tools:
+ manifest.py
+ results.py
+ renorm_spectrum.py

2 files changed, 23 insertions, 155 deletions

diff --git a/scripts/ciao_bkg_spectra.sh b/scripts/ciao_bkg_spectra.sh
index 3b2fb22..07a4bc9 100755
--- a/scripts/ciao_bkg_spectra.sh
+++ b/scripts/ciao_bkg_spectra.sh
@@ -45,7 +45,7 @@
 case "$1" in
     -[hH]*|--[hH]*)
         printf "usage:\n"
-        printf "    `basename $0` evt=<evt2_clean> reg=<reglist> blank=<blanksky_evt> nh=<nH> z=<redshift> [ grouptype=<NUM_CTS|BIN> grouptypeval=<number> binspec=<binspec> log=<log_file> ]\n"
+        printf "    `basename $0` reg=<reglist> [ evt=<evt2_clean> blank=<blanksky_evt> 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"
@@ -155,14 +155,14 @@ else
 fi
 # check given nH
 if [ -z "${nh}" ]; then
-    N_H=$(results.py get nh)
+    N_H=$(results.py -b get nh)
 else
     N_H=${nh}
 fi
 printf "## use nH: ${N_H}\n" | ${TOLOG}
 # check given redshift
 if [ -z "${z}" ]; then
-    REDSHIFT=$(results.py get z)
+    REDSHIFT=$(results.py -b get z)
 else
     REDSHIFT=${z}
 fi
diff --git a/scripts/ciao_deproj_spectra.sh b/scripts/ciao_deproj_spectra.sh
index 1ec5ed6..e134fde 100755
--- a/scripts/ciao_deproj_spectra.sh
+++ b/scripts/ciao_deproj_spectra.sh
@@ -35,12 +35,14 @@
 ## [4] CIAO 4.6 Release Notes
 ##     http://cxc.cfa.harvard.edu/ciao/releasenotes/ciao_4.6_release.html
 ##
-AUTHOR="Weitian LI <liweitianux@gmail.com>"
-CREATED="2012/07/24"
-UPDATED="2015/06/03"
-VERSION="v10.0"
+## Author: Weitian LI <liweitianux@gmail.com>
+## Created: 2012-07-24
 ##
-## ChangeLogs:
+## Change logs:
+## 2017-02-26, Weitian LI
+##   * Use 'manifest.py', 'results.py', and 'renorm_spectrum.py'
+##   * Simplify 'specextract' parameters handling
+##   * Other cleanups
 ## v10.0, 2015/06/03, Aaron LI
 ##   * Copy needed pfiles to current working directory, and
 ##     set environment variable $PFILES to use these first.
@@ -81,19 +83,14 @@ VERSION="v10.0"
 ##   * added `background renormalization'
 ##
 
-unalias -a
-export LC_COLLATE=C
-
 ## usage, help {{{
 case "$1" in
     -[hH]*|--[hH]*)
         printf "Usage:\n"
-        printf "    `basename $0` evt=<evt2_clean> reg=<radial_reg> bkgd=<blank_evt | lbkg_reg | bkg_spec> basedir=<base_dir> nh=<nH> z=<redshift> [ grouptype=<NUM_CTS|BIN> grouptypeval=<number> binspec=<binspec> log=<log_file> ]\n"
+        printf "    `basename $0` reg=<radial_reg> [ bkgd=<blank_evt | lbkg_reg | bkg_spec> evt=<evt2_clean> 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}, ${UPDATED}\n"
         exit ${ERR_USG}
         ;;
 esac
@@ -101,14 +98,9 @@ esac
 
 ## default parameters {{{
 # default `event file' which used to match `blanksky' files
-#DFT_EVT="_NOT_EXIST_"
-DFT_EVT="`\ls evt2*_clean.fits`"
+DFT_EVT=$(manifest.py -b getpath -r evt2_clean)
 # default `radial region file'
-#DFT_REG_IN="_NOT_EXIST_"
 DFT_REG_IN="rspec.reg"
-# default dir which contains `asols, asol.lis, ...' files
-# DFT_BASEDIR=".."
-DFT_BASEDIR="_NOT_EXIST_"
 # default parameters for 'dmgroup'
 DFT_GROUPTYPE="NUM_CTS"
 DFT_GROUPTYPEVAL="20"
@@ -121,18 +113,9 @@ DFT_LOGFILE="deproj_spectra_`date '+%Y%m%d'`.log"
 XSPEC_DEPROJ="xspec_deproj.xcm"
 XSPEC_PROJTD="xspec_projected.xcm"
 
-## howto find files in `basedir'
-# default `asol.lis pattern'
-DFT_ASOLIS_PAT="acis*asol?.lis"
-# default `bad pixel filename pattern'
-DFT_BPIX_PAT="acis*repro*bpix?.fits"
-# default `pbk file pattern'
-DFT_PBK_PAT="acis*pbk?.fits"
-# default `msk file pattern'
-DFT_MSK_PAT="acis*msk?.fits"
-
-## abundance standard
-ABUND="grsa"
+ASOL=$(manifest.py -b -s "," getpath -r asol)
+BPIX=$(manifest.py -b getpath -r bpix)
+MSK=$(manifest.py -b getpath -r msk)
 ## default parameters }}}
 
 ## error code {{{
@@ -164,39 +147,6 @@ getopt_keyval() {
         shift                           # shift, process next one
     done
 }
-
-## background renormalization (BACKSCAL) {{{
-# renorm background according to particle background
-# energy range: 9.5-12.0 keV (channel: 651-822)
-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 }'`
-    punlearn dmkeypar
-    EXPTIME=`dmkeypar $1 EXPOSURE echo=yes`
-    BACK=`dmkeypar $1 BACKSCAL echo=yes`
-    # fix `scientific notation' bug for `bc'
-    EXPTIME_B=`echo ${EXPTIME} | sed 's/[eE]/\*10\^/' | sed 's/+//'`
-    BACK_B=`echo "( ${BACK} )" | sed 's/[eE]/\*10\^/' | sed 's/+//'`
-    PB_FLUX=`echo "scale = 16; ${COUNTS} / ${EXPTIME_B} / ${BACK_B}" | bc -l`
-    echo ${PB_FLUX}
-}
-
-bkg_renorm() {
-    # $1: src spectrum, $2: back spectrum
-    PBFLUX_SRC=`pb_flux $1`
-    PBFLUX_BKG=`pb_flux $2`
-    BACK_OLD=`dmkeypar $2 BACKSCAL echo=yes`
-    BACK_OLD_B=`echo "( ${BACK_OLD} )" | sed 's/[eE]/\*10\^/' | sed 's/+//'`
-    BACK_NEW=`echo "scale = 16; ${BACK_OLD_B} * ${PBFLUX_BKG} / ${PBFLUX_SRC}" | bc -l`
-    printf "\`$2': BACKSCAL:\n"
-    printf "    ${BACK_OLD} --> ${BACK_NEW}\n"
-    punlearn dmhedit
-    dmhedit infile=$2 filelist=none operation=add \
-        key=BACKSCAL value=${BACK_NEW} comment="old value: ${BACK_OLD}"
-}
-## bkg renorm }}}
 ## functions end }}}
 
 ## parameters {{{
@@ -212,8 +162,6 @@ fi
 printf "## use logfile: \`${LOGFILE}'\n"
 [ -e "${LOGFILE}" ] && mv -fv ${LOGFILE} ${LOGFILE}_bak
 TOLOG="tee -a ${LOGFILE}"
-echo "process script: `basename $0`" >> ${LOGFILE}
-echo "process date: `date`" >> ${LOGFILE}
 ## log }}}
 
 # check given parameters
@@ -222,12 +170,6 @@ if [ -r "${evt}" ]; then
     EVT=${evt}
 elif [ -r "${DFT_EVT}" ]; then
     EVT=${DFT_EVT}
-else
-    read -p "clean evt2 file: " EVT
-    if [ ! -r "${EVT}" ]; then
-        printf "ERROR: cannot access given \`${EVT}' evt file\n"
-        exit ${ERR_EVT}
-    fi
 fi
 printf "## use evt file: \`${EVT}'\n" | ${TOLOG}
 # check given region file(s)
@@ -245,7 +187,7 @@ fi
 printf "## use radial reg: \`${REG_IN}'\n" | ${TOLOG}
 # check given bkgd, determine background {{{
 if [ -z "${bkgd}" ]; then
-    read -p "> background (blanksky_evt | lbkg_reg | bkg_spec): " BKGD
+    BKGD=$(manifest.py -b getpath -r bkg_spec)
 else
     BKGD=${bkgd}
 fi
@@ -309,33 +251,18 @@ fi
 # bkgd }}}
 # check given nH
 if [ -z "${nh}" ]; then
-    read -p "> value of nH: " N_H
+    N_H=$(results.py -b get nh)
 else
     N_H=${nh}
 fi
 printf "## use nH: ${N_H}\n" | ${TOLOG}
 # check given redshift
 if [ -z "${z}" ]; then
-    read -p "> value of redshift: " REDSHIFT
+    REDSHIFT=$(results.py -b get z)
 else
     REDSHIFT=${z}
 fi
 printf "## use redshift: ${REDSHIFT}\n" | ${TOLOG}
-# check given dir
-if [ -d "${basedir}" ]; then
-    BASEDIR=${basedir}
-elif [ -d "${DFT_BASEDIR}" ]; then
-    BASEDIR=${DFT_BASEDIR}
-else
-    read -p "> basedir (contains asol files): " BASEDIR
-    if [ ! -d "${BASEDIR}" ]; then
-        printf "ERROR: given \`${BASEDIR}' NOT a directory\n"
-        exit ${ERR_DIR}
-    fi
-fi
-# remove the trailing '/'
-BASEDIR=`echo ${BASEDIR} | sed 's/\/*$//'`
-printf "## use basedir: \`${BASEDIR}'\n" | ${TOLOG}
 # check given dmgroup parameters: grouptype, grouptypeval, binspec
 if [ -z "${grouptype}" ]; then
     GROUPTYPE="${DFT_GROUPTYPE}"
@@ -371,37 +298,6 @@ 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"
 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`
-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`
-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`
-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`
-if [ -z "${MSK}" ]; then
-    printf "ERROR: cannot find \"${DFT_MSK_PAT}\" in dir \`${BASEDIR}'\n"
-    exit ${ERR_MSK}
-fi
-printf "## use msk: \`${MSK}'\n" | ${TOLOG}
 ## check files }}}
 
 ## prepare parameter files (pfiles) {{{
@@ -490,39 +386,11 @@ for i in `seq ${LINES}`; do
     # NO background response files
     # NO background spectrum (generate by self)
     # 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
-    # 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_CIAO})]" \
-        outroot="r${i}_${ROOTNAME}" bkgfile="" asp="@${ASOLIS}" \
+        outroot="r${i}_${ROOTNAME}" bkgfile="" asp="${ASOL}" \
         mskfile="${MSK}" badpixfile="${BPIX}" \
-        ${P_PBKFILE} ${P_CORRECT} \
-        weight=yes ${P_WEIGHTRMF} \
+        correctpsf=no weight=yes weight_rmf=yes \
         energy="0.3:11.0:0.01" channel="1:1024:1" \
         bkgresp=no combine=no binarfwmap=2 \
         grouptype=NONE binspec=NONE \
@@ -547,7 +415,7 @@ for i in `seq ${LINES}`; do
 
     ## bkg renormalization {{{
     printf "Renormalize background ...\n"
-    bkg_renorm ${RSPEC_PI} ${RSPEC_BKG_PI}
+    renorm_spectrum.py -r ${RSPEC_PI} ${RSPEC_BKG_PI}
     ## bkg renorm }}}
 
     ## group spectrum {{{
@@ -663,7 +531,7 @@ ignore **:0.0-0.7,7.0-**
 
 method leven 1000 0.01
 # change abundance standard
-abund ${ABUND}
+abund grsa
 xsect bcmc
 cosmo 70 0 0.73
 xset delta 0.01
@@ -741,7 +609,7 @@ ignore **:0.0-0.7,7.0-**
 
 method leven 1000 0.01
 # change abundance standard
-abund ${ABUND}
+abund grsa
 xsect bcmc
 cosmo 70 0 0.73
 xset delta 0.01