1

昨日から答えを探していましたが、運がありません。そのため、各波長でのフラックス値を含む 1D スペクトル (.fits) ファイルがあります。それらを2D配列(x、y)=(波長、光束)に変換し、割り当てられた波長(x)で光束(y)を返すプログラムを書きたいと考えています。私はこれを試しました：

#modules
import scipy
import numpy as np
import pyfits as pf

#Target Global Vaiables
hdulist_tg = pf.open('cutmask1-2.0001.fits')    
hdr_tg = hdulist_tg[0].header
flux_tg = hdulist_tg[0].data
crval_tg = hdr_tg['CRVAL1']             #Starting wavelength 
cdel_tg = hdr_tg['CDELT1']              #Wavelength axis width 
wave_tg = crval_tg + np.arange(3183)*cdel_tg        #Create an x-axis
wavelist = [6207,6315,6369,6438,6490,6565,6588]

wave_flux=[]
diff = 10
for wave in wave_tg:
    for flux in flux_tg:
        wave_flux.append((wave,flux))           

for item in wave_flux:
    wave = item[0]
    flux = item[1]

#Where I got my actual wavelength that exists in wave_tg
    diffmatch = np.abs(wave - wavelist[0])
    if diffmatch < diff:
        flux_wave = flux  
        diff = diffmatch 
        wavematch = wave

print wavelist[0],flux_wave,wavematch

しかし、波長が異なっていても、プログラムは常に同じフラックス値を返します。助けてください...

4

2 に答える 2

1

私はあなたのコードにいくつかの変更を加えました:

  • numpy ot createwave_fluxを using としてndarray使用しnp.hstack()、np.repeat()np.tile()

  • ファンシーインデックスを使用して検索に一致する値を取得する

結果のコードは次のとおりです。

#modules
import scipy
import numpy as np
import pyfits as pf

#Target Global Vaiables
hdulist_tg = pf.open('cutmask1-2.0001.fits')
hdr_tg = hdulist_tg[0].header
flux_tg = hdulist_tg[0].data
crval_tg = hdr_tg['CRVAL1']             #Starting wavelength
cdel_tg = hdr_tg['CDELT1']              #Wavelength axis width
wave_tg = crval_tg + np.arange(3183)*cdel_tg        #Create an x-axis
wavelist = [6207,6315,6369,6438,6490,6565,6588]

wave_flux = np.vstack(( np.repeat(wave_tg, len(flux_tg)),
                        np.tile(flux_tg, len(wave_tg)) )).transpose()

wave_ref = wavelist[0]
diff = 10

print wave_flux[ np.abs(wave_flux[:,0]-wave_ref) < diff ]

wave_flux列の波動値と列0のフラックス値を持つサブグループを返します1：

[[ 6197.10300138   500.21020508]
 [ 6197.10300138   523.24102783]
 [ 6197.10300138   510.6390686 ]
 ...,
 [ 6216.68436446   674.94732666]
 [ 6216.68436446   684.74255371]
 [ 6216.68436446   712.20098877]]
于 2013-07-02T07:30:58.380 に答える
1

2 次元テーブルの作成を完全にスキップして、次のものを使用しますinterp。

fluxvalues = np.interp(wavelist, wave_tg, flux_tg)

投稿したファイルについては、wave_tg 配列の長さがハードコーディングされているため、投稿したコードが機能しません。したがって、むしろ使用することをお勧めします

wave_tg = crval_tg + np.arange(len(flux_tg))*cdel_tg

また、何らかの理由で、あなたが投稿したファイルは実際にあなたが調べている波長まで上がっていないようです。対応する波長を正しく計算していることを確認するか、正しい波長を検索していることを確認する必要がある場合があります。

于 2013-07-02T10:19:45.037 に答える