From 90eb3cf2b94650ce470bdf215b9b977f1e6a0a06 Mon Sep 17 00:00:00 2001 From: Aaron LI Date: Thu, 16 Feb 2017 16:30:08 +0800 Subject: Add 'renorm_spectrum.py' to replace 'chandra_bkg_rescale.sh' Also remove 'chandra_pb_flux.sh' --- scripts/chandra_bkg_rescale.sh | 80 ------------------------------------------ scripts/chandra_pb_flux.sh | 69 ------------------------------------ scripts/renorm_spectrum.py | 74 ++++++++++++++++++++++++++++++++++++++ 3 files changed, 74 insertions(+), 149 deletions(-) delete mode 100755 scripts/chandra_bkg_rescale.sh delete mode 100755 scripts/chandra_pb_flux.sh create mode 100755 scripts/renorm_spectrum.py diff --git a/scripts/chandra_bkg_rescale.sh b/scripts/chandra_bkg_rescale.sh deleted file mode 100755 index 8d4abc7..0000000 --- a/scripts/chandra_bkg_rescale.sh +++ /dev/null @@ -1,80 +0,0 @@ -#!/bin/sh -# -# background rescale, by adjusting `BACKSCAL' -# according to the photon flux values in `9.5-12.0 keV' -# -# Weitian LI -# 2012/08/14 -# -# Changelogs: -# 2015/06/03, Aaron LI -# * Copy needed pfiles to tmp directory, -# set environment variable $PFILES to use these first. -# and remove them after usage. -# - -## background rescale (BACKSCAL) {{{ -# rescale 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_rescale() { - # $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 rescale }}} - -if [ $# -ne 2 ] || [ "x$1" = "x-h" ]; then - printf "usage:\n" - printf " `basename $0` \n" - printf "\nNOTE:\n" - printf " is the spectrum to be adjusted\n" - exit 1 -fi - -## prepare parameter files (pfiles) {{{ -CIAO_TOOLS="dmstat dmkeypar dmhedit" - -PFILES_TMPDIR="/tmp/pfiles-$$" -[ -d "${PFILES_TMPDIR}" ] && rm -rf ${PFILES_TMPDIR} || mkdir ${PFILES_TMPDIR} - -# Copy necessary pfiles for localized usage -for tool in ${CIAO_TOOLS}; do - pfile=`paccess ${tool}` - [ -n "${pfile}" ] && punlearn ${tool} && cp -Lvf ${pfile} ${PFILES_TMPDIR}/ -done - -# Modify environment variable 'PFILES' to use local pfiles first -export PFILES="${PFILES_TMPDIR}:${PFILES}" -## pfiles }}} - -# perform `bkg_rescale' -bkg_rescale "$1" "$2" - -# clean pfiles -rm -rf ${PFILES_TMPDIR} - -exit 0 - diff --git a/scripts/chandra_pb_flux.sh b/scripts/chandra_pb_flux.sh deleted file mode 100755 index 59d4418..0000000 --- a/scripts/chandra_pb_flux.sh +++ /dev/null @@ -1,69 +0,0 @@ -#!/bin/sh -# -# chandra particle background -# 9.5-12.0 keV (channel: 651-822) -# PI = [ energy(eV) / 14.6 eV + 1 ] -# -# Weitian -# 2012/07/30 -# -# ChangeLog: -# v2.0, 2015/06/03, Aaron LI -# * Copy needed pfiles to tmp directory, -# set environment variable $PFILES to use these first. -# and remove them after usage. -# v1.1: August 2, 2012 -# * fix bugs with scientific notation in `bc' -# - -if [ $# -eq 0 ] || [ "x$1" = "x-h" ]; then - echo "usage:" - echo " `basename $0` ..." - exit 1 -fi - -## energy range: 9.5 -- 12.0 keV -EN_LOW="9.5" -EN_HI="12.0" -CH_LOW=651 -CH_HI=822 - -echo "Energy: $EN_LOW -- $EN_HI (keV)" -echo "Channel: $CH_LOW -- $CH_HI" - -## prepare parameter files (pfiles) {{{ -CIAO_TOOLS="dmstat dmkeypar" - -PFILES_TMPDIR="/tmp/pfiles-$$" -[ -d "${PFILES_TMPDIR}" ] && rm -rf ${PFILES_TMPDIR} || mkdir ${PFILES_TMPDIR} - -# Copy necessary pfiles for localized usage -for tool in ${CIAO_TOOLS}; do - pfile=`paccess ${tool}` - [ -n "${pfile}" ] && punlearn ${tool} && cp -Lvf ${pfile} ${PFILES_TMPDIR}/ -done - -# Modify environment variable 'PFILES' to use local pfiles first -export PFILES="${PFILES_TMPDIR}:${PFILES}" -## pfiles }}} - -while ! [ -z $1 ]; do - f=$1 - shift - echo "FILE: $f" - punlearn dmstat - COUNTS=`dmstat "$f[channel=${CH_LOW}:${CH_HI}][cols counts]" | grep 'sum:' | awk '{ print $2 }'` - punlearn dmkeypar - EXPTIME=`dmkeypar $f EXPOSURE echo=yes` - BACK=`dmkeypar $f 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 " counts / exptime / backscal: ${COUNTS} / ${EXPTIME} / ${BACK}" - echo " ${PB_FLUX}" -done - -# clean pfiles -rm -rf ${PFILES_TMPDIR} - diff --git a/scripts/renorm_spectrum.py b/scripts/renorm_spectrum.py new file mode 100755 index 0000000..a09c90d --- /dev/null +++ b/scripts/renorm_spectrum.py @@ -0,0 +1,74 @@ +#!/usr/bin/env python3 +# +# Copyright (c) 2017 Weitian LI +# MIT license + +""" +Renormalize the (background) spectrum by equaling its particle background +flux (e.g., 9.5-12.0 keV) with respect to its corresponding source spectrum. + +The ``BACKSCAL`` keyword of the background spectrum is modified to realize +the above renormalization. +""" + +import sys +import argparse +import subprocess + +from spectrum import Spectrum + + +def renorm_spectrum(specfile, specfile_ref, elow=9500, ehigh=12000): + """ + Modify the ``BACKSCAL`` of ``specfile`` in order to equal its + particle background flux to that of ``specfile_ref``. + + Parameters + ---------- + specfile : str + (background) spectrum to be renormalized/modified. + specfile_ref : str + (source/reference) spectrum + elow, ehigh : float, optional + Lower and upper energy limit of the particle background. + """ + spec = Spectrum(specfile) + spec_ref = Spectrum(specfile_ref) + flux = spec.calc_pb_flux(elow=elow, ehigh=ehigh) + flux_ref = spec_ref.calc_pb_flux(elow=elow, ehigh=ehigh) + bs_old = spec.BACKSCAL + bs_new = bs_old * flux / flux_ref + subprocess.check_call(["punlearn", "dmhedit"]) + subprocess.check_call([ + "dmhedit", "infile=%s" % specfile, + "filelist=none", "operation=add", + "key=BACKSCAL", "value=%s" % bs_new, + "comment='Old BACKSCAL: %s'" % bs_old + ]) + print("%s:BACKSCAL: %f -> %f" % (specfile, bs_old, bs_new), + file=sys.stderr) + + +def main(): + parser = argparse.ArgumentParser( + description="Renormalize background spectrum w.r.t. source spectrum") + parser.add_argument("-L", "--energy-low", dest="elow", + type=int, default=9500, + help="Lower energy limit of the particle " + + "background [eV] (default: 9500 eV)") + parser.add_argument("-H", "--energy-high", dest="ehigh", + type=int, default=12000, + help="Upper energy limit of the particle " + + "background [eV] (default: 12000 eV)") + parser.add_argument("-r", "--spec-ref", dest="spec_ref", required=True, + help="Reference (source) spectrum") + parser.add_argument("spec", dest="spec", + help="(background) spectrum to be renormalized") + args = parser.parse_args() + + renorm_spectrum(specfile=args.spec, specfile_ref=args.spec_ref, + elow=args.elow, ehigh=args.ehigh) + + +if __name__ == "__main__": + main() -- cgit v1.2.2