XSPECのメモページ
XSPECに関してのリンク一覧
本家のサイト。XSPECの(非常に分厚い)マニュアルもある。
マニュアルは
こちら。ちなみにこれはver12.僕が使っているのはver11.変更点についてはマニュアルを参照してくださいませ。
以下、解析の際に参考にするサイト
吸収やモデルを考慮した状況下での、flux変換などが、ウェブ上でできる。
XSPECに関係するもののリンク
XSPECの拡張版 (broad band SEDなどのplot, fittingをしたい人向け) として、ISISというのがあるらしい。multi-wavelength dataに対して、XSPEC likeな使い方をしたい人向け。
XSPECとは
XSPECの設定ファイルを作る
bashで言うところの.bashrcとか.bash_profile的なものを作ることができる。マックユーザーなら、
locate .xspec
として、まずは.xspecの場所を探し、.xspecに移動する。
移動したら、そこでxspec.rcを作成する
emacs -nw xspec.rc
そのxspec.rcファイルに、XSPECを呼び出すときの初期設定を書きこむ。例えば
#show PGPLOT Window
cpd /xw
#Plot parameters
setplot energy
と言った感じ。自分の好きなようにカスタマイズしていけばよい。
XSPECに入ってからSpectrumを描くまで
ここでは、XSPECに入ってから、dataを読み込み、Spectrumを描くところまでを記す。
xspecにログイン
xspec
でxspecにログインできる。
スペクトルをplotする
ターミナルで
xspec
と打ってログインし、読み込みたいデータセットを入れていく。
具体的には
XSPEC>data 1:1 xis0_3_srcmath.pha (ソースファイル)
XSPEC>back 1 xis0_3_bgd.pha (バックグラウンドファイル)
XSPEC>resp 1 xis0_3_rsp.rmf (レスポンスファイル)
XSPEC>arf 1 xis0_3_arf.arf (arfファイル)
など。ちなみに、2つめのデータセットは1:1を2:2とし、1のところを2とする。例えば、
XSPEC>data 2:2 xis1_srcmath.pha (ソースファイル)
XSPEC>back 2 xis2_bgd.pha (バックグラウンドファイル)
XSPEC>resp 2 xis2_rsp.rmf (レスポンスファイル)
XSPEC>arf 2 xis2_arf.arf (arfファイル)
など。レスポンスファイルとarfファイルをまとめてレスポンスファイルにまとめる場合もあって、その場合は(当然だが)arfファイルはかかなくてよい。
どんなデータを読み込んだかを確認したいときは
XSPEC>show files
で確認できる。
Sectrumを表示する
dataを読み込んだら、以下の方法でSpectrumを表示できる
XSPEC>setplot ene
XSPEC>plot ldata
上の"setplot ene"は横軸をenergy(kev)で表す、という意味。
channelで表したいときは、"setplot ch"とすればよい。Defaultでは"setplot"は"setplot ch"となる。
"plot ldata"はplot loading dataの意味。Spectrumが表示される。色は1,2,3の順に白・赤・緑で表示される。
以下は、フィッティングのあとに見ればよいが、
XSPEC>plot ldata del
で、モデルフィッティングしたχのグラフが下に現れる。
範囲を指定する
energy(channel)の範囲を指定したいときは、以下のようにする。
XSPEC>ignore **-2.0 10.0-**
これで、「2keV以下、10keV以上を切る」という意味になる。energyで切るときは、必ず小数点以下まで書くこと。そうしないとchannelを切ると判断されてしまう。
XSPEC>ignore 2.0-10.0
とすると、2.0~10.0keVを無視する、となる。ちなみに、
何度もignoreをすると、すべてのignoreの範囲の「かつ(and)」を取ることになるので注意。
ちなみに、これだとすべてのdatasetを上の範囲で切ることになるが、例えばdataset 1だけをきりたければ、
XSPEC>ignore 1:**-1.0 10.0-**
などとすると、dataset 1 だけが「1keV以下、10keV以上を切る」という意味になる。
範囲の指定をはずす
範囲指定を外したいときは、"notice"を使う。
XSPEC>notice 0.0-50.0
とすれば、「0keVから50keVを表示する」となる。
特定のDatasetだけに適用したければ、ignoreと同様に、
XSPEC>notice n:hogehoge
とすればよい。ここでnはdatasetのnumber。
Spectrumを保存する
plotで図を表示した後、図を保存したければ、
XSPEC>iplot
でqdpに移動し、
PLT> hardcopy filename.ps/cps
でps fileが作られる。ちなみに、
PLT>wenviron filename
としておくと、qdpとpcoファイルが出来て、後で図の編集が可能になる。
ちなみに、labelをはりたいときは、
label top hogehoge
などとすると、Topの名前がかわる。よく使うものは、
la x hogehoge (x軸の名前をhogehogeにする)
la y hogehoge (y軸の名前をhogehogeにする)
time off (右下にある時刻表示をoffにする)
fo ro (fontをromanにする)
line of/off on n (n番をline表示にする。nはline number)
lw m on n (n番の線の太さを変える。mは自然数, nはline number)
lw m on n1..n2 (n1からn2番までの線の太さをmに変える)
lw m (図の枠の太さをmに変える)
ls m on n (n番のcomponentの線のスタイルをmにかえる。m=1,2,3,4,5まで。詳細はls ?で見よ。)
color m on n (mは自然数, nはline number)
error off 2 (fitting結果をlineとして表示)
cs 1.5 (title等を大きく表示)
lwidthのmは大きいほど線が太くなる。
colorの番号mなどについては、2は赤、4は青など。
lsの各番号がどのようなスタイルかを調べたければ
ls ?
とすればよい。
また、各モデルのcomponentが何番であるのかをチェックしたいときは、iplotのあとに
info
と打てばわかる。
このラベルを替えたものを保存したければ、このあとに、hardcopyをする。
またxspecにもどりたいときは、
PLT>exit
とすれば
XSPEC>
に戻る。
iplot内で使用するコマンドをxspecで使用する
setplot command を使用すれば、iplotで使用するコマンドをxspecでも使うことができる。
setplot co color 4 on 2
setplot co r y ymin ymax
setplot co r x xmin xmax
など。
setp co r x/y xmin/ymin xmax/ymax
は特定のx軸の領域、あるいはy軸の領域のみを拡大表示したいときに使える。
読み込んだdata、設定の保存
spectrumを描くために読み込んだdataや設定(ignoreとかfitting(後で出てくる)のparameterとか)をまた読み込むのはしんどいし、必ずミスする。そんなときは、
XSPEC>save all filename.xcm
として、今の設定を保存する。また読みこみたいときは、
XSPEC>@filename.xcm
とすれば、読み込みができる。心配なら、show files,show allなどでチェックすればよい。
読み込んだデータをASCII dataとして取得する
fitting結果やスペクトルデータ等を、iplotじゃなくて、自分の好きなplotterで表示したい場合もあるかもしれない。そういう時は以下のようにすればよい。
表示したいパラメータ (xとかyとか) を選んで、
XSPEC>tclout plot ldata x
XSPEC>puts $xspec_tclout
7.84725 7.93571 8.02418 ...
などとすると、それらのパラメータが並ぶ。これらをlocalのASCII fileに落としたいときは、logコマンドなどを使ってやれば、保存できる。
XSPEC>tclout plot ldata x
XSPEC>log ascii_spec_x.tbl
XSPEC>puts $xspec_tclout
XSPEC>log none
yの値も同じようにして、後でpasteなりしてくっつければよい。
Fittingをする
ここはまだ工事中です。
Dataを読み込んで描いたSpectrumをいろいろなモデルでfittingさせる。
方法はいくつかある。
- 対象天体の論文で、fitting modelを探す
- 先生（or先輩）に聞く
- 自分で考える
モデルがまだよくわかっていない自分には上の2つの選択肢がよく使われる。
Modelの指定
上で探したfitting modelを実際に入れてやる。
そのためには、
装置依存のconstantや吸収のphabs、系外のbremsstrahlung であるzbremss、powerlawのzpowerlwなどを入れていく。
XSPEC>model constant*phabs(zbremss+ ( zpowerlw)zhighect + pexrav )
とすれば、modelを適用することができる。
このときにいろいろparameterを聞かれるが、何を聞かれるかは、xspecで事前に調べておく。(column densityであったり、zであったり・・・etc)
Modelの作成
XSPECのbuilt-inにはないモデルを使いたい場合、解析的な函数であれば、mdefineというモデル作成コマンドで、local modelを作ることができる。
使い方は
mdefine modelname func
という形でかける。例えば、(すでに存在するけど、)黒体輻射モデルを描きたいときは、
mdefine junkbb E**-3/(exp(E/(k*T))-1)
などと書けばよい。この時、Eの関数で、パラメータはk, Tである。実際のところ、kはボルツマン定数でfreezeさせ、Tをフリーとして動かせば良い。また、すごい特殊な事情で、横軸を波長λ[um]・縦軸をBv(λ)[Jy]で表現しなければいけない場合、
mdefine jbb E^(-3)/(exp(a1/(T*E))-1)
fitting modelについて
constについては、
すざくメモのAll listsから最適の値についてを探すことができる。
ちなみに、2009年12月現在は
これの8ページ目を読めばよい。これによると、xis0=1.000に対して、xis1=1.007+/- 0.011,xis3=1.009+/-0.011,PIN=1.181+/-0.016となる。
モデルを追加する
modelをした後に、modelを追加したり・削除したりとあれやこれやしたくなる。そういうときは、おしなべてeditmodというコマンドを使ってやればよい。具体的には、今modelとして、
constant*phabs(zbremss+ ( zpowerlw)zhighect + pexrav )
というのがあったとして、zphabsをzhighectのとなりにかけたい、とする。そのときは、
XSPEC>editmod constant*phabs(zbremss+ ( zpowerlw)zhighect*zphabs + pexrav )
としてやればよい。このあとは、さきほどのmodelと同様にparameterに値を入れていってやる。ちなみに、modelと違って、editmodはコマンド1回につき、1つずつしかモデルを追加・削除できない点には注意しておく。ちなみに、editmodしたあとは、
XSPEC>renorm
をしておいてやると、これはmodelがrenormalizeされて、後々のfitが速くすむ(はず)である。
modelの追加方法その2(読み飛ばしてもよい)
editmod以外にもmodelを追加する方法がある。それがaddcompというコマンド. show paraなどとすると、各modelの番号が<>で囲まれている。これをcomponent numberという。新しいmodelを追加したいだけだったら、いちいちeditmodでたくさん書くよりかは、
addcomp n gau
(ここでnはcomponent number)
などとして、新しいmodelを入れてあげればよい。
modelの分布をみる
上のモデルがそれぞれどのように効いているのかを区別して見たければ、
XSPEC>plot eeu
とすれば、data,モデル成分のそれぞれを見ることができる。
ちなみに、dataとfittingの度合いを見るときは、
XSPEC>plot ldata del
が有用。こうすれば、fittingとdataの合い具合がわかる。ここでldataにすると、y軸がlogスケールに変わる。つけなければ、y軸はlinearなまま。ちなみに、modelとdataのratioを出したいときは、
XSPEC>plot ldata ratio
とすればよい。また、modelの各componentがどの程度contributeしているかを見たい場合は、
XSPEC>setp add
と打てば良い。
有効面積を見たい
各観測装置の有効面積をみたければ、
XSPEC>plot eff
とすればよい。
SUZAKU HXDのPINは、12keV以下は有効面積と相性が悪くなるらしい。くわしくは
こちら。
また、modelだけのプロットがしたい場合は、
XSPEC>plot eem
とすればよい。
fittingさせる
上で定めたmodelをfittingさせるには、
XSPEC>fit
とすればよい。
ちなみに、fittingの途中で、
Continue fitting?
などと聞かれるので、毎回yesと答えるのが面倒臭いのであれば、
XSPEC>query yes
と答えておくと、
Querying disabled - assuming answer is yes
となって、ずっとfittingを続けてくれる状態をつくれる。そのあと、
XSPEC>fit
とすればよい。
ちなみに、query yesにしたままおいておくと、fittingのあとも勝手にyesと答えてしまう仕様になっている。元の状態に戻したければ、
XSPEC>query on
としておけばよい。
このあとは、fittingが終了したあと、plotしてfitting度合を見て、足りなければ
XSPEC>fit
を何度か繰り返す。fitを繰り返す目安としては、"Reduced chi-squared"が減りつづけている間は、fitを繰り返すとよい。
free parameterの動きをみる。
free parameterを適当な(適切な、という意味ではない)値を入れてみたものの、それが適切かどうかはわからない場合、stepparで、そのfree parameterの変動によって、どのようにxi squredが動くかがわかる。具体的には
XSPEC>steppar 17 13 20 7
などと使う。これは、「paprameter番号17を、13から20まで7等分して見たときの、それぞれのxi squredを教えてくださいまし」という意味になる。なので、一般にはparameter番号をn,動かしたい範囲を[x1,x2]として、m等分したいときは、
XSPEC>steppar n x1 x2 m
とする。
free parameterを確認する
今freeparameterが何なのかを知りたいときは
XSPEC>show free
でわかる。
parameterの変換
いろいろとparameterを入れた後、parameterを直したいときがある。この時は
XSPEC>newpar n x,+/-1,min,min,Max,Max
とする。ここで、nはparameterの番号、xは入れ直したいparameterの値、+/-1は、-1ならfixされたparameter、+1なら、free parameterになる。そして、minからMAXの範囲内に指定したいときは、",min,min,Max,Max"と入れる。ここで注意しなければならないのは、例えばDefaultの範囲が(min,Max)=(0,100)のとき、x=200などと指定して、最後のminとMaxを指定しないと、x=200は反映されず、勝手に0~100の範囲に修正される。なので、Defaultの範囲外に指定をするときは、かならず新しい範囲を指定しなければならない。
また、既存のparameterと同じにしたければ、
XSPEC>newpar n = m
とする。ここで、nが直したいparameter,mが既存のparameter。ここで注意しておくこととしては、
XSPEC>newpar n=m
などと=の間のスペースをなくすとあまりよろしくない、という点に注意しておく。
それ以外にも、各パラメータどおしの四則演算のものを入れたいとき、
具体的にはn番目のパラメータに、m番目とl番目のパラメータの和を入れたい場合は、
XSPEC>newpar n = m + l
とすればよい。また、とあるパラメータの2倍、としたい場合は
XSPEC>newpar n = 2.0 * m
というように、小数点を含めた形で書かないといけない。もし
XSPEC>newpar n = 2 * m
とすると、これは”2番めのパラメータとm番目のパラメータの積をn番目のパラメータの値とする”という意味になってしまう。
からまったparameterをほどく
今、例えばparameter番号25番が、
=par 8
で結ばれているとする。この関係を解消したいときは、
XSPEC>untie 8 25
とすればよい。
nをfix parameterにする
n番目のparameterをfixしたparameterにしたい場合は、
XSPEC>freeze n
とすればよい。
nをfree parameterにする
n(,m)番目のparameterをfreeparameterにしたい場合は、
XSPEC>thaw n m
とすればよい。
Errorの計算
Errorを計算する。非常に時間がかかる。
まず、
XSPEC>query yes
Querying disabled - assuming answer is yes
とすることで、「もし何か聞かれても、すべてyesと答えます」と宣言する。これで、ほっておいても勝手にErrorの計算をしてくれる。
Errorをかけていると、よく
Minimization may have run into a problem, check your result
とでることがある。これは、おそらく傾きの変動が多数あり、極小値が他にもあることが予測されるときにでる注意。どのparameterに対してこの注意が出ているのかを見て、各自でこの注意に耳を傾けるかどうかを決める。
XSPEC>show free
で、いまは何がfree parameterなのかを確認し、
3,5,11,17番がfree parameterなのであれば、
XSPEC>error 3 5 11 17
とする。
これで、reduced xi squred が下がれば、もう一度fitして、errorを計算して、をして、落ち着くまで繰り返す。
XSPECを出る
XSPECを出たいときは
XSPEC>quit
Do you really want to exit? yes
XSPEC: quit
とする。
Fe K line1のfitting
以下は、Feのfittingのひとつの方法である。
まず、上でfitしたdatasetのうち、
zpowerlow norm
だけを残して、他のfreeparameterをすべて、freezeさせる。
具体的には、
XSPEC>show free
をして、free parameterを確認して、zpowerlow norm以外のfree parameterの番号を
XSPEC>freeze n1 n2 n3
とfreezeしていく(n1,n2,n3はfree parameterの番号)。
ここで、Feのfitをするための準備として、
XSPEC>ignore 1:**-3.0
XSPEC>ignore 2:**-3.0
とxisの3.0keV以下をカットしておく。
これらのデータセットを保存する。具体的には、
XSPEC>save all xis_fe.xcm
などと保存しておいて、emacsなどでxis_fe.xcmを開いて、xisのデータセットだけを残す(PIN,GSOなどのデータセットは省いておく)。
ここで、いったんxspecを出て、xspecに入りなおす。
XSPEC>@xis_fe.xcm
XSPEC>fit
XSPEC>notice 5.0-7.0
として、Fe K lineのある5.0-7.0keV領域を復活させる。
ここから、Feのfitをする。Diskline成分を入れる場合と、zgaussを入れる場合を考える。
Diskline modelを入れる
editmodでdisklineを加える。もし系外の場合は、lineEには、E_{obs}を入れることに注意。
および
から、
が得られる。
[keV]などを入れる。
diskline:lineE>6.3,-1
Betor10>-3,-1
Rin>10,+1
Rout>10000,-1
などとする(上の値は例)。 Betor10はemissivity law。
renormをして、以下、fitをして、
errorをして・・・を繰り返す。
zgauss modelを入れる場合
editmodでzgaussを入れればよい。
1:zgauss:lineE>6.4,-1
1:zgauss:sigma>
1:zgauss:redshift>
などが聞かれる。
ちなみに、
WARNING*:RENORM:No variable model to allow renormalization
とか聞かれることがあるが、気にしない。これは、いまはrenormしたいのに、model のparameterとしてrenormを止めているので、出るのは当然ではある。
以下は、diskline modelを入れる場合と一緒。
番外編
ここでは、上の流れには直接関係ないものを記した。
XSPECを強制終了する。
XSPECを立ち上げているDirectoryで
ps aux | grep xspec
とすると
Name 23188 888 0.8 88888 88888 abc/8 R+ 15:22 11:23 xspec11
Name 23706 0.0 0.0 888 888 abc/9 S+ 16:10 0:00 grep xspec
などと表示される。消したい方を選んで、左から２番目の数字を入れる。
kill 23188
これで強制終了される。
XSPECでscreenを使う
screenコマンドをxspecでするときには注意が必要。
xspecの画面が出ているとscreenコマンドが上手く働かない（らしい）。
以下のような手順で進める。
screen
でscreenに入り、
export PGPLOT_TYPE=/NULL
xspec11
とすると、xspecの画面がでない。これをすると画面表示が使えないが、例えば「error計算で一晩ずっと計算を走らせたいけれど、このPC共用だしな・・・」という学部生特有の悩みをもているひともいるかもしれない。これで、そういう問題とはおさらば。
あとはターミナルを閉じて、error計算の度合を見たくなったときに、別のターミナルで
screen -r
とすれば、再読み込みがされる。
(でも、もっとうまい方法がありそう。。)
errorの計算短縮法？その1
errorの計算は時間がかかるので、いろいろな短縮方法(？)がある。その1つとして、normだけをfree parameterにして、計算してしまう方法がある。renormよりも信頼できる値になる(？)とかならないとか。
やり方は
XSPEC>show free
でfree parameterのnorm以外をすべて
XSPEC>freeze n1 n2 n3
としてやり、fittingやerrorをしてやる。ここで、n1,n2,n3はfree parameterを表す。
その後、
XSPEC>thaw n1 n2 n3
として、parameterをfreeに戻して、再度errorの計算などをする。
fとfiとfit
xspec11では、タスクを簡略化することができる。例えば、
XSPEC>show free
は
XSPEC>show f
となり、
XSPEC>show file
は
XSPEC>show fi
とできる。でも、
XSPEC>show fit
は
XSPEC>show fit
としなければならない、という不思議。
Equivalent Widthの概算を出す(工事中)
Equivalent Width(EW)が大体どれくらいなのかを出す方法。
今回は6.4keVにあるFe輝線の場合を記す。
まず、show freeなどで、Fe輝線をfittingしているmodelのnormに注目する。これは、次元が[photon/cm^2/s]で与えられている。
次に、
newpar 24 ,0
XSPEC>flux 6.3 6.5
として、6.3~6.5keVの間のcontinuumのfluxを計算する。
setplot
plotの表示を変えたいときは、
setplot hogehoge
とすればよい。XSPECのマニュアルのsetplotにもたくさんまとまっている。よく用いるのは、
setplot ene (x軸をenergy表示)
setplot xlog on/off (x軸のlog表示on/off)
plotの表示
基本は
plot ldata (del) (縦軸log表示. delを付けたら、fittingのreduced X^2も下に表示される)
plot data (縦軸linear表示)
plot eeu (モデルの各成分を見る)
など。
ログの取り方
上にも書いたように、
save all ~.xcm
と書くと、
data files
resp files
model
model paramers
fitting results
等々の情報が入ったデータが出来上がる。
ただし、これらは情報が多すぎて見づらい時がある。
そういうときは、以下のように、必要な情報だけをいれたログファイルを作成すればよい。
log hoge.log
show para
show fit
log none
などとすると、hoge.logというログファイルができる。
コマンドのみのログファイル
xspec上で使用したコマンドのみの履歴を残したいときは、
script ~.xcm
script none
で挟んであげたコマンドがすべてファイルに読み込まれる。
backgroundのcontributionを変化させる。
backgroundの影響をいじってみたいときは、corfileというのを用いる。
corfile 1 hogehoge.pha
として、
cornorm 1 0.90
などとすると、90%のbackgroundで引いてくれる(要編集)。
最終更新：2014年04月28日 14:02
