我在研究CCD阵列中的陷阱。目前我使用的是NumPy和Scipy,我已经能够将大多数的呼叫矢量化,这给了我一些加速。 目前,我的代码中的瓶颈是,我必须从代码内部循环中的大量不同插值中检索一个数字。这一步占用了97%的计算时间。在
我在这里举了一个简单的例子来说明我的问题:
import numpy as np
from scipy.interpolate import interp1d
# the CCD array containing values from 0-100
array = np.random.random(200)*100
# a number of traps at different positions in the CCD array
n_traps = 100
trap_positions = np.random.randint(0,200,n_traps)
# xvalues for the interpolations
xval = [0,10,100]
# each trap has y values corresponding to the x values
trap_yvals = [np.random.random(3)*100 for _ in range(n_traps)]
# The xval-to-yval interpolation is made for each trap
yval_interps = [interp1d(xval,yval) for yval in trap_yvals]
# moving the trap positions down over the array
for i in range(len(array)):
# calculating new trap position
new_trap_pos = trap_positions+i
# omitting traps that are outside array
trap_inside_array = new_trap_pos < len(array)
# finding the array_vals (corresponding to the xvalues in the interpolations)
array_vals = array[new_trap_pos[trap_inside_array]]
# retrieving the interpolated y-values (this is the bottleneck)
yvals = np.array([yval_interps[trap_inside_array[t]](array_vals[t])
for t in range(len(array_vals))])
# some more operations using yvals
有没有一种方法可以优化,也许使用Cython或类似的方法?在
我仔细考虑了一下,我想我找到了一个很好的解决方案,我想和大家分享,尽管这意味着我将回答我自己的问题。在
首先,我意识到,与其使用scipy.插值函数,我可以找到两个值之间的插值。这可以用这个小函数来完成
这给了我一些加速,但我想看看是否可以做得更好,所以我将函数移植到Cython并在所有陷阱上添加循环,这样我就不必在python代码中这样做了:
^{pr2}$我对不同的方法运行了一些计时(使用了一些比原问题中指出的更大的数组和更多的陷阱):
使用原始方法,即interp1d:(最佳值3)15.1秒
使用插值edunivariatespline(k=1)代替interp1d,如@M.T:(最佳值3)7.25秒
使用双值插值函数:(最佳值为3)1.34秒
使用Cython实现两个值插值(最佳值为3)0.113秒
相关问题 更多 >
编程相关推荐