我有一个巨大的文本文件，其中包含一百万颗恒星的位置(x、y、z(和速度分量(vx、vy、vz(。经过一些旋转和投影，我得到了新的恒星位置和速度分量(x'、y'、z'、vx'、vy'、vz'(。
我的最后一步是计算沿着视线的速度，这就像我必须"；平均值"；vz分量，为此，我尝试创建一个FITS文件，其中每个像素都包含vz分量的平均值。
这里有一部分我的代码：
mod = np.genfromtxt('data_bar_region.txt')
x = list(mod[:,0])
y = list(mod[:,1])
vz = mod[:,5]
x_rang_1 = np.arange(-40, 41, 1)
y_rang_1 = np.arange(-40, 41, 1)
fake_data_1 = np.array((len(x_rang_1),len(x_rang_1)))
for i in range(len(x_rang_1)-1):
for j in range(len(y_rang_1)-1):
vel_tmp = []
for index in range(len(x)):
if x_rang_1[i] <= x[index] <= x_rang_1[i+1]:
if y_rang_1[j] <= y[index] <= y_rang_1[j+1]:
vel_tmp.append(vz[index])
fake_data_1[j,i] = np.mean(vel_tmp)
hdu1 = fits.PrimaryHDU(fake_data_1)
hdu1.writeto('TEST.fits')
这个代码太慢了(在我的笔记本电脑上花了大约8个小时(，我不知道如何加快速度。
你有没有一些建议或其他方法可以更好更快地计算v_LOS？
编辑：在执行"；平均"；，我必须把图像分成几个形状和尺寸的部分(这样的部分被称为"仓"(。
这里有一张[垃圾箱的图像(在右侧面板上，有相同的垃圾箱图像，但它被放大了，以更好地证明什么是垃圾箱(]1。
因此，我有另一个FITS文件(称为bins.FITS(，其尺寸与fake_data_1相同，我只想找到这两个文件之间的对应关系，因为我想计算恒星在几个bin中分布的平均值和std。
或者，我有一个文本文件，其中包含哪个像素属于特定bin的信息，例如：
x y箱
1 1 34
12 34
2 3 34
34 56 37
34 57 37
34 58 37
bins.file的大小为(564585(，因此还有fake_data_1，改变了x和y范围的开始和停止的机会。我附上了整个脚本：
mod = np.genfromtxt('data_new_bar_scaled.txt')
# to match the correct position and size of the observation,
# I have to multiply by a factor equal to the semi-size
x = mod[:, 0]*(585-1)/200
y = mod[:, 1]*(564-1)/200
vz = mod[:,5]
A = fits.open('NGC4277_TESIkinematic.fits')
bins = A[7].data.T
start_x = -(585-1)/2
stop_x = (585-1)/2
step_x = step # step in x_rang_1
x_rang = np.arange(start_x, stop_x + step_x, step_x)
start_y = -(564-1)/2
stop_y = (564-1)/2
step_y = step # step in y_rang_1
y_rang = np.arange(start_y, stop_y + step_y, step_y)
fake_data_1 = np.empty((len(x_rang), len(y_rang)))
fake_data_1[:] = np.NaN # initialize with NaN
print(fake_data_1.shape)
print(bins.shape)
d = {}
for i in range(len(x)):
index_for_x = math.floor((x[i] - start_x) / step_x)
index_for_y = math.floor((y[i] - start_y) / step_y)
if 0 <= index_for_x < len(x_rang) and 0 <= index_for_y < len(y_rang):
key = (x_rang[index_for_x], y_rang[index_for_y])
if key in d:
d[key].append(vz[i])
else:
d[key] = [vz[i]]
bb = np.unique(bins)
print(len(bb))
for i, x in enumerate(x_rang):
for j, y in enumerate(y_rang):
key = (x, y)
for z in range(len(bb)):
j,k = np.where(bb[z]==bins)
print('index :', z)
if key in d:
fake_data_1[j,k] = np.mean(d[key])
您的代码非常慢，因为代码中的嵌套循环迭代了1600(80*80(次。您可以通过使用字典并仅对一百万颗恒星进行一次迭代来提高性能。
您可以尝试以下代码，其速度大约快1600倍：
import numpy as np
import math
mod = np.genfromtxt('data_bar_region.txt')
x = list(mod[:, 0])
y = list(mod[:, 1])
vz = mod[:, 5]
x_rang_1 = np.arange(-40, 41, 1)
y_rang_1 = np.arange(-40, 41, 1)
fake_data_1 = np.empty((len(x_rang_1), len(y_rang_1)))
fake_data_1[:] = np.NaN # initialize with NaN
d = {}
for i in range(len(x)):
key = (math.floor(x[i]), math.floor(y[i]))
if key in d:
d[key].append(vz[i])
else:
d[key] = [vz[i]]
for i, x in enumerate(x_rang_1):
for j, y in enumerate(y_rang_1):
key = (x, y)
if key in d:
fake_data_1[i, j] = np.mean(d[key])
hdu1 = fits.PrimaryHDU(fake_data_1)
hdu1.writeto('TEST.fits')
更新
对于x_rang_1(或y_rang_1(中step的通用版本，您可以尝试以下代码：
import numpy as np
import math
mod = np.genfromtxt('data_bar_region.txt')
x = list(mod[:, 0])
y = list(mod[:, 1])
vz = mod[:, 5]
start_x_rang_1 = -40
stop_x_rang_1 = 40
step_x_rang_1 = 0.5 # step in x_rang_1
x_rang_1 = np.arange(start_x_rang_1, stop_x_rang_1 + step_x_rang_1, step_x_rang_1)
start_y_rang_1 = -40
stop_y_rang_1 = 40
step_y_rang_1 = 1 # step in y_rang_1
y_rang_1 = np.arange(start_y_rang_1, stop_y_rang_1 + step_y_rang_1, step_y_rang_1)
fake_data_1 = np.empty((len(x_rang_1), len(y_rang_1)))
fake_data_1[:] = np.NaN # initialize with NaN
d = {}
for i in range(len(x)):
index_for_x_rang_1 = math.floor((x[i] - start_x_rang_1) / step_x_rang_1)
index_for_y_rang_1 = math.floor((y[i] - start_y_rang_1) / step_y_rang_1)
if 0 <= index_for_x_rang_1 < len(x_rang_1) and 0 <= index_for_y_rang_1 < len(y_rang_1):
key = (x_rang_1[index_for_x_rang_1], y_rang_1[index_for_y_rang_1])
if key in d:
d[key].append(vz[i])
else:
d[key] = [vz[i]]
for i, x in enumerate(x_rang_1):
for j, y in enumerate(y_rang_1):
key = (x, y)
if key in d:
fake_data_1[i, j] = np.mean(d[key])
hdu1 = fits.PrimaryHDU(fake_data_1)
hdu1.writeto('TEST.fits')
更新2
也许像下面这样？
当我认为输入是时
x y vz
0 0.1 10
1.8 0 4
1.2 1.9 5.2
bins = np.array(
[[34, 35, 34, 34, 36],
[37, 36, 34, 35, 36],
[34, 35, 37, 36, 34]]) # shape: (5, 3)
你想要以下代码吗？
import numpy as np
import math
x = np.array([0, 1.8, 1.2, ])
y = np.array([0.1, 0, 1.9, ])
vz = np.array([10, 4, 5.2])
start_x_rang_1 = 0
stop_x_rang_1 = 2
step_x_rang_1 = 1 # step in x_rang_1
x_rang_1 = np.arange(start_x_rang_1, stop_x_rang_1 + step_x_rang_1, step_x_rang_1)
start_y_rang_1 = 0
stop_y_rang_1 = 0.5
step_y_rang_1 = 2 # step in y_rang_1
y_rang_1 = np.arange(start_y_rang_1, stop_y_rang_1 + step_y_rang_1, step_y_rang_1)
fake_data_1 = np.empty((len(x_rang_1), len(y_rang_1))) # shape: (3, 5)
fake_data_1[:] = np.NaN # initialize with NaN
bins = np.array(
[[34, 35, 34, 34, 36],
[37, 36, 34, 35, 36],
[34, 35, 37, 36, 34]]) # shape: (3, 5)
d_bins = {}
for i in range(len(x)):
index_for_x_rang_1 = math.floor((x[i] - start_x_rang_1) / step_x_rang_1)
index_for_y_rang_1 = math.floor((y[i] - start_y_rang_1) / step_y_rang_1)
if 0 <= index_for_x_rang_1 < len(x_rang_1) and 0 <= index_for_y_rang_1 < len(y_rang_1):
key = bins[index_for_x_rang_1, index_for_y_rang_1]
if key in d_bins:
d_bins[key].append(vz[i])
else:
d_bins[key] = [vz[i]]
d_bins_mean = {}
for bin in d_bins:
d_bins_mean[bin] = np.mean(d_bins[bin])
get_corresponding_mean = np.vectorize(lambda x: d_bins_mean.get(x, np.NaN))
result = get_corresponding_mean(bins)
print(result)
它打印
[[10. nan 10. 10. nan]
[ 4.6 nan 10. nan nan]
[10. nan 4.6 nan 10. ]]