# Making a low pass filter for irregular samples

We are interested in reconstructing the red channel of a scene.

The red channel is slowly varying.

Unfortunately the camera does not have a Bayer sensor, and the red samples are irregularly spaced over the grid:

How do you design a 2D low-pass filter (convolution kernel) that:

• only depends on the samples of red that you have
• has as close to constant frequency response as possible?

I was thinking of something along the lines of solving an optimization problem relating to the NDFT. Would this idea work? :

Take the NDFT of the red components we have. Then in continuous space, form the sum of all the sinusoids from the frequency information in the NDFT. Re-sample the function along the whole grid.

• Have you considered interpolation schemes such as Kriging? You can weight the estimates with their distance from the known pixel values. – A_A Jul 8 '16 at 14:03

I realize that this is not a filter, but I do not believe any filter will give you a reasonable result given the strange sample locations. This seems to be borne out by several references on the web. This one for example and this PDF one.

I suggest something something along these lines, using projections onto convex sets (POCS).

Algorithm Descroption

Assume that the known red pixel values $R_i$ are at $(x_i,y_i)$ so that your start red image is: $$I_{\tt red}(x,y) = \sum_{i=1}^N R_i \delta(x-x_i, y-y_i)$$

Set $k=0$, and initialize your estimate of the fully sampled red image with $I_{\tt red}(x,y)$: $$\hat{I}^k_{\tt red}(x,y) = I_{\tt red}(x,y).$$

Take the 2D FFT of $\hat{I}^k_{\tt red}$: $$F_{\tt red}(n,m) = {\tt FFT}\{\hat{I}^k_{\tt red}(x,y)\}$$ Because the red channel is slowly varying, remove all frequencies above a cutoff: $$F^{lp}_{\tt red}(n,m) = F_{\tt red}(n,m) L(n,m)$$ where $L$ is a low pass mask.

Form the inverse FFT to get $\hat{I}^{k+1}_{\tt red}(x,y)$ $$\hat{I}^{k+1}_{\tt red}(x,y) = {\tt FFT}^{-1}\{F^{lp}_{\tt red}(n,m)\}$$

Set the known pixels of $\hat{I}^{k+1}_{\tt red}(x,y)$ from $I_{\tt red}(x,y)$.

If the difference between $\hat{I}^{k+1}_{\tt red}(x,y)$ and $\hat{I}^{k}_{\tt red}(x,y)$ is small enough, stop, otherwise reiterate with $k=k+1$.

Results

An example run with just a $6 \times 6$ matrix is below. The left image shows the original pixels. The right image shows the POCS version. The red dots show the original pixel locations in both images.

Different Initial Conditions

I tried setting the initial image to

Ired <- array(100000,c(ImageSize, ImageSize))


and

Ired <- array(128,c(ImageSize, ImageSize))


and the results are below, so I do not believe the initial values for the unknown pixels matter.

The plot scales the greyscale values, which is why the 100,000 example looks wrong in the left-hand plot.

R Code Below

# 32018
# Non-Bayer 2D interpolation

# Model pixels as:
# P(x,y) = (a_0 + a_1 * x + a_2 * x^2) * (b_0 + b_1 * y + b_2 * y^2)
#        = [a_0 a_1 a_2 b_0 b_1 b_2] * [1 x x^2 1 y y^2]^T
# And so we need to find the a_i and b_i from specific pixel and coordinate values.
#
# install.packages('ptw')

require(ptw)

{

return(p2)
}

Npoints <- 8
xy <- array(c(c(1,1), c(1,5), c(2,3), c(3,6), c(4,2), c(4,4), c(5,6), c(6,3)), c(2,Npoints))
xy_norm <- (xy - min(xy))/(max(xy) - min(xy))
Pxy <- array(runif(N, 0, 255), c(1,Npoints))

ImageSize <- max(xy)

Ired <- array(0,c(ImageSize, ImageSize))

for (idx in 1:Npoints)
{
Ired[xy[1,idx], xy[2,idx]] = Pxy[idx]
}

Ihatred <-array(0,c(ImageSize,ImageSize))
Ihatred2 <- Ired
count <- 0

while (sum(abs(Ihatred - Ihatred2)) > 1)
{
Ihatred <- Ihatred2

Fredsize <- dim(Fred)[1]
cutoff1 <- Fredsize / 8
cutoff2 <- 7 * Fredsize / 8

Fred[4:16,] <- 0
Fred[,4:16] <- 0

IFred <- fft(Fred, inverse=TRUE) / Fredsize / Fredsize

Ihatred2 <- Re(IFred[1:ImageSize, 1:ImageSize])

for (idx in 1:Npoints)
{
Ihatred2[xy[1,idx], xy[2,idx]] = Pxy[idx]
}

count <- count + 1
if (count > 100)
{
break
}
}

par(mfrow=c(1,2))
image(Ired, col= grey(seq(0, 1, length = 256)))
for (idx in 1:Npoints)
{
points(xy_norm[1,idx] , xy_norm[2,idx] , col="red", lwd=5)
}
title('Original with only known red pixels')
image(Ihatred2, col= grey(seq(0, 1, length = 256)))
title('POCS solution')
for (idx in 1:Npoints)
{
points(xy_norm[1,idx] , xy_norm[2,idx] , col="red", lwd=5)
}

• Won't the initial image being mostly black influence what it converges to? – enthdegree Jul 8 '16 at 4:59
• @enthdegree : I do not believe it should matter, and I have added two examples. One has the unknown pixels set to 100,000 and the other has them set to 128. Both converge to the same image. – Peter K. Jul 8 '16 at 11:50