Initial commit

Added files:
  * mass_profile: developed by Junhua GU, modified by Weitian LI
  * opt_utilities: developed by Junhua GU
  * tools/cosmo_calc: originated from 'calc_distance', modified
  * scripts: scripts used to process Chandra ACIS data
  * files: useful files used in processing
  * HOWTO_chandra_acis_process.txt
  * README.md

1 files changed, 529 insertions, 0 deletions

diff --git a/scripts/ciao_sbp_v3.1.sh b/scripts/ciao_sbp_v3.1.sh
new file mode 100755
index 0000000..d111b7f
--- /dev/null
+++ b/scripts/ciao_sbp_v3.1.sh
@@ -0,0 +1,529 @@
+#!/bin/sh
+#
+unalias -a
+export LC_COLLATE=C
+###########################################################
+## extract `surface brighness profile'                   ##
+##                                                       ##
+## NOTES:                                                ##
+## only ACIS-I (chip: 0-3) and ACIS-S (chip: 7) supported##
+## `merge_all' conflict with Heasoft's `pget' etc. tools ##
+##                                                       ##
+## LIweitiaNux <liweitianux@gmail.com>                   ##
+## August 16, 2012                                       ##
+###########################################################
+
+###########################################################
+## Changes made by zzh: (2013-02-01)
+## removes the region in ccd gap of ACIS_I 
+## removes the region in the area of point source
+## need asol file to prevent offset 
+###########################################################
+
+SCRIPT_PATH=`readlink -f $0`
+SCRIPT_DIR=`dirname ${SCRIPT_PATH}`
+CCDGAP_SCRIPT="chandra_ccdgap_rect.py"
+
+## about, used in `usage' {{{
+VERSION="v3.1"
+UPDATE="2013-02-01"
+## about }}}
+
+## error code {{{
+ERR_USG=1
+ERR_DIR=11
+ERR_EVT=12
+ERR_BKG=13
+ERR_REG=14
+ERR_CELL_REG=15
+ERR_ASOL=21
+ERR_BPIX=22
+ERR_PBK=23
+ERR_MSK=24
+ERR_BKGTY=31
+ERR_SPEC=32
+ERR_DET=41
+ERR_ENG=42
+ERR_CIAO=100
+## error code }}}
+
+## usage, help {{{
+case "$1" in
+    -[hH]*|--[hH]*)
+        printf "usage:\n"
+        printf "    `basename $0` evt_e=<evt_e_file> reg=<sbp_reg> expmap=<exp_map> cellreg=<cell_reg> aspec=<asol_file> [bkg=<bkg> log=<logfile> ]\n"
+        printf "\nversion:\n"
+        printf "${VERSION}, ${UPDATE}\n"
+        exit ${ERR_USG}
+        ;;
+esac
+## usage, help }}}
+
+## default parameters {{{
+# default energy band
+E_RANGE="700:7000"
+# default ccd edge cut pixel (20pixel)
+DFT_CCD_EDGECUT=25
+# default `event file' which used to match `blanksky' files
+#DFT_EVT="_NOT_EXIST_"
+DFT_EVT_E="`ls evt2_*_e*.fits 2> /dev/null`"
+# default expmap
+DFT_EXPMAP="`ls expmap*.fits 2> /dev/null | head -n 1`"
+# default `radial region file' to extract surface brightness
+#DFT_SBP_REG="_NOT_EXIST_"
+DFT_SBP_REG="sbprofile.reg"
+# defalut pointsource region file 
+DFT_CELL_REG="`ls celld*.reg 2> /dev/null`"
+# defalut asol file
+DFT_ASOL_FILE="`ls pcad*asol*fits 2> /dev/null`"
+
+# default `log file'
+DFT_LOGFILE="sbp_`date '+%Y%m%d'`.log"
+
+## default parameters }}}
+
+## functions {{{
+# process commandline arguments
+# cmdline arg format: `KEY=VALUE'
+getopt_keyval() {
+    until [ -z "$1" ]
+    do
+        key=${1%%=*}                    # extract key
+        val=${1#*=}                     # extract value
+        keyval="${key}=\"${val}\""
+        echo "## getopt: eval '${keyval}'"
+        eval ${keyval}
+        shift                           # shift, process next one
+    done
+}
+## functions }}}
+
+## check CIAO init {{{
+if [ -z "${ASCDS_INSTALL}" ]; then
+    printf "ERROR: CIAO NOT initialized\n"
+    exit ${ERR_CIAO}
+fi
+
+## XXX: heasoft's `pget' etc. tools conflict with some CIAO tools
+printf "set \$PATH to avoid conflicts between HEAsoft and CIAO\n"
+export PATH="${ASCDS_BIN}:${ASCDS_CONTRIB}:${PATH}"
+printf "## PATH: ${PATH}\n"
+## check CIAO }}}
+
+## parameters {{{
+# process cmdline args using `getopt_keyval'
+getopt_keyval "$@"
+
+## check log parameters {{{
+if [ ! -z "${log}" ]; then
+    LOGFILE="${log}"
+else
+    LOGFILE=${DFT_LOGFILE}
+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
+# check evt file
+if [ -r "${evt_e}" ]; then
+    EVT_E=${evt_e}
+elif [ -r "${DFT_EVT_E}" ]; then
+    EVT_E=${DFT_EVT_E}
+else
+    read -p "clean evt2 file: " EVT_E
+    if [ ! -r "${EVT_E}" ]; then
+        printf "ERROR: cannot access given \`${EVT_E}' evt file\n"
+        exit ${ERR_EVT}
+    fi
+fi
+printf "## use evt_eng file: \`${EVT_E}'\n" | ${TOLOG}
+# check evt file
+if [ -r "${expmap}" ]; then
+    EXPMAP=${expmap}
+elif [ -r "${DFT_EXPMAP}" ]; then
+    EXPMAP=${DFT_EXPMAP}
+else
+    read -p "expocure map file: " EXPMAP
+    if [ ! -r "${EXPMAP}" ]; then
+        printf "ERROR: cannot access given \`${EXPMAP}' expmap file\n"
+        exit ${ERR_EVT}
+    fi
+fi
+printf "## use expmap file: \`${EXPMAP}'\n" | ${TOLOG}
+# check given region file(s)
+if [ -r "${reg}" ]; then
+    SBP_REG=${reg}
+elif [ -r "${DFT_SBP_REG}" ]; then
+    SBP_REG=${DFT_SBP_REG}
+else
+    read -p "> surface brighness radial region file: " SBP_REG
+    if [ ! -r "${SBP_REG}" ]; then
+        printf "ERROR: cannot access given \`${SBP_REG}' region file\n"
+        exit ${ERR_REG}
+    fi
+fi
+printf "## use reg file(s): \`${SBP_REG}'\n" | ${TOLOG}
+# check bkg
+if [ ! -z "${bkg}" ]; then
+    BKG=${bkg}
+#else
+#    read -p "> bkg: " BKG
+fi
+if [ -r "${BKG}" ]; then
+    printf "## use bkg: \`${BKG}'\n" | ${TOLOG}
+else
+    BKG="NULL"
+fi
+# check cell region file
+if [ -r "${cellreg}" ]; then 
+    CELL_REG="${cellreg}"
+elif [ -r "${DFT_CELL_REG}" ] ; then
+    CELL_REG="${DFT_CELL_REG}"
+else 
+    read -p "> celldetect region file: " CELL_REG
+    if [ ! -r "${CELL_REG}" ]; then
+        printf "ERROR: cannot access given \`${CELL_REG}' region file \n"
+        exit ${ERR_CELL_REG}
+    fi
+fi
+printf "## use cell reg file(s): \`${CELL_REG}'\n" | ${TOLOG}
+# check asol file
+if [ -r "${aspec}" ]; then
+    ASOL_FILE="${aspec}"
+elif [ -r "${DFT_ASOL_FILE}" ] ; then
+    ASOL_FILE="${DFT_ASOL_FILE}"
+else
+    read -p ">asol file: " ASOL_FILE
+    if [ ! -r "${ASOL_FILE}" ] ; then
+        printf " ERROR: cannot access asol file \n"
+        exit ${ERR_ASOL}
+    fi
+fi
+printf "## use asol file(s) : \`${ASOL_FILE}'\n" | ${TOLOG}
+## parameters }}}
+
+## determine ACIS type {{{
+# consistent with `ciao_procevt'
+punlearn dmkeypar
+DETNAM=`dmkeypar ${EVT_E} DETNAM echo=yes`
+if echo ${DETNAM} | grep -q 'ACIS-0123'; then
+    printf "## \`DETNAM' (${DETNAM}) has chips 0123 -> ACIS-I\n"
+    ACIS_TYPE="I"
+    CCD="0:3"
+elif echo ${DETNAM} | grep -q 'ACIS-[0-6]*7'; then
+    printf "## \`DETNAM' (${DETNAM}) has chip 7 -> ACIS-S\n"
+    ACIS_TYPE="S"
+    CCD="7"
+else
+    printf "ERROR: unknown detector type: ${DETNAM}\n"
+    exit ${ERR_DET}
+fi
+## ACIS type }}}
+
+## check validity of pie region {{{
+INVALID=`grep -i 'pie' ${SBP_REG} | awk -F'[,()]' '$7 > 360'`
+SBP_REG_FIX="_${SBP_REG%.reg}_fix.reg"
+if [ "x${INVALID}" != "x" ]; then
+    printf "*** WARNING: some pie regions' END_ANGLE > 360\n" | ${TOLOG}
+    printf "*** script will fix ...\n"
+    cp -fv ${SBP_REG} ${SBP_REG_FIX}
+    # using `awk' to fix
+    awk -F'[,()]' '{
+        if ($7 > 360) {
+            printf "%s(%.2f,%.2f,%.2f,%.2f,%.2f,%.2f)\n", $1,$2,$3,$4,$5,$6,($7-360)
+        }
+        else {
+            print $0
+        }
+    }' ${SBP_REG} | sed '/^\#/d' > ${SBP_REG_FIX}
+else
+    cat  ${SBP_REG} | sed '/^\#/d' >${SBP_REG_FIX}
+fi
+## pie validity }}}
+
+## main process {{{
+
+## generate `skyfov'
+# XXX: omit `aspec', NOT provide `asol' file
+# otherwise the size of evt_img NOT match the `expmap'
+printf "generate skyfov ...\n"
+SKYFOV="_skyfov.fits"
+punlearn skyfov
+skyfov infile="${EVT_E}" outfile="${SKYFOV}" aspec="${ASOL_FILE}" clobber=yes
+
+## get CCD fov regions {{{
+printf "make regions for CCD ...\n"
+TMP_LIST="_tmp_list.txt"
+TMP_REC="_tmp_rec.reg"
+if [ "${ACIS_TYPE}" = "S" ]; then
+    # ACIS-S
+    punlearn dmlist
+    dmlist infile="${SKYFOV}[ccd_id=${CCD}][cols POS]" opt="data,clean" | awk '{for (i=1;i<=NF;i++) print $i }' |sed -e ':a;N;s/\n/,/;ta' | awk -F"]," '{print "polygon("$2}' | awk -F"NaN" '{print $1}' >${TMP_LIST} 
+    python ${SCRIPT_DIR}/${CCDGAP_SCRIPT} ${TMP_LIST} >${TMP_REC}
+    XC=` cat ${TMP_REC} | awk -F\( '{print $2}' |awk -F\) '{print $1}' |awk -F\, '{print $1}'`
+    YC=` cat ${TMP_REC} | awk -F\( '{print $2}' |awk -F\) '{print $1}' |awk -F\, '{print $2}'`
+    ADD_L=` cat ${TMP_REC} | awk -F\( '{print $2}' |awk -F\) '{print $1}' |awk -F\, '{print $3/2}'`
+    ANG=` cat ${TMP_REC} | awk -F\( '{print $2}' |awk -F\) '{print $1}' |awk -F\, '{print $5}'`
+    while [ 1 -eq 1 ]; do
+        if [ `echo "${ANG} < 0" |bc -l ` -eq 1 ] ; then
+            ANG=` echo " ${ANG} + 90 " | bc -l `
+        elif [ `echo "${ANG} >=90" |bc -l ` -eq 1 ] ; then
+            ANG=` echo " ${ANG} - 90 " | bc -l`
+        else
+            break
+        fi
+    done
+    ANG=`echo "${ANG}/180*3.1415926" |bc -l` 
+    CCD_1_X_RAW=` echo " ${XC}  ${ADD_L} ${ANG} "| awk '{print $1-$2*cos($3)-$2*sin($3)}' ` 
+    CCD_2_X_RAW=` echo " ${XC}  ${ADD_L} ${ANG} "| awk '{print $1+$2*cos($3)-$2*sin($3)}' `
+    CCD_3_X_RAW=` echo " ${XC}  ${ADD_L} ${ANG} "| awk '{print $1-$2*cos($3)+$2*sin($3)}' `
+    CCD_4_X_RAW=` echo " ${XC}  ${ADD_L} ${ANG} "| awk '{print $1+$2*cos($3)+$2*sin($3)}' `
+    CCD_1_Y_RAW=` echo " ${YC}  ${ADD_L} ${ANG} "| awk '{print $1+$2*cos($3)-$2*sin($3)}' `
+    CCD_2_Y_RAW=` echo " ${YC}  ${ADD_L} ${ANG} "| awk '{print $1+$2*cos($3)+$2*sin($3)}' `
+    CCD_3_Y_RAW=` echo " ${YC}  ${ADD_L} ${ANG} "| awk '{print $1-$2*cos($3)-$2*sin($3)}' `
+    CCD_4_Y_RAW=` echo " ${YC}  ${ADD_L} ${ANG} "| awk '{print $1-$2*cos($3)+$2*sin($3)}' `
+    CCD_1_RAW=` echo "${CCD_1_X_RAW},${CCD_1_Y_RAW}"`
+    CCD_2_RAW=` echo "${CCD_2_X_RAW},${CCD_2_Y_RAW}"`
+    CCD_3_RAW=` echo "${CCD_3_X_RAW},${CCD_3_Y_RAW}"`
+    CCD_4_RAW=` echo "${CCD_4_X_RAW},${CCD_4_Y_RAW}"`
+    REG_CCD_RAW="`echo "polygon(${CCD_1_RAW}, ${CCD_2_RAW}, ${CCD_4_RAW}, ${CCD_3_RAW}) " `"
+    DX_2T1=$(echo "`echo ${CCD_2_RAW} | awk -F\, '{print $1}'`-`echo  ${CCD_1_RAW} |awk -F\, '{print $1}'`" |bc -l)
+    DY_2T1=$(echo "`echo ${CCD_2_RAW} | awk -F\, '{print $2}'`-`echo  ${CCD_1_RAW} |awk -F\, '{print $2}'`" |bc -l)
+    CCD_WIDTH=`echo "sqrt(${DX_2T1}*${DX_2T1}+${DY_2T1}*${DY_2T1})" | bc -l`
+    CCD_2T1_MOV_X=`echo "${DFT_CCD_EDGECUT}*${DX_2T1}/${CCD_WIDTH}" | bc -l`
+    CCD_2T1_MOV_Y=`echo "${DFT_CCD_EDGECUT}*${DY_2T1}/${CCD_WIDTH}" | bc -l`
+    DX_3T1=$(echo "`echo ${CCD_3_RAW} | awk -F\, '{print $1}'`-`echo  ${CCD_1_RAW} |awk -F\, '{print $1}'`" |bc -l)
+    DY_3T1=$(echo "`echo ${CCD_3_RAW} | awk -F\, '{print $2}'`-`echo  ${CCD_1_RAW} |awk -F\, '{print $2}'`" |bc -l)
+    CCD_3T1_MOV_X=`echo "${DFT_CCD_EDGECUT}*${DX_3T1}/${CCD_WIDTH}" | bc -l`
+    CCD_3T1_MOV_Y=`echo "${DFT_CCD_EDGECUT}*${DY_3T1}/${CCD_WIDTH}" | bc -l`
+    CCD_1_X=$(echo "` echo ${CCD_1_RAW} |awk -F\, '{print $1}' ` + `echo ${CCD_2T1_MOV_X}` +`echo ${CCD_3T1_MOV_X}` "| bc -l)
+    CCD_1_Y=$(echo "` echo ${CCD_1_RAW} |awk -F\, '{print $2}' ` + `echo ${CCD_2T1_MOV_Y}` +`echo ${CCD_3T1_MOV_Y}` "| bc -l)
+    DX_1T2=$(echo "`echo ${CCD_1_RAW} | awk -F\, '{print $1}'`-`echo  ${CCD_2_RAW} |awk -F\, '{print $1}'`" |bc -l)
+    DY_1T2=$(echo "`echo ${CCD_1_RAW} | awk -F\, '{print $2}'`-`echo  ${CCD_2_RAW} |awk -F\, '{print $2}'`" |bc -l)
+    CCD_1T2_MOV_X=`echo "${DFT_CCD_EDGECUT}*${DX_1T2}/${CCD_WIDTH}" | bc -l`
+    CCD_1T2_MOV_Y=`echo "${DFT_CCD_EDGECUT}*${DY_1T2}/${CCD_WIDTH}" | bc -l`
+    DX_4T2=$(echo "`echo ${CCD_4_RAW} | awk -F\, '{print $1}'`-`echo  ${CCD_2_RAW} |awk -F\, '{print $1}'`" |bc -l)
+    DY_4T2=$(echo "`echo ${CCD_4_RAW} | awk -F\, '{print $2}'`-`echo  ${CCD_2_RAW} |awk -F\, '{print $2}'`" |bc -l)
+    CCD_4T2_MOV_X=`echo "${DFT_CCD_EDGECUT}*${DX_4T2}/${CCD_WIDTH}" | bc -l`
+    CCD_4T2_MOV_Y=`echo "${DFT_CCD_EDGECUT}*${DY_4T2}/${CCD_WIDTH}" | bc -l`
+    CCD_2_X=$(echo "` echo ${CCD_2_RAW} |awk -F\, '{print $1}' ` + `echo ${CCD_1T2_MOV_X}` +`echo ${CCD_4T2_MOV_X}` "| bc -l)
+    CCD_2_Y=$(echo "` echo ${CCD_2_RAW} |awk -F\, '{print $2}' ` + `echo ${CCD_1T2_MOV_Y}` +`echo ${CCD_4T2_MOV_Y}` "| bc -l)
+    DX_1T3=$(echo "`echo ${CCD_1_RAW} | awk -F\, '{print $1}'`-`echo  ${CCD_3_RAW} |awk -F\, '{print $1}'`" |bc -l)
+    DY_1T3=$(echo "`echo ${CCD_1_RAW} | awk -F\, '{print $2}'`-`echo  ${CCD_3_RAW} |awk -F\, '{print $2}'`" |bc -l)
+    CCD_1T3_MOV_X=`echo "${DFT_CCD_EDGECUT}*${DX_1T3}/${CCD_WIDTH}" | bc -l`
+    CCD_1T3_MOV_Y=`echo "${DFT_CCD_EDGECUT}*${DY_1T3}/${CCD_WIDTH}" | bc -l`
+    DX_4T3=$(echo "`echo ${CCD_4_RAW} | awk -F\, '{print $1}'`-`echo  ${CCD_3_RAW} |awk -F\, '{print $1}'`" |bc -l)
+    DY_4T3=$(echo "`echo ${CCD_4_RAW} | awk -F\, '{print $2}'`-`echo  ${CCD_3_RAW} |awk -F\, '{print $2}'`" |bc -l)
+    CCD_4T3_MOV_X=`echo "${DFT_CCD_EDGECUT}*${DX_4T3}/${CCD_WIDTH}" | bc -l`
+    CCD_4T3_MOV_Y=`echo "${DFT_CCD_EDGECUT}*${DY_4T3}/${CCD_WIDTH}" | bc -l`
+    CCD_3_X=$(echo "` echo ${CCD_3_RAW} |awk -F\, '{print $1}' ` + `echo ${CCD_1T3_MOV_X}` +`echo ${CCD_4T3_MOV_X}` "| bc -l)
+    CCD_3_Y=$(echo "` echo ${CCD_3_RAW} |awk -F\, '{print $2}' ` + `echo ${CCD_1T3_MOV_Y}` +`echo ${CCD_4T3_MOV_Y}` "| bc -l)
+    DX_2T4=$(echo "`echo ${CCD_2_RAW} | awk -F\, '{print $1}'`-`echo  ${CCD_4_RAW} |awk -F\, '{print $1}'`" |bc -l)
+    DY_2T4=$(echo "`echo ${CCD_2_RAW} | awk -F\, '{print $2}'`-`echo  ${CCD_4_RAW} |awk -F\, '{print $2}'`" |bc -l)
+    CCD_2T4_MOV_X=`echo "${DFT_CCD_EDGECUT}*${DX_2T4}/${CCD_WIDTH}" | bc -l`
+    CCD_2T4_MOV_Y=`echo "${DFT_CCD_EDGECUT}*${DY_2T4}/${CCD_WIDTH}" | bc -l`
+    DX_3T4=$(echo "`echo ${CCD_3_RAW} | awk -F\, '{print $1}'`-`echo  ${CCD_4_RAW} |awk -F\, '{print $1}'`" |bc -l)
+    DY_3T4=$(echo "`echo ${CCD_3_RAW} | awk -F\, '{print $2}'`-`echo  ${CCD_4_RAW} |awk -F\, '{print $2}'`" |bc -l)
+    CCD_3T4_MOV_X=`echo "${DFT_CCD_EDGECUT}*${DX_3T4}/${CCD_WIDTH}" | bc -l`
+    CCD_3T4_MOV_Y=`echo "${DFT_CCD_EDGECUT}*${DY_3T4}/${CCD_WIDTH}" | bc -l`
+    CCD_4_X=$(echo "` echo ${CCD_4_RAW} |awk -F\, '{print $1}' ` + `echo ${CCD_2T4_MOV_X}` +`echo ${CCD_3T4_MOV_X}` "| bc -l)
+    CCD_4_Y=$(echo "` echo ${CCD_4_RAW} |awk -F\, '{print $2}' ` + `echo ${CCD_2T4_MOV_Y}` +`echo ${CCD_3T4_MOV_Y}` "| bc -l)
+    REG_CCD_CUT=`echo "polygon(${CCD_1_X},${CCD_1_Y},${CCD_2_X},${CCD_2_Y},${CCD_4_X},${CCD_4_Y},${CCD_3_X},${CCD_3_Y})"`
+    REG_FILE_CCD="_ccd.reg"
+    [ -e "${REG_FILE_CCD}" ] && mv -f ${REG_FILE_CCD} ${REG_FILE_CCD}_bak
+    echo "${REG_CCD_CUT}" >>${REG_FILE_CCD}
+
+elif [ "${ACIS_TYPE}" = "I" ]; then
+    # ACIS-I
+    TMP_REG_FILE_CCD="_ccd_tmp.reg"
+    [ -e "${TMP_REG_FILE_CCD}" ] && mv -f ${TMP_REG_FILE_CCD} ${TMP_REG_FILE_CCD}_bak
+    for i in `seq 0 3` ; do 
+    punlearn dmlist
+    dmlist infile="${SKYFOV}[ccd_id=${i}][cols POS]" opt="data,clean" | awk '{for (i=1;i<=NF;i++) print $i }' |sed -e ':a;N;s/\n/,/;ta' | awk -F"]," '{print "polygon("$2}' | awk -F"NaN" '{print $1}' >${TMP_LIST} 
+    python ${SCRIPT_DIR}/${CCDGAP_SCRIPT} ${TMP_LIST} >${TMP_REC}
+    XC=` cat ${TMP_REC} | awk -F\( '{print $2}' |awk -F\) '{print $1}' |awk -F\, '{print $1}'`
+    YC=` cat ${TMP_REC} | awk -F\( '{print $2}' |awk -F\) '{print $1}' |awk -F\, '{print $2}'`
+    ADD_L=` cat ${TMP_REC} | awk -F\( '{print $2}' |awk -F\) '{print $1}' |awk -F\, '{print $3/2}'`
+    ANG=` cat ${TMP_REC} | awk -F\( '{print $2}' |awk -F\) '{print $1}' |awk -F\, '{print $5}'`
+    while [ 1 -eq 1 ]; do
+        if [ `echo "${ANG} < 0" |bc -l ` -eq 1 ] ; then
+            ANG=` echo " ${ANG} + 90 " | bc -l `
+        elif [ `echo "${ANG} >=90" |bc -l ` -eq 1 ] ; then
+            ANG=` echo " ${ANG} - 90 " | bc -l`
+        else
+            break
+        fi
+    done
+    ANG=`echo "${ANG}/180*3.1415926" |bc -l` 
+    CCD_1_X_RAW=` echo " ${XC}  ${ADD_L} ${ANG} "| awk '{print $1-$2*cos($3)-$2*sin($3)}' ` 
+    CCD_2_X_RAW=` echo " ${XC}  ${ADD_L} ${ANG} "| awk '{print $1+$2*cos($3)-$2*sin($3)}' `
+    CCD_3_X_RAW=` echo " ${XC}  ${ADD_L} ${ANG} "| awk '{print $1-$2*cos($3)+$2*sin($3)}' `
+    CCD_4_X_RAW=` echo " ${XC}  ${ADD_L} ${ANG} "| awk '{print $1+$2*cos($3)+$2*sin($3)}' `
+    CCD_1_Y_RAW=` echo " ${YC}  ${ADD_L} ${ANG} "| awk '{print $1+$2*cos($3)-$2*sin($3)}' `
+    CCD_2_Y_RAW=` echo " ${YC}  ${ADD_L} ${ANG} "| awk '{print $1+$2*cos($3)+$2*sin($3)}' `
+    CCD_3_Y_RAW=` echo " ${YC}  ${ADD_L} ${ANG} "| awk '{print $1-$2*cos($3)-$2*sin($3)}' `
+    CCD_4_Y_RAW=` echo " ${YC}  ${ADD_L} ${ANG} "| awk '{print $1-$2*cos($3)+$2*sin($3)}' `
+    CCD_1_RAW=` echo "${CCD_1_X_RAW},${CCD_1_Y_RAW}"`
+    CCD_2_RAW=` echo "${CCD_2_X_RAW},${CCD_2_Y_RAW}"`
+    CCD_3_RAW=` echo "${CCD_3_X_RAW},${CCD_3_Y_RAW}"`
+    CCD_4_RAW=` echo "${CCD_4_X_RAW},${CCD_4_Y_RAW}"`
+    REG_CCD_RAW=`echo "polygon(${CCD_1_RAW}, ${CCD_2_RAW}, ${CCD_4_RAW}, ${CCD_3_RAW}) " `
+    DX_2T1=$(echo "`echo ${CCD_2_RAW} | awk -F\, '{print $1}'`-`echo  ${CCD_1_RAW} |awk -F\, '{print $1}'`" |bc -l)
+    DY_2T1=$(echo "`echo ${CCD_2_RAW} | awk -F\, '{print $2}'`-`echo  ${CCD_1_RAW} |awk -F\, '{print $2}'`" |bc -l)
+    CCD_WIDTH=`echo "sqrt(${DX_2T1}*${DX_2T1}+${DY_2T1}*${DY_2T1})" | bc -l`
+    CCD_2T1_MOV_X=`echo "${DFT_CCD_EDGECUT}*${DX_2T1}/${CCD_WIDTH}" | bc -l`
+    CCD_2T1_MOV_Y=`echo "${DFT_CCD_EDGECUT}*${DY_2T1}/${CCD_WIDTH}" | bc -l`
+    DX_3T1=$(echo "`echo ${CCD_3_RAW} | awk -F\, '{print $1}'`-`echo  ${CCD_1_RAW} |awk -F\, '{print $1}'`" |bc -l)
+    DY_3T1=$(echo "`echo ${CCD_3_RAW} | awk -F\, '{print $2}'`-`echo  ${CCD_1_RAW} |awk -F\, '{print $2}'`" |bc -l)
+    CCD_3T1_MOV_X=`echo "${DFT_CCD_EDGECUT}*${DX_3T1}/${CCD_WIDTH}" | bc -l`
+    CCD_3T1_MOV_Y=`echo "${DFT_CCD_EDGECUT}*${DY_3T1}/${CCD_WIDTH}" | bc -l`
+    CCD_1_X=$(echo "` echo ${CCD_1_RAW} |awk -F\, '{print $1}' ` + `echo ${CCD_2T1_MOV_X}` +`echo ${CCD_3T1_MOV_X}` "| bc -l)
+    CCD_1_Y=$(echo "` echo ${CCD_1_RAW} |awk -F\, '{print $2}' ` + `echo ${CCD_2T1_MOV_Y}` +`echo ${CCD_3T1_MOV_Y}` "| bc -l)
+    DX_1T2=$(echo "`echo ${CCD_1_RAW} | awk -F\, '{print $1}'`-`echo  ${CCD_2_RAW} |awk -F\, '{print $1}'`" |bc -l)
+    DY_1T2=$(echo "`echo ${CCD_1_RAW} | awk -F\, '{print $2}'`-`echo  ${CCD_2_RAW} |awk -F\, '{print $2}'`" |bc -l)
+    CCD_1T2_MOV_X=`echo "${DFT_CCD_EDGECUT}*${DX_1T2}/${CCD_WIDTH}" | bc -l`
+    CCD_1T2_MOV_Y=`echo "${DFT_CCD_EDGECUT}*${DY_1T2}/${CCD_WIDTH}" | bc -l`
+    DX_4T2=$(echo "`echo ${CCD_4_RAW} | awk -F\, '{print $1}'`-`echo  ${CCD_2_RAW} |awk -F\, '{print $1}'`" |bc -l)
+    DY_4T2=$(echo "`echo ${CCD_4_RAW} | awk -F\, '{print $2}'`-`echo  ${CCD_2_RAW} |awk -F\, '{print $2}'`" |bc -l)
+    CCD_4T2_MOV_X=`echo "${DFT_CCD_EDGECUT}*${DX_4T2}/${CCD_WIDTH}" | bc -l`
+    CCD_4T2_MOV_Y=`echo "${DFT_CCD_EDGECUT}*${DY_4T2}/${CCD_WIDTH}" | bc -l`
+    CCD_2_X=$(echo "` echo ${CCD_2_RAW} |awk -F\, '{print $1}' ` + `echo ${CCD_1T2_MOV_X}` +`echo ${CCD_4T2_MOV_X}` "| bc -l)
+    CCD_2_Y=$(echo "` echo ${CCD_2_RAW} |awk -F\, '{print $2}' ` + `echo ${CCD_1T2_MOV_Y}` +`echo ${CCD_4T2_MOV_Y}` "| bc -l)
+    DX_1T3=$(echo "`echo ${CCD_1_RAW} | awk -F\, '{print $1}'`-`echo  ${CCD_3_RAW} |awk -F\, '{print $1}'`" |bc -l)
+    DY_1T3=$(echo "`echo ${CCD_1_RAW} | awk -F\, '{print $2}'`-`echo  ${CCD_3_RAW} |awk -F\, '{print $2}'`" |bc -l)
+    CCD_1T3_MOV_X=`echo "${DFT_CCD_EDGECUT}*${DX_1T3}/${CCD_WIDTH}" | bc -l`
+    CCD_1T3_MOV_Y=`echo "${DFT_CCD_EDGECUT}*${DY_1T3}/${CCD_WIDTH}" | bc -l`
+    DX_4T3=$(echo "`echo ${CCD_4_RAW} | awk -F\, '{print $1}'`-`echo  ${CCD_3_RAW} |awk -F\, '{print $1}'`" |bc -l)
+    DY_4T3=$(echo "`echo ${CCD_4_RAW} | awk -F\, '{print $2}'`-`echo  ${CCD_3_RAW} |awk -F\, '{print $2}'`" |bc -l)
+    CCD_4T3_MOV_X=`echo "${DFT_CCD_EDGECUT}*${DX_4T3}/${CCD_WIDTH}" | bc -l`
+    CCD_4T3_MOV_Y=`echo "${DFT_CCD_EDGECUT}*${DY_4T3}/${CCD_WIDTH}" | bc -l`
+    CCD_3_X=$(echo "` echo ${CCD_3_RAW} |awk -F\, '{print $1}' ` + `echo ${CCD_1T3_MOV_X}` +`echo ${CCD_4T3_MOV_X}` "| bc -l)
+    CCD_3_Y=$(echo "` echo ${CCD_3_RAW} |awk -F\, '{print $2}' ` + `echo ${CCD_1T3_MOV_Y}` +`echo ${CCD_4T3_MOV_Y}` "| bc -l)
+    DX_2T4=$(echo "`echo ${CCD_2_RAW} | awk -F\, '{print $1}'`-`echo  ${CCD_4_RAW} |awk -F\, '{print $1}'`" |bc -l)
+    DY_2T4=$(echo "`echo ${CCD_2_RAW} | awk -F\, '{print $2}'`-`echo  ${CCD_4_RAW} |awk -F\, '{print $2}'`" |bc -l)
+    CCD_2T4_MOV_X=`echo "${DFT_CCD_EDGECUT}*${DX_2T4}/${CCD_WIDTH}" | bc -l`
+    CCD_2T4_MOV_Y=`echo "${DFT_CCD_EDGECUT}*${DY_2T4}/${CCD_WIDTH}" | bc -l`
+    DX_3T4=$(echo "`echo ${CCD_3_RAW} | awk -F\, '{print $1}'`-`echo  ${CCD_4_RAW} |awk -F\, '{print $1}'`" |bc -l)
+    DY_3T4=$(echo "`echo ${CCD_3_RAW} | awk -F\, '{print $2}'`-`echo  ${CCD_4_RAW} |awk -F\, '{print $2}'`" |bc -l)
+    CCD_3T4_MOV_X=`echo "${DFT_CCD_EDGECUT}*${DX_3T4}/${CCD_WIDTH}" | bc -l`
+    CCD_3T4_MOV_Y=`echo "${DFT_CCD_EDGECUT}*${DY_3T4}/${CCD_WIDTH}" | bc -l`
+    CCD_4_X=$(echo "` echo ${CCD_4_RAW} |awk -F\, '{print $1}' ` + `echo ${CCD_2T4_MOV_X}` +`echo ${CCD_3T4_MOV_X}` "| bc -l)
+    CCD_4_Y=$(echo "` echo ${CCD_4_RAW} |awk -F\, '{print $2}' ` + `echo ${CCD_2T4_MOV_Y}` +`echo ${CCD_3T4_MOV_Y}` "| bc -l)
+    REG_CCD_CUT=`echo "polygon(${CCD_1_X},${CCD_1_Y},${CCD_2_X},${CCD_2_Y},${CCD_4_X},${CCD_4_Y},${CCD_3_X},${CCD_3_Y})"`
+    echo ${REG_CCD_CUT} >>${TMP_REG_FILE_CCD}
+done
+    REG_FILE_CCD="_ccd.reg"
+    [ -e "${REG_FILE_CCD}" ] && mv -fv ${REG_FILE_CCD} ${REG_FILE_CCD}_bak
+ #   echo "` cat ${TMP_REG_FILE_CCD} | head -n 1 | tail -n 1` + ` cat ${TMP_REG_FILE_CCD} | head -n 2 | tail -n 1` +`cat ${TMP_REG_FILE_CCD} | head -n 3 | tail -n 1`+`cat ${TMP_REG_FILE_CCD} | head -n 4 | tail -n 1 `" >${REG_FILE_CCD}
+    cat "${TMP_REG_FILE_CCD}">${REG_FILE_CCD}
+else
+    #
+    printf "*** ERROR ACIS_TYPE ***\n"
+    exit 255
+fi
+## }}}
+
+## cut ccd region edge
+echo "${REG_CCD_RAW}" >_ccd_raw.reg
+
+# exit 233
+
+## generate new regions within CCD for dmextract
+SBP_REG_INCCD="_${SBP_REG%.reg}_inccd.reg"
+[ -e "${SBP_REG_INCCD}" ] && mv -fv ${SBP_REG_INCCD} ${SBP_REG_INCCD}_bak
+
+echo "CMD: cat ${CELL_REG} | grep \( | sed -e ':a;N;s/\n/-/;ta'"
+CELL_REG_USE=`cat ${CELL_REG} | grep \( | sed -e ':a;N;s/\n/-/;ta'`
+# exit 233
+
+if [ "${ACIS_TYPE}" = "S" ]; then
+    grep -iE '^(pie|annulus)' ${SBP_REG_FIX} | sed "s/$/\ \&\ `cat ${REG_FILE_CCD}`/" | sed "s/$/\ \-\ ${CELL_REG_USE}/" > ${SBP_REG_INCCD}
+    else 
+    L=`cat ${SBP_REG_FIX} | wc -l `
+
+    for i in `seq 1 $L` ; do 
+        echo "`cat ${SBP_REG_FIX} |head -n $i | tail -n 1 ` & `cat ${REG_FILE_CCD} | head -n 1 `- ${CELL_REG_USE} | `cat ${SBP_REG_FIX} |head -n $i | tail -n 1` & `cat ${REG_FILE_CCD} | head -n 2| tail -n 1 `- ${CELL_REG_USE} |`cat ${SBP_REG_FIX} |head -n $i | tail -n 1 ` & `cat ${REG_FILE_CCD} | head -n 3 | tail -n 1 `- ${CELL_REG_USE} |`cat ${SBP_REG_FIX} |head -n $i | tail -n 1 ` & `cat ${REG_FILE_CCD} | tail -n 1 `- ${CELL_REG_USE}  " >>${SBP_REG_INCCD}
+    done
+fi
+# ds9 ${EVT_E} -region ${SBP_REG_INCCD}
+
+## `surface brightness profile' related data {{{
+## extract sbp
+printf "extract surface brightness profile ...\n"
+SBP_DAT="${SBP_REG%.reg}.fits"
+[ -e "${SBP_DAT}" ] && mv -fv ${SBP_DAT} ${SBP_DAT}_bak
+if [ -r "${BKG}" ]; then
+    EXPO_EVT=`dmkeypar ${EVT_E} EXPOSURE echo=yes`
+    EXPO_BKG=`dmkeypar ${BKG} EXPOSURE echo=yes`
+    BKG_NORM=`echo "${EXPO_EVT} ${EXPO_BKG}" | awk '{ printf("%g", $1/$2) }'`
+    printf "   == (BKG subtracted; bkgnorm=${BKG_NORM}, energy:${E_RANGE}) ==\n"
+    punlearn dmextract
+    dmextract infile="${EVT_E}[bin sky=@${SBP_REG_INCCD}]" outfile="${SBP_DAT}" \
+        exp="${EXPMAP}" bkg="${BKG}[energy=${E_RANGE}][bin sky=@${SBP_REG_INCCD}]" \
+        bkgexp=")exp" bkgnorm=${BKG_NORM} opt=generic clobber=yes
+else
+    punlearn dmextract
+    dmextract infile="${EVT_E}[bin sky=@${SBP_REG_INCCD}]" outfile="${SBP_DAT}" \
+        exp="${EXPMAP}" opt=generic clobber=yes
+fi
+
+## add `rmid' column
+printf "add \`RMID' & \`R_ERR' column ...\n"
+SBP_RMID="${SBP_DAT%.fits}_rmid.fits"
+[ -e "${SBP_RMID}" ] && mv -fv ${SBP_RMID} ${SBP_RMID}_bak
+punlearn dmtcalc
+dmtcalc infile="${SBP_DAT}" outfile="${SBP_RMID}" \
+    expression="RMID=(R[0]+R[1])/2,R_ERR=(R[1]-R[0])/2" \
+    clobber=yes
+
+## output needed sbp data to files
+printf "output needed sbp data ...\n"
+SBP_TXT="${SBP_DAT%.fits}.txt"
+SBP_QDP="${SBP_DAT%.fits}.qdp"
+[ -e "${SBP_TXT}" ] && mv -fv ${SBP_TXT} ${SBP_TXT}_bak
+[ -e "${SBP_QDP}" ] && mv -fv ${SBP_QDP} ${SBP_QDP}_bak
+punlearn dmlist
+dmlist infile="${SBP_RMID}[cols RMID,R_ERR,SUR_FLUX,SUR_FLUX_ERR]" \
+    outfile="${SBP_TXT}" opt="data,clean"
+
+## QDP for sbp {{{
+printf "generate a handy QDP file for sbp ...\n"
+cp -fv ${SBP_TXT} ${SBP_QDP}
+# change comment sign
+sed -i'' 's/#/!/g' ${SBP_QDP}
+# add QDP commands
+sed -i'' '1 i\
+READ SERR 1 2' ${SBP_QDP}
+sed -i'' '2 i\
+LABEL Y "Surface Flux (photons/cm\\u2\\d/pixel\\u2\\d/s)"' ${SBP_QDP}
+sed -i'' '2 i\
+LABEL X "Radius (pixel)"' ${SBP_QDP}
+sed -i'' '2 i\
+LABEL T "Surface Brightness Profile"' ${SBP_QDP}
+## QDP }}}
+
+printf "generate sbp fitting needed files ...\n"
+SBP_RADIUS="radius_sbp.txt"
+SBP_FLUX="flux_sbp.txt"
+[ -e "${SBP_RADIUS}" ] && mv -fv ${SBP_RADIUS} ${SBP_RADIUS}_bak
+[ -e "${SBP_FLUX}" ] && mv -fv ${SBP_FLUX} ${SBP_FLUX}_bak
+punlearn dmlist
+dmlist infile="${SBP_RMID}[cols R]" \
+    opt="data,clean" | awk '{ print $2 }' > ${SBP_RADIUS}
+# change the first line `R[2]' to `0.0'
+sed -i'' 's/R.*/0\.0/' ${SBP_RADIUS}
+dmlist infile="${SBP_RMID}[cols SUR_FLUX,SUR_FLUX_ERR]" \
+    opt="data,clean" > ${SBP_FLUX}
+# remove the first comment line
+sed -i'' '/#.*/d' ${SBP_FLUX}
+
+## sbp data }}}
+
+## main }}}
+
+exit 0
+