Advice on vectorizing block-wise operations in Numpy
By using apply_along_axis
, you can get rid of cauchy_schwartz
. Since you are not overly concerned with the pre-processing time, assume you have obtained the array params
which contains the flattened patches
params = np.random.rand(3,2,100)
as you can see the shape of params
is (3,2,100)
, the three numbers 3, 2, and 100 are just randomly chosen to create an auxiliary array to demonstrate the logic of using apply_along_axis
. 3 corresponds to the number of patches you have (determined by the patch shape and the image size), 2 corresponds to the two images, and 100 corresponds to the flattened patches. Therefore, the axes of params
is (idx of patches, idx of images, idx of entries of a flattened patch)
, this exactly matches the list params
created by your code
params = []
for i in range(0,patch1.shape[0],1):
for j in range(0,patch1.shape[1],1):
window1 = np.copy(imga[i:i+N,j:j+N]).flatten()
window2 = np.copy(imgb[i:i+N,j:j+N]).flatten()
params.append((window1, window2))
With the auxiliary array params
, here is my solution:
hist = np.apply_along_axis(lambda x: np.histogram(x,bins=11)[0],2,params)
hist = hist / np.sum(hist,axis=2)[...,None]
n_d = np.sum(np.product(hist,axis=1),axis=1)
d_d = np.sum(np.product(np.power(hist,2),axis=1),axis=1)
res = -1.0 * np.log10(n_d, d_d)