# XSPEC 使い方
黄色は藤原パソコン内情報
## フィッティングの流れ
①xspecの立ち上げ
②データの読み込み
③モデル読み込み→パラメータ書き込み
⑤フィッティング
## データ読み込み
### データ読み込み
観測データ、バックグラウンド、rmf、arfをxspecで読み込む。吉武さんのメモ通りに行っている場合、これらはひとまとめのデータとして紐づけされているはずなので以下のようにして読み込む。
```bash=
data 1:1 データファイル
```
2つ目のデータを読み込むのは下のようにする。
```bash=
data 2:2 データファイル
```
:::info
紐づけがうまくできていないときは4つのデータをひとつづつ読み込む必要がある。
```bash=
data dataファイル
resp rmfファイル
arf arfファイル
back bkgファイル
```
:::
### ignore
データはうまい具合に自分で切り取って使う。例えば一番低いエネルギーの点は必ず除去しないといけないし、高エネルギー側も統計に外れるような点は除去する必要がある。ignoreの使い方は以下の通り。
例：データ1の3eV以下と80eV以上をignoreするとき。
```bash=
ignore 1: **-3.0 80.0-**
```
.0までつけないとエネルギーではなく点の数でignoreされてしまうので注意する。
逆にnotice でignoreされた点を復活できる。
show daで読み込んだデータの確認ができる。
:::warning
吉武さんのメモで作ったデータは
nustar :rspがarfとrmfのまとめ
newton :rmfがarfとrmfのまとめ？
→mos1だけrmfがないけど大丈夫？
→nongrpとpiの違いは？
suzaku :まだ途中だけど多分BIはarfとrmfはまとめられていない。FIはaddascaspecするからまとめられる？
:::
## モデル読み込み
### モデル
```bash=
model モデルの形
```
を打つ。「モデルの形」には
```bash=
const*phabs*(phabs*zcutoffpl+const*zcutoffpl+pexmon)
```
のような、いろいろなモデルを四則演算等で結んだモデルを入れる。どんなモデルを入れるかは物理的状況を考えるか、論文を再現する場合はその論文の最初の方を見るかする。
show mo で読み込んだモデルを見ることができる。
また、モデルの理論スペクトルを見たいときは
```bash=
dummyrsp (energy min) (energy max) (bin number)
```
とし、pl emやpl eemを打つ。(bin number)の後にlogを入力するとbinがログスケールになる。
### パラメータ
モデルを読み込んだらすぐにパラメータを設定するよう求められる。パラメータは固定するパラメータとフィッティングによって動かされるパラメータの2つがある。前者は考える物理状況に合わせて書き込み、パラメータの後にスペースを一つ開けて-1を書くと固定される。後者は適当な数字を打つか、空欄のままでエンターを打つ。同じパラメータを二回設定しないといけないときは、二回目は「=p(数字)」で一回目のパラメータと紐づけることができる。
:::info
赤方偏移
vではなくzを代入する。おすすめのサイトはNASAのNED(論文に出展を書きやすいから)。
https://ned.ipac.caltech.edu/
:::
:::info
銀河水素柱密度
地球から目標天体を見た時にその視線方向の天の川銀河の水素柱密度を出してくれる神サイトがある。
https://www.swift.ac.uk/analysis/nhtot/
天体名を入れると柱密度が出てくるが、そのうち一番右の値を使うらしい(中谷さん)。
ただしxspecでは柱密度の単位が$[{10}^{22}atoms{cm}^{-2}]$なので例えば上の画像のような柱密度の時にはxspecには0.155を設定する。
:::
:::info
ノルム(norm)
ノルムの設定でモデルの形が上下するけどなんでかは忘れてしまった。あとで聞く。
:::
:::warning
念のために231201のモデルとパラメータ設定の様子を貼っておく。
```bash=
XSPEC12>model const*phabs*(phabs*zcutoffpl+const*zcutoffpl+pexmon)
Input parameter value, delta, min, bot, top, and max values for ...
1 0.01( 0.01) 0 0 1e+10 1e+10
1:data group 1::constant:factor>1 -1
1 0.001( 0.01) 0 0 100000 1e+06
2:data group 1::phabs:nH>0.155 -1
1 0.001( 0.01) 0 0 100000 1e+06
3:data group 1::phabs:nH>10
1 0.01( 0.01) -3 -2 9 10
4:data group 1::zcutoffpl:PhoIndex>1.9
15 0.01( 0.15) 0.01 1 500 500
5:data group 1::zcutoffpl:HighECut>370 -1
0 -0.01( 0.01) -0.999 -0.999 10 10
6:data group 1::zcutoffpl:Redshift>0013 -1
1 0.01( 0.01) 0 0 1e+20 1e+24
7:data group 1::zcutoffpl:norm>
1 0.01( 0.01) 0 0 1e+10 1e+10
8:data group 1::constant:factor>0.01
1 0.01( 0.01) -3 -2 9 10
9:data group 1::zcutoffpl:PhoIndex>=p4
15 0.01( 0.15) 0.01 1 500 500
10:data group 1::zcutoffpl:HighECut>=p5
0 -0.01( 0.01) -0.999 -0.999 10 10
11:data group 1::zcutoffpl:Redshift>=p6
1 0.01( 0.01) 0 0 1e+20 1e+24
12:data group 1::zcutoffpl:norm>=p7
2 0.01( 0.02) 1.1 1.1 2.5 2.5
13:data group 1::pexmon:PhoIndex>=p4
1000 -1( 10) 1 1 1e+06 1e+06
14:data group 1::pexmon:foldE>=p5
-1 -1( 0.01) -1e+06 -1e+06 1e+06 1e+06
15:data group 1::pexmon:rel_refl>-0.001 0.001 -1 -1 0 0
0 -0.01( 0.01) 0 0 4 4
16:data group 1::pexmon:redshift>=p6
1 -0.01( 0.01) 0 0 1e+06 1e+06
17:data group 1::pexmon:abund>1 -1
1 -0.01( 0.01) 0 0 100 100
18:data group 1::pexmon:Fe_abund>1 -1
60 1( 0.6) 0 0 85 85
19:data group 1::pexmon:Incl>60
1 0.01( 0.01) 0 0 1e+20 1e+24
20:data group 1::pexmon:norm>=p7
***Error: Desired Value 13 is outside hard range -0.999 - 10
Re-enter parameter value, delta, min, bot, top, and max values for ...
0 -1( 1) -0.999 -0.999 10 10
data group 1::zcutoffpl:Redshift>0.013 -1
Input parameter value, delta, min, bot, top, and max values for ...
1 -1( 0.01) 0 0 1e+10 1e+10
21:data group 2::constant:factor>1 1
0.155 -1( 0.00155) 0 0 100000 1e+06
22:data group 2::phabs:nH>
10 0.001( 0.1) 0 0 100000 1e+06
23:data group 2::phabs:nH>
1.9 0.01( 0.019) -3 -2 9 10
24:data group 2::zcutoffpl:PhoIndex>
370 -1( 3.7) 0.01 1 500 500
25:data group 2::zcutoffpl:HighECut>
0.013 -1( 0.00013) -0.999 -0.999 10 10
26:data group 2::zcutoffpl:Redshift>
1 0.01( 0.01) 0 0 1e+20 1e+24
27:data group 2::zcutoffpl:norm>
0.01 0.01( 0.0001) 0 0 1e+10 1e+10
28:data group 2::constant:factor>
= p4
29:data group 2::zcutoffpl:PhoIndex>
= p5
30:data group 2::zcutoffpl:HighECut>
= p6
31:data group 2::zcutoffpl:Redshift>
= p7
32:data group 2::zcutoffpl:norm>
= p4
33:data group 2::pexmon:PhoIndex>
= p5
34:data group 2::pexmon:foldE>
-0.001 0.001( 1e-05) -1 -1 0 0
35:data group 2::pexmon:rel_refl>
= p6
36:data group 2::pexmon:redshift>
1 -1( 0.01) 0 0 1e+06 1e+06
37:data group 2::pexmon:abund>
1 -1( 0.01) 0 0 100 100
38:data group 2::pexmon:Fe_abund>
60 1( 0.6) 0 0 85 85
39:data group 2::pexmon:Incl>
= p7
40:data group 2::pexmon:norm>
========================================================================
Model constant<1>*phabs<2>(phabs<3>*zcutoffpl<4> + constant<5>*zcutoffpl<6> + pexmon<7>) Source No.: 1 Active/On
Model Model Component Parameter Unit Value
par comp
Data group: 1
1 1 constant factor 1.00000 frozen
2 2 phabs nH 10^22 0.155000 frozen
3 3 phabs nH 10^22 10.0000 +/- 0.0
4 4 zcutoffpl PhoIndex 1.90000 +/- 0.0
5 4 zcutoffpl HighECut keV 370.000 frozen
6 4 zcutoffpl Redshift 1.30000E-02 frozen
7 4 zcutoffpl norm 1.00000 +/- 0.0
8 5 constant factor 1.00000E-02 +/- 0.0
9 6 zcutoffpl PhoIndex 1.90000 = p4
10 6 zcutoffpl HighECut keV 370.000 = p5
11 6 zcutoffpl Redshift 1.30000E-02 = p6
12 6 zcutoffpl norm 1.00000 = p7
13 7 pexmon PhoIndex 1.90000 = p4
14 7 pexmon foldE keV 370.000 = p5
15 7 pexmon rel_refl -1.00000E-03 +/- 0.0
16 7 pexmon redshift 1.30000E-02 = p6
17 7 pexmon abund 1.00000 frozen
18 7 pexmon Fe_abund 1.00000 frozen
19 7 pexmon Incl deg 60.0000 +/- 0.0
20 7 pexmon norm 1.00000 = p7
Data group: 2
21 1 constant factor 1.00000 +/- 0.0
22 2 phabs nH 10^22 0.155000 = p2
23 3 phabs nH 10^22 10.0000 = p3
24 4 zcutoffpl PhoIndex 1.90000 = p4
25 4 zcutoffpl HighECut keV 370.000 = p5
26 4 zcutoffpl Redshift 1.30000E-02 = p6
27 4 zcutoffpl norm 1.00000 = p7
28 5 constant factor 1.00000E-02 = p8
29 6 zcutoffpl PhoIndex 1.90000 = p9
30 6 zcutoffpl HighECut keV 370.000 = p10
31 6 zcutoffpl Redshift 1.30000E-02 = p11
32 6 zcutoffpl norm 1.00000 = p12
33 7 pexmon PhoIndex 1.90000 = p13
34 7 pexmon foldE keV 370.000 = p14
35 7 pexmon rel_refl -1.00000E-03 = p15
36 7 pexmon redshift 1.30000E-02 = p16
37 7 pexmon abund 1.00000 = p17
38 7 pexmon Fe_abund 1.00000 = p18
39 7 pexmon Incl deg 60.0000 = p19
40 7 pexmon norm 1.00000 = p20
________________________________________________________________________
Compton reflection from neutral medium.
See help for details.
If you use results of this model in a paper,
please refer to Magdziarz & Zdziarski 1995 MNRAS, 273, 837
Fit statistic : Chi-Squared 3.460244e+07 using 92 bins, spectrum 1, group 1.
Chi-Squared 1.105473e+09 using 112 bins, spectrum 2, group 2.
Total fit statistic 1.140076e+09 with 197 d.o.f.
Test statistic : Chi-Squared 1.140076e+09 using 204 bins.
Null hypothesis probability of 0.000000e+00 with 197 degrees of freedom
Current data and model not fit yet.
```
:::
一度設定したパラメータを書き換える場合は、
```bash=
newp (パラメータの番号) (新しいパラメータ)
```
を書き込む。
### 輝線
ガウシアンで輝線のモデルをいれる方法。最初に
```bash=
addc (輝線を入れる場所の番号) zgauss
```
を打つ。するとパラメータを入れるよう指示されるので、打つ。特にエネルギーEは直接輝線の位置に対応するので注意する。
下は輝線を入れてみた具体例。
```bash=
XSPEC12>show mo
Current model list:
Model constant<1>*phabs<2>(phabs<3>*zcutoffpl<4> + constant<5>*zcutoffpl<6> + pexmon<7>) Source No.: 1 Active/On
For Data Group(s): 1 2
Using energies from responses.
XSPEC12>addc 7 zgauss
Input parameter value, delta, min, bot, top, and max values for ...
6.5 0.05( 0.065) 0 0 1e+06 1e+06
13:data group 1::zgauss:LineE>2.3
0.1 0.05( 0.001) 0 0 10 20
14:data group 1::zgauss:Sigma>0.01 -1
0 -0.01( 0.01) -0.999 -0.999 10 10
15:data group 1::zgauss:Redshift>0.013 -1
1 0.01( 0.01) 0 0 1e+20 1e+24
16:data group 1::zgauss:norm>
Input parameter value, delta, min, bot, top, and max values for ...
6.5 0.05( 0.065) 0 0 1e+06 1e+06
37:data group 2::zgauss:LineE>
0.1 0.05( 0.001) 0 0 10 20
38:data group 2::zgauss:Sigma>
0 -0.01( 0.01) -0.999 -0.999 10 10
39:data group 2::zgauss:Redshift>
1 0.01( 0.01) 0 0 1e+20 1e+24
40:data group 2::zgauss:norm>
Fit statistic : Chi-Squared 152.07 using 92 bins, spectrum 1, group 1.
Chi-Squared 3.648008e+11 using 112 bins, spectrum 2, group 2.
Total fit statistic 3.648008e+11 with 195 d.o.f.
Test statistic : Chi-Squared 3.648008e+11 using 204 bins.
Null hypothesis probability of 0.000000e+00 with 195 degrees of freedom
Current data and model not fit yet.
```
zgaussのnormで輝線の大体の大きさ(高さ)を設定できるのでモデルの形をうまく実際のスペクトルの形に近づけてからフィッティングを行う。
delcomp モデル番号 でモデルを消せる
untie パラメータ番号1 パラメータ番号2 でパラメータの紐づけを解除できる
# errorの出し方
基本は
```bash=
error (free parameter number)
```
で出るがすぐにエラーが起こってしまうので、出てきたエラーと解決方法をまとめる。まとめたい。
:::info
```bash=
***Warning: Zero alpha-matrix diagonal element for parameter 14
Parameter 14 is pegged at 0.000466135 due to zero or negative pivot element, likely
caused by the fit being insensitive to the parameter.
```
fitを押してこれがでてきたときはフィットはできる？けど気になるときはパラメータ14があるモデルを消すとこのWarningが消える。
:::
:::info
```bash=
Parameter pegged at hard limit: 10
```
σは10から90に設定されているのでフィット中に10以下になろうとすると強制的に値が10で止まってしまう。
どうやって直したらいい？
:::
## pdp 出力
http://www.heal.phy.saitama-u.ac.jp/~terada/web_archive/cosmic.riken.jp/suzaku/help/guide/fstep_web/node11.html
```bash=
#@xspec
XSPEC12>iplot
PLT>we (filename)
```
とすると(filename).qdp (filename).pcoというファイルができる。このうちqdpの方には多分スペクトルの情報が数値で書いてある。weの前にdummyrspを実行しておけばqdpにはモデルのみのスペクトルが数値で出力される。
これを使えば、xspec上でfluxコマンドを実行して標準出力をとって...としなくてもよくなる。
```bash=
#@シェル
~/..../fit> qdp (filename).qdp
To produce plot, please enter
PGPLOT file/type: /xs
PLT> pl
```
とすると図が出力される。
# background 解析
cstatで解析を行うとき、sourceからbkgを引く通常の解析だとポワソン統計でなくなる。
ので特別なbackground解析をしなければならない。
https://heasarc.gsfc.nasa.gov/xanadu/xspec/manual/node27.html
※rsl_nxbgenのbackscaleを1にしないと、backgroundが引かれるどころか足されてしまうしまう。
# iplot
iplotコマンドをできるだけまとめる
以下全てxspecで```iplot```と打った後に実行する
```bash=
r x AAA BBB #x軸をAAA keV~BBB keVに変更
r y AAA BBB #y軸をAAA keV~BBB keVに変更
co AAA on BBB #BBB番の線をAAA色に変更 AAA=1:黒, 2:赤, 3:緑, 4:青,...
#example
co 2 on 2 #2番の線を赤色に変更
co 3 on 2..5 #2~5番の色を緑に変更
co 1 on 1..100 #1~100番の色(=全ての線の色)を黒に変更
co off AAA # AAA番の線を表示しない
log x on #X軸をlog表示、xをyにするとy軸
log x off #X軸を線形表示、xをyにするとy軸
wi 2 #window2に移動、pl ld delの際はmodel+dataがwindoiw1, residualがwindow2
ls AAA on BBB #BBB番の線の表示をAAAに変更
#example
ls 1 on 1..100 #全ての線の表示を実線に変更
lw AAA on BBB #BBB番の線の太さをAAAに変更
#example
lw 5 on 1 #1番の線の太さを５(若干太め)に変更
```
# memo
::: info
fits event file
HDU 0 (Primary HDU) → ヘッダーのみ（データ無し）
HDU 1 (Binary Table) → イベントリスト（光子ごとの位置・エネルギー・時間）
HDU 2 (Image) → 2D イメージ（スカイ座標での画像）
HDU 3 (Good Time Intervals) → 観測が有効だった時間の表
:::
# phabsとcabs
phabsは束縛電子を仮定していて(でないと吸収が起きない)、
のexpの中がΣ各元素の電子数＊吸収の断面積、みたいな形になっている
なのでアバンダンス依存する
一方cabsは自由電子によるトムソン散乱、コンプトン散乱を仮定しているので、アバンダンスには依存しない
# fakeit
https://hackmd.io/@tenoto/H1Q70UrGP
fakeit
# binnnng
```bash=
grppha input output
chkey respfile ***
chkey ancrfile ***
chkey backfile ***
```
Sign in with Wallet
Connect another wallet