aboutsummaryrefslogtreecommitdiffstats
diff options
context:
space:
mode:
-rw-r--r--xspec/xspec_avg_tz.tcl109
-rw-r--r--xspec/xspec_tprofile.tcl14
2 files changed, 115 insertions, 8 deletions
diff --git a/xspec/xspec_avg_tz.tcl b/xspec/xspec_avg_tz.tcl
new file mode 100644
index 0000000..7720f64
--- /dev/null
+++ b/xspec/xspec_avg_tz.tcl
@@ -0,0 +1,109 @@
+# Copyright (c) 2017 Weitian LI <liweitianux@live.com>
+# MIT license
+#
+# XSPEC Tcl script to calculate the average temperature and abundance
+# by ting the temperatures and abundances of all regions.
+#
+# NOTE:
+# * Use XSPEC model `projct*wabs*apec`
+# * if error < 0.01, then assume error = 0.01
+# * if fabs(error) < 1e-7, then set error as `NULL' in output
+# (this may be caused by frozen or tied parameters)
+#
+# Output:
+# * tz_average.txt : average temperature and abundance (and errors)
+#
+# Change logs:
+# 2017-02-06, Weitian LI
+# * Initial version based on `xspec_tprofile.tcl`
+
+
+set NAME "xspec_avg_tz.tcl"
+set DATE [ exec date ]
+# Set the error flag
+set ERR_FLAG "FALSE"
+
+# Return TCL results for XSPEC commands.
+set xs_return_result 1
+
+# Keep going until fit converges.
+query yes
+
+# Get the value of the specified parameter
+proc get_value {pn} {
+ tclout param $pn
+ scan $xspec_tclout "%f" value
+ return $value
+}
+
+# Calculate the error of the specified parameter
+proc get_error {pn} {
+ global ERR_FLAG
+ error 1.0 $pn
+ tclout error $pn
+ scan $xspec_tclout "%f %f" val_l val_u
+ set err [ expr {($val_u - $val_l) / 2.0} ]
+ # error treatment
+ if {abs($err) < 1.0e-7} {
+ puts "*** WARNING: invalid error value ***"
+ set err "NULL"
+ set ERR_FLAG "TRUE"
+ } elseif {$err < 0.01} {
+ set err 0.01
+ }
+ return $err
+}
+
+# Output file
+set tz_fn "tz_average.txt"
+if {[ file exists $tz_fn ]} {
+ exec mv -fv $tz_fn ${tz_fn}_bak
+}
+set tz_fd [ open $tz_fn w ]
+
+# Header
+puts $tz_fd "# Average temperature and abundance"
+puts $tz_fd "#"
+puts $tz_fd "# Generated by: ${NAME}"
+puts $tz_fd "# Created: ${DATE}"
+puts $tz_fd "#"
+puts $tz_fd "# Item Value Error"
+
+# Number of data group
+tclout datasets
+set datasets $xspec_tclout
+
+# Untie and thaw the temperature and abundance of the first region
+untie 5 6
+thaw 5 6
+
+# Tie all other temperature and abundance to that of the first region
+for {set i 2} {$i <= ${datasets}} {incr i} {
+ # Parameter number of temperature and abundance
+ set temp_pn [ expr {8 * $i - 3} ]
+ set abund_pn [ expr {8 * $i - 2} ]
+ newpar $temp_pn = 5
+ newpar $abund_pn = 6
+}
+
+# Fit the spectra again
+fit
+
+# Get the values and errors of average temperature and abundance
+set temp_val [ get_value 5 ]
+set temp_err [ get_error 5 ]
+set abund_val [ get_value 6 ]
+set abund_err [ get_error 6 ]
+
+# Output the average temperature and abundance
+puts "Average temperature: $temp_val, $temp_err"
+puts "Average abundance: $abund_val, $abund_err"
+puts $tz_fd "temp $temp_val $temp_err"
+puts $tz_fd "abund $abund_val $abund_err"
+
+close $tz_fd
+
+# Print a WARNING if there are any errors
+if {[ string equal $ERR_FLAG "TRUE" ]} {
+ puts "*** WARNING: there are invalid error values ***"
+}
diff --git a/xspec/xspec_tprofile.tcl b/xspec/xspec_tprofile.tcl
index e792d3c..ec7263c 100644
--- a/xspec/xspec_tprofile.tcl
+++ b/xspec/xspec_tprofile.tcl
@@ -2,20 +2,18 @@
## XSPEC Tcl script to extract the temperature profile from the
## fitted *deprojection spectra analysis* `projct*wabs*apec` results.
##
-## Requirements:
+## NOTE:
## * Spectra loaded *in order* (in order to correctly calculate radii)
## * Use XSPEC model `projct*wabs*apec`
-##
-## NOTES:
## * if error < 0.01, then assume error = 0.01
## * if fabs(error) < 1e-7, then set error as `NULL' in output
## (this may be caused by frozen or tied parameters)
##
## Output:
-## * xspec_saveall.xcm (save all for the current state)
-## * fitted_params.log (log the fitted parameters)
-## * tprofile.qdp (temp profile QDP file for graphing)
-## * tprofile.txt (temp profile data file for later fit)
+## * xspec_saveall.xcm : save all for the current state
+## * fitted_params.log : log the fitted parameters
+## * tprofile.qdp : temp profile QDP file for graphing
+## * tprofile.txt : temp profile data file for later fit
##
## Weitian LI <liweitianux@live.com>
## Created: 2012-08-12
@@ -35,7 +33,7 @@
## basic variables {{{
-## record process date
+set NAME "xspec_avg_tz.tcl"
set DATE [ exec date ]
set ERR_FLAG "FALSE"
## basic vars }}}