その輝度の見積もりがむずかしくbackground estimationで苦労することもしばしば。
すざくの良質なスペクトルを元にH. Uchiyamaさんがその表面輝度をモデル化しています。
H. Uchiyama PhD (2010) PDFはこちら
Sgr A*を中心とする２つの2-D exponentialの和を用いたモデル (eq.4.2) は特に便利なので、
銀河座標を放り込むと表面輝度を表示するかんたんなスクリプトを書きました。
座標を２セットいれると２領域間での輝度の比を表示します。
無駄にnumpyつかってるのはご愛嬌。。
※初期版にひどいバグ:
l0, b0で絶対値とるの忘れてました！
12/5に参照した方は修正した下記をご覧下さい。
※初期版にひどいバグ:
l0, b0で絶対値とるの忘れてました！
12/5に参照した方は修正した下記をご覧下さい。
#!/usr/bin/env python
'''
The model expression and the parameters are respectively
taken from eq.4.2 and table 4.11 of H. Uchiyama 2010, PhD Thesis
http://repository.tksc.jaxa.jp/pl/dr/IS8000028000/en
The parameters were measured with the K-alpha lines from FeXXV
and FeXXVI and the 5-10 keV continuum bands.
The unit of surface brightness is 1e-7 ph/s/cm2/am2 for
Fe Ka lines and 1e-13 erg/s/cm2/am2 for 5-10keV continuum.
'''
import numpy as np
gpos_SgrAs = np.array([-0.056, -0.046])
a1 = (26.3, 10.0, 3.4)
a2 = (1.6, 0.25, 0.20)
hl1 = (0.47, 0.51, 0.48)
hb1 = (0.22, 0.22, 0.19)
hl2 = (30, 47, 55)
hb2 = (1.7, 17, 2.5)
band2colnum = {
'FeXXV':0,
'FeXXVI':1,
'hard':2
}
def calc_sb(l, b, band):
i = band2colnum[band]
l0, b0 = np.abs(np.array([l, b]) - gpos_SgrAs)
sb1 = a1[i] * np.exp(-l0/hl1[i]) * np.exp(-b0/hb1[i])
sb2 = a2[i] * np.exp(-l0/hl2[i]) * np.exp(-b0/hb2[i])
sb = sb1 + sb2
return(sb)
def calc_ratio(l1,b1, l2,b2, band):
sb1 = calc_sb(l1, b1, band)
sb2 = calc_sb(l2, b2, band)
ratio = sb1 / sb2
return(ratio)
if __name__ == '__main__':
import sys
argvs = sys.argv
try:
l, b, band = argvs[1:4]
cnum = band2colnum[band]
sb = calc_sb(float(l), float(b), band)
print(sb)
except:
try:
l1, b1, l2, b2, band = argvs[1:6]
cnum = band2colnum[band]
rt = calc_ratio(float(l1), float(b1), float(l2), float(b2), band)
print(rt)
except:
print('Usage1 calc sb : calcgdxe.py l b band')
print('Usage2 calc sb ratio : calcgdxe.py l1 b1 l2 b2 band')
print('band should be selected from FeXXV, FeXXVI, or hard')
exit()
追記。
標準で入ってないないnumpyというライブラリを使っているので、
import numpyでこける場合は以下の３点を修正すれば動くと思います。
めんどくさいので動作みかくにん。
1)
import numpy as np
gpos_SgrAs = np.array([-0.056, -0.046])
↓
import math
l_SgrAs = -0.056
b_SgrAs = -0.046
2)
l0, b0 = np.abs(np.array([l, b]) - gpos_SgrAs)
↓
l0 = abs(l - l_SgrAs)
b0 = abs(b - b_SgrAs)
3)
np.exp()
↓
math.exp()
