The other day, I asked this question [3x3x3 Laplace Kernel?](https://mathoverflow.net/questions/31594/3x3x3-laplace-kernel), regarding what the 3x3x3 kernel was for applying a Laplacian convolution. 


On that page, it mentions the kernels were "deduced by using discrete differential quotients." Does anyone know how this is done? For example, what if I need a 5x5x5 kernel? How do I "solve" for it? I tried looking in Matlab but didn't find anything. 


The reason I want to do this is because if I have a 100x100x10 pixel image I want to apply a Laplacian convolution to, I can convolve with my 3x3x3 kernel, OR I can apply a Fourier transform to both, and then use point-wise multiplication to get the same solution, then use the inverse Fourier to get my result back, which should match the result I got using the direct convolution method. As I mentioned, to do the Fourier approach, I need two images - my original, and the "Laplacian" one. I figured this Laplacian one needs to be generated in the same way as the 3x3x3 one my previous question identified, but needs to also be the same dimensions as my real image (100x100x10 in this example) in order for the Fourier transform approach to work. So my question above is - how do I 'fill in' the 3-dimensional matrix that is my discrete Laplacian operator?


Sorry if this question doesn't make sense. In a very dumbed down sense, I just want to know how to solve for the 3x3x3 kernel myself, so that I can solve for 'larger kernels' like a 100x100x10 one.



How to "fill in" 3-dimensional Laplacian kernels
If you want the "5x5x5" kernel, then it is no longer the Laplacian in the usual sense. 


Let me quickly describe where that kernel you saw on the Wikipedia page comes from. The Laplacian, as a differential operator, is $\sum\_i (\partial\_i)^2$. Now if we discretize the space into a grid, we can approximate the partial derivative in the $x$ direction by the finite difference of the function at a point $f(x\_0,y\_0,z\_0)$ and $f(x\_0 + 1, y\_0,z\_0)$. Here we use the usual formula
$$ \frac{\partial}{\partial x}f(x\_0,y\_0,z\_0) = \lim\_{h\to 0} \frac{f(x\_0 + h, y\_0, z\_0) - f(x\_0,y\_0,z\_0)}{h} $$
And since that the space has discrete steps, the smallest non-vanishing interval you can take is of size step 1. So morally speaking $f(x\_0 + 1 , y\_0,z\_0) - f(x\_0,y\_0,z\_0)$ is the derivative in the $+x$ direction of the function $f$ defined on the grid. 


The second derivative can be expressed as
$$ \frac{\partial^2}{\partial x^2} f(x\_0,y\_0,z\_0) = \lim\_{h\to 0} \frac{\partial\_xf(x\_0 + h, y\_0, z\_0) - \partial\_xf(x\_0,y\_0,z\_0)}{h}$$
$$ = \lim\_{h\to 0} \frac{ f(x\_0 + h, y\_0, z\_0) - f(x\_0,y\_0,z\_0) - f(x\_0,y\_0,z\_0) + f(x\_0 - h, y\_0, z\_0)}{h^2} $$
So the discrete version is simply the difference between the discrete derivative in the $+x$ direction and that in the $-x$ direction. (Notice that since we discretize, the forward and backward derivatives are not necessary the same.) In other words, the discrete second derivative is 
$$ f(x\_0 + 1, y\_0, z\_0) + f(x\_0-1,y\_0,z\_0) - 2 f(x\_0,y\_0,z\_0) $$


Now you just need to add it over the three spatial directions. 


So naturally the discrete Laplace kernel has support with in a cube that is 3 on a side. So what you want to do, I think, is probably just trivially extend the filter: make the filter 100 on a side and put the 3x3x3 kernel smack in the center (or as close to it as you can), and pad the rest with zero. Then you can take the Fourier transform and multiply against the Fourier transform of your original image. Note that I suspect whatever software solution you have will take the Fourier transform by assuming the data is periodic, so if you do this the edges will wrap around and bleed into the opposite sides a little. So the two results probably will not match on a two pixels wide band around the boundary. But they should match in the interior. 


