def Cal_Correlogram(data):
global START_POINT
global TIME_FRAME
global NUM_CHANNEL
global MAX_DELAY
global NEW_WINDOW_HALF
numTimeFrame=int(math.ceil(float(data.shape[1]-START_POINT)/float(TIME_FRAME)))
corr=np.zeros((numTimeFrame,NUM_CHANNEL,MAX_DELAY),dtype=np.float64)
for chan in xrange(NUM_CHANNEL):#range(NUM_CHANNEL):
#blip()
for timeFrame in xrange(numTimeFrame):
tmp_corr1=0
for timeStep in xrange(-NEW_WINDOW_HALF,NEW_WINDOW_HALF):
tmp_corr1=tmp_corr1+data[chan,START_POINT+timeFrame*TIME_FRAME+timeStep]**2
tmp_corr1=math.sqrt(tmp_corr1)
for delay in xrange(MAX_DELAY):
tmp_corr=0.
tmp_corr2=0.
for timeStep in xrange(-NEW_WINDOW_HALF,NEW_WINDOW_HALF):
tmp_corr +=data[chan,START_POINT+timeFrame*TIME_FRAME+timeStep]*data[chan,START_POINT+timeFrame*TIME_FRAME-delay+timeStep]
tmp_corr2 +=data[chan,START_POINT+timeFrame*TIME_FRAME-delay+timeStep]**2
if(math_util.REALZEROP(tmp_corr1) or math_util.REALZEROP(tmp_corr2)):
corr[timeFrame,chan,delay]=0.0
else:
corr[timeFrame,chan,delay] = tmp_corr/math.sqrt(tmp_corr2)/tmp_corr1;
return corr
通过编写如下的build.py 文件
from distutils.core import setup
from distutils.extension import Extension
from Cython.Distutils import build_ext
setup(cmdclass = {'build_ext': build_ext},ext_modules = [Extension("Cal_Correlogram", ["corr_c.pyx"])])
然后运行python setup.py build_ext --inplace
得到Cal_Correlogram.pyd 和libCal_Correlogram.a等文件。运行测试。
发现其实速度并不能提升多少,跑了一次700*100图像的自相关,很慢2-3分钟(我的笔记本)。
那么像使用了numpy等库的程序如何提升速度,后面发现原来其实上面的代码主要是慢在for循环,不知道怎么python的for循环这么慢。因此把矩阵单位元素的计算,用numpy的内部函数代替,换成矩阵与矩阵或行列之间的运算,这样就极大程度的减少了for循环。代码修改如下:
import numpy as np
from common import *
import math_util
import math
def Cal_Correlogram(data):
global START_POINT
global TIME_FRAME
global NUM_CHANNEL
global MAX_DELAY
global NEW_WINDOW_HALF
numTimeFrame=int(math.ceil(float(data.shape[1]-START_POINT)/float(TIME_FRAME)))
corr=np.zeros((numTimeFrame,NUM_CHANNEL,MAX_DELAY),dtype=np.float64)
cdef timeFrame=0
for timeFrame in xrange(numTimeFrame):
tmp_corr1=0.
tmp_corr1+=np.sum(data[:,START_POINT+timeFrame*TIME_FRAME-NEW_WINDOW_HALF:START_POINT+timeFrame*TIME_FRAME+NEW_WINDOW_HALF]**2,axis=1)
tmp_corr1=np.sqrt(tmp_corr1)
for delay in xrange(MAX_DELAY):
tmp_corr=0.
tmp_corr2=0.
tmp_corr +=np.sum(data[:,START_POINT+timeFrame*TIME_FRAME-NEW_WINDOW_HALF:START_POINT+timeFrame*TIME_FRAME+NEW_WINDOW_HALF]*data[:,START_POINT+timeFrame*TIME_FRAME-delay-NEW_WINDOW_HALF:START_POINT+timeFrame*TIME_FRAME-delay+NEW_WINDOW_HALF],axis=1)
tmp_corr2 +=np.sum(data[:,START_POINT+timeFrame*TIME_FRAME-delay-NEW_WINDOW_HALF:START_POINT+timeFrame*TIME_FRAME-delay+NEW_WINDOW_HALF]**2,axis=1)
corr[timeFrame,:,delay] = tmp_corr/np.sqrt(tmp_corr2)/tmp_corr1;
return corr
速度提升不少,通过time()计算,运行一次100*700图像的自相关,用时21秒。