diff options
| -rw-r--r-- | xspec/xspec_tprofile.tcl | 160 |
1 files changed, 34 insertions, 126 deletions
diff --git a/xspec/xspec_tprofile.tcl b/xspec/xspec_tprofile.tcl index 0bfcbf5..e792d3c 100644 --- a/xspec/xspec_tprofile.tcl +++ b/xspec/xspec_tprofile.tcl @@ -1,50 +1,38 @@ -##################################################################### -## XSPEC Tcl script ## -## Task: -## To generate `Temperature profile' -## and `cooling function data file' -## according to `deprojection spectra analysis' results +## XSPEC Tcl script to extract the temperature profile from the +## fitted *deprojection spectra analysis* `projct*wabs*apec` results. ## ## Requirements: -## deprojection spectra analysis, -## spectra loaded *in order* (in order to correctly calc *radius*) -## use XSPEC model `projct*wabs*apec' (AS IT IS) -## because this script cannot handle other *parameter numbers* +## * Spectra loaded *in order* (in order to correctly calculate radii) +## * Use XSPEC model `projct*wabs*apec` ## ## NOTES: -## *solar abundance std* assumed to be `GRSA' (1998) -## 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 params) +## * 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: -## tcl_xspec_saveall.xcm (save all for the current state) -## tcl_fitted_params.log (log the fitted parameters) -## tcl_temp_profile.qdp (temp profile QDP file for graphing) -## tcl_temp_profile.txt (temp profile data file for later fit) -## tcl_coolfunc_dat.txt (coolfunc 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) ## -## LIweitiaNux <liweitianux@gmail.com> -## August 12, 2012 +## Weitian LI <liweitianux@live.com> +## Created: 2012-08-12 ## -## v1.1, 2012/08/12 -## improve description -## modify error treatment -## v2.0, 2012/08/14, LIweitiaNux -## add `output a QDP file for temperature profile' -## cancel `calc errors of abundance' -## improve error treatment (add warning) -##################################################################### +## Change logs: +## 2017-02-05, Weitian LI +## * Drop the cooling function calculation +## * Rename to 'xspec_tprofile.tcl' +## * Cleanups +## v2.0, 2012-08-14, Weitian LI +## * add `output a QDP file for temperature profile' +## * cancel `calc errors of abundance' +## * improve error treatment (add warning) +## v1.1, 2012-08-12 +## * improve description +## * modify error treatment -## about {{{ -set NAME "xspec_coolfunc_v2.tcl" -set VERSION "v2, 2012-08-14" -## about }}} - -## ask `normalization' used to calculate cooling function -puts -nonewline "LY> Input the NORMALIZATION used to calc cooling function: " -set norm [ gets stdin ] ## basic variables {{{ ## record process date @@ -53,10 +41,6 @@ set ERR_FLAG "FALSE" ## basic vars }}} ## xspec settings {{{ -# solar abundance std -# use `grsa' (1998) -set abund_std "grsa" - # Return TCL results for XSPEC commands. set xs_return_result 1 @@ -65,7 +49,7 @@ query yes ## xspec settings }}} ## xspec, save current results {{{ -set fn_save "tcl_xspec_saveall.xcm" +set fn_save "xspec_saveall.xcm" if {[ file exists $fn_save ]} { exec mv -fv $fn_save ${fn_save}_bak } @@ -73,7 +57,7 @@ save all $fn_save ## xspec save all }}} ## xspec log, log current fitted params {{{ -set fn_log "tcl_fitted_params.log" +set fn_log "fitted_params.log" if {[ file exists $fn_log ]} { exec mv -fv $fn_log ${fn_log}_bak } @@ -83,12 +67,8 @@ log none ## xspec log }}} ## set output file {{{ -set tpro_qdp_fn "tcl_temp_profile.qdp" -set tpro_fn "tcl_temp_profile.txt" -set cool_fn "tcl_coolfunc_dat.txt" -# tmp file used to calc `cooling function' -set tmpxspec_fn "_tcl_xspec_tmp.xcm" -set tmpres_fn "_tcl_coolres_tmp.txt" +set tpro_qdp_fn "tprofile.qdp" +set tpro_fn "tprofile.txt" # check file status if {[ file exists $tpro_qdp_fn ]} { exec mv -fv $tpro_qdp_fn ${tpro_qdp_fn}_bak @@ -96,24 +76,14 @@ if {[ file exists $tpro_qdp_fn ]} { if {[ file exists $tpro_fn ]} { exec mv -fv $tpro_fn ${tpro_fn}_bak } -if {[ file exists $cool_fn ]} { - exec mv -fv $cool_fn ${cool_fn}_bak -} set tpro_qdp_fd [ open $tpro_qdp_fn w ] set tpro_fd [ open $tpro_fn w ] -set cool_fd [ open $cool_fn w ] ## output files }}} ## datasets, number of data group tclout datasets set datasets $xspec_tclout -## get `nH' and `redshift' at first -tclout param 4 -scan $xspec_tclout "%f" nh -tclout param 7 -scan $xspec_tclout "%f" redshift - ## QDP file header {{{ puts $tpro_qdp_fd "! generated by: ${NAME}" puts $tpro_qdp_fd "! created date: ${DATE}" @@ -144,7 +114,7 @@ for {set i 1} {$i <= ${datasets}} {incr i} { set r [ expr {($r_in + $r_out) / 2.0} ] set r_err [ expr {($r_out - $r_in) / 2.0} ] puts "radius: $r, radius_err: $r_err ($r_in, $r_out)" - # output values: temperature and abundance + # output fitted results # determine the param number of temperature and abundance # temperature set temp_pn [ expr {8 * $i - 3} ] @@ -163,84 +133,22 @@ for {set i 1} {$i <= ${datasets}} {incr i} { set temp_err 0.01 } puts "temperature (p${temp_pn}): $temp_val, $temp_err" - # abundance - set abund_pn [ expr {8 * $i - 2} ] - tclout param $abund_pn - scan $xspec_tclout "%f" abund_val - #error 1.0 $abund_pn - #tclout error $abund_pn - #scan $xspec_tclout "%f %f" abund_val_l abund_val_u - #set abund_err [ expr {($abund_val_u - $abund_val_l) / 2.0} ] - #if {abs($abund_err) < 1.0e-7} { - # puts "*** WARNING: WRONG error values" - # set abund_err "NULL" - # set ERR_FLAG "TRUE" - #} elseif {$abund_err < 0.01} { - # set abund_err 0.01 - #} - #puts "abundance (p${abund_pn}): $abund_val $abund_err" + ## abundance + #set abund_pn [ expr {8 * $i - 2} ] + #tclout param $abund_pn + #scan $xspec_tclout "%f" abund_val # output `temp_profile' related results puts $tpro_qdp_fd "$r $r_err $temp_val $temp_err" puts $tpro_fd "$r $r_err $temp_val $temp_err" - - ## generate a xspec script to calc cooling function {{{ - # remove previous files - if {[ file exists $tmpxspec_fn ]} { - exec rm -fv $tmpxspec_fn - } - if {[ file exists $tmpres_fn ]} { - exec rm -fv $tmpres_fn - } - # open script file - set tmpxspec_fd [ open $tmpxspec_fn w ] - puts $tmpxspec_fd "set tmpres_fn \"$tmpres_fn\"" - puts $tmpxspec_fd "set tmpres_fd \[ open \$tmpres_fn w \]" - puts $tmpxspec_fd "query yes" - puts $tmpxspec_fd "statistic chi" - puts $tmpxspec_fd "method leven 10 0.01" - puts $tmpxspec_fd "abund ${abund_std}" - puts $tmpxspec_fd "xsect bcmc" - puts $tmpxspec_fd "cosmo 70 0 0.73" - puts $tmpxspec_fd "xset delta 0.01" - puts $tmpxspec_fd "systematic 0" - puts $tmpxspec_fd "dummyrsp 0.3 11.0 1024" - puts $tmpxspec_fd "model wabs*apec & ${nh} & ${temp_val} & ${abund_val} & ${redshift} & $norm & /*" - puts $tmpxspec_fd "flux 0.7 7.0" - puts $tmpxspec_fd "tclout flux 1" - puts $tmpxspec_fd "scan \$xspec_tclout \"%f %f %f %f\" holder holder holder res_flux" - puts $tmpxspec_fd "puts \$res_flux" - puts $tmpxspec_fd "puts \$tmpres_fd \$res_flux" - puts $tmpxspec_fd "close \$tmpres_fd" - puts $tmpxspec_fd "tclexit" - close $tmpxspec_fd - ## xspec script }}} - ## calc cooling function and output results {{{ - # exec XSPEC to calculate - exec xspec - $tmpxspec_fn - # read results from files - set tmpres_fd [ open $tmpres_fn r ] - set tmpline [ gets $tmpres_fd ] - scan $tmpline "%f" coolf_val - puts "cooling function: $coolf_val" - close $tmpres_fd - puts $cool_fd "$r $coolf_val" - ## cooling function }}} -} - -## remove tmp files -if {[ file exists $tmpres_fn ]} { - exec rm -fv $tmpres_fn } ## close file close $tpro_qdp_fd close $tpro_fd -close $cool_fd ## check `ERR_FLAG', print WARNING if {[ string equal $ERR_FLAG "TRUE" ]} { puts "*** WARNING: there are WRONG error values" } -# EOF # vim: set ts=8 sw=4 tw=0 fenc= ft=tcl: # |