**1.4 Methods of matrix data-to-image conversion**

Let's consider used algorithmic primitives of interpolation the colors applied to form the color image in digital photographic cameras.

Let light filters of primary colors are allocated in Bayer's grid according to a picture 1.

The algorithms used for recovery of missing color components, are represent "know-how" of vendors and, as a rule, vary depending on model of the camera and type of a photosensitive matrix. However most often they are constructed on the basis of linear and median filtrations primitives, threshold gradients and persistence of color tone.


Fig. 1. Color filter array in the Bayer structure

Digital Camera Identification Based on Original Images 305

Color interpolation can also be performed by median filtration. Application of the median filter offered in [8] is carried out in two stages. On the first by bilinear interpolation components *\* R ( ),* x *\* G( )* x and *\* B( )* x are calculated, and then the difference between channels with the subsequent median filtration. Let *Mrg( )* x , *Mgb( )* x , *Mrb( )* x designate differences after a median filtration. For each pixel sampling of missing colors is estimated as a difference between current value of a component and an appropriate difference after a

Recovery of colors can be performed also by the gradient method offered in [9] and for the first time used in photocamera of Kodak DCS 200. The method is based on three-stage process which saves boundaries at interpolation led in a direction, perpendicular their orientation. In the beginning the G-channel along boundaries is interpolated. For example, in case of interpolation of "green" pixel in a position (4,4) horizontal and vertical gradients

> *4,4 4,2 4,6 4,4 4,4 2,4 6.4 4,4*

If horizontal gradient *ǻ4,4* greater than vertical gradient *V4,4*, it specifies to possible boundary in a horizontal direction and then interpolation of value of "green" pixel is performed only

*G( 4,4) ( 3,4 5,4 g g )/2* .

And on the contrary. If horizontal and vertical gradients are equal, values of pixels of the

*G( 4,4) ( 3,4 4,3 4,5 5,4 gggg )/4* .

Missing *R( )* x and *B( )* x channels are recovered by interpolation on the basis of constant color tone. For example, the missing blue component of pixels with coordinates (3,4) and (4,4)

*B( 3,4) b G( 3,3) b G( 3,5 )) / 2 G( 3,4) 3,3 3,5* ,

*B( 4,4) (b G( 3,3) b G( 3,5 ) b G( 5,3) b G( 5,5 )) / 4 G( 4,4). 3,3 3,5 5,3 5,5*

The method on the basis of the variable gradients activated on a threshold (Threshold Based Variable Number of Gradient) is based on a variable amount of the gradients which usage is controlled by exceeding of threshold values. In the given primitive possibility to use gradients on all eight directions, namely in two horizontal, two vertical (N, S, E and W

In each direction on a matrix of pixels the gradient for the selected point on the basis of an array of 5x5 adjacent pixels is calculated. The choice of a configuration of a neighborhood is

"green" channel calculated by averaging four adjacent pixels:

according to [10] is interpolated by following expressions:

**1.8 Interpolation based on variable threshold gradients** 

accordingly) and four diagonal NW, SW, NE and SE are added.

*H (b b ) / 2 b , V b b /2 b .* 

**1.7 Interpolation based on median filtering** 

median filtration.

for "blue" are calculated:

in a vertical direction:

#### **1.5 Interpolation based on linear filtering**

The elementary primitive of color interpolation is the algorithm of a bilinear filtration which is applied to each channel independently. For channel G ("green") the filter kernel represents:

$$k = \frac{1}{4} \begin{bmatrix} 0 & 1 & 0 \\ 1 & 4 & 1 \\ 0 & 1 & 0 \end{bmatrix}' $$

And for channels "red" and "blue" accordingly:

$$k = \frac{1}{4} \begin{bmatrix} 1 & 2 & 1 \\ 2 & 4 & 2 \\ 1 & 2 & 1 \end{bmatrix}\_+$$

Other algorithm of the general application is bicubic interpolation, at which kernels for channels of the primary colors are the following:

$$k\_{\mathcal{G}} = \frac{1}{256} \begin{bmatrix} 0 & 0 & 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & -9 & 0 & -9 & 0 & 0 \\ 0 & -9 & 0 & 81 & 0 & -9 & 0 \\ 1 & 0 & 81 & 256 & 81 & 0 & 1 \\ 0 & -9 & 0 & 81 & 0 & -9 & 0 \\ 0 & 0 & -9 & 0 & -9 & 0 & 0 \\ 0 & 0 & 0 & 1 & 0 & 0 & 0 \end{bmatrix},$$

$$k\_{\mathcal{R},\mathcal{B}} = \frac{1}{256} \begin{bmatrix} 1 & 0 & -9 & -16 & -9 & 0 & 1 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ -9 & 0 & 81 & 144 & 81 & 0 & -9 \\ -16 & 0 & 144 & 256 & 81 & 0 & -16 \\ -9 & 0 & 81 & 144 & 81 & 0 & -9 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 1 & 0 & -9 & -16 & 9 & 0 & 1 \end{bmatrix}$$

.

#### **1.6 Interpolation based on color hue constance**

Color interpolation can be led also on the basis of assumptions of persistence of color tone in localized areas. Generally, selection of a color tone constant is possible considering property of orderliness of colors within a color circle. Interpolation of a constant of the color tone, offered in [7], is one of the most widespread methods used up to professional cameras. The constant of color tone is defined as a difference between the main color components. At the first stage the algorithm interpolates green channel *G*, using the bilinear method considered earlier. For an estimation of an error of "red" pixels bilinear interpolation of a difference *\* R ( ) G( )* x x , which then incremented by *G( )* x . The channel "blue" is recovered similarly.

### **1.7 Interpolation based on median filtering**

304 Applications of Digital Signal Processing

The elementary primitive of color interpolation is the algorithm of a bilinear filtration which is applied to each channel independently. For channel G ("green") the filter kernel

> *010 <sup>1</sup> k 141 4*

,

*010* <sup>ª</sup> <sup>º</sup> « » « » « » <sup>¬</sup> <sup>¼</sup>

> 121 <sup>1</sup> <sup>242</sup>

> > 121

*00 0 1 0 00 00 9 0 900 0 9 0 81 0 9 0*

*00 9 0 900 00 0 1 0 00*

*1 0 9 16 9 0 1 0 0 0 0 00 0 9 0 81 144 81 0 9*

*0 0 0 0 00 0 1 0 9 16 9 0 1*

,

.

<sup>ª</sup> <sup>º</sup> « » « » « » <sup>¬</sup> <sup>¼</sup>.

Other algorithm of the general application is bicubic interpolation, at which kernels for

*<sup>1</sup> <sup>k</sup> 1 0 81 256 81 0 1 <sup>256</sup> 0 9 0 81 0 9 0*

*<sup>1</sup> <sup>k</sup> 16 0 144 256 81 0 16 <sup>256</sup> 9 0 81 144 81 0 9*

Color interpolation can be led also on the basis of assumptions of persistence of color tone in localized areas. Generally, selection of a color tone constant is possible considering property of orderliness of colors within a color circle. Interpolation of a constant of the color tone, offered in [7], is one of the most widespread methods used up to professional cameras. The constant of color tone is defined as a difference between the main color components. At the first stage the algorithm interpolates green channel *G*, using the bilinear method considered earlier. For an estimation of an error of "red" pixels bilinear interpolation of a difference *\* R ( ) G( )* x x , which then incremented by *G( )* x . The channel "blue" is recovered

<sup>ª</sup> <sup>º</sup> « » « » « » « » « » « » « » « » « » <sup>¬</sup> <sup>¼</sup>

<sup>ª</sup> <sup>º</sup> « » « » « » « » « » « » « » « » « » <sup>¬</sup> <sup>¼</sup>

4

*k*

**1.5 Interpolation based on linear filtering** 

And for channels "red" and "blue" accordingly:

channels of the primary colors are the following:

*G*

*R,B*

**1.6 Interpolation based on color hue constance** 

represents:

similarly.

Color interpolation can also be performed by median filtration. Application of the median filter offered in [8] is carried out in two stages. On the first by bilinear interpolation components *\* R ( ),* x *\* G( )* x and *\* B( )* x are calculated, and then the difference between channels with the subsequent median filtration. Let *Mrg( )* x , *Mgb( )* x , *Mrb( )* x designate differences after a median filtration. For each pixel sampling of missing colors is estimated as a difference between current value of a component and an appropriate difference after a median filtration.

Recovery of colors can be performed also by the gradient method offered in [9] and for the first time used in photocamera of Kodak DCS 200. The method is based on three-stage process which saves boundaries at interpolation led in a direction, perpendicular their orientation. In the beginning the G-channel along boundaries is interpolated. For example, in case of interpolation of "green" pixel in a position (4,4) horizontal and vertical gradients for "blue" are calculated:

$$\begin{aligned} H\_{4,4} &= \left| (b\_{4,2} + b\_{4,6}) \right\rangle / \left| 2 - b\_{4,4} \right\rangle. \\ V\_{4,4} &= \left| b\_{2,4} + b\_{6,4} \right| / \left| 2 - b\_{4,4} \right\rangle. \end{aligned}$$

If horizontal gradient *ǻ4,4* greater than vertical gradient *V4,4*, it specifies to possible boundary in a horizontal direction and then interpolation of value of "green" pixel is performed only in a vertical direction:

$$G(4,4) = (g\_{3,4} + g\_{5,4}) / \, ^\prime 2 \dots$$

And on the contrary. If horizontal and vertical gradients are equal, values of pixels of the "green" channel calculated by averaging four adjacent pixels:

$$G(4,4) = (\mathcal{g}\_{3,4} + \mathcal{g}\_{4,3} + \mathcal{g}\_{4,5} + \mathcal{g}\_{5,4}) / \, ^\dagger \text{A} \dots$$

Missing *R( )* x and *B( )* x channels are recovered by interpolation on the basis of constant color tone. For example, the missing blue component of pixels with coordinates (3,4) and (4,4) according to [10] is interpolated by following expressions:

*B( 3,4) b G( 3,3) b G( 3,5 )) / 2 G( 3,4) 3,3 3,5* ,

$$B(4,4) = (b\_{3,3} - G(3,3) + b\_{3,5} - G(3,5) + b\_{5,3} - G(5,3) + b\_{5,5} - G(5,5)) / \, / \, 4 + G(4,4).$$

#### **1.8 Interpolation based on variable threshold gradients**

The method on the basis of the variable gradients activated on a threshold (Threshold Based Variable Number of Gradient) is based on a variable amount of the gradients which usage is controlled by exceeding of threshold values. In the given primitive possibility to use gradients on all eight directions, namely in two horizontal, two vertical (N, S, E and W accordingly) and four diagonal NW, SW, NE and SE are added.

In each direction on a matrix of pixels the gradient for the selected point on the basis of an array of 5x5 adjacent pixels is calculated. The choice of a configuration of a neighborhood is

Digital Camera Identification Based on Original Images 307

allocated between *1iN* d d is considered as a set of the pixels correlated with adjacent pixels. Considering periodic layout of filters (lattice filters - color filter array (CFA) the given correlations will show periodicity. Being based on it in considered article the assumption about identity of scales of pixel sets with different *x* and *y* that a set of the correlated pixels,

Let's consider the right member of equation (1.1) as *filter F* applied to *I( )* x , designating operation of a filtration *F(I( ))* x as well as summed averaged square errors from both sides

*<sup>1</sup> MSE(F(I( ))) I(x x ,y y ) I(x,y)| WH*

Where H and W - height and width of an image accordingly. Adding the virtual correlated

*x 1y 1 i 1 <sup>1</sup> MSE(F(I( ))) I(x x ,y y ) WH*

*M <sup>T</sup> SE(F(I( )),I( )) X AX* x x ,

The coefficient of a matrix *A* contains the full information for determination of variable vector *ȃ*, however, in article obtaining ȃ optionally and for the further analysis enough

It was empirically revealed that the correlated pixels mask shown in a figure 2, yields good

*i i j j*

¦ ¦ ' ' ' ' , *1 i,j N 1* <sup>d</sup> d .

The extension of the equation (1.3) gives the square form rather ȃ = {ǂ 1, ǂ2, …, ǂN+1} T:

*i ii*

*<sup>2</sup> WH N*

x ¦¦¦D ' ' (2)

*i ii*

x ¦¦¦D ' ' . (3)

*2*

*WH N*

*x 1y 1 i 1*

pixel ǂN+1 =-1, ƦxN+1 = Ʀ yN+1=0, the equation (1.2) assumes more arranged air:

*<sup>1</sup> A(i, j) I(x x ,y y ) I(x x ,y y ) WH*

Fig. 2. Correlated pixels mask, where — is a center and *Ʀ* — correlated pixels

*(1 i, j N 1),*

d d

 *\* A i, j A(i, j) A / A,*

ª º ¬ ¼

and according to their weight for each pixel in *I( )* x are identical.

*W H*

*x 1y 1*

On a following step the analysis of principal components is done.

Numerical values of elements *A* after obtaining are normalized:

from *I( )* x , we receive:

Where

matrix affirms that *A*.

result (N=12).

done by empirically detected feeble dependence of a difference of the calculated gradient from colors and considered pixels.

For example, vertical, horizontal and diagonal gradients for the "red" pixel allocated in a point (3,3) will be equal accordingly:

$$\begin{split} \mathcal{N} &= \left| \mathcal{g}\_{2,3} - \mathcal{g}\_{4,3} \right| + \left| r\_{1,3} - r\_{3,3} \right| + \left| b\_{2,2} - b\_{4,2} \right| / 2 + \left| b\_{2,4} - b\_{4,4} \right| / 2 + \left| \mathcal{g}\_{1,2} - \mathcal{g}\_{3,2} \right| / 2 + \left| \mathcal{g}\_{1,4} - \mathcal{g}\_{3,4} \right| / 2 \\ \mathcal{E} &= \left| \mathcal{g}\_{3,2} - \mathcal{g}\_{3,4} \right| + \left| r\_{3,3} - r\_{3,5} \right| + \left| b\_{2,2} - b\_{4,2} \right| / 2 + \left| b\_{4,2} - b\_{3,4} \right| / 2 + \left| \mathcal{g}\_{2,3} - \mathcal{g}\_{2,5} \right| / 2 + \left| \mathcal{g}\_{4,3} - \mathcal{g}\_{4,5} \right| / 2 \\ \mathcal{S} &= \left| b\_{2,4} - b\_{4,2} \right| + \left| r\_{5,1} - r\_{3,3} \right| + \left| \mathcal{g}\_{2,3} - \mathcal{g}\_{3,2} \right| / 2 + \left| \mathcal{g}\_{3,4} - \mathcal{g}\_{4,3} \right| / 2 + \left| \mathcal{g}\_{3,2} - \mathcal{g}\_{4,1} \right| / 2 + \left| \mathcal{g}\_{4,3} + \mathcal{g}\_{5,2} \right| / 2. \end{split}$$

On the basis of a set containing 8 gradients, threshold *T* is calculated, allowing to define, what directions were used. *T* it is defined as*T k min k (max min) 1 2* , where *min* and *max* are the minimum and maximum gradients accordingly, and *k1* and *k2* constants. Author's values are *k1=1,5* and *k2=0,5*. Those directions which gradient is less than a threshold are selected, and for each selected direction mean values for "blue", "red" and "green" are calculated. For example, at coordinates (3,3) mean values for directions N, E, SW are the following:

*<sup>N</sup> R (r r ) / 2 1,3 3,3* , *<sup>N</sup> G g 2,3* , *NB (b b ) / 2 2,2 2,4* , *<sup>E</sup> R (r r ) / 2 3,3 3,5* , *<sup>E</sup> G g 3,4* , *<sup>E</sup> B (b b ) / 2 2,4 4,4* , *SW R (r r ) / 2 3,3 5,1* , *SW G (g g g g ) / 4 3,2 4,1 4,3 5,2* , *SW B b 4,2* .

Let's designate mean values red, blue and green as *Ravg*, *Gavg Bavg* accordingly. Then for the selected pixel mean averaging values for red, dark blue and green in the selected directions will be: *Ravg = (Rs+RE+Rse)*/3*, Gavg = (Gs+GE+GSE)/3, Bavg = (BS+BE+BSE)* (for pixel pixel (3,3) and directions *S, E, SE*). A final estimation of missing color components levels are: *G* (3,3) *=r3,3* + *(Gavg-Ravg)* and *B* (3,3) *=r3,3* + *(Bavg+Ravg)* [11].

#### **2. Cameras identification techniques**

#### **2.1 Camera identification based on artifacts of color interpolation**

There are several approaches to the implementation of identification systems for digital cameras based on the above characteristics.

In [12] cameras identification is done based on color interpolation features. The recognition process involves the following steps.

Designating *I( )* as one of *R( ) G( )* , *B( )* channels provided that the pixel in coordinates *(x,y)* is correlated linearly with other pixels, it is possible to express value of brightness of a color component as the weighed total of brightness of components of adjacent pixels:

$$I(\mathbf{x}, \mathbf{y}) = \sum\_{i=1}^{N} \alpha\_i I(\mathbf{x} + \Delta \mathbf{x}\_i, \mathbf{y} + \Delta \mathbf{y}\_i) \, , \tag{1}$$

Where *N* is a number of correlated pixels, ǂi, Ʀxj, Ʀyj - weight and offset on an axis *x* and an axis *y* of the pixel correlated from *i* th pixel accordingly. The set of such coordinates Ʀȣi, Ʀyi

allocated between *1iN* d d is considered as a set of the pixels correlated with adjacent pixels. Considering periodic layout of filters (lattice filters - color filter array (CFA) the given correlations will show periodicity. Being based on it in considered article the assumption about identity of scales of pixel sets with different *x* and *y* that a set of the correlated pixels, and according to their weight for each pixel in *I( )* x are identical.

Let's consider the right member of equation (1.1) as *filter F* applied to *I( )* x , designating operation of a filtration *F(I( ))* x as well as summed averaged square errors from both sides from *I( )* x , we receive:

$$MSE(F(I(\bullet))) = \frac{1}{WH} \sum\_{x=1}^{W} \sum\_{y=1}^{H} \left| \sum\_{i=1}^{N} \alpha\_i I(x + \Delta x\_i, y + \Delta y\_i) - I(x, y) \right|^2 \tag{2}$$

Where H and W - height and width of an image accordingly. Adding the virtual correlated pixel ǂN+1 =-1, ƦxN+1 = Ʀ yN+1=0, the equation (1.2) assumes more arranged air:

$$MSE(F(I(\bullet))) = \frac{1}{V \text{NH}} \sum\_{x=1}^{\mathcal{W}} \sum\_{y=1}^{H} \left| \sum\_{i=1}^{N} \alpha\_i I(\mathbf{x} + \Delta \mathbf{x}\_i, y + \Delta y\_i) \right|^2 \tag{3}$$

The extension of the equation (1.3) gives the square form rather ȃ = {ǂ 1, ǂ2, …, ǂN+1} T:

$$MSE(F(I(\bullet)), I(\bullet)) = X^T A X \ , \ \ \ $$

Where

306 Applications of Digital Signal Processing

done by empirically detected feeble dependence of a difference of the calculated gradient

For example, vertical, horizontal and diagonal gradients for the "red" pixel allocated in a

*2,3 4,3 1,3 3,3 2,2 4,2 2,4 4,4 1,2 3,2 1,4 3,4 3,2 3,4 3,3 3,5 2,2 4,2 4,2 3,4 2,3 2,5 4,3 4,5 2,4 4,2 5 ,1 3,3 2,3 3,2 3,4 4,3 3,2 4,1 4,*

*3 5,2 g / 2.*

On the basis of a set containing 8 gradients, threshold *T* is calculated, allowing to define, what directions were used. *T* it is defined as*T k min k (max min) 1 2* , where *min* and *max* are the minimum and maximum gradients accordingly, and *k1* and *k2* constants. Author's values are *k1=1,5* and *k2=0,5*. Those directions which gradient is less than a threshold are selected, and for each selected direction mean values for "blue", "red" and "green" are calculated. For example, at coordinates (3,3) mean values for directions N, E, SW are the

*<sup>N</sup> R (r r ) / 2 1,3 3,3* , *<sup>N</sup> G g 2,3* , *NB (b b ) / 2 2,2 2,4* ,

*<sup>E</sup> R (r r ) / 2 3,3 3,5* , *<sup>E</sup> G g 3,4* , *<sup>E</sup> B (b b ) / 2 2,4 4,4* ,

*SW R (r r ) / 2 3,3 5,1* , *SW G (g g g g ) / 4 3,2 4,1 4,3 5,2* , *SW B b 4,2* .

Let's designate mean values red, blue and green as *Ravg*, *Gavg Bavg* accordingly. Then for the selected pixel mean averaging values for red, dark blue and green in the selected directions will be: *Ravg = (Rs+RE+Rse)*/3*, Gavg = (Gs+GE+GSE)/3, Bavg = (BS+BE+BSE)* (for pixel pixel (3,3) and directions *S, E, SE*). A final estimation of missing color components levels are: *G* (3,3)

There are several approaches to the implementation of identification systems for digital

In [12] cameras identification is done based on color interpolation features. The recognition

Designating *I( )* as one of *R( ) G( )* , *B( )* channels provided that the pixel in coordinates *(x,y)* is correlated linearly with other pixels, it is possible to express value of brightness of a

*I(x,y) I(x x ,y y )*

Where *N* is a number of correlated pixels, ǂi, Ʀxj, Ʀyj - weight and offset on an axis *x* and an axis *y* of the pixel correlated from *i* th pixel accordingly. The set of such coordinates Ʀȣi, Ʀyi

*i ii*

D ' ' ¦ , (1)

color component as the weighed total of brightness of components of adjacent pixels:

*N*

*i 1*

*=r3,3* + *(Gavg-Ravg)* and *B* (3,3) *=r3,3* + *(Bavg+Ravg)* [11].

**2.1 Camera identification based on artifacts of color interpolation** 

**2. Cameras identification techniques** 

cameras based on the above characteristics.

process involves the following steps.

*N g g r r b b /2 b b /2 g g /2 g g /2 E g g r r b b /2 b b /2 g g /2 g g /2*

 

*SW b b r r g g / 2 g g / 2 g g / 2 g*

from colors and considered pixels.

point (3,3) will be equal accordingly:

following:

$$A(\mathbf{i}, \mathbf{j}) = \frac{1}{\text{VHI}} \sum\_{x=1}^{\text{IV}} \sum\_{y=1}^{H} I(\mathbf{x} + \Delta \mathbf{x}\_i, y + \Delta y\_i) \cdot I(\mathbf{x} + \Delta \mathbf{x}\_j, y + \Delta y\_j) \,, \ 1 \le \mathbf{i}, \mathbf{j} \le \mathbf{N} + \mathbf{1} \ldots$$

The coefficient of a matrix *A* contains the full information for determination of variable vector *ȃ*, however, in article obtaining ȃ optionally and for the further analysis enough matrix affirms that *A*.

It was empirically revealed that the correlated pixels mask shown in a figure 2, yields good result (N=12).

On a following step the analysis of principal components is done.

Fig. 2. Correlated pixels mask, where — is a center and *Ʀ* — correlated pixels

Numerical values of elements *A* after obtaining are normalized:

$$\begin{aligned} A^\* \left( i, j \right) &= \left\lceil A(i, j) - \overline{A} \right\rceil / \left\lceil \overline{A} \right\rceil, \\ \left( 1 \le i, j \le N + 1 \right), \end{aligned}$$

Digital Camera Identification Based on Original Images 309

photocamera. The method also can be used for check on authenticity of pictures since in those areas at which there are signs of editing the responses of neural network increase by 2- 3 times. However, usage of replaced fragments of the image from the same photocamera as a background makes detection impossible [15]. Mehdi Kharrazi, etc. in [16] proposed the method of photocameras identification on the basis of image features. The task of determination of the camera with which help the analyzable picture has been received was thus considered. Proceeding from known sequence of information processing from a photosensitive matrix, it is possible to select two stages, importing the most essential distortions: a stage of debayerization, i.e. full-color image restoration and a postprocessing


1 2 *G*


Along with enumerated features metrics of image quality proposed in [16] has been used.

To classify vectors the SVM-based classifier has been used. At learning stage 120 of 300 images were used, with 180 at test stage. An average accuracy of camera identification in "1

In [17] an identification method based on proprietary interpolation algorithms used in camera. The basis of algorithm is pixel correlation estimation listed in [18] with two estimations: estimation of pixel value by adjacent pixels' values and demosaic kernel used for raw data processing. As precise configuration of area used for interpolation is uknown, several different configurations were used, with additional assumbtion about different interpolation algorithms used in gradient and texturized areas. Camera identification experiments were done on a basis of two cameras: Sony DSC-P51 Ȗ Nikon E-2100. It has been acknowledged that filter kernel increase leads to accuracy increase (for kernels 3x3,

Camera identification based on postprocessing algorithms features posess several disadvantages, the most fundamental is impossibility of practical use for one-model camera

In [19] camera identification method based on defective ("hot" and "dead" pixels) are presented but its effectiveness is limited for cameras without build-in pixel defects

*R* ,

*E*

2

*E*

2 1 2 *B*

*R* .

stage. Totally authors select 34 signs of classification, among them: - Cross-channel correlation R-G, R-B, B-G (3 scalar features).


*E*

and averaging each sub-band) (9 features).

4x4, 5x5, accuracies were from 89.3 to 95.7%).

identification, even in "1 ot of 2" case.

correction and "dark frame" subtraction.

**2.2 Camera identification based on matrix defects** 

All used metrics can be divided into following groups: - pixelwise difference metrics (MSE, AMSE).

2

1 2 *G*


out of 2" were 98,73% with 88,02% when images were regular photos.

*B* ,

values (3 scalar features).


$$\text{Where } \qquad \qquad \overline{A} = \frac{1}{\left(N+1\right)^2} \sum\_{1 \le i,j \le N+1} A(i,j) \dots$$

Let *\* A* it is *<sup>2</sup> N* a-dimensional vector of signs E . Accepting total number of vectors for training of neural network as *L 12 L { , ,..., }* E E E and their average according:

$$
\overline{\mathbb{B}} = \frac{1}{L} \sum\_{i=1}^{L} \mathbb{B}\_{i} \, \prime \, \prime
$$

$$
\mathbb{B}\_{i}^{\*} = \mathbb{B}\_{i} - \overline{\mathbb{B}}\_{i} \, \prime
$$

$$
i = 1, 2, \dots L.
$$

The covariace matrix will be:

$$\mathbf{C} = \{ \mathtt{\mathfrak{F}}\_1, \mathtt{\mathfrak{F}}\_2, \ldots, \mathtt{\mathfrak{F}}\_L \} \{ \mathtt{\mathfrak{F}}\_1, \mathtt{\mathfrak{F}}\_2, \ldots, \mathtt{\mathfrak{F}}\_L \}^T \ne (\mathtt{L} - \mathtt{1}) \dots$$

Let eigenvalues and eigenvectors C - {*nj1, nj2, … <sup>nj</sup>L*} and {*Ǐ1, Ǐ2, … ǏL*}, *1 2 L1 L { ... }* ] d] d ] d] accordingly.

Eigenvectors corresponding *M* greatest eigenvalues, form a vector of features *V* = [*nj1, nj2, … <sup>nj</sup>M*] T. The experiments led by authors, shown that *M 15* is enough. *\** E*<sup>i</sup>* as a result transforms to *\** \* E *i i V* with dimensionality reduction. Recognition of the image belonging to the specific camera was carried out by trained neural network of direct propagation with 15 input neurons, 50 neurons in the hidden layer (with tangential activation function) and one output neuron (sigmoidal activation function). If we denote a set of color interpolation algorithms D: *\* D D{}* where is the empty set corresponding initial *I( )* . Identification by color interpolation consists in defining conversion *\* d D* which with the greatest probability has been fulfilled over *I( )* , i.e. the purpose is to reference available *I( )* to one of classes *\* D* of the learning images set nearest to *I( )* , in space of conversion characteristics of debayerization. Thus, to each class of images one neural network should be set in correspondence. To select total number of neural networks the following conditions according to authors is necessary to consider. For the different *<sup>i</sup> d* it is necessary to use different neural networs ( *\* d ,d D ,d d 12 1 2* z ) considering essential difference of debayerization operators, applied to the channel of "green" and channels "red"/"blue" should use different neural networks for each color component. By authors of a method the best result were shown when three networks was used, one for each channel. Thus the total of neural networks makes *\* D D1* for each channel and *3( D 1)* totally [12]. The given method does not depend on color channels used for identification, on each channel the independent decision which is pseudo-independent, as channels are mutually correlated as shown in [13-14]. Accuracy of identification has been checked by authors on learning and test samplings on 100 images [13]. Accuracy of recognition of 7 algorithms of the interpolation was 100 % (errors of the first and second type are equal to zero). Accuracy of classification by offered methods oȠ real photocameras made 95-100 % depending on a 308 Applications of Digital Signal Processing

Let *\* A* it is *<sup>2</sup> N* a-dimensional vector of signs E . Accepting total number of vectors for

*L i i 1 1 L* E E ¦ ,

*i i , i 1,2,...L.* E E E 

*\* \* \* \* \* \*T C [ , ,..., ][ , ,..., ] / (L 1)* EE E EE E *12 L 12 L* .

Let eigenvalues and eigenvectors C - {*nj1, nj2, … <sup>nj</sup>L*} and {*Ǐ1, Ǐ2, … ǏL*}, *1 2 L1 L { ... }* ] d] d ] d]

Eigenvectors corresponding *M* greatest eigenvalues, form a vector of features *V* = [*nj1, nj2, … <sup>nj</sup>M*] T. The experiments led by authors, shown that *M 15* is enough. *\** E*<sup>i</sup>* as a result transforms to *\** \* E *i i V* with dimensionality reduction. Recognition of the image belonging to the specific camera was carried out by trained neural network of direct propagation with 15 input neurons, 50 neurons in the hidden layer (with tangential activation function) and one output neuron (sigmoidal activation function). If we denote a set of color interpolation algorithms D: *\* D D{}* where is the empty set corresponding initial *I( )* . Identification by color interpolation consists in defining conversion *\* d D* which with the greatest probability has been fulfilled over *I( )* , i.e. the purpose is to reference available *I( )* to one of classes *\* D* of the learning images set nearest to *I( )* , in space of conversion characteristics of debayerization. Thus, to each class of images one neural network should be set in correspondence. To select total number of neural networks the following conditions according to authors is necessary to consider. For the different *<sup>i</sup> d* it is necessary to use different neural networs ( *\* d ,d D ,d d 12 1 2* z ) considering essential difference of debayerization operators, applied to the channel of "green" and channels "red"/"blue" should use different neural networks for each color component. By authors of a method the best result were shown when three networks was used, one for each channel. Thus the total of neural networks makes *\* D D1* for each channel and *3( D 1)* totally [12]. The given method does not depend on color channels used for identification, on each channel the independent decision which is pseudo-independent, as channels are mutually correlated as shown in [13-14]. Accuracy of identification has been checked by authors on learning and test samplings on 100 images [13]. Accuracy of recognition of 7 algorithms of the interpolation was 100 % (errors of the first and second type are equal to zero). Accuracy of classification by offered methods oȠ real photocameras made 95-100 % depending on a

*\**

training of neural network as *L 12 L { , ,..., }* E E E and their average according:

*1 i,j N 1 <sup>1</sup> <sup>A</sup> A(i, j) (N 1)* d d ¦ .

Where *<sup>2</sup>*

The covariace matrix will be:

accordingly.

photocamera. The method also can be used for check on authenticity of pictures since in those areas at which there are signs of editing the responses of neural network increase by 2- 3 times. However, usage of replaced fragments of the image from the same photocamera as a background makes detection impossible [15]. Mehdi Kharrazi, etc. in [16] proposed the method of photocameras identification on the basis of image features. The task of determination of the camera with which help the analyzable picture has been received was thus considered. Proceeding from known sequence of information processing from a photosensitive matrix, it is possible to select two stages, importing the most essential distortions: a stage of debayerization, i.e. full-color image restoration and a postprocessing stage. Totally authors select 34 signs of classification, among them:


$$E\_1 = \frac{\left\|G\right\|^2}{\left\|B\right\|^2} \; \; \; E\_1 = \frac{\left\|G\right\|^2}{\left\|R\right\|^2} \; \; \; \; E\_1 = \frac{\left\|B\right\|^2}{\left\|R\right\|^2} \; \; \; \; \; $$


Along with enumerated features metrics of image quality proposed in [16] has been used. All used metrics can be divided into following groups:


To classify vectors the SVM-based classifier has been used. At learning stage 120 of 300 images were used, with 180 at test stage. An average accuracy of camera identification in "1 out of 2" were 98,73% with 88,02% when images were regular photos.

In [17] an identification method based on proprietary interpolation algorithms used in camera. The basis of algorithm is pixel correlation estimation listed in [18] with two estimations: estimation of pixel value by adjacent pixels' values and demosaic kernel used for raw data processing. As precise configuration of area used for interpolation is uknown, several different configurations were used, with additional assumbtion about different interpolation algorithms used in gradient and texturized areas. Camera identification experiments were done on a basis of two cameras: Sony DSC-P51 Ȗ Nikon E-2100. It has been acknowledged that filter kernel increase leads to accuracy increase (for kernels 3x3, 4x4, 5x5, accuracies were from 89.3 to 95.7%).

#### **2.2 Camera identification based on matrix defects**

Camera identification based on postprocessing algorithms features posess several disadvantages, the most fundamental is impossibility of practical use for one-model camera identification, even in "1 ot of 2" case.

In [19] camera identification method based on defective ("hot" and "dead" pixels) are presented but its effectiveness is limited for cameras without build-in pixel defects correction and "dark frame" subtraction.

Digital Camera Identification Based on Original Images 311

where *P* is a non-linear function of pixel's value, its coordinates and its neighborhood

Structure noise can be suppressed by subtracting additive noise *ij c* , then dividing pixels

*ij ij ij ij x' ( y c )/ f '*,

where *ij x '* is a corrected pixels value, *ij f '* is an approximation of *ij f* by averaging multiple

*(k) ij k ij (k) ij i,j,k*

¦ .

*f f ' mn f*

¦

This operation cannot be done on *ij p* , only over raw data from photosensitive matrix *ij y*

To get better understanding influence of structural noise onto resulted images and

Using ambient light 118 images were made on Canon camera with automatic exposure and

All obtained images possessed pronounced brightness gradient (vignetting). To eliminate that low-frequency distortion the HF-filter with cutoff frequency at (150/1136) S . Then images were averaged thus random noise was suppressed and structural noise summed. Spectrum of the signal resembles white spectrum with decrease of HF-components area, which is explainable as consequences of color interpolation over pixel neighborhood. PNUnoises are not presented in saturated and completely dark areas where FPN prevails. Owing to noise-like of the PNU-components of matrix noise, it is natural to use correlation method

In the absence of access in consumer-grade cameras to sensors output *ij y* , usually it is impossible to extract PNU from gray-frame. However it is possible to approximate noise by averaging multiple images *p(k) k= 1,…,Np.* Process speed-up is performed by filtering and

)( *kk* )()( *<sup>k</sup>* )( *pFpn* .

Other advantage of operation with residual noise that low-frequency component of PRNU is automatically suppressed. It is obvious that, the more the number of images (N> 50), the less influence of the single source image will take place. Originally, the filter based on wavelet

*N( ij y )* .

value by normalized frame's values:

prior successive image processing. Properties of pixel non-uniformity noise

*ij f* , *k 1..K* :

determine its characteristics in the following experiments were done:

**2.3 Identification based on non-uniformity of pixels sensitivity** 

transform was used. So advantages of this method are: - No access to internals of camera is required;


focused on infinity. White balance was set to create neutral gray images.

flat-exposure frames *(k)*

for its detection [16].

averaging of residual noise *n(k)*:

Camera identification based on dark frame correction along with obvious advantage of identification of concrete camera sample inherent critical disadvantage namely requirement of dark frames to identify cameras, which makes this method nearly completely useless in practical sense.

In [20] the camera identification method based on non-uniformity of pixel matrix namely different photosensitivity of pixels.

There are many sources of defects and noises which are generated at different image processing stages. Even if sensors form several images of absolutely static scene, the resulted digital representations may posess insignificant alterations of intensity between "same pixel" of image. It appears partly from shot noise [14,15] which is random, and partially because of structure non-uniformity, which is deterministic and slowly changed across even very large sets of image for similar conditions.

Structural non-uniformity presented in every image and can be used for camera identification. Due to similarity of non-uniformity's nature and random noise it is frequently named structural noise.

By averaging multiple images context impact is reduced and structural noises are separated structural matrix noise can be viewed as two components — fixed pattern noise (FPN) and photo-response non-uniform noise (PRNU). Fixed pattern noise is induced by dark currents and defined primarily by pixels non-uniformity in absence of light on sensitive sensor area. Due to additive nature of FPN, modern digital cameras suppress it automatically by subtracting the dark frame from every image [14]. FPN depends on matrix temperature and time of exposure. Natural images primary structural noise component is PRNU. It is caused by pixels non-uniformity (PNU), primarily non-uniform photosensitivity due to nonhomogeneity of silicon wafers and random fluctuations in sensor manufacturing process. Source and character of noise induced by pixels non-uniformities make correlation of noise extracted from two even identical matrixes small. Also temperature and humidity don't render influence to PNU-noise. Light refraction on dust particles and optical system also also induced its contribution to PRNU-noise, but these effects are not stable (dust can migrate over the matrix surface, vignette type changes with focal length or lens change) hence, can't be used for reliable identification.

The model of image obtaining process is the following. Let absolute photon number on pixel's area with coordinates *(i,j)* corresponds *ij x* , where *i 1..m, j 1..n* , *m n* u photosensitive matrix resolution. If we designate shooting noise as K*ij* , additive noise due to reading and other noises as *ij* H , dark currents as *ij c* . Then sensor's output *ij y* is:

$$\{y\_{i\dot{j}} = f\_{i\dot{j}}(\mathfrak{x}\_{i\dot{j}} + \mathfrak{η}\_{i\dot{j}}) + \mathfrak{c}\_{i\dot{j}} + \mathfrak{c}\_{i\dot{j}}\,. \tag{\*}$$

Here *ij f* is almost 1 and is multiplicative PRNU-noise.

Final image pixels *ij p* are completely formed after multiple-stage processing of *ij y* including, interpolation over adjacent pixels, color correction and image filtering. Many of that operations are non-linear like gamma correction white balance estimation, adaptive Bayer structure interpolation based on strategies for missing color recoveries. So:

$$p\_{\vec{\imath}\vec{\jmath}} = P(y\_{\vec{\imath}\vec{\jmath}\prime}N(y\_{\vec{\imath}\vec{\jmath}}), i, j) \; .$$

310 Applications of Digital Signal Processing

Camera identification based on dark frame correction along with obvious advantage of identification of concrete camera sample inherent critical disadvantage namely requirement of dark frames to identify cameras, which makes this method nearly completely useless in

In [20] the camera identification method based on non-uniformity of pixel matrix namely

There are many sources of defects and noises which are generated at different image processing stages. Even if sensors form several images of absolutely static scene, the resulted digital representations may posess insignificant alterations of intensity between "same pixel" of image. It appears partly from shot noise [14,15] which is random, and partially because of structure non-uniformity, which is deterministic and slowly changed across even

Structural non-uniformity presented in every image and can be used for camera identification. Due to similarity of non-uniformity's nature and random noise it is frequently

By averaging multiple images context impact is reduced and structural noises are separated structural matrix noise can be viewed as two components — fixed pattern noise (FPN) and photo-response non-uniform noise (PRNU). Fixed pattern noise is induced by dark currents and defined primarily by pixels non-uniformity in absence of light on sensitive sensor area. Due to additive nature of FPN, modern digital cameras suppress it automatically by subtracting the dark frame from every image [14]. FPN depends on matrix temperature and time of exposure. Natural images primary structural noise component is PRNU. It is caused by pixels non-uniformity (PNU), primarily non-uniform photosensitivity due to nonhomogeneity of silicon wafers and random fluctuations in sensor manufacturing process. Source and character of noise induced by pixels non-uniformities make correlation of noise extracted from two even identical matrixes small. Also temperature and humidity don't render influence to PNU-noise. Light refraction on dust particles and optical system also also induced its contribution to PRNU-noise, but these effects are not stable (dust can migrate over the matrix surface, vignette type changes with focal length or lens change)

The model of image obtaining process is the following. Let absolute photon number on pixel's area with coordinates *(i,j)* corresponds *ij x* , where *i 1..m, j 1..n* , *m n* u photosensitive matrix resolution. If we designate shooting noise as K*ij* , additive noise due to

Final image pixels *ij p* are completely formed after multiple-stage processing of *ij y* including, interpolation over adjacent pixels, color correction and image filtering. Many of that operations are non-linear like gamma correction white balance estimation, adaptive

*ij ij ij p P( y ,N(y ),i, j)* ,

*ij ij ij ij ij ij y f (x ) c* K H . (\*)

reading and other noises as *ij* H , dark currents as *ij c* . Then sensor's output *ij y* is:

Bayer structure interpolation based on strategies for missing color recoveries. So:

practical sense.

named structural noise.

different photosensitivity of pixels.

very large sets of image for similar conditions.

hence, can't be used for reliable identification.

Here *ij f* is almost 1 and is multiplicative PRNU-noise.

where *P* is a non-linear function of pixel's value, its coordinates and its neighborhood *N( ij y )* .

Structure noise can be suppressed by subtracting additive noise *ij c* , then dividing pixels value by normalized frame's values:

$$\propto\_{\vec{\eta}} = (y\_{\vec{\eta}} - c\_{\vec{\eta}}) / \langle f\_{\vec{\eta}} \rangle\_{\vec{\nu}}$$

where *ij x '* is a corrected pixels value, *ij f '* is an approximation of *ij f* by averaging multiple flat-exposure frames *(k) ij f* , *k 1..K* :

$$f\_{\vec{\imath}\vec{\jmath}} = \frac{\sum f\_{\vec{\imath}\vec{\jmath}}^{(k)}}{\sum f\_{\vec{\imath}\vec{\jmath}}^{(k)}} m \cdot n \; .$$

This operation cannot be done on *ij p* , only over raw data from photosensitive matrix *ij y*

prior successive image processing.

Properties of pixel non-uniformity noise

To get better understanding influence of structural noise onto resulted images and determine its characteristics in the following experiments were done:

Using ambient light 118 images were made on Canon camera with automatic exposure and focused on infinity. White balance was set to create neutral gray images.

All obtained images possessed pronounced brightness gradient (vignetting). To eliminate that low-frequency distortion the HF-filter with cutoff frequency at (150/1136) S . Then images were averaged thus random noise was suppressed and structural noise summed. Spectrum of the signal resembles white spectrum with decrease of HF-components area, which is explainable as consequences of color interpolation over pixel neighborhood. PNUnoises are not presented in saturated and completely dark areas where FPN prevails. Owing to noise-like of the PNU-components of matrix noise, it is natural to use correlation method for its detection [16].

#### **2.3 Identification based on non-uniformity of pixels sensitivity**

In the absence of access in consumer-grade cameras to sensors output *ij y* , usually it is impossible to extract PNU from gray-frame. However it is possible to approximate noise by averaging multiple images *p(k) k= 1,…,Np.* Process speed-up is performed by filtering and averaging of residual noise *n(k)*:

$$m^{(k)} = p^{(k)} - F(p^{(k)})\_{\cdot}$$

Other advantage of operation with residual noise that low-frequency component of PRNU is automatically suppressed. It is obvious that, the more the number of images (N> 50), the less influence of the single source image will take place. Originally, the filter based on wavelet transform was used. So advantages of this method are:


Digital Camera Identification Based on Original Images 313

*N*

The best results were achieved by *5 5* u mask. It has been shown that Wiener filter provides better separation, comparing with wavelet transform filter in [21]. The offered identification technique has been researched for identification possibility of 13 cameras [22-27], each with 100 images. Images from every camera were divided into 2 sets – training set used for camera fingerprinting and the test set, used for identity check [21]. The central crop of an

To process an image *I* for fingerprint creation or identification the color-to-grayscale conversion has been applied. Fingerprint is an averaged sum of all HF-components, forming *Wapprx* value. To check identity of an image *qI* , the correlation coefficient is evaluated:

( ( ), ) *q apprx p cc F I W* ,

**1 2 3 4 5 6 7 8 9 10 11 12 13 1 0.0908** 0.0010 0.0042 0.0004 0.0037 0.0026 0.0028 0.0020 0.0040 0.0033 0.0031 0.0032 0.0035 0.0006 **0.1494** -0.0001 0.0015 0.0005 0.0006 -0.0001 0.0000 0.0005 0.0003 0.0004 0.0001 0.0007 0.0028 -0.0001 **0.1364** 0.0003 0.0018 0.0004 0.0017 0.0020 0.0030 0.0023 0.0017 0.0024 0.0016 -0.0001 0.0009 -0.0004 **0.1889** 0.0002 -0.0007 -0.0004 -0.0009 -0.0000 -0.0003 -0.0007 -0.0001 -0.0005 0.0054 0.0019 0.0035 0.0000 **0.0727** 0.0025 0.0042 0.0024 0.0044 0.0022 0.0033 0.0029 0.0038 0.0022 0.0004 0.0004 0.0000 0.0015 **0.1423** 0.0006 0.0001 0.0016 0.0009 0.0013 0.0007 0.0036 0.0010 -0.0007 0.0005 -0.0001 0.0010 0.0001 **0.2645** -0.0012 0.0010 0.0008 0.0012 0.0017 0.0014 0.0003 0.0008 0.0009 -0.0001 0.0004 -0.0001 -0.0005 **0.7079** 0.0000 0.0012 -0.0001 0.0009 -0.0005 0.0049 0.0004 0.0046 0.0001 0.0031 0.0018 0.0027 0.0029 **0.1038** 0.0017 0.0033 0.0036 0.0019 0.0027 0.0013 0.0031 0.0002 0.0023 0.0015 0.0017 0.0032 0.0025 **0.1005** 0.0090 0.0030 0.0014 0.0011 0.0007 0.0020 -0.0010 0.0013 0.0004 0.0012 -0.0002 0.0006 0.0013 **0.3776** 0.0006 0.0007 0.0025 0.0002 0.0023 -0.0006 0.0018 0.0004 0.0023 0.0016 0.0028 0.0016 0.0020 **0.1294** 0.0016 0.0015 0.0001 0.0010 -0.0006 0.0017 0.0022 0.0009 -0.0006 0.0008 0.0009 0.0009 0.0003 **0.2747** 

On intersection of columns and lines with identical indexes there are correlation coefficients of images and a fingerprint, received by the same camera. Thus, at matrix coincidence,

Photosensitive matrix of a modern digital camera naturally possesses non-uniformity of its elements, both photosensitive and signal amplifiers. As the charge is transferred by columns, the well-known phenomena called banding occurs, resulting high-frequency noise. After image reconstruction process [3] and subjective quality improvements completing, the resulted image is compressed, usually according to JPEG standard, which introduces

*apprx*

*W*

image with 1024x1024 pixels size was used for identification purposes.

where *p* — is a correlation coefficient, and *cc* — cross correlation.

Table 1. An averaged correlation coefficients for 13 cameras.

**2.6 Image rotation detection based on Radon transform** 

correlation value is 0.1 - 0.7 and for incoincident cameras is 0.001 - 0.054.

blocking effect, and contributes regular pattern to rows and columns as well.

where *F(I )* is a filter operation.

( ( ))

*I FI*

¦ ,

*N*

#### **2.4 Detection based on correlation coefficient**

To detect image *Ȟ* belonging to specific camera *ǿ* it is possible to calculate correlation ǒǿ between residual and structural noise *n=p-F(p)* for the camera:

$$\rho\_C = \operatorname{corr}(n, P\_C) = \frac{(n - \overline{n})(P\_C - \overline{P}\_C)}{\|n - \overline{n}\| \|P\_C - \overline{P}\_C\|} \cdot \frac{1}{n}$$

Now it is possible to define distribution *ǒǿ (q)* for different images *q* made by the camera *C* and distribution *ǒC (q ')* for images *q '* made not by the camera *C. Based on* Neumann-Pirson approach and minimizing error rate the reached accuracy of classification made from 78 % to 95 % on 320 images from 9 digital cameras.

#### **2.5 Identification technique of digital recording devices based on correlation of digital images**

For development of a technique of identification of photocameras under images it is necessary to consider architecture of prospective system of identification. The system includes units:


An input format for identification system should be lossless format like full-color BMP to which all images and video streams are convertible. Typical output formats of modern cameras are JPEG and TIFF. In the feature vector former, digital image is converted to the feature vector represents an image for identification an storage purposes.

In the unit of device identification the estimation of likeness of two or more vectors is estimated allowing to accept or reject device similarity hypothesis.

#### **2.5.1 Feature vector forming for digital cameras identification**

Feature vector former is based on photosensitive matrix identification techniques, namely PRNU-features. As there will always be both signal and noise (PRNU-components and image context and (or) other noises) it is preferable to use filters to increase signal-noise ratio. To select HF-components, which represent PRNU can be done by Wiener filtering:

$$\begin{aligned} \mu &= \frac{1}{NM} \sum\_{n\_1, n\_2 \in \eta} a(n\_1, n\_2) \\ \sigma^2 &= \frac{1}{NM} \sum\_{n\_1, n\_2 \in \eta} a^2(n\_1, n\_2) - \mu^2 \\ b(n\_1, n\_2) &= \mu + \frac{\sigma^2 - \nu^2}{\sigma^2} (a(n\_1, n\_2) - \mu)\_{\mu} \end{aligned}$$

where *N* and *M* are number of pixels of neighborhood by y and x axis respectively. *1 2 a(n ,n )* — is a value of pixel with *1 2 (n ,n )* coordinates. Thus averaged values for specific matrix is:

$$W\_{\text{approx}} = \frac{\sum\_{N} (I - F(I))}{N} \prime$$

where *F(I )* is a filter operation.

312 Applications of Digital Signal Processing

To detect image *Ȟ* belonging to specific camera *ǿ* it is possible to calculate correlation ǒǿ

*<sup>C</sup> <sup>C</sup> PPnn PPnn Pncorr* ))((

Now it is possible to define distribution *ǒǿ (q)* for different images *q* made by the camera *C* and distribution *ǒC (q ')* for images *q '* made not by the camera *C. Based on* Neumann-Pirson approach and minimizing error rate the reached accuracy of classification made from 78 %

**2.5 Identification technique of digital recording devices based on correlation of digital** 

For development of a technique of identification of photocameras under images it is necessary to consider architecture of prospective system of identification. The system

An input format for identification system should be lossless format like full-color BMP to which all images and video streams are convertible. Typical output formats of modern cameras are JPEG and TIFF. In the feature vector former, digital image is converted to the

In the unit of device identification the estimation of likeness of two or more vectors is

Feature vector former is based on photosensitive matrix identification techniques, namely PRNU-features. As there will always be both signal and noise (PRNU-components and image context and (or) other noises) it is preferable to use filters to increase signal-noise ratio. To select HF-components, which represent PRNU can be done by Wiener filtering:

1 2

*an n*

2 2 1 2 2 1 2

*a nn*

( , ) (( , ) )

V

1 2

 P

P

,

feature vector represents an image for identification an storage purposes.

1 2

¦

,

*n n*

*NM*

P

V

*1 2 a(n ,n )* — is a value of pixel with *1 2 (n ,n )* coordinates.

Thus averaged values for specific matrix is:

*NM*

1 2

,

*n n*

P

2 22

*bn n an n*

where *N* and *M* are number of pixels of neighborhood by y and x axis respectively.

V Q

K

<sup>1</sup> (, )

<sup>1</sup> (, )

K

¦

estimated allowing to accept or reject device similarity hypothesis.

**2.5.1 Feature vector forming for digital cameras identification** 

),( .

*CC CC*

**2.4 Detection based on correlation coefficient** 

to 95 % on 320 images from 9 digital cameras.

**images** 

includes units:






between residual and structural noise *n=p-F(p)* for the camera:

U

The best results were achieved by *5 5* u mask. It has been shown that Wiener filter provides better separation, comparing with wavelet transform filter in [21]. The offered identification technique has been researched for identification possibility of 13 cameras [22-27], each with 100 images. Images from every camera were divided into 2 sets – training set used for camera fingerprinting and the test set, used for identity check [21]. The central crop of an image with 1024x1024 pixels size was used for identification purposes.

To process an image *I* for fingerprint creation or identification the color-to-grayscale conversion has been applied. Fingerprint is an averaged sum of all HF-components, forming *Wapprx* value. To check identity of an image *qI* , the correlation coefficient is evaluated:

$$p = \text{cc}(F(I\_q), W\_{\text{approx}})\_{\text{rev}}$$

where *p* — is a correlation coefficient, and *cc* — cross correlation.


Table 1. An averaged correlation coefficients for 13 cameras.

On intersection of columns and lines with identical indexes there are correlation coefficients of images and a fingerprint, received by the same camera. Thus, at matrix coincidence, correlation value is 0.1 - 0.7 and for incoincident cameras is 0.001 - 0.054.

#### **2.6 Image rotation detection based on Radon transform**

Photosensitive matrix of a modern digital camera naturally possesses non-uniformity of its elements, both photosensitive and signal amplifiers. As the charge is transferred by columns, the well-known phenomena called banding occurs, resulting high-frequency noise. After image reconstruction process [3] and subjective quality improvements completing, the resulted image is compressed, usually according to JPEG standard, which introduces blocking effect, and contributes regular pattern to rows and columns as well.

Digital Camera Identification Based on Original Images 315

In drawing there are maxima at values of rotating angle with added 90º multiples.. At transition from a corner 89º to a corner 90º occurrence of maxima in a peak spectrum is observed. Similar change of character of a peak spectrum gives a possibility to establish

(Ȏ) (b) Fig. 4. Spectrums of the projections corresponding to angles 89º and 90º (a-b) for 1024x1024

Result (an average of a spectrum for the Radon-transformed image for projection angles from 0**º** to 360**º** with 10**º** step) is presented in figure 5 and a dissection of it — an average of Radon projection spectrograms for image fragment 1024x1024 pixels in figure 5. Local maximums at 10**º**, 100**º**, 190**º**, 280**º** in figure 6 correspond diagonally-placed maximums in figure 6. To determine the influence of image size change (resize operation), defined as the relation of the linear sizes of the initial image to resulted one on possibility of rotation detection by Radon transform, different scales of original image has also been investigated.

Fig. 5. An average of a spectrum of capacity of transformation of Radon (corners from 0º to

360º) at corners of turn from 0º to 360º with step 10º

value of a image rotation degree.

pixels image fragment

This phenomena can be used to detect angle of rotation. To detect an angle of image rotation the Radon transform of its fragment could be performed with successive analysis of Radon projections with Fourier transform. Radon transform can be defined as follows: Let function *f(x, y)*, is defined in *D*. We will consider some straight line *L* on a plane *ȣȡ*, crossing area *D*. Then, integrating function *f (x, y)* along line *L*, we receive a projection or linear integral of function *f*. Integration along all possible lines *L* on a plane allows to define Radon transform:

$$f^\* = Rf = \int\_L f(x, y)ds,$$

where *ds* - an increment of length along *L*.

For minimization of edge effects impact of analyzed area on high-frequency part of an image it is advisable to apply Radon transform over circular fragment with smoothed borders. Selection of a fragment from the image and smoothing of its borders were done [4] by normalized two-dimensional Gaussian window

$$h(t) = \frac{\exp\left(\frac{-t^2}{2\delta^2}\right)}{\sqrt{2\pi}\cdot\delta},$$

$$\delta = \frac{\sqrt{\ln(2)}}{2\pi BT}$$

shown in figure 3.

Refinement to an angle which 90º degrees multiple is possible to make due to uncompensated banding traces, which are consequences of non-uniformity of image brightness component obtained from CCD or CMOS matrixes [2] and traces of compression artifacts. A consequence of the given phenomenon will be unequal level of maxima of a Fourier spectrum obtained from result of Radon transform that allows to select only 2 or (in some cases) 4 angles. Examples of columns spectrograms for a matrix of Radon transformed image fragment 1024x1024 pixel size are represented in figure 4.

Fig. 3. Two-dimensional normalized Gaussian window used to select an image fragment

314 Applications of Digital Signal Processing

This phenomena can be used to detect angle of rotation. To detect an angle of image rotation the Radon transform of its fragment could be performed with successive analysis of Radon projections with Fourier transform. Radon transform can be defined as follows: Let function *f(x, y)*, is defined in *D*. We will consider some straight line *L* on a plane *ȣȡ*, crossing area *D*. Then, integrating function *f (x, y)* along line *L*, we receive a projection or linear integral of function *f*. Integration along all possible lines *L* on a plane allows to define Radon transform:

> *L <sup>f</sup> Rf f(x,y)ds,* ³

For minimization of edge effects impact of analyzed area on high-frequency part of an image it is advisable to apply Radon transform over circular fragment with smoothed borders. Selection of a fragment from the image and smoothing of its borders were done [4]

<sup>2</sup> exp <sup>2</sup> ( ) <sup>2</sup>

*h t*

G

image fragment 1024x1024 pixel size are represented in figure 4.

2

*t*

§ · ¨ ¸ © ¹

G

S G

*ln 2 2 BT*

S

Refinement to an angle which 90º degrees multiple is possible to make due to uncompensated banding traces, which are consequences of non-uniformity of image brightness component obtained from CCD or CMOS matrixes [2] and traces of compression artifacts. A consequence of the given phenomenon will be unequal level of maxima of a Fourier spectrum obtained from result of Radon transform that allows to select only 2 or (in some cases) 4 angles. Examples of columns spectrograms for a matrix of Radon transformed

Fig. 3. Two-dimensional normalized Gaussian window used to select an image fragment

*\**

where *ds* - an increment of length along *L*.

shown in figure 3.

by normalized two-dimensional Gaussian window

In drawing there are maxima at values of rotating angle with added 90º multiples.. At transition from a corner 89º to a corner 90º occurrence of maxima in a peak spectrum is observed. Similar change of character of a peak spectrum gives a possibility to establish value of a image rotation degree.

Fig. 4. Spectrums of the projections corresponding to angles 89º and 90º (a-b) for 1024x1024 pixels image fragment

Result (an average of a spectrum for the Radon-transformed image for projection angles from 0**º** to 360**º** with 10**º** step) is presented in figure 5 and a dissection of it — an average of Radon projection spectrograms for image fragment 1024x1024 pixels in figure 5. Local maximums at 10**º**, 100**º**, 190**º**, 280**º** in figure 6 correspond diagonally-placed maximums in figure 6. To determine the influence of image size change (resize operation), defined as the relation of the linear sizes of the initial image to resulted one on possibility of rotation detection by Radon transform, different scales of original image has also been investigated.

Fig. 5. An average of a spectrum of capacity of transformation of Radon (corners from 0º to 360º) at corners of turn from 0º to 360º with step 10º

Digital Camera Identification Based on Original Images 317

The methods of digital cameras identification allow defining the fact of origin of digital images from the specific camera. In comparison with artificial digital watermarks embedded by either the special software, or device modification, identification based on innate difference of every single camera allows to identify cameras by the analysis of statistical regularities in digital media. Explicit advantages of the given identification methods are their applicability to images of consumer cameras without necessity of internals access or camera firmware modification. Methods of cameras identification on the basis of processing differences allow to identify cameras vendors. Methods of identification based on nonuniformities of record path allow to identify separate copies of one model of camera. Essential hindrance for correct identification of cameras is scaling and rotation of the images which are exposed to identification process. To ascertain the fact of rotation and its reversing the Radon transform with the subsequent projections processing by Fourier

[1] Grubunin V. G., Okov I. N., Turintsev I. V. Tsifrovaya steganografiya. – M.: Solon-Press,

[2] Osborne C., van Schyndel R., Tirkel A. A digital watermark. – IEEE International

[3] Ramkumar M. Data Hiding in Multimedia. PhD Thesis. – New Jersey Institute of

[4] Hartung F., Su J., Girod B. Spread Spectrum Watermarking: Malicious Attacks and

[5] Petitcolas F., Anderson R., Kuhn M. Attacks on Copyright Marking Systems. – Lecture

[6] Langelaar G., Lagendijk R., Biemond J. Removing spatial spread spectrum watermarks

[8] D. R. Cok, Signal processing method and apparatus for producing interpolated

[9] W. T. Freeman, Median filter for reconstructing missing color samples, U.S. Patent, No.

[10] C. A. Laroche and M. A. Prescott, Apparatus and method for adaptively interpolating a full color image utilizing chrominance gradients, U.S. Patent, No. 5,373,322, 1994. [11] Yangjing Long Yizhen Huang Image Based Source Camera Identification using

[12] B. K. Gunturk, Y. Altunbasak and R. M. Mersereau, Color plane interpolation using

[13] L. Lam and C.Y. Suen, Application of majority voting to pattern recognition: An

alternating projections, IEEE Transactions on Image Processing, 11(9):997-1013,

analysis of its behavior and performance, IEEE Transactions on Systems, Man and

chrominance values in a sampled color image signal, U.S.Patent, No. 4,642,678,

Conference on Image Processing, 1994. – P. 86-90

Notes in Computer Science, 1998. – P. 218-238.

by non-linear filtering. – Proceedings EUSIPCO-98, 1998. [7] B. E. Bayer, Color imaging array, U.S. Patent, No. 3,971,065, 1976

**3. Conclusion** 

transform can be used.

2002. – 272 p.

Counterattacks.

4,724,395, 1988.

Demosaicking

Cybernetics, 27(5):553–568, 1997.

1986.

2002.

Technology,1999. – 72 p.

**4. References** 

Fig. 6. An average of the spectrogram of projections of Radon of a fragment of the image dimension 1024x1024 pixel for corners from 0 to 360º

In figure 7 the two-dimensional dependence diagram of a magnitude average calculated on a set of image Radon transforms is shown where peak of normalized averaged spectrum located at 10**º** of Radon transform for angles from [5**º**..15**º**], applied to the image with an initial rotation angle of 10º and image scaling factor varied from 1 to 0,1 by 0,1 step. Even at 0.2 scale coefficient the maximum, which corresponds to correct rotation angle, is visible, so rotation operation can be undone. In figure 6 values of spectrogram average of Radon transform (rotation angles from (0**º**..20**º**)) for 80 images obtained from one camera and rotated by 10**º** (a) with histograms (b) are shown.

Fig. 7. Dependence of a normalized mean value of averaged spectrum on scaling factor and angle of a Radon projection
