# Registration and Sequencing of Vessels Section Images at Macroscopic Levels

Aneta Górniak, Ewa Skubalska-Rafajlowicz

► To cite this version:

Aneta Górniak, Ewa Skubalska-Rafajlowicz. Registration and Sequencing of Vessels Section Images at Macroscopic Levels. 14th Computer Information Systems and Industrial Management (CISIM), Sep 2015, Warsaw, Poland. pp.399-410, 10.1007/978-3-319-24369-6_33. hal-01444483

HAL Id: hal-01444483

https://inria.hal.science/hal-01444483

Submitted on 24 Jan 2017

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Registration and sequencing of vessels section images at macroscopic levels

Aneta Górniak¹ and Ewa Skubalska-Rafajłowicz¹

Wroclaw University of Technology
Faculty of Electronics, Department of Computer Engineering
aneta.gorniak@pwr.edu.pl, ewa.rafajlowicz@pwr.edu.pl

**Abstract.** In this paper we present a new approach to registration and sequencing of microscope images obtained from serial sections of large blood vessels at macroscopic levels of magnification. It is assumed that subsequent section images may be located inadequately in the image series. Translations and rotations of the object of interest can occur. Some images can be also reflected vertically or horizontally. The proposed algorithm is based on the center of gravity estimation and the phase-only correlation (POC) and uses standard image normalization as a preprocessing procedure. The method is fully-automatic and robust to common image distortions. The quality of registration is measured by the mean value of the sum of absolute difference between images. This criterion can be used also for slide images sequencing, when the image acquisition is performed independently for each section. A set of experiments, carried out using sampled microscopic images of a vein section, proves experimentally the effectiveness of the proposed approach.

**Keywords:** image processing, microscopic image, rigid image registration, image normalization, object matching, phase-only correlation

## 1 Introduction

Medical images are widely used within healthcare for diagnosis, planning treatment, guiding treatment and monitoring disease progression. In many cases, multiple images are acquired from subjects at different times and with different imaging modalities. There is also a notion of comparing images obtained from patient cohorts in place of single subjects captured multiple times [6].

A medical image can refer to various types of images that possess very different underlying physical principles and very different applications. It covers images varying from microscopic images of histological sections to video images used for remote consultation. It becomes crucial to find ways of accurately aligning the information in the different images and providing tools for visualizing the combined images.

Image registration is a very important step in processing sequences of images [6], [11], [16]. A large number of algorithms have been developed to perform registration of medical and biological images. Image registration is the process
of aligning multiple images representing the same scene that were captured at
a different time, at or by a different set of modalities [6]. The process involves the
transformation of the coordinate systems of the reference image and the input
image into the joint coordinate system. It must ensure that the correspond-
ing structures are precisely located at the same positions in the images to be
matched [3]. One type of approach in image registration consists of estimating
the geometric transformations of translation, rotation, scaling or perspective be-
tween the reference and input images, based on the pixel intensity. The process
is iterative and aims to optimize the elected measure of similarity between the
images. Other approaches may involve the identification of corresponding points
or areas in the images [5], [11], [16].

The most common problems of medical image registration cover the varying
deformations and distortions of the imagined tissues, the alignment of images
with different dimensionality, the alignment of images from different subjects or
the geometric distortion formed in the imaging process.

The other set of problems is finding the correct ways of assessing the reg-
istration accuracy. The required accuracy may vary between applications. The
most promising approach in ensuring acceptable accuracy is visual assessment
of the registered images before they are used [6], [9].

Image registration has application in many fields. It is used in remote sens-
ing for multi-spectral classification, environment monitoring, change detection
or weather forecasting. In medicine, it is used in signal fusion from multiple
sources to acquire more accurate information about the patient, to monitor
tumor growth or to compare patient data with anatomical atlases [6]. In [4],
a method of MRI image registration with sub-pixel accuracy and its application
in fMRI is presented. Other applications in biomedical image registration are
presented in [2], [10].

In this paper we present a new approach to the registration and sequencing
of microscope images [9], [13], [16] obtained from serial sections of large blood
vessels at macroscopic levels of magnification. It is assumed that subsequent
section images may be located inadequately in the image series. Furthermore,
not only translations and rotations of the object of interest can occur. Some
images can be also reflected vertically or horizontally. In general, the process
of acquiring the data is imperfect and leads to geometric and pixel intensity-
based distortions in the images. Due to the macroscopic levels of magnification
(less than 100) we can restrict the rigid registration of images, since object
deformation obtained during sectioning is here negligible. The small differences
between geometry of the subsequent images should not be corrected, because
these differences are caused by the natural morphological structure changes and
their viability depends on section thickness. The proposed algorithm is based on
the center of gravity estimation [1] and the phase-only correlation (POC) [3], [8],
[15], [12] and it uses a standard image normalization as preprocessing procedure
[7]. The method is fully-automatic and robust to common image distortions.
The quality of registration is measured by the mean value of the sum of absolute
difference between images (23). This criterion can be used also for slide image sequencing, when image acquisition is performed independently for each section.

In the next Section the problem of microscopic section images registration is formulated. Section 3 provides a concise description of specialized methods used in the proposed registration procedure. Section 4 presents the detailed registration procedure adapted to the blood vessel section images obtained at macroscopic levels of magnification. Section 5 describes a set of experiments for evaluating the performance of the proposed algorithm. The experimental data consists of two sets of microscopic images containing vein sections. Finally, in Section 6 some brief conclusions are presented.

## 2 Image Registration

An image is a two-dimensional sampled function $f(x, y)$ with a discrete set of coordinates $x = 1, ..., N_1$ and $y = 1, ..., N_2$, where $N_1$ and $N_2$ denote the width and the height of the image. Color intensity of the pixels in the image takes up values from the interval $[0, 1]$. Images are presented in greyscale.

Images $f$ and $g$ represent the same type of an object, captured with the same or different modalities. The acquisition can result in the appearance of differences between acquired objects. We may classify the differences into the following categories: naturally occurring differences, noise-induced differences, geometric distortions and intensity distortions.

The first type of differences stems from the two acquired objects being naturally different from each other. It translates into the lack of correspondence between the pixels of the two objects and the difference in intensity values of those pixels. The difference is described with the value $\beta(x, y)$

$$ \beta(x, y) = |f(x, y) - g(x, y)|. \qquad (1) $$

If the objects in the images are non-identical, the sum of $\beta(x, y)$ values will be grater than zero. These differences should be retained, because they contain true information about subsequent changes of the morphological structure under consideration.

The other type of naturally occurring differences comes in the form of noise induced during the acquisition process. This type of difference is always present in the images.

The second group of differences, called distortions, corresponds to geometric displacement of the captured scene and disparity of intensity levels in the image.

The difference in intensity levels of the images results in the objects appearing in the images with varying intensity. Let us consider that each image's intensity changes can be described with parameters $\alpha$ and $\gamma$, which are called scaling factor and offset, respectively. The intensity transformation is denoted with $\mathcal{F}$ and describes the following dependency

$$ \mathcal{F}\{f(x, y)\} = \alpha f(x, y) + \gamma. \qquad (2) $$
The result of the intensity transformation $\mathcal{F}$ of the image $f(x, y)$ is the image $\alpha f(x, y) + \gamma$. These image differences can be easily removed by the image normalization procedure [7].

The geometric distortion in the acquired objects requires the geometric transformation to map an object from one image to another, because the corresponding points in the images do not match the corresponding points on the objects in these images. The adequate geometric transformation is denoted with $\mathcal{T}$ and maps both position and the intensity at that position from one image to another. This dependency is given by

$$f(x, y) = \mathcal{T}\{g(x, y)\}. \qquad (3)$$

We assume that the image $g(x, y)$ is a translated, rotated and possibly reflected (vertically or horizontally) replica of the image $f(x, y)$.

The difference between these two images can be described by an unknown displacement function $(d_1(x), d_2(y))$. So, we obtain

$$g(x, y) = f(x + d_1(x), y + d_2(y)) \qquad (4)$$

or

$$f(x, y) = g(x - d_1(x), y - d_2(y)).$$

The objective of the image registration is to estimate the affine displacement function $(d_1(x), d_2(y))$. We assume that the displacement is a result (sum) of geometric transformations of translation, rotation and reflection in a prescribed order. We also assume that the values of the displacement $(d_1(x), d_2(y))$ can be obtained by estimating the parameters of the respective composite geometric transformations.

1. Rotation of the image by the angle $\theta$ around the point $(x_0, y_0)$

$$f(x, y) = g(x_{\theta}, y_{\theta}) \qquad (5)$$

where $x_{\theta}$ and $y_{\theta}$ are given by

$$\begin{aligned} x_{\theta} &= (x - x_0) \cdot \cos \theta + (y - y_0) \cdot \sin \theta + x_0 \\ y_{\theta} &= -(x - x_0) \cdot \sin \theta + (y - y_0) \cdot \cos \theta + y_0. \end{aligned} \qquad (6)$$

2. Translation of the image by the translation vector $(\tau_1, \tau_2)$

$$f(x, y) = g(x - \tau_1, y - \tau_2). \qquad (7)$$

3. Reflection of the image. We consider two types of reflections: vertical and horizontal. Reflection about a vertical axis of abscissa $x_0$ is

$$f(x, y) = g(-x + 2x_0, y) \qquad (8)$$

and the reflection about a horizontal axis of ordinate $y_0$ is

$$f(x, y) = g(x, -y + 2y_0). \qquad (9)$$

The displacement $(d_1(x), d_2(y))$ is a sum of displacements generated by each of the transformations. The order by which the transformations are performed is significant and it may influence the accuracy of the results of displacement function estimation.
# 3 Used methods

We have used normalization of the intensity levels in the image as a preprocessing step.

## 3.1 Normalization of the intensity levels in the image

The purpose of the normalization of the intensity levels in the image is elimination of the intensity-based distortion so the images that differ only in intensity become identical [7]. Let $f_1$ and $f_2$ be two images which differ in the pixel intensity space (see (2)). For arbitrary values of $\alpha$ and $\gamma$, the following dependency holds

$$ \mathcal{N}\{\alpha_1 f + \gamma_1\} = \mathcal{N}\{\alpha_2 f + \gamma_2\} \qquad (10) $$

where $\alpha_1$ and $\alpha_2$ are the intensity levels of the images $f_1$ and $f_2$, and $\gamma_1$ and $\gamma_2$ are offsets of pixel intensities in both images, respectively. $\mathcal{N}$ denotes the normalization function of the intensity in the image.

To normalize the intensity level of the image, the influence of $\alpha$ and $\gamma$ coefficients needs to be eliminated. The standard normalization function [7] is given by

$$ \mathcal{N}\{f(x, y)\} = \frac{\alpha f(x, y) - \alpha\mu}{\alpha f_{max} - \alpha f_{min}} = \frac{\alpha(f(x, y) - \mu)}{\alpha(f_{max} - f_{min})} = \frac{f(x, y) - \mu}{f_{max} - f_{min}}, \quad (11) $$

where $f_{max}$ is a maximum intensity value and $f_{min}$ is a minimum intensity value in the image. The $\mu$ parameter is the mean value of intensity in the image

$$ \mu = \frac{1}{N_1 N_2} \sum_{x,y} f(x, y) \qquad (12) $$

where $\sum_{x,y}$ denotes $\sum_{x=1}^{N_1} \sum_{y=1}^{N_2}$.

## 3.2 The center of gravity of the object in the image

Calculating the position of the center of gravity in an image allows us for detection of the object and estimation of its position in the image. In practice it allows us to reduce the size of the processed image to the area where the object is localized.

We assume that a set of coordinates $(x, y)$ in the image is the center of the object to be detected. A detection algorithm, defined with a probability-based function $c(x, y)$, which results in a local maximum of $c(x, y)$ around the center of each detected object [1].

For an object in a greyscale image the center of gravity is weighted and is defined as

$$ (x_c, y_c) = \left( \frac{\sum_{x,y} x w(x, y)}{\sum_{x,y} w(x, y)}, \frac{\sum_{x,y} y w(x, y)}{\sum_{x,y} w(x, y)} \right) \qquad (13) $$
where $\sum_{x,y}$ describes $\sum_{x=1}^{N_1} \sum_{y=1}^{N_2}$. The weight function $w(x, y)$ is defined as

$$w(x, y) = a(f(x, y) - m) \qquad (14)$$

where $m < \min_{x,y}(f(x,y))$, when $a > 0$ and $m > \max_{x,y}(f(x,y))$, when $a < 0$. Parameter $a$ decides if $(x_c, y_c)$ drifts toward bright pixels ($a > 0$) or dark pixels ($a < 0$) in the image.

The results give a sub-pixel estimate of the central coordinates of an object. However, since we are not working with sub-pixel precision the results are rounded to pixels. Of course, sub-pixel accuracy can be used, when more subtle precision of image registration is necessary.

### 3.3 Rotation

We propose a simple, but relatively accurate, method of rotation estimation. The image under consideration is rotated around the center of gravity $(x_c, y_c)$ for different $\theta$ values. New pixel values are calculated

$$\begin{aligned} x_{\theta} &= (x - x_c) \cdot \cos \theta - (y - y_c) \cdot \sin \theta + x_c \\ y_{\theta} &= (x - x_c) \cdot \sin \theta + (y - y_c) \cdot \cos \theta + y_c \end{aligned} \qquad (15)$$

and used subsequently in the algorithm presented in detail in the next Section.

### 3.4 The Phase-Only Correlation (POC) function

Correlation between phase-only versions of the two images to be aligned is used for image registration. This matching technique uses the phase component of two-dimensional discrete Fourier Transform (2D-DFT) of the two images [3], [8]. Furthermore, the POC method can be also used for rotation estimation [12].

We assume that $f(x, y)$ is the reference image and $g(x, y)$ is the translated image of $f(x, y)$ by $(\tau_1, \tau_2)$, both of the size $N_1 \times N_2$. For mathematical simplicity we assume that the index ranges are $x = -M_1 \cdots M_1$ ($M_1 > 0$) and $y = -M_2 \cdots M_2$ ($M_2 > 0$) with $N_1 = 2M_1 + 1$ and $N_2 = 2M_2 + 1$ [15].

The 2D Discrete Fourier Transforms (2D DFTs) of the images $f(x, y)$ and $g(x, y)$, denoted $F(k, l)$ and $G(k, l)$, are defined as follows

$$F(k, l) = \sum_{x,y} f(x, y) W_{N_1}^{kx} W_{N_2}^{ly}, \qquad (16)$$

$$G(k, l) = \sum_{x,y} g(x, y) W_{N_1}^{kx} W_{N_2}^{ly}, \qquad (17)$$

where $k = -M_1 \cdots M_1$, $l = -M_2 \cdots M_2$, $W_{N_1} = e^{-j \frac{2\pi}{N_1}}$, $W_{N_2} = e^{-j \frac{2\pi}{N_2}}$, and operator $\sum_{x,y}$ describes $\sum_{x=-M_1}^{M_1} \sum_{y=-M_2}^{M_2}$. For the 2D DFTs of the images, the cross-spectrum $R(k, l)$ is given by

$$R(k, l) = F(k, l) \overline{G(k, l)}, \qquad (18)$$
where $\overline{G(k, l)}$ is the complex conjugate of $G(k, l)$. The cross-phase spectrum $\hat{R}(k, l)$ is a normalized cross spectrum $R(k, l)$ defined as

$$ \hat{R}(k, l) = \frac{F(k, l) \overline{G(k, l)}}{|F(k, l) \overline{G(k, l)}|} \qquad (19) $$

The Phase-Only Correlation function $\hat{r}(x, y)$ is the 2D Inverse Discrete Fourier Transform (2D IDFT) of $\hat{R}(k, l)$ and is given by

$$ \hat{r}(x, y) = \frac{1}{N_1 N_2} \sum_{k,l} \hat{R}(k, l) W_{N_1}^{-kx} W_{N_2}^{-ly} \qquad (20) $$

If the images $f(x, y)$ and $g(x, y)$ are the same, i.e., $f(x, y) = g(x, y)$, the POC function is defined as

$$ \hat{r}(x, y) = \frac{1}{N_1 N_2} \sum_{k,l} W_{N_1}^{-kx} W_{N_2}^{-ly} = \delta(x, y) = \begin{cases} 1 & \text{for } x = y = 0 \\ 0 & \text{for others} \end{cases} \qquad (21) $$

where $\delta(x, y)$ is the Kronecker's delta function. If the images are translated by some values $(\tau_1, \tau_2)$ as in (7) then $\hat{r}(x, y)$ is a translated Kronecker's delta function [17] given by

$$ \hat{r}(x, y) = \begin{cases} 1 & \text{for } (x, y) = (-\tau_1, -\tau_2) \\ 0 & \text{for others} \end{cases} \qquad (22) $$

The position of the maximal value of POC function $\hat{r}(x, y)$ identifies the translation vector $(-\tau_1, -\tau_2)$. It has an opposite sign to $(\tau_1, \tau_2)$. Of course, there are known many other methods of translation estimation (see for example [14] and references cited therein). We have used the POC method, since it is relatively robust to the image differences.

# 4 Algorithm Implementation

This section presents the proposed algorithm of object matching in the series of vessel section images with the use of the Phase-Only Correlation function, the center of gravity calibration and intensity levels correction. The algorithm matches the input image to the reference image by estimating the values of the parameters of the geometric transformations (5)-(9) forming the displacement between the images.

In the proposed algorithm we can list the following steps: estimation of the object position in the image, normalization of intensity in the image, rotation alignment, displacement alignment, reflection alignment and assessment of the match. The algorithm has a modular structure that allows for easy inclusion of new procedures or modification and removal of already applied procedures.
The main loop of the algorithm consists of the following steps:

1. Input the series of images

2. Normalize the images with the function (11).

3. **for each** pair of normalized images $f(x, y)$ and $g(x, y)$ in the series **do**:
   Calculate the new position of the image center $(x_c, y_c)$ for $f(x, y)$ and $g(x, y)$ using the center of gravity function (13).

Extract the $W_1 \times W_2$ window with the center at $(x_c, y_c)$ from $f(x, y)$ and $g(x, y)$.

4. **loop for** $\theta$ on range $[-\theta_{max}, \theta_{max}]$ **do**:
   Rotate image $g^W(x, y)$ around its $(x_c, y_c)$ coordinates by the angle $\theta$. The rotated image is denoted $g_{\theta}^W(x, y)$.

5. **for** the non-reflected image $g^W(x, y)$ **do**:
   Estimate the displacement $(\tau_1, \tau_2)$ with the POC function for the images $f^W(x, y)$ and $g_{\theta}^W(x, y)$.

Translate image $g_{\theta}^W(x, y)$ by the displacement $(\tau_1, \tau_2)$.

Evaluate the match of $f^W(x, y)$ and $g_{\theta}^W(x - \tau_1, y - \tau_2)$ per criterion (23).

Record the results.

6. **repeat for** the vertically and horizontally reflected image $g_{\theta}^W(x, y)$.

7. **repeat loop for** the next value of $\theta$.

8. Return the best result.

Estimating the center of gravity of an image allows for the initial detection of the object in the image. The center of gravity $(x_c, y_c)$ drifts toward the center of the object. The estimation allows us to reduce the size of the processing images $f(x, y)$ and $g(x, y)$ to the area around $(x_c, y_c)$. The new set of images $f^W(x, y)$ and $g^W(x, y)$ of size $W_1 \times W_2$ (where $W_1 < N_1$ and $W_2 < N_2$) is focused on the mass of the object in the image and has most of the peripheral values removed, making it better suited for the POC function. For $f^W(x, y)$ and $g^W(x, y)$ further adjustments of the image center can be made.

Estimation of the center of gravity of the object in the image consists of the following steps:

1. Calculate the preliminary coordinates of the center of gravity points $(x_c^f, y_c^f)$ and $(x_c^g, y_c^g)$ of images $f(x, y)$ and $g(x, y)$.

2. Assign values to indexes $x_W$ and $y_W$ of window $W_1 \times W_2$ with the center at $(x_c, y_c)$: $x_W = x_c - \frac{W_1}{2}, \dots, x_c + \frac{W_1}{2}$ and $y_W = y_c - \frac{W_2}{2}, \dots, y_c + \frac{W_2}{2}$

3. Extract the window area from images $f(x, y)$ and $g(x, y)$. The new set of images, named $f^W(x, y)$ and $g^W(x, y)$, is the size of $W_1 \times W_2$ with centers in $(x_c^f, y_c^f)$ and $(x_c^g, y_c^g)$ respectively; $x = x_W$ and $y = y_W$ for $f^W(x, y)$ and $g^W(x, y)$.

The rotation procedure is performed on $g^W(x, y)$. The new image is rotated by the angle $\theta$ around the center point of $g^W(x, y)$. For the rotated image $g_{\theta}^W(x, y)$ a new space of coordinates is assigned (see (5)). The new values are interpolated from $g(x, y)$ using a linear method. The angular range of $\theta$ is $-\theta_{max} \le \theta \le \theta_{max}$. In practical application, we used $\theta_{max} = 180^{\circ}$ with $1^{\circ}$ spacing.

Estimation of object translation with the POC function provides the translation vector $(\tau_1, \tau_2)$ for the rotated by the angle $\theta$ image. This part consists of the following steps:
1. Calculate DFT for images $f^W(x, y)$ and $g_\theta^W(x, y)$. $F(k, l) = \text{DFT}[f^W(x, y)]$, $G(k, l) = \text{DFT}[g_\theta^W(x, y)]$.

2. Calculate $\hat{R}(k, l)$ for $F(k, l)$ and $G(k, l)$.

3. Calculate $\hat{r}(x, y)$ of $\hat{R}(k, l)$.

4. Find coordinates of maximum value of $\hat{r}(x, y)$. It will be the vector of translation $(\tau_1, \tau_2)$.

The value of translation vector $(\tau_1, \tau_2)$ is obtained from the POC function $\hat{r}(x, y)$ of $f^W(x, y)$ and $g_\theta^W(x, y)$. The image $g_\theta^W(x, y)$ is translated by $(-\tau_1, -\tau_2)$. The estimation with the POC function is repeated for the image $g_\theta^W(x, y)$ reflected horizontally and vertically.

The evaluation of the match is performed on the images $f^W(x, y)$ and $g_\theta^W(x - \tau_1, y - \tau_2)$. The criterion is the mean value of the sum of absolute difference between $f^W(x, y)$ and $g_\theta^W(x - \tau_1, y - \tau_2)$ given by

$$ \frac{1}{W_1 W_2} \sum_{x,y} \left| f^W(x, y) - g^W(x - \tau_1, y - \tau_2) \right| \quad (23) $$

The smaller the value of (23), the better the match is.

The final result consists of the value (23) and the set of parameters describing the transformed image $g^W(x, y)$. These are the translation vector $(\tau_1, \tau_2)$, the rotation angle $\theta$, the center of gravity point $(x_c, y_c)$ and the information on the type of reflection. The parameters allow for the geometric transformation of the image $g$ into the image $f$ by using the set of equations (5)-(9).

## 5 Experiments and Discussion

This section describes a set of experiments for evaluating the performance of the proposed algorithm. The experimental data consists of two sets of microscopic images containing vein sections. We can distinguish two types of shapes for these sections. The sets consist of six and five images of sizes 690 × 615 pixels and 746 × 1144 pixels respectively. The sample images of each set are shown in Fig. 1 and Fig. 2. The testing of the algorithm is performed on each image sequence. The acquired results consist of the estimated values of the parameters

<figure><img src="image_2.png" /><figcaption>Fig. 1. The image set I consisting of 6 vein section images (I-1,I-2,I-3,I-4,I-5,I-6) of the size 690 × 615 pixels.</figcaption></figure>
<figure><img src="image_3.png" /><figcaption>**Fig. 2.** The image set II (II-1,II-2,II-3,II-4,II-5) consisting of 5 vein section images of the size 746 × 1144 pixels.</figcaption></figure>

needed to match the input image to the reference image using the basic geometric transformations. The parameters include: the center of gravity point ($x_c, y_c$), the rotation angle $\theta$, the translation vector ($\tau_1, \tau_2$) and the type of reflection applied. The final results were chosen based on the value of the error criterion (23) used

<figure><img src="image_4.png" /><figcaption>**Fig. 3.** The results of the section images registration using the proposed method. I-5 registered with respect to I-1 (on the left) and I-2 registered with respect to I-3 (on the right).</figcaption></figure>

<table><caption>**Table 1.** The value of criterion (23) for the best match of the pair of images in set I.</caption><tr><td>I-1</td><td>0.0708</td><td>0.0781</td><td>0.0716</td><td>0.0569</td><td>0.0611</td></tr><tr><td>0.0651</td><td>I-2</td><td>0.0633</td><td>0.0708</td><td>0.0528</td><td>0.0511</td></tr><tr><td>0.0825</td><td>0.0690</td><td>I-3</td><td>0.0853</td><td>0.0787</td><td>0.0706</td></tr><tr><td>0.0669</td><td>0.0720</td><td>0.0776</td><td>I-4</td><td>0.0688</td><td>0.0616</td></tr><tr><td>0.0594</td><td>0.0526</td><td>0.0714</td><td>0.0717</td><td>I-5</td><td>0.0546</td></tr><tr><td>0.0650</td><td>0.0459</td><td>0.0619</td><td>0.0665</td><td>0.0557</td><td>I-6</td></tr></table>

to measure the accuracy of the match for each pair of images in the sequence. The results are shown in Tables 1 and 2. A row in the table denotes the reference image and a column denotes the input image to be matched. From these tables, we can see that the error value varies depending on the reference image and the input image used. The cause of this disparity lies in the initial differences between the objects in the images. Therefore, assessing the value of the error criterion (23), we have to consider evaluating each result based on which pair of images it comes from. There is no constant threshold that applies to all results. The best considered course of action is to choose the best match at the cross-section of the
<table><caption>**Table 2.** The value of criterion (23) for the best match of the pair of images in set II.</caption><tr><td>II-1</td><td>0.1175</td><td>0.0895</td><td>0.0639</td><td>0.0676</td></tr><tr><td>0.1225</td><td>II-2</td><td>0.0979</td><td>0.1238</td><td>0.1265</td></tr><tr><td>0.0732</td><td>0.0953</td><td>II-3</td><td>0.0768</td><td>0.0733</td></tr><tr><td>0.0576</td><td>0.1062</td><td>0.0751</td><td>II-4</td><td>0.0627</td></tr><tr><td>0.0612</td><td>0.1248</td><td>0.0832</td><td>0.0605</td><td>II-5</td></tr></table>

analyzed image, i.e., the best result in both the row and the column containing the image. The most reliable assessment remains the visual confirmation of the acquired results.

The sample application of the acquired results to match the input image to reference image is presented in Fig. 3. The matched image was obtained by geometrically transforming the input image using the required parameters.

It should be noted that the method allows for the retrieval (building) of the match for a pair of images using the composition of geometric transformations from the corresponding pairs. Considering equation (4), it is also possible to acquire a reverse match by transforming the reference image to match the input image.

## 6 Conclusion

The paper presents a proposed algorithm of object matching for a series of images. The approach makes use of the distinctive features of the objects to estimate the match. The algorithm is shown to perform on the objects that are not identical, different on intensity level and geometrically transformed. It estimates the parameters of the geometric transformation differentiating the images without the influence of natural differences between the objects and the intensity of the images. It should be indicated, that the proposed approach is rather simple. It works in the image space domain and in the frequency domain consecutively, reducing distortion introduced by images displacements.

The method allows for further series analysis of the images. There is a possibility of applying the method to image scheduling and sequence construction.

*Acknowledgments* The authors express their special thanks to Dr. Agnieszka Malińska for the section images and to professor Maciej Zabel from University School of Medicine, Poznań, Poland, for indicating problems of microscopic section images registration.

## References

1. van Assen, H.C., Egmont-Petersen, M., Reiber, J.H.C.: Accurate object localization in gray level images using the center of gravity measure; accuracy versus precision. (2011)
2. Berberidis, K., Evangelidis, G. D., Karybali, I. G., Psarakis, E. Z.: An efficient spatial domain technique for subpixel image registration. Signal Processing: Image Communication, 23, 711–724 (2008)

3. Brunelli, R.: Template matching techniques in computer vision. Theory and practice. Southern Gate, Chichester: Wiley (2009)

4. DeLaPaz, R.L., Ma, Q.Y., Perera, G.M., Tang, H., Wu, E.X.: FFT-based Subpixel MRI Image Registration and Its Application in FMRI. Proc. Intl. Sot. Mag. Reson. Med., 8, 1751 (2000)

5. Fitzpatrick, J. M., Sonka, M.: Handbook of Medical Imaging, Volume 2. Medical Image Processing and Analysis (SPIE Press Monograph Vol. PM80/SC). SPIE Publications, (2009)

6. Flusser, J., Zitova, B.: Image registration methods: a survey. Image and Vision Computing, 21, 977–1000 (2003)

7. Gonzalez R.C., Woods R.E.: Digital Image Processing. Prentice Hall (2007)

8. Hwang, J.J., Kim, D.N., Rao, K.R.: Fast Fourier Transform: Algorithms And Applications. Springer (2010)

9. Ourselin, S., Roche, A., Subsol, G., Pennec, X., and Ayache, N. : Reconstructing a 3D structure from serial histological sections. Image and vision computing. 19(1), 25–31 (2001)

10. Lippolis, G., Edsjo, A., Helczynski, L., Bjartell, A., Overgaard, N.: Automatic registration of multi-modal microscopy images for integrative analysis of prostate tissue sections. MBC Cancer 13:408 (2013)

11. Modersitzki, J.: Numerical Methods for Image Registration (Numerical Mathematics and Scientific Computation). Oxford University Press, Oxford, New York (2004)

12. Reddy, B. and Chatterji, B.: An FFT-based technique for translation, rotation, and scale-invariant image registration. IEEE Transactions on Image Processing. 5, 1266–1271 (1996)

13. Schwier, M., Böhler, T., Hahn, H.K., Dahmen, U., Dirsch, O.: Registration of histological whole slide images guided by vessel structures. Journal of Pathology Informatics. 4:2(Suppl):S10 (2013)

14. Skubalska-Rafajlowicz, E.: Estimation of horizontal and vertical translations of large images based on columns and rows mean energy matching. Multidimensional Systems and Signal Processing. 25:2, 273–294 (2014)

15. Takita, K., Aoki, T., Sasaki, Y, Higuchi, T. and Kobayashi, K.: High-Accuracy Subpixel Image Registration Based on Phase-Only Correlation. IEICE Trans. Fundamentals of Electronics, Communication and Computer Sciences E86-A (8), 1925–1934 (2003)

16. Wang, C.W., Ka, S.M. and Chen, A.: Robust image registration of biological microscopic images. Scientific Reports 4: 6050 (2014)

17. Zhang, X., Abe, M., Kawamata, M.: An Efficient Subpixel Image Registration Based on the Phase-Only Correlations of Image Projections. The 2010 International Symposium on Communications and Information Technologies, 997–1001 (2010)