Polynomial Approximation of 2D image -Part 2

Pi19404
February 18, 2014

Contents

Contents
Polynomial Approximation of 2D image -Part 2 3

Introduction . . . 0.1.1 Poynomial 0.1.2 Code . . . . 0.2 code . . . . . . . . . References . . . . . . .

0.1

. . . . . . . . . . . Interpolation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

. . . . .

3 3 8 13 13

2 | 14

Polynomial Approximation of 2D image -Part 2

Polynomial Approximation of 2D image -Part 2

0.1 Introduction
In this article we will look at polynomial approximation of image neighborhood and efficient implementation in C/C++ OpenCV platform.

    

In the earlier articles we looked at how to represent a image region using quadratic polynomials Here we will look at efficient way to compute the polynomial basis representation of an image. The concept of seperable symmetric/anit symmetric convolution filters described in the earlier articles will be also useful. The first step is the computation of seperable filter coefficients The basis functions are as shown in figure 1a

0.1.1 Poynomial Interpolation

    

Let us assume that the neighborhood is approximated by a quadratic polynomial f (x; y ) = a + bx + cy + dx2 + ey 2 + f xy We known the value of in the image
f (x; y )

ie the image values at all points

2N + 1

Let us consider a rectangular neighborhood of size a point. The aim is to find the coefficients borhood. The x and y value ranges from
a; b; c; d; e; f

about

for this neigh-

N toN

3 | 14

Polynomial Approximation of 2D image -Part 2

(a) 2D basis function

4 | 14

Polynomial Approximation of 2D image -Part 2

For each (x; y ) eg (0; 0); (0; 1) : : : (N; N ) we can equate it with general form of f (x; y ) to get a equation.
f (0; f (1;

0) =

(1) (2) (3) (4)

1) =

+b+c+d+e+f
f (0;

1) =

+c+e
:::

These set of equations can be written in form of the matrix

2 3 2 61 6 41
:::

0 1 0

0 1 1

0 1 0

6b7 2 3 6 7 f (0; 0) 7 6 7 1 17 6 c 7 = 4 f (0; 1) 5 1 05 6 d7 6 7 f (0; 1)::: 4e5

0 0
f

(5)

In general it can be written as

2 3 2 61 6 41
:::

x1 x2 x3

y1 y2 y3

x1 x2

2 2

y1 y2 y3

2 2 2

2 x3

6b7 2 3 6 7 f1 7 6 7 x2 y2 7 6c7 = 4 f2 5 x3 y3 5 6 d 7 6 7 f 3::: 4e5

x1 y1 f

(6)

This can be expressed as

Ax

=
A
T T

b b

(7) (8) (9)

(A
x

= (A

A)1(A
:::

A)x

b)

The matrix

(AT b)

can be expressed as
1
x1 y1

2 6 6 6 6 6 6 4

1
x2 y2

1
x3 y3

32

2 x1 2 y1

2 x2 2 y2
x2 y2

x3 y3

6b7 : : :7 6f2 7 76 7 6 7 6 7 : : :7 6f3 7 7 6 7 = 6c7 7 6 7 6d7 ::: 76 : 7 6 7 4e5 : : :5 4 : 5

::: : f

f1

2 3
a

(10)

x1 y1

x3 y3

This matrix simply compute the projection of basis onto image

5 | 14

Polynomial Approximation of 2D image -Part 2 neighborhood being considered

a

XX
fi
x y

(11)

= = = =

XX XX
x y

xi

fi

(12)

yi

fi

(13)

XX XX
x y x y

xi

fi

(14)

yi

fi

(15)

XX
x y

yi xi

fi

(16) (17)

  

This can be efficiently computed by performining 6 convolution operations on the image The basis matrix for convolution are [1; x; y; x2 ; y 2 ; xy ] are N xN matrices The matries are displayed below of a neighhborhood of 7 which

6 | 14

Polynomial Approximation of 2D image -Part 2 are

7x7

matrices.
2 61 6 61 6 =6 61 61 6 41
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1

3 7 7 7 7 7X 7 7 7 5

3 63 6 63 6 =6 63 63 6 43 3
3

2 2 2 2 2 2 2
2
9 4 4 4 4 4 4 4

1 1 1 1 1 1 1
1 1 1 1 1 1 1

0 0 0 0 0 0 0 0 0 0 0 0 0 0

1 1 1 1 1 1 1 1 1 1 1 1 1 1

2 2 2 2 2 2 2 4 4 4 4 4 4 4

3 3 3 3 3 3 3 9 9 9 9 9 9 9

3 7 7 7 7 7 (18) 7 7 7 5
(19)

2 6 6 6 6 =6 6 6 6 4

3 3 3 3 3 3 3 2 2 2 2 2 2 2 1 1 1 1 1 1 1
0 1 2 3 0 1 2 3 9 4 1 0 1 4 9 9 4 1 0 1 4 9 0 1 2 3 9 4 1 0 1 4 9 9 4 1 0 1 4 9 0 1 2 3 9 4 1 0 1 4 9 0 1 2 3 9 4 1 0 1 4 9 0 1 2 3 0 1 2 3

7 69 7 6 7 69 7 6 7 X 2 = 69 7 6 7 69 7 6 5 49
9

7 7 7 7 7 (20) 7 7 7 5
(21)

64 6 61 6 2 Y =6 60 61 6 44
9

6 6 7 6 7 6 3 7 6 7 7 XY = 6 0 6 7 6 3 7 6 7 4 6 5

6 4

3 2

0 2 4 6

0 1 2 3

0 3 6 9 0 2 4 67 7 0 1 2 37 7
0 0 0 0 0 0 1 2 3 2 4 6 6 5 9 3 7 7

0 7 7(22)

Each of these filters can be expressed as product of 2 matrices

7 | 14

Polynomial Approximation of 2D image -Part 2 R,C where R is row vector and C is a column vector
2 3 617 6 7 617 6 7 7 =6 617 617 6 7 415
1 1

2 3 617 6 7 617 6 7 7 =6 617 617 6 7 415

1 1

;X

3 2 1
2 3
1

3 (23)

7 =6 6 0 7 6 1 7 6 7 4 2 5
3

3 627 6 7 617 6 7

617 6 7 617 6 7 2 7 1 1 ;X = 6 617 617 6 7 415 2 6 6 6 6 =6 6 6 6 4

9 (24)

2 3 647 6 7 617 6 7 2 7 Y =6 607 617 6 7 445

9 9

; XY

3 27 7 17 7
0 7 7 1 7 7 2 5 3

3 2 1

3 (25)

All the filters are weighted by gaussian kernel that gives more importance to pixels near the origin,and can be considered as weighted interpolation.

0.1.2 Code
The source code can be found at git repository

To realize the basis function ,seperable filter for a size 7 are shown below

//kernel is 1D gaussian kernel kernel=[0.004433048; 0.054005582; 0.24203624; 0.39905027; 0.24203624; 0.054005582; 0.004433048] //separable basis 1*kernel and 1&kernel -- 1 row=[0.004433048; 0.054005582; 0.24203624; 0.39905027; 0.24203624; 0.054005582; 0.004433048] col=[0.004433048; 0.054005582; 0.24203624; 0.39905027; 0.24203624;

8 | 14

Polynomial Approximation of 2D image -Part 2

0.054005582; 0.004433048] //separable filters for 1*kernel and x*kernel --x row=[0.004433048; 0.054005582; 0.24203624; 0.39905027; 0.24203624; 0.054005582; 0.004433048] col=[-0.013299144, -0.10801116, -0.24203624, 0, 0.24203624, 0.10801116, 0.013299144] //separable filters for y*kernel and 1*kernel --y row=[-0.013299144, -0.10801116, -0.24203624, 0, 0.24203624, 0.10801116, 0.013299144] col=[0.004433048; 0.054005582; 0.24203624; 0.39905027; 0.24203624; 0.054005582; 0.004433048] //separable filters for 1*kernel and x^2*kernel --x^2 row=[0.004433048; 0.054005582; 0.24203624; 0.39905027; 0.24203624; 0.054005582; 0.004433048] col=[0.039897431, 0.21602233, 0.24203624, 0, 0.24203624, 0.21602233, 0.039897431] //separable filter for y^2*kernel and 1*kernel --y^2 row=[0.039897431, 0.21602233, 0.24203624, 0, 0.24203624, 0.21602233, 0.039897431] col=[0.004433048; 0.054005582; 0.24203624; 0.39905027; 0.24203624; 0.054005582; 0.004433048] //separable filters for x*kernel and y*kernel --xy row=[-0.013299144; -0.10801116; -0.24203624; 0; 0.24203624; 0.10801116; 0.013299144] col=[-0.013299144, -0.10801116, -0.24203624, 0, 0.24203624, 0.10801116, 0.013299144]

  

a dot product of row col kernels specified above will provided a 2D kernel corresponding to basis function specfied in the figure 1a based on the desired size of the neighborhood being approximated.
T

It can be noted that the filter coefficients are symmetric and are either symmetic or anti symmetric and are all of the same size. The function to generate the filter coeffients is implemented by the function

9 | 14

Polynomial Approximation of 2D image -Part 2

//input the function is the size, //size of neighborhood over which the basis is to be computed void computeSeparableFilters(int size);

    

The result of the various filters can be computed by making a single pass over the entire image The class
S eperableS C onvolution

is used for achieving this.

The details of seperable convolution using symmetric/anti-symmetric filters can be found in the article The polyBasis class implements the polynomial basis computation The method voidcomputeS eparableF ilters(intsize) implements the computation of seperable filter coefficients.

SeperableSConvolution filter; void computeSeparableFilters(int size) { ..... //get the 1D gaussian weighing function Mat k=cv::getGaussianKernel(size,1,CV_32F); //c is row and k is column matrix filter.setKernel(c,k,true,true); //x and y kernels filter.setKernel(x,k,true,true); filter.setKernel(k,y,true,true); //xx and yy and xy kernels filter.setKernel(xx,k,true,true); filter.setKernel(yy,k,true,true); // Mat yk=multiply(y,k); filter.setKernel(x,yk,true,true); ..... } PolyBasis poly; poly.computeFilters(7); poly.filter.apply(i1,res[prev_index]); poly.filter.apply(i2,res[cur_index]);

 

The last step in estimation required is the matrix The matrix

A
T

(AT A)1

is 6xN matrix

10 | 14

Polynomial Approximation of 2D image -Part 2

Thus the matrix

2 6 6 6 6 6 6 4

(AT A) 1
x1 y1

is

6x6 1
x3 y3

matrix which is expressed as

: : :7 7 : : :7 7 : : :7 7 : : :5 ::: :::

1
x2 y2

2 x1 2 y1

2 x2 2 y2
x2 y2

2 x3 2 y3
x3 y3

61 6 61 6 4: : :
:::

x1 x2 x3

y1 y2 y3

x1

2 2

y1

2 2

x1 y1

3 7 7 5

x2

2 x3

y2

2 y3

x2 y2 7

x3 y3 7 (26)

x1 y1

2 P 1 P x;y 6 x 6 Px;y i 6 yi P x;y 2 =6 6 6 Px;y xi 2 4 y P x;y i

x;y

i P x;y 2 x P x;y i x;y P

P P P P
x;y x;y P

yi

xi yi yi

xi yi xi

2 3

Px;y i 3 x P x;y i 2
x;y P

P P P P
x;y x;y P

yi

2 2 3

xi yi yi

x;y

x;y

xi yi

i Px;y i 2 x y x;y i

x y

2
i

x;y P

yi xi yi

yi xi xi

x;y

x;y

x;y

x;y

yi xi

Px;y
x;y

xi yi

2 2 3

x;y P

xi yi yi

2 2 4

x;y x;y

xi yi

yi xi

i i P x;y 2 x yi 7 7 Px;y i 2 xi yi 7 x;y 7(28) P 3 xi yi 7 x;y 7 P 35 xi yi P x;y 2 2 x;y

x y

(27)

xi yi

 

The functions tions

xi ; yi

are odd function while

xi ; yi

are even func-

This implies that

X
xi

=0 =0 =0 =0 =0 =0 =0 =0

(29)

X
yi

x;y

(30)

X X
x;y x;y

x;y

xi yi

(31)

xi yi

(32)

X X
x;y

xi

(33)

yi

(34)

X
x;y

x;y

xi yi

(35)

X
x;y

xi yi

(36) (37)

11 | 14

Polynomial Approximation of 2D image -Part 2

this implies
2
1

(A

6 0 6 6 0 A) = 6P 2 6 6Px;y xi 2 4 y
x;y i

0
x;y

xi

P
x;y

xi

P
x;y

yi

0 0 0 0 0

3 7 7 7 7 (38) 7 7 5
2 2

0 0 0 0

0
x;y

yi

0
xi

0 0 0

P P

0
x;y

x;y

xi yi

2 2

x;y P

xi yi yi

2 2 4

x;y

x;y

xi yi

for a 7x7 neighborhood size we get(39)

0 0 0:9959 0 0 2:9304 0:9918 0 0:9959 0 0 0:9918 2:9304 0 0 0 0 0 0 0:9918

6 0 0:9959 0 6 6 0 0 0:9959 6 60:9959 0 0 6 40:9959 0 0

0 0 0

1:0000

7 7 7 7(40) 7 7 5

 

This computation is performed by the function

void computeBasisMatrix(int size,Mat &A);

The matrix
(AT

A)

is called as Vandermonde matrix .

//main function to be called to precompute the filter coefficients //and vandermore matrix void computeFilters(int size);

   

The last step remaining is normalization of the filter coefficients as a result of Least sequares estimation process that require application of the matrix G 1 on each pixel of multi channel image. It can be noted that result of applying the kernel on vector formed by the channel of length N results in a matrix of length N. It can also be noted that the transformation is linear. This is achieved using the class in the article
C hannelF ilter

that can be found

PolyBasis poly; Mat dst1,dst2; poly.computeFilters(7);

12 | 14

Polynomial Approximation of 2D image -Part 2

//cfilter is the channel filter poly.cfilter.apply(res[prev_index],dst1); poly.cfilter.apply(res[cur_index],dst2);

Thus a modular and efficient method to compulte the local polynomial basis of image is provided in the article.

0.2 code
The code for the same can be found in the git repository https: //github.com/pi19404/OpenVision in files ImgBasis/polybasis.cpp and ImgBasis/polybasis.hpp files.

13 | 14

Bibliography

Bibliography
[1] Kenneth Andersson and Hans Knutsson.  Continuous normalized convolution. In: ICME (1). IEEE, 2002, pp. 725728.

dblp.uni-trier.de/db/conf/icmcs/icme2002-1.html#AnderssonK02.
[2]

isbn:

0-7803-7304-9.

url: http://

Kenneth Andersson, Carl-Fredrik Westin, and Hans Knutsson.  Prediction from 87.3 (Mar. 22, 2007), pp. 353365. url: http://dblp.uni- trier.de/db/ journals/sigpro/sigpro87.html#AnderssonWK07. o-grid samples using continuous normalized convolution. In: Signal Processing

[3]

Gunnar Farnebck.  Motion-based Segmentation of Image Sequences. LiTH-ISYEX-1596. MA thesis. SE-581 83 Linkping, Sweden: Linkping University, 1996.

[4]

Gunnar Farnebck.  Polynomial Expansion for Orientation and Motion Estimation. Dissertation No 790, ISBN 91-7373-475-6. PhD thesis. SE-581 83 Linkping, Sweden: Linkping University, Sweden, 2002.

[5]

Gunnar Farneback.  Two-Frame Motion Estimation Based on Polynomial Expansion. In: SCIA. LNCS 2749. Gothenburg, Sweden, 2003, pp. 363370.

14 | 14