Python中的离散拉普拉斯算子（del2等价）

我需要Matlab（Octave）的Python/Numpy等价离散拉普拉斯算子（函数）del2（）。我尝试了几个Python解决方案，其中没有一个与del2的输出相匹配。在八音我Python中的离散拉普拉斯算子（del2等价）

image = [3 4 6 7; 8 9 10 11; 12 13 14 15;16 17 18 19] 
del2(image)

这给

0.25000 -0.25000 -0.25000 -0.75000 
    -0.25000 -0.25000 0.00000 0.00000 
    0.00000 0.00000 0.00000 0.00000 
    0.25000 0.25000 0.00000 0.00000

Python的我试过

import numpy as np 
from scipy import ndimage 
import scipy.ndimage.filters 

image = np.array([[3, 4, 6, 7],[8, 9, 10, 11],[12, 13, 14, 15],[16, 17, 18, 19]]) 
stencil = np.array([[0, 1, 0],[1, -4, 1], [0, 1, 0]]) 
print ndimage.convolve(image, stencil, mode='wrap')

其给出结果

[[ 23 19 15 11] 
[ 3 -1 0 -4] 
[ 4 0 0 -4] 
[-13 -17 -16 -20]]

我也试过

结果

scipy.ndimage.filters.laplace(image)

这给出结果

[[ 6 6 3 3] 
[ 0 -1 0 -1] 
[ 1 0 0 -1] 
[-3 -4 -4 -5]]

所以任何输出的似乎与海誓山盟。 Octave代码del2.m表明它是拉普拉斯算子。我错过了什么吗？

来源

2011-01-14 user423805

在内部，操作员都具有相同的实现（matlab显然除以4，其中Python不会）。在边界上，通过向`laplace（）`提供'mode ='wrap'`，你可以使两个Python版本相同。但通过查看Matlab结果，我不知道Matlab在边界上做了什么。 – 2011-01-15 00:32:27

实际上它在边缘做立方外推： http://www.mathworks.it/it/help/matlab/ref/del2.html 所以，如果你尝试使用laplace（）的最后一个例子，也可以在边界上获得正确的结果。 – astrojuanlu 2013-04-08 17:33:41

也许您在寻找scipy.ndimage.filters.laplace()。基于

来源

2011-01-14 14:44:30

的代码在这里

http://cns.bu.edu/~tanc/pub/matlab_octave_compliance/datafun/del2.m

我试图写一个Python等同。它似乎工作，任何反馈将不胜感激。

import numpy as np 

def del2(M): 
    dx = 1 
    dy = 1 
    rows, cols = M.shape 
    dx = dx * np.ones ((1, cols - 1)) 
    dy = dy * np.ones ((rows-1, 1)) 

    mr, mc = M.shape 
    D = np.zeros ((mr, mc)) 

    if (mr >= 3): 
     ## x direction 
     ## left and right boundary 
     D[:, 0] = (M[:, 0] - 2 * M[:, 1] + M[:, 2])/(dx[:,0] * dx[:,1]) 
     D[:, mc-1] = (M[:, mc - 3] - 2 * M[:, mc - 2] + M[:, mc-1]) \ 
      /(dx[:,mc - 3] * dx[:,mc - 2]) 

     ## interior points 
     tmp1 = D[:, 1:mc - 1] 
     tmp2 = (M[:, 2:mc] - 2 * M[:, 1:mc - 1] + M[:, 0:mc - 2]) 
     tmp3 = np.kron (dx[:,0:mc -2] * dx[:,1:mc - 1], np.ones ((mr, 1))) 
     D[:, 1:mc - 1] = tmp1 + tmp2/tmp3 

    if (mr >= 3): 
     ## y direction 
     ## top and bottom boundary 
     D[0, :] = D[0,:] + \ 
      (M[0, :] - 2 * M[1, :] + M[2, :])/(dy[0,:] * dy[1,:]) 

     D[mr-1, :] = D[mr-1, :] \ 
      + (M[mr-3,:] - 2 * M[mr-2, :] + M[mr-1, :]) \ 
      /(dy[mr-3,:] * dx[:,mr-2]) 

     ## interior points 
     tmp1 = D[1:mr-1, :] 
     tmp2 = (M[2:mr, :] - 2 * M[1:mr - 1, :] + M[0:mr-2, :]) 
     tmp3 = np.kron (dy[0:mr-2,:] * dy[1:mr-1,:], np.ones ((1, mc))) 
     D[1:mr-1, :] = tmp1 + tmp2/tmp3 

    return D/4

来源

2011-01-15 14:05:59 user423805

可以使用卷积通过卷积用适当stencil阵列来计算拉普拉斯算子：

from scipy.ndimage import convolve 
stencil= (1.0/(12.0*dL*dL))*np.array(
     [[0,0,-1,0,0], 
     [0,0,16,0,0], 
     [-1,16,-60,16,-1], 
     [0,0,16,0,0], 
     [0,0,-1,0,0]]) 
convolve(e2, stencil, mode='wrap')

来源

2011-01-18 14:25:11 2daaa

Python中的离散拉普拉斯算子（del2等价）

回答

相关问题