We have the following equation (weighted phase Laplacian):
\$c_{i,j} = U(i,j)\Delta^{^x}_{i,j} - U(i-l, j)\Delta^{x}_{i-1,j} + V(i,j)\Delta^{y}_{i,j} - V(i,j-1)\Delta^{y}_{i,j-1}\$
Where
\$U(i,j) = min(w^2_{i+1}, w^2_{i+1})V(i,j) = min(w^2_{j+1}, w^2_{j+1})\$
with \$U, V, \Delta \$ (N,N)
matrices
I coded the following code:
def compute_Weights(w):
U = np.zeros_like(w)
V = np.zeros_like(w)
U = np.minimum(w[1:, :]**2, w[:-1, :]**2)
V = np.minimum(w[:, 1:]**2, w[:, :-1]**2)
"""add the last ones rows and cols
"""
U = np.vstack((U, U[-1]))
V = np.vstack((V.T, V.T[-1])).transpose()
return U, V
def Gradient(p1, p2):
r = p1 - p2
if r > np.pi:
return r - 2 * np.pi
if r < -np.pi:
return r + 2 * np.pi
return r
"""
phase is not flattened
wts are flattened
"""
def ComputePhaseLaplacian(phase, U, V):
rows, cols = phase.shape
phi = phase.flatten()
rho = np.zeros_like(phi)
for j in range(rows):
for i in range(cols):
k = j * rows + i
k1 = k + 1 if i < cols - 1 else k - 1
k2 = k - 1 if i > 0 else k + 1
k3 = k + cols if j < rows - 1 else k - cols
k4 = k - cols if j > 0 else k + cols
w1 = U[k]
w2 = U[k-1] if i > 0 else U[k]
w3 = V[k]
w4 = V[k-cols] if j > 0 else V[k]
rho[k] = w1 * Gradient(phi[k], phi[k1]) \
+ w2 * Gradient(phi[k], phi[k2]) \
+ w3 * Gradient(phi[k], phi[k3]) \
+ w4 * Gradient(phi[k], phi[k4])
return rho
I want to know if there is a "pythonic way" to write this. I'm a C/C++ programmer so this was my perspective. Maybe something to avoid to work with flattened arrays, and to work with 2D arrays.
I.e:
phase = np.load('phase.npy') #phase 512x512 float array
U,V = ComputeWeigths(W) #W 512x512 float array with values in (0,1)
uphase = ComputePhaseLaplacian(phase, U, V) #uphase 512x512 float array