Move shell/python scripts from 'mass_profile/' to 'bin/'

12 files changed, 0 insertions, 1406 deletions

diff --git a/mass_profile/analyze_lxfx.py b/mass_profile/analyze_lxfx.py
deleted file mode 100755
index 62859ff..0000000
--- a/mass_profile/analyze_lxfx.py
+++ /dev/null
@@ -1,57 +0,0 @@
-#!/usr/bin/env python3
-#
-# Calculate the mean and standard deviation of the (Monte Carlo)
-# Lx and Fx data
-#
-# Junhua GU
-# Weitian LI
-# 2016-06-07
-#
-
-import sys
-import argparse
-import numpy as np
-
-
-def read_bands(bands):
-    """
-    Read energy bands list, each band per line.
-    """
-    bands = map(str.split, open(bands).readlines())
-    bands = ["-".join(b) for b in bands]
-    return bands
-
-
-def output(name, bands, means, sigmas, outfile=sys.stdout):
-    if outfile is not sys.stdout:
-        outfile = open(outfile, "w")
-    for b, m, s in zip(bands, means, sigmas):
-        print("%s(%s)= %4.2E +/- %4.2E erg/s" % (name, b, m, s),
-              file=outfile)
-    if outfile is not sys.stdout:
-        outfile.close()
-
-
-def main():
-    parser = argparse.ArgumentParser(
-            description="Analyze Lx/Fx results")
-    parser.add_argument("name", help="Lx or Fx")
-    parser.add_argument("infile", help="input data file")
-    parser.add_argument("outfile", help="output results file")
-    parser.add_argument("bands", help="energy bands of the input data columns")
-    args = parser.parse_args()
-
-    data = np.loadtxt(args.infile)
-    bands = read_bands(args.bands)
-    if len(bands) != data.shape[1]:
-        raise ValueError("number of energy bands != number of data columns")
-
-    means = np.mean(data, axis=0)
-    sigmas = np.std(data, axis=0)
-    output(name=args.name, bands=bands, means=means, sigmas=sigmas)
-    output(name=args.name, bands=bands, means=means, sigmas=sigmas,
-           outfile=args.outfile)
-
-
-if __name__ == "__main__":
-    main()
diff --git a/mass_profile/analyze_mass_profile.py b/mass_profile/analyze_mass_profile.py
deleted file mode 100755
index 3ee128c..0000000
--- a/mass_profile/analyze_mass_profile.py
+++ /dev/null
@@ -1,195 +0,0 @@
-#!/usr/bin/env python
-
-import sys
-import numpy
-import scipy.interpolate
-
-confidence_level=.68
-def read_file(param):
-    delta=float(param[0])
-
-    file_mass_center=open("mass_int_center.qdp").readlines();
-    file_delta_center=open("overdensity_center.qdp").readlines();
-    
-    center_r=0
-    center_m=0
-    center_gm=0
-    center_gf=0
-
-    
-    for i in range(0,len(file_mass_center)):
-        lm=file_mass_center[i].strip();
-        ld=file_delta_center[i].strip();
-        r,m=lm.split()
-        r,d=ld.split()
-        r=float(r)
-        d=float(d)
-        m=float(m)
-        if m<1e11:
-            continue
-        if d<delta:
-            center_r=r
-            center_m=m
-            for j in open("gas_mass_int_center.qdp"):
-                rgm,gm=j.strip().split()
-                rgm=float(rgm)
-                gm=float(gm)
-                if rgm>r:
-
-                    center_gm=gm
-                    center_gf=gm/m
-                    break
-            break
-    if len(param)>1 and param[1]=='c':
-        #print("%s(<r%d)=%E solar mass"%("mass",delta,center_m))
-        #print("%s%d=%E kpc"%("r",delta,center_r))
-        #print("%s(<r%d)=%E solar mass"%("gas mass",delta,center_gm))
-        #print("%s(<r%d)=%E"%("gas fraction",delta,center_gf))
-        return center_m,center_r,center_gm,center_gf,None,None,None,None
-    
-
-#print(center_gm,center_gf)
-    file_mass=open('summary_mass_profile.qdp').readlines()
-    file_delta=open('summary_overdensity.qdp').readlines()
-    file_gm=open('summary_gas_mass_profile.qdp')
-
-
-    flag=True
-    rlist=[]
-    mlist=[]
-    gmlist=[]
-    gflist=[]
-    old_m=0
-    invalid_count=0
-    for i in range(0,len(file_mass)):
-        lm=file_mass[i].strip()
-        ld=file_delta[i].strip()
-        if lm[0]=='n':
-            flag=True
-            old_m=0
-            continue
-        if not flag:
-            continue
-        r,m=lm.split()
-        m=float(m)
-        if m<1e12:
-            continue
-        if m<old_m:
-            invalid_count+=1
-            flag=False
-            continue
-        r,d=ld.split()
-        d=float(d)
-        r=float(r)
-
-        if d<delta:
-            #print("%s %e"%(d,m))
-            mlist.append(m)
-            rlist.append(r)
-            flag1=True
-            while True:
-                lgm=file_gm.readline().strip()
-                if lgm[0]=='n':
-                    break
-                rgm,gm=lgm.split()
-                rgm=float(rgm)
-                gm=float(gm)
-                if rgm>r and flag1:
-                    gmlist.append(gm)
-
-                    flag1=False
-                    gflist.append(gm/mlist[-1])
-                #print(gm,gflist[-1])
-            flag=False
-        old_m=m
-    print("%d abnormal data dropped"%(invalid_count))
-
-
-    return center_m,center_r,center_gm,center_gf,mlist,rlist,gmlist,gflist
-#center_m=numpy.mean(mlist)
-#center_r=numpy.mean(rlist)
-
-center_m,center_r,center_gm,center_gf,mlist,rlist,gmlist,gflist=read_file(sys.argv[1:])
-delta=float(sys.argv[1])
-
-if len(sys.argv)>2 and sys.argv[2]=='c':
-    print("%s(<r%d)=%E solar mass"%("mass",delta,center_m))
-    print("%s%d=%E kpc"%("r",delta,center_r))
-    print("%s(<r%d)=%E solar mass"%("gas mass",delta,center_gm))
-    print("%s(<r%d)=%E"%("gas fraction",delta,center_gf))
-    sys.exit(0)
-
-
-mlist.sort()
-rlist.sort()
-gflist.sort()
-gmlist.sort()
-
-m_idx=-1
-r_idx=-1
-gm_idx=-1
-gf_idx=-1
-delta=float(sys.argv[1])
-for i in range(len(mlist)-1):
-    if (center_m-mlist[i])*(center_m-mlist[i+1])<=0:
-        m_idx=i
-        break
-
-for i in range(len(rlist)-1):
-    if (center_r-rlist[i])*(center_r-rlist[i+1])<=0:
-        r_idx=i
-        break
-
-for i in range(len(gmlist)-1):
-    if (center_gm-gmlist[i])*(center_gm-gmlist[i+1])<=0:
-        gm_idx=i
-        break
-
-for i in range(len(gflist)-1):
-    if (center_gf-gflist[i])*(center_gf-gflist[i+1])<=0:
-        gf_idx=i
-        break
-
-
-if m_idx==-1 or r_idx==-1 or gf_idx==-1 or gm_idx==-1:
-    print("Error, the center value is not enclosed by the Monte-Carlo realizations, please check the result!")
-    print("m:%E %E %E"%(center_m,mlist[0],mlist[-1]))
-    print("gm:%E %E %E"%(center_gm,gmlist[0],gmlist[-1]))
-    print("gf:%E %E %E"%(center_gf,gflist[0],gflist[-1]))
-    print("r:%E %E %E"%(center_r,rlist[0],rlist[-1]))
-    sys.exit(1)
-
-
-mlidx=int(m_idx*(1-confidence_level))
-muidx=m_idx-1+int((len(mlist)-m_idx)*confidence_level)
-
-
-rlidx=int(r_idx*(1-confidence_level))
-ruidx=r_idx-1+int((len(rlist)-r_idx)*confidence_level)
-
-gmlidx=int(gm_idx*(1-confidence_level))
-gmuidx=gm_idx-1+int((len(gmlist)-gm_idx)*confidence_level)
-
-gflidx=int(gf_idx*(1-confidence_level))
-gfuidx=gf_idx-1+int((len(gflist)-gf_idx)*confidence_level)
-
-
-merr1=mlist[mlidx]-center_m
-merr2=mlist[muidx]-center_m
-
-rerr1=rlist[rlidx]-center_r
-rerr2=rlist[ruidx]-center_r
-
-gmerr1=gmlist[gmlidx]-center_gm
-gmerr2=gmlist[gmuidx]-center_gm
-
-gferr1=gflist[gflidx]-center_gf
-gferr2=gflist[gfuidx]-center_gf
-
-#print("%d %d %d"%(mlidx,m_idx,muidx))
-#print("%d %d %d"%(rlidx,r_idx,ruidx))
-
-print("m%d=\t%e\t %e/+%e solar mass (1 sigma)"%(delta,center_m,merr1,merr2))
-print("gas_m%d=\t%e\t %e/+%e solar mass (1 sigma)"%(delta,center_gm,gmerr1,gmerr2))
-print("gas_fraction%d=\t%e\t %e/+%e (1 sigma)"%(delta,center_gf,gferr1,gferr2))
-print("r%d=\t%d\t %d/+%d kpc (1 sigma)"%(delta,center_r,rerr1,rerr2))
diff --git a/mass_profile/calc_coolfunc.sh b/mass_profile/calc_coolfunc.sh
deleted file mode 100755
index d1d0b35..0000000
--- a/mass_profile/calc_coolfunc.sh
+++ /dev/null
@@ -1,113 +0,0 @@
-#!/bin/sh
-##
-## Calculate the 'cooling function' profile with respect to the
-## given 'temperature profile' and the average abundance, redshift,
-## and column density nH, using the XSPEC model 'wabs*apec'.
-##
-## NOTE:
-## The output cooling function values should be the 'flux' values
-## with unit 'photon/s/cm^2' (different to 'calc_coolfunc_bands.sh').
-## These results will be used by 'fit_{beta,dbeta}_sbp' to derive the
-## (3D) gas density profile from (2D) surface brightness profile,
-## whose values have unit 'photon/cm^2/pixel/s'.
-##
-## Weitian LI
-## Created: 2012-08-17
-##
-## Change logs:
-## 2017-02-17, Weitian LI
-##   * Rename from 'coolfunc_calc.sh' to 'calc_coolfunc.sh'
-##   * Clean up that do not calculate and output <coolfunc_bolo>
-## 2017-02-16, Weitian LI
-##   * Do not calculate and output 'flux_cnt_ratio.txt'
-##
-
-## cmdline arguments {{{
-if [ $# -ne 5 ]; then
-    printf "usage:\n"
-    printf "    `basename $0` <tprofile> <avg_abund> <nH> <redshift> <coolfunc_outfile>\n"
-    exit 1
-fi
-TPROFILE=$1
-ABUNDANCE=$2
-N_H=$3
-REDSHIFT=$4
-COOLFUNC_DAT=$5
-NORM=`cosmo_calc ${REDSHIFT} | grep 'norm.*cooling_function' | awk -F':' '{ print $2 }'`
-
-if [ ! -r "${TPROFILE}" ]; then
-    printf "ERROR: given tprofile '${TPROFILE}' NOT accessiable\n"
-    exit 2
-fi
-[ -e "${COOLFUNC_DAT}" ] && rm -f ${COOLFUNC_DAT}
-## arguments }}}
-
-XSPEC_CF_XCM="_calc_coolfunc.xcm"
-[ -e "${XSPEC_CF_XCM}" ] && rm -f ${XSPEC_CF_XCM}
-
-## generate xspec script {{{
-cat > ${XSPEC_CF_XCM} << _EOF_
-## XSPEC Tcl script
-## Calculate the cooling function profile w.r.t the temperature profile.
-##
-## Generated by: `basename $0`
-## Date: `date`
-
-set xs_return_results 1
-set xs_echo_script 0
-# set tcl_precision 12
-## set basic data {{{
-set nh ${N_H}
-set redshift ${REDSHIFT}
-set abundance ${ABUNDANCE}
-set norm ${NORM}
-## basic }}}
-
-## xspec related {{{
-# debug settings {{{
-chatter 0
-# debug }}}
-query yes
-abund grsa
-dummyrsp 0.01 100.0 4096 linear
-# load model 'wabs*apec' to calc cooling function
-model wabs*apec & \${nh} & 1.0 & \${abundance} & \${redshift} & \${norm} & /*
-## xspec }}}
-
-## set input and output filename & open files
-set tpro_fn "${TPROFILE}"
-set cf_fn "${COOLFUNC_DAT}"
-if { [ file exists \${cf_fn} ] } {
-    exec rm -fv \${cf_fn}
-}
-
-## open files
-set tpro_fd [ open \${tpro_fn} r ]
-set cf_fd [ open \${cf_fn} w ]
-
-## read data from tprofile line by line
-while { [ gets \${tpro_fd} tpro_line ] != -1 } {
-    scan \${tpro_line} "%f %f" radius temperature
-    #puts "radius: \${radius}, temperature: \${temperature}"
-    # set temperature value
-    newpar 2 \${temperature}
-    # calc flux & tclout
-    flux 0.7 7.0
-    tclout flux 1
-    scan \${xspec_tclout} "%f %f %f %f" _ _ _ cf_photon
-    #puts "cf: \${cf_photon}"
-    puts \${cf_fd} "\${radius}    \${cf_photon}"
-}
-
-## close opened files
-close \${tpro_fd}
-close \${cf_fd}
-
-## exit
-tclexit
-_EOF_
-## xcm generation }}}
-
-## invoke xspec to calc
-printf "invoking XSPEC to calculate cooling function profile ...\n"
-xspec - ${XSPEC_CF_XCM} > /dev/null
diff --git a/mass_profile/calc_coolfunc_bands.sh b/mass_profile/calc_coolfunc_bands.sh
deleted file mode 100755
index bebdce2..0000000
--- a/mass_profile/calc_coolfunc_bands.sh
+++ /dev/null
@@ -1,119 +0,0 @@
-#!/bin/sh
-##
-## Calculate the 'cooling function' profile for each of the energy band
-## specified in a bands file, with respect to the given 'temperature profile'
-## and the average abundance, redshift, and column density nH, using the
-## XSPEC model 'wabs*apec'.
-##
-## NOTE:
-## To calculate the luminosity and flux from the source using the
-## 'calc_lx_{beta,dbeta}', set 'nH=0'.
-## Also the output cooling function values should be the 'flux' values
-## with unit 'erg/s/cm^2'.
-##
-## Weitian LI
-## Updated: 2017-02-17
-##
-
-## cmdline arguments {{{
-if [ $# -ne 6 ]; then
-    printf "usage:\n"
-    printf "    `basename $0` <tprofile> <avg_abund> <nH> <redshift> <coolfunc_prefix> <band_list>\n"
-    exit 1
-fi
-TPROFILE=$1
-ABUNDANCE=$2
-N_H=$3
-REDSHIFT=$4
-COOLFUNC_PREFIX=$5
-BLIST=$6
-NORM=`cosmo_calc ${REDSHIFT} | grep 'norm.*cooling_function' | awk -F':' '{ print $2 }'`
-
-if [ ! -r "${TPROFILE}" ]; then
-    printf "ERROR: given tprofile '${TPROFILE}' NOT accessiable\n"
-    exit 2
-fi
-## arguments }}}
-
-XSPEC_CF_XCM="_calc_coolfunc_bands.xcm"
-[ -e "${XSPEC_CF_XCM}" ] && rm -f ${XSPEC_CF_XCM}
-
-## generate xspec script {{{
-cat > ${XSPEC_CF_XCM} << _EOF_
-## XSPEC Tcl script
-## Calculate the cooling function profile w.r.t the temperature profile,
-## for each specified energy band.
-##
-## Generated by: `basename $0`
-## Date: `date`
-
-set xs_return_results 1
-set xs_echo_script 0
-# set tcl_precision 12
-## set basic data {{{
-set nh ${N_H}
-set redshift ${REDSHIFT}
-set abundance ${ABUNDANCE}
-set norm ${NORM}
-## basic }}}
-
-## xspec related {{{
-# debug settings {{{
-chatter 0
-# debug }}}
-query yes
-abund grsa
-dummyrsp 0.01 100.0 4096 linear
-# load model 'wabs*apec' to calc cooling function
-model wabs*apec & \${nh} & 1.0 & \${abundance} & \${redshift} & \${norm} &
-## xspec }}}
-
-## set input and output filename
-set tpro_fn "${TPROFILE}"
-set blist_fn "${BLIST}"
-set cf_prefix "${COOLFUNC_PREFIX}"
-set blist_fd [ open \${blist_fn} r ]
-
-## loop over each energy band
-while { [ gets \${blist_fd} blist_line ] != -1 } {
-    if { "\${blist_line}" == "bolo" } {
-        set e1 0.01
-        set e2 100.0
-        set name_suffix bolo
-    } else {
-        set e1 [ lindex \${blist_line} 0 ]
-        set e2 [ lindex \${blist_line} 1 ]
-        set name_suffix "\${e1}-\${e2}"
-    }
-    set cf_fn "\${cf_prefix}\${name_suffix}.dat"
-    if { [ file exists \${cf_fn} ] } {
-        exec rm -fv \${cf_fn}
-    }
-    set cf_fd [ open \${cf_fn} w ]
-    set tpro_fd [ open \${tpro_fn} r ]
-
-    ## read data from tprofile line by line
-    while { [ gets \${tpro_fd} tpro_line ] != -1 } {
-        scan \${tpro_line} "%f %f" radius temperature
-        #puts "radius: \${radius}, temperature: \${temperature}"
-        # set temperature value
-        newpar 2 \${temperature}
-        # calc flux & tclout
-        flux \${e1} \${e2}
-        tclout flux 1
-        scan \${xspec_tclout} "%f" cf_erg
-        #puts "cf: \${cf_erg}"
-        puts \${cf_fd} "\${radius}    \${cf_erg}"
-    }
-    close \${tpro_fd}
-    close \${cf_fd}
-}
-
-## exit
-tclexit
-_EOF_
-## xcm generation }}}
-
-## invoke xspec to calc
-printf "invoking XSPEC to calculate cooling function profile ...\n"
-xspec - ${XSPEC_CF_XCM} > /dev/null
diff --git a/mass_profile/calc_entropy.py b/mass_profile/calc_entropy.py
deleted file mode 100755
index 736fe7e..0000000
--- a/mass_profile/calc_entropy.py
+++ /dev/null
@@ -1,104 +0,0 @@
-#!/usr/bin/env python3
-#
-# Calculate the entropy within the specified radius.
-#
-# Junhua GU
-# Weitian LI
-# 2016-06-07
-#
-
-import argparse
-import re
-from itertools import groupby
-import numpy as np
-
-
-def isplit(iterable, splitters):
-    """
-    Credit: https://stackoverflow.com/a/4322780/4856091
-    """
-    return [list(g) for k, g in groupby(iterable,
-                                        lambda x:x in splitters) if not k]
-
-
-def get_entropy(data, r):
-    """
-    Get the entropy *at* the specified radius.
-
-    XXX: whether to interpolate first?
-    """
-    radius = data[:, 0]
-    entropy = data[:, 1]
-    s = np.min(entropy[radius > r])
-    return s
-
-
-def read_merged_qdp(infile):
-    """
-    Read merged QDP with multiple group of data separated by "no no no".
-    """
-    lines = map(lambda line: re.sub(r"^\s*no\s+no\s+no.*$", "X",
-                                    line.strip(), flags=re.I),
-                open(infile).readlines())
-    lines = isplit(lines, ("X",))
-    data_groups = []
-    for block in lines:
-        data = [list(map(float, l.split())) for l in block]
-        data.append(np.row_stack(data))
-    return data_groups
-
-
-def calc_error(center_value, mc_values, ci=0.683):
-    """
-    Calculate the uncertainties/errors.
-    """
-    data = np.concatenate([[center_value], mc_values])
-    median, q_lower, q_upper = np.percentile(data, q=(50, 50-50*ci, 50+50*ci))
-    mean = np.mean(data)
-    std = np.std(data)
-    return {
-        "mean":    mean,
-        "std":     std,
-        "median":  median,
-        "q_lower": q_lower,
-        "q_upper": q_upper,
-    }
-
-
-def main():
-    parser = argparse.ArgumentParser(
-            description="Calculate the entropy within the given radius")
-    parser.add_argument("-C", "--confidence-level", dest="ci",
-                        type=float, default=0.683,
-                        help="confidence level to estimate the errors")
-    parser.add_argument("center_data",
-                        help="calculate central entropy profile " +
-                             "(e.g., entropy_center.qdp)")
-    parser.add_argument("mc_data",
-                        help="Merged QDP file of all the Monte Carlo " +
-                             "simulated entropy profiles " +
-                             "(e.g., summary_entropy.qdp)")
-    parser.add_argument("rout", type=float, help="outer radius (kpc)")
-    args = parser.parse_args()
-
-    center_data = np.loadtxt(args.center_data)
-    center_s = get_entropy(center_data, r=args.rout)
-
-    data_groups = read_merged_qdp(args.mc_data)
-    entropy_list = []
-    for dg in data_groups:
-        s = get_entropy(dg, r=args.rout)
-        entropy_list.append(s)
-    results = calc_error(center_s, entropy_list, ci=args.ci)
-    s_err_lower = results["q_lower"] - center_s
-    s_err_upper = results["q_upper"] - center_s
-
-    print("entropy= %e %+e/%+e keV cm^2 (ci=%.1f%%)" %
-          (center_s, s_err_lower, s_err_upper, args.ci * 100))
-    print("entropy(mean)= %e" % results["mean"])
-    print("entropy(median)= %e" % results["median"])
-    print("entropy(std)= %e" % results["std"])
-
-
-if __name__ == "__main__":
-    main()
diff --git a/mass_profile/calc_lxfx.sh b/mass_profile/calc_lxfx.sh
deleted file mode 100755
index f1954db..0000000
--- a/mass_profile/calc_lxfx.sh
+++ /dev/null
@@ -1,162 +0,0 @@
-#!/bin/sh
-#
-# Wrapper script used to calculate the luminosity (Lx) and flux (Fx) data,
-# which invokes the programming 'calc_lx_beta' (single-beta SBP) or
-# 'calc_lx_dbeta' (double-beta SBP).
-#
-# Output:
-#   * lx_result.txt
-#   * fx_result.txt
-#   * summary_lx.dat
-#   * summary_fx.dat
-#   * lx_beta_param.txt / lx_dbeta_param.txt
-#
-# Author: Junhua GU
-# Created: 2013-06-24
-#
-# Weitian LI
-# 2016-06-07
-#
-
-if [ $# -eq 2 ] || [ $# -eq 3 ]; then
-    :
-else
-    echo "usage:"
-    echo "    `basename $0` <mass.conf> <rout_kpc> [c]"
-    echo ""
-    echo "arguments:"
-    echo "    <mass.conf>: config file for mass profile calculation"
-    echo "    <rout_kpc>: outer/cut radius within which to calculate Lx & Fx"
-    echo "                e.g., r500, r200 (unit: kpc)"
-    echo "    [c]: optional; if specified, do not calculate the errors"
-    exit 1
-fi
-export PATH="/usr/local/bin:/usr/bin:/bin:$PATH"
-
-mass_cfg="$1"
-rout="$2"
-case "$3" in
-    [cC])
-        F_C="YES"
-        ;;
-    *)
-        F_C="NO"
-        ;;
-esac
-
-base_path=$(dirname $(realpath $0))
-
-## Extract settings/values from the config file
-abund=`grep         '^abund'         ${mass_cfg} | awk '{ print $2 }'`
-tprofile_data=`grep '^tprofile_data' ${mass_cfg} | awk '{ print $2 }'`
-tprofile_cfg=`grep  '^tprofile_cfg'  ${mass_cfg} | awk '{ print $2 }'`
-sbp_cfg=`grep       '^sbp_cfg'       ${mass_cfg} | awk '{ print $2 }'`
-
-sbp_data=`grep      '^sbp_data'      ${sbp_cfg} | awk '{ print $2 }'`
-tprofile=`grep      '^tprofile'      ${sbp_cfg} | awk '{ print $2 }'`
-z=`grep             '^z'             ${sbp_cfg} | awk '{ print $2 }'`
-cm_per_pixel=`grep  '^cm_per_pixel'  ${sbp_cfg} | awk '{ print $2 }'`
-
-if grep -q '^beta2' $sbp_cfg; then
-    MODEL="dbeta"
-else
-    MODEL="beta"
-fi
-
-PROG_TPROFILE="fit_wang2012_model"
-tprofile_dump="wang2012_dump.qdp"
-${base_path}/${PROG_TPROFILE} ${tprofile_data} ${tprofile_cfg} \
-            ${cm_per_pixel} 2> /dev/null
-mv -fv ${tprofile_dump} ${tprofile}
-
-# energy bands for which the cooling function data will be calculated
-BLIST="blist.txt"
-[ -e "${BLIST}" ] && mv -f ${BLIST} ${BLIST}_bak
-cat > ${BLIST} << _EOF_
-bolo
-0.7 7
-0.1 2.4
-_EOF_
-
-# NOTE:
-# Set 'nh=0' when calculating the cooling function values, and use the
-# value given by 'flux' with unit 'erg/s/cm^2'.
-${base_path}/calc_coolfunc_bands.sh ${tprofile} ${abund} \
-            0 ${z} "cfunc_" ${BLIST}
-
-PROG="calc_lx_${MODEL}"
-LXF_RES="lx_${MODEL}_param.txt"
-${base_path}/${PROG} ${sbp_cfg} ${rout} \
-            cfunc_bolo.dat \
-            cfunc_0.7-7.dat \
-            cfunc_0.1-2.4.dat 2> /dev/null
-LX1=`grep '^Lx1' ${LX_RES} | awk '{ print $2 }'`
-LX2=`grep '^Lx2' ${LX_RES} | awk '{ print $2 }'`
-LX3=`grep '^Lx3' ${LX_RES} | awk '{ print $2 }'`
-FX1=`grep '^Fx1' ${LX_RES} | awk '{ print $2 }'`
-FX2=`grep '^Fx2' ${LX_RES} | awk '{ print $2 }'`
-FX3=`grep '^Fx3' ${LX_RES} | awk '{ print $2 }'`
-
-echo "${LX1} ${LX2} ${LX3}" >summary_lx.dat
-echo "${FX1} ${FX2} ${FX3}" >summary_fx.dat
-
-# save the calculated central values
-mv ${LX_RES} ${LX_RES%.txt}_center.txt
-mv lx_sbp_fit.qdp lx_sbp_fit_center.qdp
-mv lx_rho_fit.dat lx_rho_fit_center.dat
-
-# only calculate the central values
-if [ "${F_C}" = "YES" ]; then
-    echo "Calculate the central values only ..."
-    ${base_path}/analyze_lxfx.py "Lx" summary_lx.dat lx_result.txt ${BLIST}
-    ${base_path}/analyze_lxfx.py "Fx" summary_fx.dat fx_result.txt ${BLIST}
-    exit 0
-fi
-
-
-###########################################################
-# Estimate the errors of Lx and Fx by Monte Carlo simulation
-MC_TIMES=100
-for i in `seq 1 ${MC_TIMES}`; do
-    ${base_path}/shuffle_profile.py ${tprofile_data} tmp_tprofile.txt
-    ${base_path}/shuffle_profile.py ${sbp_data} tmp_sbprofile.txt
-
-    # temperature profile
-    ${base_path}/${PROG_TPROFILE} tmp_tprofile.txt ${tprofile_cfg} \
-                ${cm_per_pixel} 2> /dev/null
-    mv -f ${tprofile_dump} ${tprofile}
-
-    TMP_SBP_CFG="tmp_sbp.cfg"
-    [ -e "${TMP_SBP_CFG}" ] && rm -f ${TMP_SBP_CFG}
-    cat ${sbp_cfg} | while read l; do
-        if echo "${l}" | grep -q '^sbp_data' >/dev/null; then
-            echo "sbp_data  tmp_sbprofile.txt" >> ${TMP_SBP_CFG}
-        elif echo "${l}" | grep -q '^tprofile' >/dev/null; then
-            echo "tprofile  ${tprofile}" >> ${TMP_SBP_CFG}
-        else
-            echo "${l}" >> ${TMP_SBP_CFG}
-        fi
-    done
-
-    echo "### `pwd -P`"
-    echo "### ${i} / ${MC_TIMES} ###"
-    ${base_path}/calc_coolfunc_bands.sh ${tprofile} ${abund} \
-                0 ${z} "cfunc_" ${BLIST}
-    ${base_path}/${PROG} ${TMP_SBP_CFG} ${rout} \
-                cfunc_bolo.dat \
-                cfunc_0.7-7.dat \
-                cfunc_0.1-2.4.dat 2> /dev/null
-    LX1=`grep '^Lx1' ${LX_RES} | awk '{ print $2 }'`
-    LX2=`grep '^Lx2' ${LX_RES} | awk '{ print $2 }'`
-    LX3=`grep '^Lx3' ${LX_RES} | awk '{ print $2 }'`
-    FX1=`grep '^Fx1' ${LX_RES} | awk '{ print $2 }'`
-    FX2=`grep '^Fx2' ${LX_RES} | awk '{ print $2 }'`
-    FX3=`grep '^Fx3' ${LX_RES} | awk '{ print $2 }'`
-
-    echo "${LX1} ${LX2} ${LX3}" >>summary_lx.dat
-    echo "${FX1} ${FX2} ${FX3}" >>summary_fx.dat
-done # end of 'for'
-
-# analyze Lx & Fx Monte Carlo results
-${base_path}/analyze_lxfx.py "Lx" summary_lx.dat lx_result.txt ${BLIST}
-${base_path}/analyze_lxfx.py "Fx" summary_fx.dat fx_result.txt ${BLIST}
diff --git a/mass_profile/calc_lxfx_wrapper.sh b/mass_profile/calc_lxfx_wrapper.sh
deleted file mode 100755
index 7dea0d0..0000000
--- a/mass_profile/calc_lxfx_wrapper.sh
+++ /dev/null
@@ -1,76 +0,0 @@
-#!/bin/sh
-#
-# Calculate the Lx & Fx data.
-#
-# Based on 'loop_lx.sh', but only process one source
-#
-# Weitian LI
-# 2013-10-30
-#
-
-
-full_path=`readlink -f $0`
-base_dir=`dirname $full_path`
-
-if [ $# -lt 2 ]; then
-    printf "usage:\n"
-    printf "    `basename $0` <global.cfg> [c] < 500 | 200 > ...\n"
-    exit 1
-fi
-
-cfg_file="$1"
-pre_results="final_result.txt"
-
-case "$2" in
-    [cC]*)
-        F_C="YES"
-        shift
-        ;;
-    *)
-        F_C="NO"
-        ;;
-esac
-
-shift
-echo "delta: $@"    # 'printf' not work
-
-if [ "${F_C}" = "YES" ]; then
-    printf "MODE: center\n"
-fi
-
-# exit
-
-if [ ! -r "${cfg_file}" ]; then
-    printf "ERROR: global cfg not accessible\n"
-    exit 11
-elif [ ! -r "${pre_results}" ]; then
-    printf "ERROR: previous '${pre_results}' not accessible\n"
-    exit 12
-else
-    sbp_cfg=`grep '^sbp_cfg' $cfg_file | awk '{ print $2 }'`
-    ##
-    for delta in $@; do
-        if grep -q '^beta2' $sbp_cfg; then
-            MODEL="dbeta"
-        else
-            MODEL="beta"
-        fi
-        rout=`grep "^r${delta}" ${pre_results} | sed -e 's/=/ /' | awk '{ print $2 }'`
-        if [ "${F_C}" = "YES" ]; then
-            lx_res="lx_result_${delta}_c.txt"
-            fx_res="fx_result_${delta}_c.txt"
-            CMD="$base_dir/calc_lxfx.sh $cfg_file $rout c"
-        else
-            lx_res="lx_result_${delta}.txt"
-            fx_res="fx_result_${delta}.txt"
-            CMD="$base_dir/calc_lxfx.sh $cfg_file $rout"
-        fi
-        [ -e "${lx_res}" ] && mv -f ${lx_res} ${lx_res}_bak
-        [ -e "${fx_res}" ] && mv -f ${fx_res} ${fx_res}_bak
-        ${CMD}
-        mv -f lx_result.txt ${lx_res}
-        mv -f fx_result.txt ${fx_res}
-    done
-fi
-
-exit 0
diff --git a/mass_profile/fg_2500_500.py b/mass_profile/fg_2500_500.py
deleted file mode 100755
index 67a1a11..0000000
--- a/mass_profile/fg_2500_500.py
+++ /dev/null
@@ -1,153 +0,0 @@
-#!/usr/bin/env python
-
-import sys
-import numpy
-import scipy.interpolate
-
-confidence_level=.68
-def read_file(param):
-    delta=float(param[0])
-
-    file_mass_center=open("mass_int_center.qdp").readlines();
-    file_delta_center=open("overdensity_center.qdp").readlines();
-    
-    center_r=0
-    center_m=0
-    center_gm=0
-    center_gf=0
-
-    
-    for i in range(0,len(file_mass_center)):
-        lm=file_mass_center[i].strip();
-        ld=file_delta_center[i].strip();
-        r,m=lm.split()
-        r,d=ld.split()
-        r=float(r)
-        d=float(d)
-        m=float(m)
-        if m<1e11:
-            continue
-        if d<delta:
-            center_r=r
-            center_m=m
-            for j in open("gas_mass_int_center.qdp"):
-                rgm,gm=j.strip().split()
-                rgm=float(rgm)
-                gm=float(gm)
-                if rgm>r:
-
-                    center_gm=gm
-                    center_gf=gm/m
-                    break
-            break
-    if len(param)>1 and param[1]=='c':
-        #print("%s(<r%d)=%E solar mass"%("mass",delta,center_m))
-        #print("%s%d=%E kpc"%("r",delta,center_r))
-        #print("%s(<r%d)=%E solar mass"%("gas mass",delta,center_gm))
-        #print("%s(<r%d)=%E"%("gas fraction",delta,center_gf))
-        return center_m,center_r,center_gm,center_gf,None,None,None,None
-    
-
-#print(center_gm,center_gf)
-    file_mass=open('summary_mass_profile.qdp').readlines()
-    file_delta=open('summary_overdensity.qdp').readlines()
-    file_gm=open('summary_gas_mass_profile.qdp')
-
-
-    flag=True
-    rlist=[]
-    mlist=[]
-    gmlist=[]
-    gflist=[]
-    old_m=0
-    invalid_count=0
-    for i in range(0,len(file_mass)):
-        lm=file_mass[i].strip()
-        ld=file_delta[i].strip()
-        if lm[0]=='n':
-            flag=True
-            old_m=0
-            continue
-        if not flag:
-            continue
-        r,m=lm.split()
-        m=float(m)
-        if m<1e12:
-            continue
-        if m<old_m:
-            invalid_count+=1
-            flag=False
-            continue
-        r,d=ld.split()
-        d=float(d)
-        r=float(r)
-
-        if d<delta:
-            #print("%s %e"%(d,m))
-            mlist.append(m)
-            rlist.append(r)
-            flag1=True
-            while True:
-                lgm=file_gm.readline().strip()
-                if lgm[0]=='n':
-                    break
-                rgm,gm=lgm.split()
-                rgm=float(rgm)
-                gm=float(gm)
-                if rgm>r and flag1:
-                    gmlist.append(gm)
-
-                    flag1=False
-                    gflist.append(gm/mlist[-1])
-                #print(gm,gflist[-1])
-            flag=False
-        old_m=m
-    print("%d abnormal data dropped"%(invalid_count))
-
-
-    return center_m,center_r,center_gm,center_gf,mlist,rlist,gmlist,gflist
-#center_m=numpy.mean(mlist)
-#center_r=numpy.mean(rlist)
-
-if len(sys.argv)>1:
-    center_m2500,center_r2500,center_gm2500,center_gf2500,mlist2500,rlist2500,gmlist2500,gflist2500=read_file([2500,sys.argv[1]])
-    center_m500,center_r500,center_gm500,center_gf500,mlist500,rlist500,gmlist500,gflist500=read_file([500,sys.argv[1]])
-else:
-    center_m2500,center_r2500,center_gm2500,center_gf2500,mlist2500,rlist2500,gmlist2500,gflist2500=read_file([2500])
-    center_m500,center_r500,center_gm500,center_gf500,mlist500,rlist500,gmlist500,gflist500=read_file([500])
-
-if mlist2500!=None and len(mlist2500)!=len(mlist500):
-    raise Exception("Something wrong, the number of 2500 and 500 data are different")
-
-
-if mlist2500==None:
-    print("gas fraction between r2500 and r500 is %E"%((center_gm500-center_gm2500)/(center_m500-center_m2500)))
-    sys.exit(0)
-
-gf_2500_500=[]
-
-for i in range(0,len(mlist500)):
-    if mlist500[i]-mlist2500[i]<=0:
-        continue
-    gf_2500_500.append((gmlist500[i]-gmlist2500[i])/(mlist500[i]-mlist2500[i]))
-
-gf_2500_500.sort();
-
-
-center_gf_2500_500=(center_gm500-center_gm2500)/(center_m500-center_m2500)
-gf_idx=-1
-
-for i in range(len(gf_2500_500)-1):
-    if (center_gf_2500_500-gf_2500_500[i])*(center_gf_2500_500-gf_2500_500[i+1])<=0:
-        gf_idx=i
-        break
-if gf_idx==-1:
-    raise Exception("Something wrong!")
-    
-gflidx=int(gf_idx*(1-confidence_level))
-gfuidx=gf_idx-1+int((len(gf_2500_500)-gf_idx)*confidence_level)
-
-gferr1=gf_2500_500[gflidx]-center_gf_2500_500
-gferr2=gf_2500_500[gfuidx]-center_gf_2500_500
-
-print("gas_fraction between r2500 and r500=\t%e\t %e/+%e (1 sigma)"%(center_gf_2500_500,gferr1,gferr2))
diff --git a/mass_profile/fit_mass.sh b/mass_profile/fit_mass.sh
deleted file mode 100755
index f349f24..0000000
--- a/mass_profile/fit_mass.sh
+++ /dev/null
@@ -1,240 +0,0 @@
-#!/bin/sh
-#
-# Front-end script used to calculate the mass and related values.
-#
-# Output:
-#   * final_result.txt / center_only_results.txt
-#   * beta_param.txt / dbeta_param_center.txt
-#   * gas_mass_int_center.qdp
-#   * mass_int_center.qdp
-#   * nfw_fit_center.qdp
-#   * nfw_param_center.txt
-#   * overdensity_center.qdp
-#   * rho_fit_center.dat
-#   * rho_fit_center.qdp
-#   * sbp_fit_center.qdp
-#   * entropy_center.qdp
-#   * wang2012_param_center.txt
-#   * tprofile_dump_center.qdp
-#   * tprofile_fit_center.qdp
-#   * summary_mass_profile.qdp
-#   * summary_overdensity.qdp
-#   * summary_gas_mass_profile.qdp
-#   * summary_entropy.qdp
-#
-# Junhua Gu
-# Weitian LI
-# 2016-06-07
-#
-
-if [ $# -eq 1 ] || [ $# -eq 2 ]; then
-    :
-else
-    echo "usage:"
-    echo "    `basename $0` <mass.conf> [c]"
-    echo ""
-    echo "arguments:"
-    echo "    <mass.conf>: config file for mass profile calculation"
-    echo "    [c]: optional; if specified, do not calculate the errors"
-    exit 1
-fi
-
-if ! which xspec > /dev/null 2>&1; then
-    printf "*** ERROR: please initialize HEASOFT first\n"
-    exit 2
-fi
-
-export PATH="/usr/local/bin:/usr/bin:/bin:$PATH"
-
-base_path=$(dirname $(realpath $0))
-printf "## base_path: \`${base_path}'\n"
-mass_cfg="$1"
-printf "## use configuration file: \`${mass_cfg}'\n"
-case "$2" in
-    [cC])
-        F_C="YES"
-        ;;
-    *)
-        F_C="NO"
-        ;;
-esac
-
-nh=`grep            '^nh'            ${mass_cfg} | awk '{ print $2 }'`
-abund=`grep         '^abund'         ${mass_cfg} | awk '{ print $2 }'`
-nfw_rmin_kpc=`grep  '^nfw_rmin_kpc'  ${mass_cfg} | awk '{ print $2 }'`
-tprofile_data=`grep '^tprofile_data' ${mass_cfg} | awk '{ print $2 }'`
-tprofile_cfg=`grep  '^tprofile_cfg'  ${mass_cfg} | awk '{ print $2 }'`
-sbp_cfg=`grep       '^sbp_cfg'       ${mass_cfg} | awk '{ print $2 }'`
-
-# sbp config file
-sbp_data=`grep      '^sbp_data'      ${sbp_cfg} | awk '{ print $2 }'`
-tprofile=`grep      '^tprofile'      ${sbp_cfg} | awk '{ print $2 }'`
-cfunc_profile=`grep '^cfunc_profile' ${sbp_cfg} | awk '{ print $2 }'`
-z=`grep             '^z'             ${sbp_cfg} | awk '{ print $2 }'`
-cm_per_pixel=`cosmo_calc ${z} | grep 'cm/pixel' | awk -F':' '{ print $2 }'`
-sed -i'' "s/^cm_per_pixel.*$/cm_per_pixel   ${cm_per_pixel}/" ${sbp_cfg}
-
-if grep -q '^beta2' $sbp_cfg; then
-    MODEL="dbeta"
-    MODEL_NAME="double-beta"
-else
-    MODEL="beta"
-    MODEL_NAME="single-beta"
-fi
-
-PROG_TPROFILE="fit_wang2012_model"
-tprofile_dump="wang2012_dump.qdp"
-tprofile_param_center="wang2012_param_center.txt"
-tprofile_fit_center="tprofile_fit_center.qdp"
-tprofile_center="tprofile_dump_center.qdp"
-
-${base_path}/${PROG_TPROFILE} ${tprofile_data} ${tprofile_cfg} \
-            ${cm_per_pixel} 2> /dev/null | tee ${tprofile_param_center}
-cp -fv ${tprofile_dump} ${tprofile}
-mv -fv ${tprofile_dump} ${tprofile_center}
-mv -fv fit_result.qdp ${tprofile_fit_center}
-
-${base_path}/calc_coolfunc.sh ${tprofile_center} \
-            ${abund} ${nh} ${z} ${cfunc_profile}
-cfunc_profile_center="coolfunc_profile_center.txt"
-cp -f ${cfunc_profile} ${cfunc_profile_center}
-
-PROG_SBPFIT="fit_${MODEL}_sbp"
-RES_SBPFIT="${MODEL}_param.txt"
-RES_SBPFIT_CENTER="${MODEL}_param_center.txt"
-${base_path}/${PROG_SBPFIT} ${sbp_cfg} 2> /dev/null
-mv -fv ${RES_SBPFIT} ${RES_SBPFIT_CENTER}
-cat ${RES_SBPFIT_CENTER}
-mv -fv sbp_fit.qdp sbp_fit_center.qdp
-mv -fv rho_fit.qdp rho_fit_center.qdp
-mv -fv rho_fit.dat rho_fit_center.dat
-mv -fv entropy.qdp entropy_center.qdp
-${base_path}/fit_nfw_mass mass_int.dat ${z} ${nfw_rmin_kpc} 2> /dev/null
-mv -fv nfw_param.txt      nfw_param_center.txt
-mv -fv nfw_fit_result.qdp nfw_fit_center.qdp
-mv -fv nfw_dump.qdp       mass_int_center.qdp
-mv -fv overdensity.qdp    overdensity_center.qdp
-mv -fv gas_mass_int.qdp   gas_mass_int_center.qdp
-
-## only calculate central value {{{
-if [ "${F_C}" = "YES" ]; then
-    RES_CENTER="center_only_results.txt"
-    [ -e "${RES_CENTER}" ] && mv -f ${RES_CENTER} ${RES_CENTER}_bak
-    ${base_path}/analyze_mass_profile.py  200 c | tee -a ${RES_CENTER}
-    ${base_path}/analyze_mass_profile.py  500 c | tee -a ${RES_CENTER}
-    ${base_path}/analyze_mass_profile.py 1500 c | tee -a ${RES_CENTER}
-    ${base_path}/analyze_mass_profile.py 2500 c | tee -a ${RES_CENTER}
-    ${base_path}/fg_2500_500.py c               | tee -a ${RES_CENTER}
-    exit 0
-fi
-## central value }}}
-
-## ------------------------------------------------------------------
-
-# clean previous files
-rm -f summary_overdensity.qdp
-rm -f summary_mass_profile.qdp
-rm -f summary_gas_mass_profile.qdp
-rm -f summary_entropy.qdp
-
-# Estimate the errors of Lx and Fx by Monte Carlo simulation
-printf "\n+++++++++++++++++++ Monte Carlo +++++++++++++++++++++\n"
-MC_TIMES=100
-for i in `seq 1 ${MC_TIMES}`; do
-    ${base_path}/shuffle_profile.py ${tprofile_data} tmp_tprofile.txt
-    ${base_path}/shuffle_profile.py ${sbp_data} tmp_sbprofile.txt
-
-    # temperature profile
-    ${base_path}/${PROG_TPROFILE} tmp_tprofile.txt ${tprofile_cfg} \
-                ${cm_per_pixel} 2> /dev/null
-    mv -f ${tprofile_dump} ${tprofile}
-
-    TMP_SBP_CFG="tmp_sbp.cfg"
-    [ -e "${TMP_SBP_CFG}" ] && rm -f ${TMP_SBP_CFG}
-    cat ${sbp_cfg} | while read l; do
-        if echo "${l}" | grep -q '^sbp_data' >/dev/null; then
-            echo "sbp_data  tmp_sbprofile.txt" >> ${TMP_SBP_CFG}
-        elif echo "${l}" | grep -q '^tprofile' >/dev/null; then
-            echo "tprofile  ${tprofile}" >> ${TMP_SBP_CFG}
-        else
-            echo "${l}" >> ${TMP_SBP_CFG}
-        fi
-    done
-
-    printf "## ${i} / ${MC_TIMES} ##\n"
-    printf "## `pwd -P` ##\n"
-    ${base_path}/calc_coolfunc.sh ${tprofile} ${abund} ${nh} ${z} ${cfunc_profile}
-    ${base_path}/${PROG_SBPFIT} ${TMP_SBP_CFG} 2> /dev/null
-    cat ${RES_SBPFIT}
-    ${base_path}/fit_nfw_mass mass_int.dat ${z} ${nfw_rmin_kpc} 2> /dev/null
-    cat nfw_dump.qdp     >> summary_mass_profile.qdp
-    echo "no no no"      >> summary_mass_profile.qdp
-    cat overdensity.qdp  >> summary_overdensity.qdp
-    echo "no no no"      >> summary_overdensity.qdp
-    cat gas_mass_int.qdp >> summary_gas_mass_profile.qdp
-    echo "no no no"      >> summary_gas_mass_profile.qdp
-    cat entropy.qdp      >> summary_entropy.qdp
-    echo "no no no"      >> summary_entropy.qdp
-done  # end of `for'
-
-# recover the files of original center values
-cp -f ${cfunc_profile_center} ${cfunc_profile}
-cp -f ${tprofile_center} ${tprofile}
-printf "\n+++++++++++++++++ MONTE CARLO END +++++++++++++++++++\n"
-
-## analyze results
-RES_TMP="_tmp_result_mrl.txt"
-RES_FINAL="final_result.txt"
-[ -e "${RES_TMP}" ]   && mv -fv ${RES_TMP}   ${RES_TMP}_bak
-[ -e "${RES_FINAL}" ] && mv -fv ${RES_FINAL} ${RES_FINAL}_bak
-
-${base_path}/analyze_mass_profile.py  200 | tee -a ${RES_TMP}
-${base_path}/analyze_mass_profile.py  500 | tee -a ${RES_TMP}
-${base_path}/analyze_mass_profile.py 1500 | tee -a ${RES_TMP}
-${base_path}/analyze_mass_profile.py 2500 | tee -a ${RES_TMP}
-
-R200_VAL=`grep  '^r200'  ${RES_TMP} | awk '{ print $2 }'`
-R500_VAL=`grep  '^r500'  ${RES_TMP} | awk '{ print $2 }'`
-R1500_VAL=`grep '^r1500' ${RES_TMP} | awk '{ print $2 }'`
-R2500_VAL=`grep '^r2500' ${RES_TMP} | awk '{ print $2 }'`
-
-R200E=`grep   '^r200'             ${RES_TMP} | tail -n 1 | awk '{ print $2,$3 }'`
-R500E=`grep   '^r500'             ${RES_TMP} | tail -n 1 | awk '{ print $2,$3 }'`
-R1500E=`grep  '^r1500'            ${RES_TMP} | tail -n 1 | awk '{ print $2,$3 }'`
-R2500E=`grep  '^r2500'            ${RES_TMP} | tail -n 1 | awk '{ print $2,$3 }'`
-M200E=`grep   '^m200'             ${RES_TMP} | tail -n 1 | awk '{ print $2,$3 }'`
-M500E=`grep   '^m500'             ${RES_TMP} | tail -n 1 | awk '{ print $2,$3 }'`
-M1500E=`grep  '^m1500'            ${RES_TMP} | tail -n 1 | awk '{ print $2,$3 }'`
-M2500E=`grep  '^m2500'            ${RES_TMP} | tail -n 1 | awk '{ print $2,$3 }'`
-MG200E=`grep  '^gas_m200'         ${RES_TMP} | tail -n 1 | awk '{ print $2,$3 }'`
-MG500E=`grep  '^gas_m500'         ${RES_TMP} | tail -n 1 | awk '{ print $2,$3 }'`
-MG1500E=`grep '^gas_m1500'        ${RES_TMP} | tail -n 1 | awk '{ print $2,$3 }'`
-MG2500E=`grep '^gas_m2500'        ${RES_TMP} | tail -n 1 | awk '{ print $2,$3 }'`
-FG200E=`grep  '^gas_fraction200'  ${RES_TMP} | tail -n 1 | awk '{ print $2,$3 }'`
-FG500E=`grep  '^gas_fraction500'  ${RES_TMP} | tail -n 1 | awk '{ print $2,$3 }'`
-FG1500E=`grep '^gas_fraction1500' ${RES_TMP} | tail -n 1 | awk '{ print $2,$3 }'`
-FG2500E=`grep '^gas_fraction2500' ${RES_TMP} | tail -n 1 | awk '{ print $2,$3 }'`
-
-printf "\n+++++++++++++++ RESULTS (${MODEL_NAME}) +++++++++++++++\n"
-printf "model: ${MODEL_NAME}\n"                | tee -a ${RES_FINAL}
-cat ${RES_SBPFIT_CENTER}                       | tee -a ${RES_FINAL}
-printf "\n"                                    | tee -a ${RES_FINAL}
-printf "r200= ${R200E} kpc\n"                  | tee -a ${RES_FINAL}
-printf "m200= ${M200E} Msun\n"                 | tee -a ${RES_FINAL}
-printf "gas_m200= ${MG200E} Msun\n"            | tee -a ${RES_FINAL}
-printf "gas_fraction200= ${FG200E} x100%%\n"   | tee -a ${RES_FINAL}
-printf "r500= ${R500E} kpc\n"                  | tee -a ${RES_FINAL}
-printf "m500= ${M500E} Msun\n"                 | tee -a ${RES_FINAL}
-printf "gas_m500= ${MG500E} Msun\n"            | tee -a ${RES_FINAL}
-printf "gas_fraction500= ${FG500E} x100%%\n"   | tee -a ${RES_FINAL}
-printf "r1500= ${R1500E} kpc\n"                | tee -a ${RES_FINAL}
-printf "m1500= ${M1500E} Msun\n"               | tee -a ${RES_FINAL}
-printf "gas_m1500= ${MG1500E} Msun\n"          | tee -a ${RES_FINAL}
-printf "gas_fraction1500= ${FG1500E} x100%%\n" | tee -a ${RES_FINAL}
-printf "r2500= ${R2500E} kpc\n"                | tee -a ${RES_FINAL}
-printf "m2500= ${M2500E} Msun\n"               | tee -a ${RES_FINAL}
-printf "gas_m2500= ${MG2500E} Msun\n"          | tee -a ${RES_FINAL}
-printf "gas_fraction2500= ${FG2500E} x100%%\n" | tee -a ${RES_FINAL}
-printf "\n"                                    | tee -a ${RES_FINAL}
-${base_path}/fg_2500_500.py                    | tee -a ${RES_FINAL}
-printf "\n+++++++++++++++++++++++++++++++++++++++++++++++++++++\n"
diff --git a/mass_profile/fit_sbp.sh b/mass_profile/fit_sbp.sh
deleted file mode 100755
index 5989ff1..0000000
--- a/mass_profile/fit_sbp.sh
+++ /dev/null
@@ -1,61 +0,0 @@
-#!/bin/sh
-#
-# Handy script for SBP fitting.
-# This script wraps the 'fit_beta_sbp' and 'fit_dbeta_sbp',
-# and automatically determine the sbp model according to the config file.
-#
-# Weitian LI
-# 2013-02-20
-#
-# Change logs:
-# 2017-02-07, Weitian LI
-#   * Use `sbp_cfg` from command line argument, instead of the one specified
-#     in the `mass.conf`
-#   * Update the variable names according to the updated config files
-#   * Some cleanups
-#
-
-if [ $# -ne 2 ]; then
-    echo "usage: $0 <sbp.conf> <mass.conf>"
-    exit 1
-fi
-
-sbp_cfg="$1"
-mass_cfg="$2"
-
-if [ "$0" = `basename $0` ]; then
-    script_path=`which $0`
-    base_path=`dirname ${script_path}`
-else
-    base_path=`dirname $0`
-fi
-
-nh=`grep '^nh' ${mass_cfg} | awk '{ print $2 }'`
-abund=`grep '^abund' ${mass_cfg} | awk '{ print $2 }'`
-tprofile_data=`grep '^tprofile_data' ${mass_cfg} | awk '{ print $2 }'`
-tprofile_cfg=`grep '^tprofile_cfg' ${mass_cfg} | awk '{ print $2 }'`
-
-z=`grep '^z' ${sbp_cfg} | awk '{ print $2 }'`
-cm_per_pixel=`cosmo_calc ${z} | grep 'cm/pixel' | awk -F':' '{ print $2 }'`
-sed -i'' "s/^cm_per_pixel.*$/cm_per_pixel   ${cm_per_pixel}/" ${sbp_cfg}
-cfunc_profile=`grep '^cfunc_profile' ${sbp_cfg} | awk '{ print $2 }'`
-tprofile=`grep '^tprofile' ${sbp_cfg} | awk '{ print $2 }'`
-
-if grep -q '^beta2' ${sbp_cfg}; then
-    MODEL="double-beta"
-    PROG=fit_dbeta_sbp
-else
-    MODEL="single-beta"
-    PROG=fit_beta_sbp
-fi
-
-${base_path}/fit_wang2012_model ${tprofile_data} ${tprofile_cfg} \
-            ${cm_per_pixel} 2> /dev/null
-cp wang2012_dump.qdp ${tprofile}
-if [ ! -f ${cfunc_profile} ]; then
-    ${base_path}/calc_coolfunc.sh ${tprofile} ${abund} ${nh} ${z} \
-                ${cfunc_profile}
-fi
-${base_path}/${PROG} ${sbp_cfg}
-echo "## MODEL: ${MODEL}"
-echo "## z: ${z}"
diff --git a/mass_profile/get_lxfx_data.sh b/mass_profile/get_lxfx_data.sh
deleted file mode 100755
index dbb07dd..0000000
--- a/mass_profile/get_lxfx_data.sh
+++ /dev/null
@@ -1,89 +0,0 @@
-#!/bin/sh
-#
-# collect data for 'loop_lx.sh' & 'calc_lxfx_simple.sh'
-#
-
-if [ $# -lt 2 ]; then
-    printf "usage:\n"
-    printf "    `basename $0` <dir> [c] <delta ...>\n"
-    exit 1
-fi
-
-DIR="$1"
-shift
-case "$1" in
-    [cC]*)
-        F_C="YES"
-        shift
-        ;;
-    *)
-        F_C="NO"
-        ;;
-esac
-printf "CENTER_MODE: $F_C\n"
-echo "DELTA: $@"
-
-cd $DIR
-pwd -P
-
-INFO=`ls ../*_INFO.json 2> /dev/null`
-if [ ! -z "$INFO" ]; then
-    OI=`grep '"Obs\.\ ID' ${INFO} | sed 's/.*"Obs.*":\ //' | sed 's/\ *,$//'`
-    NAME=`grep '"Source\ Name' ${INFO} | sed 's/.*"Source.*":\ //' | sed 's/^"//' | sed 's/"\ *,$//'`
-    UNAME=`grep '"Unified\ Name' ${INFO} | sed 's/.*"Unified.*":\ //' | sed 's/^"//' | sed 's/"\ *,$//'`
-    Z=`grep '"redshift' ${INFO} | sed 's/.*"redshift.*":\ //' | sed 's/\ *,$//'`
-fi
-
-printf "# OI,NAME,UNAME,Z"
-if [ "${F_C}" = "YES" ]; then
-    for DELTA in $@; do
-        printf ",L${DELTA}(bolo),L${DELTA}(0.7-7),L${DELTA}(0.1-2.4),F${DELTA}(bolo),F${DELTA}(0.7-7),F${DELTA}(0.1-2.4)"
-    done
-    printf "\n"
-else
-    for DELTA in $@; do
-        printf ",L${DELTA}(bolo),L${DELTA}ERR(bolo),L${DELTA}(0.7-7),L${DELTA}ERR(0.7-7),L${DELTA}(0.1-2.4),L${DELTA}ERR(0.1-2.4),F${DELTA}(bolo),F${DELTA}ERR(bolo),F${DELTA}(0.7-7),F${DELTA}ERR(0.7-7),F${DELTA}(0.1-2.4),F${DELTA}ERR(0.1-2.4)"
-    done
-    printf "\n"
-fi
-
-printf "# $OI,$NAME,$UNAME,$Z"
-
-if [ "${F_C}" = "YES" ]; then
-    for DELTA in $@; do
-        LX_RES="lx_result_${DELTA}_c.txt"
-        FX_RES="fx_result_${DELTA}_c.txt"
-        if [ -r ${LX_RES} ] && [ -r ${FX_RES} ]; then
-            Lbolo=`grep '^Lx(bolo' ${LX_RES} | awk '{ print $2 }'`
-            L077=`grep '^Lx(0\.7-7' ${LX_RES} | awk '{ print $2 }'`
-            L0124=`grep '^Lx(0\.1-2\.4' ${LX_RES} | awk '{ print $2 }'`
-            Fbolo=`grep '^Fx(bolo' ${FX_RES} | awk '{ print $2 }'`
-            F077=`grep '^Fx(0\.7-7' ${FX_RES} | awk '{ print $2 }'`
-            F0124=`grep '^Fx(0\.1-2\.4' ${FX_RES} | awk '{ print $2 }'`
-            printf ",$Lbolo,$L077,$L0124,$Fbolo,$F077,$F0124"
-        fi
-    done
-    printf "\n"
-else
-    for DELTA in $@; do
-        LX_RES="lx_result_${DELTA}.txt"
-        FX_RES="fx_result_${DELTA}.txt"
-        if [ -r ${LX_RES} ] && [ -r ${FX_RES} ]; then
-            Lbolo=`grep '^Lx(bolo' ${LX_RES} | awk '{ print $2 }'`
-            LboloERR=`grep '^Lx(bolo' ${LX_RES} | awk '{ print $4 }'`
-            L077=`grep '^Lx(0\.7-7' ${LX_RES} | awk '{ print $2 }'`
-            L077ERR=`grep '^Lx(0\.7-7' ${LX_RES} | awk '{ print $4 }'`
-            L0124=`grep '^Lx(0\.1-2\.4' ${LX_RES} | awk '{ print $2 }'`
-            L0124ERR=`grep '^Lx(0\.1-2\.4' ${LX_RES} | awk '{ print $4 }'`
-            Fbolo=`grep '^Fx(bolo' ${FX_RES} | awk '{ print $2 }'`
-            FboloERR=`grep '^Fx(bolo' ${FX_RES} | awk '{ print $4 }'`
-            F077=`grep '^Fx(0\.7-7' ${FX_RES} | awk '{ print $2 }'`
-            F077ERR=`grep '^Fx(0\.7-7' ${FX_RES} | awk '{ print $4 }'`
-            F0124=`grep '^Fx(0\.1-2\.4' ${FX_RES} | awk '{ print $2 }'`
-            F0124ERR=`grep '^Fx(0\.1-2\.4' ${FX_RES} | awk '{ print $4 }'`
-            printf ",$Lbolo,$LboloERR,$L077,$L077ERR,$L0124,$L0124ERR,$Fbolo,$FboloERR,$F077,$F077ERR,$F0124,$F0124ERR"
-        fi
-    done
-    printf "\n"
-fi
-
diff --git a/mass_profile/shuffle_profile.py b/mass_profile/shuffle_profile.py
deleted file mode 100755
index 124f505..0000000
--- a/mass_profile/shuffle_profile.py
+++ /dev/null
@@ -1,37 +0,0 @@
-#!/usr/bin/env python3
-#
-# Shuffle the profile data point values according to their errors.
-#
-# Weitian LI
-# 2017-02-07
-
-import sys
-import numpy as np
-
-
-if len(sys.argv) != 3:
-    print("Usage: %s <input_profile> <shuffled_profile>")
-    sys.exit(1)
-
-
-# 4-column data: radius, err, temperature/brightness, err
-data = np.loadtxt(sys.argv[1])
-
-x1 = data[:, 2]
-xe = data[:, 3]
-x2 = np.zeros(shape=x1.shape)
-
-for i in range(len(x2)):
-    if x1[i] <= 0 or xe[i] <= 0:
-        # Skip shuffle
-        x2[i] = x1[i]
-
-    v = -1.0
-    while v <= 0:
-        v = np.random.normal(0.0, 1.0) * xe[i] + x1[i]
-    x2[i] = v
-
-# Replace original values
-data[:, 2] = x2
-
-np.savetxt(sys.argv[2], data)