我正在尝试使用 memoryview 将以前的 python 代码转换为 cython,并且我不断从第一 bolo 行(靠近底部)得到编译器崩溃:
from __future__ import print_function
from builtins import range
from builtins import object
import cython
from cython.view cimport array as cvarray
from cpython cimport bool
from libc.math cimport round
from libc.stdlib cimport malloc, free
import numpy as np
cimport numpy as np
class LegCache(object):
def __init__(self):
self.d = {}
pass
def prep_legendre(self, n, polyorder):
p = (n, polyorder)
if p not in self.d:
self.d[p] = prep_legendre(n, polyorder)
return self.d[p]
@cython.boundscheck(False)
@cython.wraparound(False)
cdef prep_legendre(int n, int polyorder):
'''make array of legendre's'''
assert type(n) == int and type(polyorder) == int
cdef int[:,:] legendres = np.empty([n, polyorder + 1], dtype=int)
cdef int l = 0
legendres[:, 0] = np.ones(n)
if polyorder > 0:
legendres[:, 1] = np.linspace(-1, 1, n)
for i in range(polyorder - 1):
l = i + 1
np.multiply(l/(l+1), legendres[:,l-1])
cdef double[:,:] q = np.empty([polyorder + 1, polyorder + 1], dtype=double)
cdef double[:,:] r = np.empty([n, polyorder + 1], dtype=double)
cdef double[:,:] qt = np.empty([polyorder+1, polyorder+1], dtype=double)
cdef double[:,:] rinv = np.empty([polyorder + 1, n], dtype=double)
q, r = np.linalg.qr(legendres)
rinv = np.linalg.inv(r)
qt = q.T.copy()
return legendres, rinv, qt
def filter_slice_legendre_qr_mask_precalc(bolo,mask,legendres):
m=legendres.shape[1]
n=legendres.shape[0]
l2 = legendres*np.tile(mask.reshape(n,1),[1,m])
q,r=np.linalg.qr(l2)
rinv = np.linalg.inv(r)
p = np.dot(q.T,bolo)
coeff=np.dot(rinv,p)
out=bolo-np.dot(legendres,coeff)
return out,coeff
@cython.boundscheck(False)
@cython.wraparound(False)
cdef poly_filter_array(
double[:,:] array,
np.ndarray[DTYPE3_t, cast=True, ndim=2] mask_remove, # I think this casting should still work like this
np.ndarray[DTYPE3_t, cast=True, ndim=2] mask,
int[:] scan_list,
int ibegin,
int polyorder,
double minfrac=.75):
""" writes over input array
"""
cdef double nold = -1
# do nothing
if polyorder < 0:
return array
#damn, work
cdef int nch = array.shape[0]
cdef int nt = array.shape[1]
cdef int ns = len(scan_list)
cdef double[:,:,:] coeff_out = np.empty([nch, ns, nt], dtype = double)
legcache = LegCache()
cdef int istart = 0
cdef int n = 0
cdef int start = 0
cdef double mean = 0.0
cdef int[:,:] legendres = np.empty([n, polyorder + 1], dtype=int)
cdef double[:,:] qt = np.empty([polyorder+1, polyorder+1], dtype=double)
cdef double[:,:] rinv = np.empty([polyorder + 1, n], dtype=double)
#cdef double[:,:] bolo # I think you can get away without giving it a value. bolo changes in size throughout the loop
cdef int[:] goodhits = np.empty(np.shape(mask)[1], dtype = int)
# I'm not sure about the size of these
cdef double[:,:] coeff = np.empty([]) # I can't remember how dot product work right now but this should be easy to sort out
# remove mean
if polyorder == 0:
for s in range(len(scan_list)):
istart, n = scan_list[s]
start = istart - ibegin
for i in range(nch):
if np.any(mask[i, start:start + n]):
mean = np.average(
array[i, start:start + n], weights=mask[i, start:start + n])
array[i, start:start + n] -= mean
coeff_out[i, s, 0] = mean
# other cases
if polyorder > 0:
for s in range(len(scan_list)):
istart, n = scan_list[s]
start = istart - ibegin
if n <= polyorder: # otherwise cannot compute legendre polynomials
for i in range(nch):
mask[i, start:start + n] = 0 # flag it
# remove this region from actual data as well
mask_remove[i, start:start + n] = 0
print('Not enough points (%d) to build legendre of order (%d)' % (n, polyorder))
continue
goodhits = np.sum(mask[:, start:start + n], axis=1)
if n != nold:
legendres, rinv, qt = legcache.prep_legendre(n, polyorder)
rinvqt = np.dot(rinv, qt)
nold = n
# handle no masked ones
for i in range(nch):
if goodhits[i] != n:
continue # skip for now
bolo[i, :] = array[i, start:start + n] #where problem starts
coeff = np.dot(rinvqt, bolo)
coeff_out[i, s, :] = coeff
bolo -= np.dot(legendres, coeff)
for i in range(nch):
if goodhits[i] == n:
continue # skip since dealt with above
if goodhits[i] < minfrac * n: # not enough points
mask[i, start:start + n] = 0 # flag it
# remove this region from actual data as well
mask_remove[i, start:start + n] = 0
continue
bolo, coeff = filter_slice_legendre_qr_mask_precalc(
array[i, start:start + n], mask[i, start:start + n], legendres)
array[i, start:start + n] = bolo
coeff_out[i, s, :] = coeff
return coeff_out
当我尝试编译代码时,它会引发非特定错误“ExpandInplaceOperators 中的编译器崩溃”。我完全迷路了。