I am working on simulating traps in CCD arrays. Currently I am using NumPy and Scipy, and I have been able to vectorize most of the calls which have given me some speed-up. At the moment the bottleneck in my code is that I have to retrieve a number from a large number of different interpolations in the inner loop of my code. This particular step takes up ~97% of the computing time.
I have made a simple example of my problem here:
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 Is there a way this can be optimized, maybe using Cython or similar?
