2D-массив Python - Как подключить х и получить значение у?

Question

Я искал ответ со вчерашнего дня, но не повезло. Поэтому у меня есть файл одномерного спектра (.fits) со значением потока на каждой длине волны. Я преобразовал их в двумерный массив (x, y) = (длина волны, поток) и хочу написать программу, которая будет возвращать поток (y) на некоторых назначенных длинах волн (x). Я попробовал это:

#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

но программа всегда возвращает одно и то же значение потока, даже если длина волны различна. Пожалуйста помоги...

Answer 1

Я бы вообще пропустил создание двумерной таблицы и просто использовал interp:

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

Для файла, который вы разместили, код, который вы опубликовали, не работает из-за жестко заданной длины массива wave_tg. Поэтому я бы рекомендовал вам использовать

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

Кроме того, по какой-то причине кажется, что файл, который вы разместили, на самом деле не соответствует длине волны, которую вы ищете. Возможно, вам придется проверить, правильно ли вы рассчитываете соответствующие длины волн, или проверить, что вы ищете правильные длины волн.

Answer 2

Я сделал некоторые изменения в вашем коде:

• используя NumPy не создать wave_flux как 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]]