Python-RADEX
项目描述
RADEX（www.sron.rug.nl/~vdtak/radex/）在Python中的包装器。
截至v0.2版本，创建于2013年10月26日，此包包括命令行程序的Python包装器和用f2py创建的直接包装Fortran代码。
f2py包装版本的安装过程
您需要在您的路径上安装gfortran和f2py。如果您已成功从源代码构建numpy，那么您应该都有。
您需要做的只是
$ python setup.py install_radex install这将调用一个名为install_radex的过程，该过程从radex主页下载RADEX的最新版本，修补源代码，并构建一个名为radex.so的文件，这是一个可以导入的Python共享对象。
有关更多详细信息，请参阅安装页面。
如果您想让pyradex在启动Python之前查找特定目录中的分子数据文件，您可以在启动Python之前指定环境变量RADEX_DATAPATH。也可以通过datapath关键字交互式指定。
使用f2py包装版本
Fortran代码的直接包装器使用一个名为Radex的类作为其基础结构。此类对于直接操作RADEX输入和直接访问其输出非常有用。
示例
import pyradex
import numpy as np
R = pyradex.Radex()
Tlvg = R(collider_densities={'oH2':900,'pH2':100}, column=1e16, species='co',method='lvg')
Tslab = R(collider_densities={'oH2':900,'pH2':100}, column=1e16, species='co',method='slab')
Tsphere = R(collider_densities={'oH2':900,'pH2':100}, column=1e16, species='co',method='sphere')
Tlvg[:3].pprint()
Tslab[:3].pprint()
Tsphere[:3].pprint()结果
Tex tau frequency upperstateenergy upperlevel lowerlevel upperlevelpop lowerlevelpop flux
------------- -------------- ----------- ---------------- ---------- ---------- ---------------- --------------- -----------------
15.2747101724 0.937692338925 115.2712018 5.53 2 1 0.273140336953 0.453621905471 2.93964536078e-14
10.8673211326 2.74275175782 230.538 16.6 3 2 0.0518618367484 0.273140336953 9.26125039465e-14
8.30670325364 2.01021823976 345.7959899 33.19 4 3 0.00379591658449 0.0518618367484 8.16324298598e-14
Tex tau frequency upperstateenergy upperlevel lowerlevel upperlevelpop lowerlevelpop flux
------------- -------------- ----------- ---------------- ---------- ---------- ---------------- -------------- -----------------
17.8076937528 0.681341951256 115.2712018 5.53 2 1 0.312979158313 0.394862780876 2.89304678735e-14
14.8865118666 1.96024230849 230.538 16.6 3 2 0.102821702575 0.312979158313 1.38012283784e-13
11.448407058 2.03949857132 345.7959899 33.19 4 3 0.00920322307626 0.102821702575 1.6139902821e-13
Tex tau frequency upperstateenergy upperlevel lowerlevel upperlevelpop lowerlevelpop flux
------------- ------------- ----------- ---------------- ---------- ---------- ---------------- -------------- -----------------
14.38256087 1.06765591906 115.2712018 5.53 2 1 0.243400727834 0.480559204909 2.93394133644e-14
9.28920337666 3.1666639484 230.538 16.6 3 2 0.037299201561 0.243400727834 7.24810556601e-14
7.50189023571 1.84556901411 345.7959899 33.19 4 3 0.00307839203073 0.037299201561 6.19215196139e-14
请注意，由于RADEX的编写方式，即使用公共块，因此这些对象中存储的值是相同的！您不能在任何时候拥有RADEX类的两个独立副本。
命令行版本的推荐安装过程
像往常一样使用make radex，但通过每次注释掉以下三行中的一行来构建两个可执行文件：radex_sphere、radex_lvg和radex_slab。
c parameter (method = 1) ! uniform sphere parameter (method = 2) ! expanding sphere (LVG) c parameter (method = 3) ! plane parallel slab (shock)将它们复制到您的系统路径
使用python setup.py install安装pyradex
简单示例
使用一些平凡的默认值
In [1]: import pyradex
In [2]: T = pyradex.radex(collider_densities={'H2':1000})
WARNING: Assumed thermal o/p ratio since only H2 was given but collider file has o- and p- H2 [pyradex.core]
In [3]: T.pprint(show_units=True)
J_up J_low E_UP FREQ WAVE T_EX TAU T_R POP_UP POP_LOW FLUX_Kkms FLUX_Inu
K GHz um K K K km / s erg / (cm2 s)
---- ----- ---- -------- --------- ----- --------- ------- ------ ------- --------- -------------
1 0 5.5 115.2712 2600.7576 5.044 0.0004447 0.00086 0.4709 0.47 0.0009155 1.806e-11
In [4]: T.meta
Out[4]:
{'Column density [cm-2]': '1.000E+12',
'Density of H2 [cm-3]': '1.000E+03',
'Density of oH2 [cm-3]': '3.509E-04',
'Density of pH2 [cm-3]': '1.000E+03',
'Geometry': 'Uniform sphere',
'Line width [km/s]': '1.000',
'Molecular data file': '/Users/adam/repos/Radex/data/co.dat',
'Radex version': '20nov08',
'T(background) [K]': '2.730',
'T(kin) [K]': '10.000'}
时间信息
即，有多快？
%timeit T = pyradex.radex(collider_densities={'H2':1000})
1 loops, best of 3: 149 ms per loop
for n in 10**np.arange(6):
%timeit T = pyradex.radex(collider_densities={'H2':n})
10 loops, best of 3: 149 ms per loop
10 loops, best of 3: 150 ms per loop
10 loops, best of 3: 149 ms per loop
10 loops, best of 3: 151 ms per loop
10 loops, best of 3: 150 ms per loop
10 loops, best of 3: 149 ms per loop
for n in 10**np.arange(12,18):
....: %timeit T = pyradex.radex(collider_densities={'H2':1000}, column_density=n)
10 loops, best of 3: 149 ms per loop
10 loops, best of 3: 149 ms per loop
10 loops, best of 3: 149 ms per loop
10 loops, best of 3: 150 ms per loop
10 loops, best of 3: 152 ms per loop
10 loops, best of 3: 157 ms per loop
这些结果表明，即使在需要更多迭代的非常光厚的案例中，执行时间主要由Python开销主导。
如果您重新进行这些测试，比较Fortran包装器与“原始”版本，差异将非常巨大。以下测试可以在timing.py中看到
Python: 0.892609834671 Fortran: 0.0151958465576 py/fortran: 58.7403822016 Python: 0.902825832367 Fortran: 0.0102920532227 py/fortran: 87.7206727205 Python: 0.876524925232 Fortran: 0.0730140209198 py/fortran: 12.0048850096 Python: 0.836034059525 Fortran: 0.0925290584564 py/fortran: 9.03536762906 Python: 0.880390882492 Fortran: 0.0725519657135 py/fortran: 12.1346248008 Python: 0.96048283577 Fortran: 0.0753719806671 py/fortran: 12.7432346512
制作网格
使用脚本更高效，但您仍然可以这样做...
for n in 10**np.arange(12,18):
T = pyradex.radex(collider_densities={'H2':1000}, column_density=n)
T.pprint()
Row# Line# E_UP FREQ WAVE T_EX TAU T_R POP_UP POP_LOW FLUX_Kkms FLUX_Inu
---- ----- ---- -------- --------- ----- --------- ------- ------ ------- --------- ---------
1 0 5.5 115.2712 2600.7576 5.044 0.0004447 0.00086 0.4709 0.47 0.0009155 1.806e-11
Row# Line# E_UP FREQ WAVE T_EX TAU T_R POP_UP POP_LOW FLUX_Kkms FLUX_Inu
---- ----- ---- -------- --------- ----- -------- -------- ------ ------- --------- ---------
1 0 5.5 115.2712 2600.7576 5.047 0.004444 0.008589 0.471 0.4698 0.009143 1.803e-10
Row# Line# E_UP FREQ WAVE T_EX TAU T_R POP_UP POP_LOW FLUX_Kkms FLUX_Inu
---- ----- ---- -------- --------- ----- ------- ------- ------ ------- --------- ---------
1 0 5.5 115.2712 2600.7576 5.075 0.04415 0.08473 0.4721 0.4681 0.0902 1.779e-09
Row# Line# E_UP FREQ WAVE T_EX TAU T_R POP_UP POP_LOW FLUX_Kkms FLUX_Inu
---- ----- ---- -------- --------- ----- ------ ------ ------ ------- --------- ---------
1 0 5.5 115.2712 2600.7576 5.336 0.4152 0.7475 0.4817 0.4527 0.7957 1.569e-08
Row# Line# E_UP FREQ WAVE T_EX TAU T_R POP_UP POP_LOW FLUX_Kkms FLUX_Inu
---- ----- ---- -------- --------- ----- ----- ---- ------ ------- --------- ---------
1 0 5.5 115.2712 2600.7576 6.929 2.927 3.49 0.5057 0.3745 3.715 7.327e-08
Row# Line# E_UP FREQ WAVE T_EX TAU T_R POP_UP POP_LOW FLUX_Kkms FLUX_Inu
---- ----- ---- -------- --------- ----- ----- ---- ------ ------- --------- ---------
1 0 5.5 115.2712 2600.7576 9.294 18.09 5.96 0.4696 0.2839 6.345 1.252e-07
如果您想使用直接包装的版本创建网格，请使用恒定温度的循环：每次加载新的温度时，RADEX都必须读取分子数据文件并跨碰撞率值进行插值，这可能是一个相当大的开销。
如果您想构建网格，请不要每次都创建一个astropy表！这似乎在每个迭代中主导了开销。
关于LVG计算的自我一致性说明
LVG计算的单位很奇怪。线的吸收率仅取决于视线方向上的速度相干柱，即每千米每秒的柱。
LVG Sobolev近似的关键假设是每个“单元”可以独立处理，这样就没有非局部辐射效应。
这种独立性意味着局部体积密度和总视线柱密度之间有一个分离。
然而，典型代码（如RADEX、DESPOTIC）报告的量是积分的视线值。柱密度、丰度和局部体积密度并不独立。
为了有一个自洽的云（或视线），您必须假设某些长度尺度。通常，人们指定每长度尺度的速度梯度而不是绝对长度尺度，但长度尺度很重要。
如果指定了氢的总柱密度 N(H) 和密度 n(H)，则长度尺度很简单：N(H)/n(H) = L。如果您增加密度，这个长度尺度会减小——到目前为止一切都很好。
在RADEX中，标准自由变量是感兴趣分子的柱。如果您更改分子的柱，这是可能的，并保持RADEX中所有其他内容（n(H)、dV）固定，则这种变化可以解释为尺寸尺度或柱的变化。
可以考虑将长度尺度作为自由参数的替代可能性，但这种方法包含改变涉及过程解释的危险：如果长度尺度减小，而固定delta-V，则速度梯度 dv/dl 必须更大。应避免这种解释，因为它承担着打破LVG假设的风险。速度梯度也经常通过观察到的线宽作为强制约束，而长度尺度在大多数情况下只受到弱约束。
在DESPOTIC中，自由变量是总柱密度、密度、丰度和速度梯度。因此，长度被留作依赖变量，与上述内容一致。
类（Despotic & Radex）被构造得长度是依赖变量，而所有其他变量都可以更改。由于丰度不是RADEX的显式输入，这是通过一些后台属性机制实现的。
