The Problem
When applying high-pass filters to images, such as height maps, a common issue is edge artifacts. These are caused by having to extrapolate pixels and most common methods cannot mitigate these artifacts.
Even Photoshop suffers from the same issue as shown below:
Ideally, the high-pass filter should produce a uniform output. However, there is a gradient near the edge. This is caused by having to extrapolate pixels outside the image when convolving it, which is often achieved using replicate or reflection padding. Here is a plot of the pixel intensities for the same example:
The Solution
Extrapolating pixels based on a local linear regression can help to solve this issue. A naive implementation for 1 dimension would look like this:
1. Calculate the slope $b$ and intercept $a$ of the neighbourhood of each pixel weighted by the convolution kernel using linear regression
2. Apply the kernel to the image and extrapolate missing samples by $a + b*d_{ij}$. Here $d_{ij}$ is the distance from the current pixel to the pixel being extrapolated.
This means that pixels outside of the boundary are extrapolated using the slope of the local neighbourhood of each pixel. The only challenge remaining is to find an algorithm which can do this in a single pass, otherwise the memory accesses would need to be doubled or tripled, which would significantly affect the performance of a GPU implementation.
Finding the local slope requires calculating a weighted linear regression for the neighbourhood of each pixel and this is indeed possible in a single pass which can be combined with the convolution itself:
$ S_w = \sum{w_i} $
$ S_x = \sum{w_i \cdot x_i} $
$ S_y = \sum{w_i \cdot y_i} $
$ S_{xx} = \sum{w_i \cdot x_i^2} $
$ S_{xy} = \sum{w_i \cdot x_i \cdot y_i} $
Then we can calculate the slope $ b $ and intercept $ a $ using:
$ b = \frac{S_w \cdot S_{xy} – S_x \cdot S_y}{S_w \cdot S_{xx} – S_x^2} $
$ a = \frac{S_y – b \cdot S_x}{S_w} $
Results
With the extrapolation applied to the same linear gradient, the blurred version perfectly matches the original and the difference is zero:
Here is a more complex example with random noise added to a sine wave:
The same example without extrapolation yields the following result:
This improvement is particularly noticeable in high-pass filtered data but also makes a difference for low-pass filtering.
The Matlab implementation of this method in one dimension can be downloaded here (rename to .m, WordPress doesn’t like the Matlab extension and I’m too lazy to fix it):