1. Introduction

Nowadays, many researches have been done in medical image processing to present systems that can detect and analyze various diseases automatically. Minimizing expenditures, increasing detection speed, reducing subjective mistakes and enhancing accuracy are the chief goals of these systems. Particularly in ophthalmology, a lot of researches have been conducted to propose and apply systems for diagnosis and analysis of several retinal diseases including Diabetic Retinopathy (DR) and Age-related Macular Degeneration (AMD) [1].

Image registration is the procedure of aligning two or more images into a unique common coordinate system. Applying registration on multimodal images leads to integration of information, thus benefitting from complementary nature of images. There are various ocular imaging modalities such as Fluorescein Angiography (FA) and Optical Coherence Tomography (OCT), each of which shows some features of retina better than the others. FA is an invasive imaging technique that clearly illustrates retina vasculature through injection of fluorescein dye. FA is known as the gold standard ocular imaging for presenting neurovascular structure and assessing vasculature diseases such as DR and Diabetic Macular Edema (DME) [2]. On the other hand, OCT is non-invasive, prevalent ocular imaging that works based on the Michelson interferometry rule using near-infrared light instead of sound in ultrasound and presents a cross-sectional view of retinal structure. OCT can be used to identify macular edema, macular cysts, vitreomacular traction, sub-retinal fluid, pigment epithelial detachment, choroidal neovascularization and measuring retinal thickness [3].

As mentioned before, FA is used to diagnose DR as a common retinal disease with hyperfluorescent spots as hallmarks of Microaneurysm (MA) [4]. The greater number of MA spots indicates the DR progression. Albeit FA depicts neurovascular structure of retina exactly, it utilizes intravenous injection that is harmful for human body and requires an expert photographer with the patient remaining strictly still. Unlike FA, OCT is non-invasive with retinal in-depth data. Therefore, clinicians usually use information of both imaging modalities for better diagnosis and evaluation of retinal diseases. Great number of papers have evaluated MA visual properties via intuitive comparison of FA and OCT images simultaneously [4–8] whereas information of both images can be added together using automatic and precise retinal image registration. Furthermore, due to reduction in subject-dependent errors, this can be more helpful for clinical purposes.

So far, many retinal image registration techniques have been proposed which align images of various modalities into a unique coordinate system. From a technical point of view, these techniques can be classified into two categories: segmentation-based methods and non-segmentation-based ones [9]. In the former, the retinal image is partitioned into meaningful segments such as vessels and optic disc making the analysis of the image simpler [10]. In the latter, the full image contents are processed using image processing techniques. A comprehensive review, analysis, and classification of the retinal vessel segmentation algorithms and methodologies is provided by [11], also a survey of state-of-the-art techniques for extracting retinal vessels is presented in [12]. In [13], a fast and accurate registration algorithm is proposed based on Salient Feature Regions (SFR) followed by local rigid transformation. The alignment process is then augmented with a global second order polynomial transformation. This method is only applicable to mono-modal fundus images. In [14], a registration algorithm for poor quality multimodal retinal images using Harris corner detector and Partially Intensity Invariant Feature Descriptor (PIIFD) is presented. A multimodal registration of Spectral Domain OCT (SD-OCT) volumes and fundus photographs is suggested in [15] which uses Features from Accelerated Segment Test (FAST) feature detector and Histogram of Oriented Gradients (HOG) descriptor. In this method, the affine transformation is computed using RANdom SAmple Consensus (RANSAC) [16]. A multi-resolution difference of Gaussian pyramid with Saddle detector (D-Saddle) is proposed in [17], to detect feature points on the low-quality region that consists of vessels with varying contrasts and sizes. In [18] a two-step framework is proposed to register multimodal retinal images. In the first step, HOG descriptor matching is applied to the mean phase image for transformation computation and in the second step, the deformation registration method based on Modality Independent Neighborhood Descriptor (MIND) is used to improve the registration accuracy. In [19], reliable feature point sets are extracted based on SURF detector and PIIFD descriptor from model and target retinal images. Then Mixture Feature and Structure Preservation (MFSP) is used to map feature point sets for finding exact point correspondences. An Adaptive RANSAC (A-RANSAC) method is proposed in [20] for multimodal retinal image registration. In this method, the threshold value is chosen so that the Root Mean Square Error (RMSE) and the number of removed matches are optimized simultaneously. It should be noted that all of the aforementioned methods and techniques are included in none segmentation-based category. The following methods are from the segmentation-based one.

In [21] a registration of OCT projection images with color fundus photographs is proposed based on a hierarchical local feature matching that utilizes the search of local maximization of the similarity function. Reference [22] registers pairs of FA frames for the sake of segmenting fluorescein leakage in patients with DME where two-step registration method performs global and local registration in turn. In global registration step, after vessel extraction using exploratory dijkstra forest algorithm, Speeded Up Robust Features (SURF) features are detected, extracted, and matched, then M-estimator Sample Consensus (MSAC) [23] algorithm is applied for transformation computation. Thereafter in local registration step, intensity multi-resolution registration is performed on local patches to obtain optimal results. In [24], a joint vessel segmentation and deformable registration model is proposed based on Convolutional Neural Network (CNN). In vessel segmentation, a style loss guides the model to generate segmentation maps and helps transform images of different modalities into consistent representations. In deformable registration, a content loss helps finding dense correspondence for multi-modal images. An image registration algorithm is proposed in [25] to trace changes in the retina structure across modalities using vessel segmentation and automatic landmark detection. The segmentation of the vessels is done using a U-Net and the detection of the vessel junctions is achieved with Mask R-CNN. In [26] a new method is presented to register SD-OCT and FA images. For this purpose SLO images captured by Heidelberg Spectrally HRA2/OCT device are used as intermediate images. To the best of our knowledge, only [26] registers OCT images with FA photographs via corresponding SLO images.

As a result of differences in multimodal retinal images, when common feature-based registration methods are applied in RANSAC or MSAC, the outlier matches are not removed completely. State-of-the-art mismatch removal techniques are proposed in [27–30]. An efficient approach is designed in [27] based on maintaining the local neighborhood structures of the potential true matches. The problem is formulated into a mathematical model, and a closed-form solution with linearithmic time and linear space complexities is derived. In [28], the mismatch removal is converted to a two-class classification problem, and learning a general classifier to determine the correctness of an arbitrary putative match. Method in [29] adaptively clusters the putative matches into several motion consistent clusters together with an outlier/mismatch cluster. The classic Density Based Spatial Clustering method of Applications with Noise (DBSCAN) is customized in the context of feature matching to implement spatial clustering. An iterative clustering strategy is designed to promote the matching performance in case of severely degraded data. In [30], a feature guided Gaussian Mixture Model (GMM) is proposed for the non-rigid registration of retinal images. The problem is formulated as an estimation of a feature guided mixture of densities: a GMM is fitted to one point set in which the centers of the Gaussian densities characterized by spatial positions associated with local appearance descriptors are constrained to coincide with the other point set.

Due to differences in contrast, resolution and brightness of the multimodal retinal images, prevalent feature detection, extraction and matching schemes do not result in perfect matches. Also, the images resulted from vessel extraction of image pairs are not exactly the same and the accuracy of retinal image registration methods that use vessel segmentation is largely dependent on the accuracy of vessel segmentation. In addition, the relationships between retinal image pairs are usually modeled by affine transformation, which cannot generate accurate alignments due to the non-planar retina surface.

In this paper, a novel method is proposed for registering FA and SD-OCT images using SLO photographs. For doing so, after applying some preprocessing techniques on FA and SLO images, the vessels are extracted in both images and the gray level images are converted to new black and white vessel-map ones. Next, a new feature detector-descriptor method is applied on these vessel-map photographs to remove outlier matches completely. For this goal, first, SURF features are detected and extracted from both images. Likewise, FAST-HOG framework is applied on both images to detect and extract new feature set. Applying only one of these feature detector-descriptors separately in RANSAC or MSAC, does not result in perfect matches. So, SURF feature descriptors from two images are matched using approximate nearest neighbors [31]. The same procedure is iterated with HOG feature descriptors. Since this matching method may include some erroneous matches, the MSAC method is applied separately on matched SURF and HOG features to refine them, but because of vessel segmentation errors and different modalities, some outliers still remained. After that the two vectors of matched featuresF resulted from the previous steps are concatenated to enhance matched features. After applying MSAC again on concatenated features and obtaining perfect matches, global rigid transformation is calculated and applied to the vessel-map and the original images. Combining these features and applying the MSAC for three times, guarantee complete elimination of outlier matches. In addition, to consider non-planar surface of retina, the transformed images are globally aligned again applying a Gaussian model for the surface of retina which increases the accuracy of the previous step. This reduces the amount of deformation effects resulted from the subsequent local registration step. High image deformations in multimodal image registration, could affect the primary image, especially in parts of image that contain no vessels. Eventually displacement field between transformed and reference original images are estimated using method presented in [32] to register two images perfectly. Our dataset contains 36 image pairs of 21 subjects with diabetic retinopathy. Using accurate registration of FA and SD-OCT images, enables us to exactly detect location and morphological appearance of diabetic retinopathy symptoms such as MAs and leakages in OCT B-scans.

The proposed methods are explained in section 2. section 3 presents experimental results and finally, this paper is concluded in section 4.

2. Methods

Because of the different nature of two modalities, direct registration of OCT images and FA photographs is not possible; one presenting a depth view of retina and the other showing the surface of it. On the other hand, Heidelberg Spectralis HRA2/OCT device provides the capability of simultaneous SLO/OCT imaging by means of a single light source that guarantees exact correspondence between SLO image of retina surface and cross sectional B-scans of OCT. In this paper precise automatic registration of FA images with OCT photographs is proposed using registration of FA images with SLO photographs. Comprehensive block diagram of the proposed four-step method is demonstrated in Fig. 1. As shown, overall process is divided into four main sections: data acquisition, preprocessing, global registration, and local registration. In the following, the steps are explained in detail.

 

Fig. 1. Block diagram of overall process of registration.

Download Full Size | PDF

2.1 Data acquisition

In this research, dataset is comprised of 36 pairs of FA and SLO images of 21 subjects with diabetic retinopathy, where SLO image pixels are perfectly in correspondence with OCT B-scans (Fig. 2 and Fig. 3). The FA, OCT, and SLO images are captured via Heidelberg Spectralis HRA2/OCT device in Didavaran Eye Clinic, Isfahan. Moreover the FA and SLO images are the same size as 768×768 pixels and FA images were captured with two different fields of views (30 and 55 degrees). Here, the FA and SLO images are considered as moving (source) and fixed (target) images, respectively.

 

Fig. 2. (a) SLO image. (b) the OCT B-scan corresponds to the yellow line in the SLO image.

Download Full Size | PDF

 

Fig. 3. Our data set includes 36 pairs of FA and SLO images. Odd columns illustrate FA images of different field of view (30 and 55 degrees), even columns depict the SLO images corresponding to the left side column.

Download Full Size | PDF

2.2 Preprocessing

In this step, first the contrast of FA and SLO images are enhanced utilizing the Contrast Limited Adaptive Histogram Equalization (CLAHE) method [33]. Then FA image pixels intensities are complemented to become similar to SLO image. Next, the vessels of both moving and fixed photographs are extracted. In order to do so, the method presented in [34] is applied to precisely segment the main vessel structure of the smoothed image from original one, then binarization using Otsu's threshold selection method [35] and elimination of small objects from binary image are performed. Due to differences in contrast, resolution and brightness of the multimodal retinal images, the binary images resulted from vessel extraction of FA and SLO images are not exactly the same (Fig. 4).

 

Fig. 4. Vessel Extraction. (a) FA image. (b) vessels extracted from FA. (c) SLO image. (d) vessels extracted from SLO.

Download Full Size | PDF

2.3 Global registration

Global registration process consists of bringing FA and SLO images into spatial alignment. Here, FA image is mapped into the space of SLO photograph. For the registration scheme, it should be taken into account that the retina surface is not flat [1,36,37]. The registration schemes in [36,37] use the quadratic model for the surface of retina. Since these methods consider only the small shifts in eye and then model them as translations and rotations, they cannot be applied directly to multimodal and multi-scale images. Also, quadratic model does not adequately describe the retina surface for field of view of our image data set and leads to inaccurate results. On the other hand, studies in [1] describes the retina as a Gaussian surface, but they did not use this theory in related applications. Therefore, a new two-step global registration is proposed applying Gaussian model to improve the global registration accuracy. In our proposed method, the first step aligns two images globally using a new feature-based method and making them of the same scale, and then Gaussian model is applied to the system to calculate displacements to decrease the effects of deformations in the consequent local registration step. Rigid global registration is obtained using the following steps:

2.3.1 Control points detection

From another point of view, image registration techniques can be classified into two categories: intensity based techniques and feature based ones. In the former, whole image or sub- images are registered whereas the latter uses distinct image features as control points and matches corresponding features of two images for registration. Control points are well-defined image pixels that can be detected robustly and are distributed all over the image. In the proposed method, features are detected, extracted and matched to be used as input of transformation computation stage. Therefore, after vessel extraction, the control points should be detected in vessel-map images. As mentioned before, due to inevitable errors in vessel segmentation algorithms, using a feature detector-descriptor and MSAC alone does not lead to perfect matches. To address this problem, two feature detection algorithms are utilized: Features from Accelerated Segment Test (FAST) corner points [38] and Speeded Up Robust Features (SURF) [39].

2.3.1.1 FAST corner detector

FAST corner detector utilizes a Bresenham circle of radius 3 containing 16 pixels around the candidate pixel (Fig. 6(a)). If the pixel satisfies either of the two following conditions, it can be determined as the corner point:

2.3.1.2 SURF (Fast-Hessian) blob detector

Blobs are image regions which differ from surrounding regions in some image properties like color, brightness, etc. The principles behind the SURF blob detector are based on that of Scale-Invariant Feature Transform (SIFT) [40], however they differ in details. SURF works on integral images instead of original ones and uses box filters to increase the speed of filtering. Utilizing integral image, the summation of values within a rectangle of a different size in image can be computed just through evaluating four corners of that rectangle. SURF employs determinant of hessian matrix Eq. (1) to specify the location and scale of feature and scale of feature points.

 

Fig. 5. (a) Fast corner detector and Bresenham circle. (b) FAST corners in FA vessel map (green dots). (c) FAST corners in SLO vessel map (green dots).

Download Full Size | PDF

 

Fig. 6. SURF features (a) 150 strongest SURF features on SLO vessel map (green circles). (b) 150 strongest SURF features on FA vessel map (green circles).

Download Full Size | PDF

In Eq. (2), $\det ({H_{approx}})$, ${D_{xx}}$, ${D_{yy}}$ and ${D_{xy}}$ are approximations for $\det (H)$, ${L_{xx}}$, ${L_{yy}}$ and ${L_{xy}}$, respectively [38]. Figure 6 indicates the 150 strongest SURF features in both SLO and FA vessel maps.

2.3.2 Feature extraction

When features have been detected, the local patch around each feature can be described as a compact representation named feature vector. In proposed method, the SURF (fast-hessian) and FAST detected features are extracted using SURF and HOG feature descriptors.

2.3.2.1 SURF feature descriptor

To achieve rotation-invariant descriptor, main orientation is assigned to the previously-detected SURF features. For doing so, Haar wavelet responses are calculated in vertical and horizontal directions in circular neighborhood of radius 6×s around detected feature, where s denotes the image scale of the feature. After that, appropriate Gaussian weights are applied. The main orientation is appraised via computing the sum of all wavelet responses within a sliding orientation window of 60 degrees angle. Then for a $(4 \times 4)$ sub-region aligned to the estimated orientation around detected feature, horizontal and vertical wavelet responses are taken and a vector $v = (\sum {{d_x},\sum {{d_y}} ,\sum {|{{d_x}} |} } ,\sum {|{{d_y}} |} )$ is formed for each pixel leading to $64$-dimension $(4 \times 4 \times 4)$ feature vector [39].

2.3.2.2 HOG feature descriptor

Here HOG descriptor is used to describe FAST features [41]. Individual gradients may be affected by noise. This is why histogram of gradients is utilized to represent image patches that makes descriptor robust against noise. For this aim, first image is divided to blocks and each block is partitioned into cells. Thereafter gradient value and orientation are calculated for each pixel within the cell. $N$ bins are considered, each of which represents an angle starting from 0 degree, ending by $180 - ({\raise0.7ex\hbox{${180}$} \!\mathord{\left/ {\vphantom {{180} N}} \right.}\!\lower0.7ex\hbox{$N$}})$ degree with steps of ${\raise0.7ex\hbox{${180}$} \!\mathord{\left/ {\vphantom {{180} N}} \right.}\!\lower0.7ex\hbox{$N$}}$ degree for unsigned gradients. Within the cell, the gradient values corresponding to the same orientation are summed leading to $N$-dimension vector where each dimension contains aggregate of corresponding gradient values. Taking into account that if the gradient orientation is between two bins, the corresponding value should be split proportionally into both adjacent bins. Finally, feature vector is normalized to be not affected by lighting variations. In the proposed method, $N$=9, cell size = [8,8] pixels, and block size = [2,2] cells, so HOG vector has 36 dimensions (9×2×2). Figure 7 visualizes HOG descriptors on vessel map images.

 

Fig. 7. HOG visualization. (a) SLO vessel map. (b) HOG visualization of (a). (c) FA vessel map. (d) HOG visualization of (c).

Download Full Size | PDF

2.3.3 Feature matching

After FAST and SURF features are detected and extracted utilizing HOG and SURF descriptors, HOG descriptors of both binary images containing vessel maps should be matched. This procedure is applied for SURF features in the same way. To do this, the method in [31] is used which presents an automatic algorithm for searching approximate nearest neighbor considering data sets and desired precision to speed up the matching process. According to [31], depending on the data set and required precision, two methods namely searching hierarchical k-means trees with a priority search order [31] and multiple randomized kd-trees [42] can perform perfectly. Afterwards, since the abovementioned matching methods may include some erroneous matches, the MSAC method is applied separately on matched SURF and HOG features to refine them and obtain better matches. This results in more precise inlier set of matched features, but as can be seen in Fig. 8(a) some outliers still remained. MSAC method is a variant of RANSAC. For more details, the thresholding in RANSAC can be understood as a top-hat cost function whereas MSAC uses truncated quadratic cost function [43]. The steps of MSAC are as follows:

 

Fig. 8. Final MSAC step refines the previously-matched features. (a) matched features after joining separately refined SURF and HOG matches (still contains few outliers). (b) obtaining perfect matches after final MSAC step.

Download Full Size | PDF

2.3.4 Transformation computation

To compute the similarity transformation which is registration objective, the matched FAST and SURF features of fixed and moving images are joined. Thereafter to calculate transformation, the MSAC method is utilized again and final transformation is re-computed utilizing all inliers associated with the largest number of ${S_i}$. Using SURF and FAST features together in this manner has two main advantages: first, it benefits from properties of both feature detection-description frameworks and secondly, refining matched features using the MSAC for three times, causes the erroneous matches to be removed completely. Figure 8 shows how final MSAC step refines the previously-matched features and perfects them. Then the computed transformation is applied to the original images. The proposed feature detection, description and matching scheme is flexible and can be generalized. In fact, other feature detection, extraction and matching techniques can be added to the scheme to increase the number of inlier matches under the penalty of longer execution time. In general, for $k$ number of separate feature detectors, m-sac should be applied for $k + 1$ times.

2.3.5 Applying the Gaussian model for retina surface

As mentioned before, the retina surface is not flat and in the proposed method, Gaussian model is applied to improve the global registration accuracy. Figure 9 illustrates a point on retina surface that is expressed in two different 3D coordinate systems related to two different cameras namely ${C_p}$ and ${C_q}$. These points are represented as ${P^T} = (X,Y,Z)$ and ${Q^T} = (X^{\prime},Y^{\prime},Z^{\prime})$ in coordinate systems of ${C_p}$ and ${C_q}$. Assume the transformed FA and SLO images are shown by ${I_p}$ and ${I_q}$ corresponding to $P$ and $Q$ coordinate systems. The Gaussian equation which models the retina surface and relates $X,Y$ and $Z$ in the coordinate system of ${C_p}$ is:

 

Fig. 9. Two different camera viewpoints of the retina.

Download Full Size | PDF

Here, weak perspective projection model adequately projects the 3D retina surface to 2D image plane [36,37]. The projections of P and $Q$ in images ${I_p}$ and ${I_q}$ considering homogenous coordinates is as follows:

From Eq. (6):

Now from Eq. (5), Eq. (8), and Eq. (9), the $X^{\prime}$ and $Y^{\prime}$, are calculated in terms of $x$ and $y$. Finally $x^{\prime}$ and $y^{\prime}$ can be obtained from Eq. (8) as follows:

 

Fig. 10. Qualitative evaluation of three global registration schemes (Yellow rectangles highlight the differences). (a) global registration using first step(false color). (b) global registration exploiting quadratic model for retina surface(false color). (c) global registration using Gaussian model for retina surface(false color). In false color method, first and second images are shown in green and magenta, and the overlapped regions are depicted in white.

Download Full Size | PDF

2.4 Local registration

Since global registration aligns images approximately (Fig. 11(a), Fig. 11(b), Fig. 11(d) and Fig. 11(e)), local registration is required to estimate displacement field and register images accurately (Fig. 11(c) and Fig. 11(f)). For local registration, the method proposed in [32] is applied. To achieve this, one image called diffusing model should be locally deformed to become similar to the other fixed one, named scene. In this approach, diffusing model (transformed FA image) is assumed as a grid model that must be deformed. This model contains vertices that are particles. On the other hand, contours of the scene objects (SLO image objects) are assumed to be membrane with effectors named demons. These demons are placed throughout the membrane and the presumptive vectors which are perpendicular to this contour, oriented from inside of image objects to the outside. In grid model, vertices polarity are labeled as either inside or outside and accordingly in the other image, demons locally apply forces to push the model inside the other image objects or vice versa [44]. For example, the demon corresponding to vertex with inside polarity, pushes that vertex of model toward inside of the image object. This method applies transformation iteratively. Magnitudes of forces applied from demons are constant within an iteration, but decreases between iterations for convergence. On the other words, local registration using diffusing model can be divided into two steps:

 

Fig. 11. Vessel maps and images after different steps of the proposed registration. (a) vessel maps after first step of global registration (false color). (b) vessel maps after applying Gaussian model to the first step of global registration (false color). (c) vessel maps after local registration step (false color) (arrows point out the differences in registration accuracy and yellow circles show the key points). (d-f) images correspond to the vessel maps (checker board). (g-i) zoomed images correspond to the ROI specified via red rectangle in second row (checker board). (j-l) zoomed images correspond to the ROI specified via magenta rectangle in second row (checker board). (m-o) zoomed images correspond to the ROI specified via yellow rectangle in second row (checker board). In checker board method, image is composed of alternating rectangular regions from two images.

Download Full Size | PDF

2.4.1 Pre-computation of demons set

The demons should be extracted from the scene (SLO). They can be elicited from all image pixels, image contours or edges. Each demon includes information of its spatial position, direction from the inside to the outside, current displacement from the corresponding point of the model.

2.4.2 Iterative process

After applying an initial transformation to the model, the iterative process is started. Each iteration is composed of two steps:

 

Fig. 12. All data set images (zoomed ROI) after local registration step (checker board). Since registration of pair of images 36 is failed, this image is not shown.

Download Full Size | PDF

3. Results

In this paper, three standard and prevalent parameters have been considered to assess the performance of proposed registration method: Root-Mean-Square (RMS) error, success rate and speed of registration.

3.1 RMS error

In mono-modal image registration, the pixel intensity values of two images can be utilized to determine distance between corresponding pixels. However in this case, for multimodal images, distances are calculated according to Eq. (11) where n and $||{{p_i}} ||$ stand for the number of landmarks and Euclidean distances of $i$-th corresponding landmarks.

Since vessel branches and crossovers are easy to locate even in multimodal retinal images, these points have been determined as landmarks in each image. Here $n$ equals 15. The mean and standard deviation of RMS error of proposed method after global and local registration steps is indicated in Table 1, also the box-plot of RMS errors of global and local registration steps are demonstrated in Fig. 13. The process of measuring RMS errors is done by two observers, and the differences between mean and standard deviation of RMS errors measured by them is negligible. Nevertheless, the maximum values are reported in Table 1. Table 2 compares mean RMS error of proposed method with other methods. It can be seen from this table that due to combining SURF and FAST features, applying Gaussian model for retina surface, and final local registration, in the proposed method the RMS error have been dramatically reduced.

Table 1. Mean and maximum RMS errors of global and local registration

View Table | View all tables in this article

 

Fig. 13. The box-plot of RMS errors of global and local registration steps.

Download Full Size | PDF

Table 2. Quantitative evaluation of mean and maximum RMS errors of different methods

View Table | View all tables in this article

3.2 Success rate

Success rate determines what percentage of images are aligned correctly. Registration success or failure, is directly observable. The success rate of different methods is compared in Table 3. Figure 14 indicates qualitative evaluation of the proposed method.

 

Fig. 14. Qualitative evaluation of proposed method. (a) matching features using SURF on image case 3. (b) matching features using proposed method. (c) registration failure utilizing SURF features.(in color) (d) registration using proposed method. (e) matching features using FAST-HOG on image case 8. (f) matching features using proposed method. (g) registration failure utilizing FAST-HOG. (in color) (h) registration using proposed method.

Download Full Size | PDF

Table 3. Success rate of different methods on the dataset

View Table | View all tables in this article

3.3 Speed

Running time of the proposed method is compared with other methods in Table 4. It is worth mentioning that, our experiments were done using a PC with Windows 7 64-bit operating system, 4 GB RAM, and Intel Core i5 CPU 2.5 GHz whereas other reported results were obtained utilizing a PC with Windows 7 64-bit operating system, 64 GB RAM, and Intel Xeon i5 CPU 3.7 GHz.

Table 4. Running time of different methods

View Table | View all tables in this article

4. Discussions and conclusions

In this paper, a novel method is proposed to register OCT images with FA photographs in an accurate manner. This method uses SLO images as intermediate, and includes four main steps namely image acquisition, preprocessing, global registration, and local registration. Because of using SURF and FAST features together and refining the matches for three times, the proposed method has inherited the advantages of both methods resulting in enhanced success rate. On the other hand, applying Gaussian model for retina surface and local registration after global registration leads to dramatic reduction in RMS error. Precise registration of FA and OCT images via SLO photographs enables us to exactly detect location and morphological appearance of diabetic retinopathy symptoms such as MAs and leakages in OCT B-scans. In this research, after registration of FA and OCT images using proposed method, morphological appearance and location of MAs in the OCT b-scans are estimated exploiting the same properties in the corresponding FA images. Using the fact that MAs in the FA images are indicated by hyper-fluorescent dots, and utilizing the position of these dots in FA photographs, we analyzed location and morphological appearance of MAs in OCT b-scans (Fig. 15). For this purpose, 38 MAs which are located in 20 images are analyzed. In 8 out of 38 cases, MAs had complete or incomplete ring sign, whereas remaining MAs have no ring sign. Some leakages in OCT B-scans were observed as blob signs whereas other leakages did not show any visual sign. Figure 16 illustrates blob sign of a leakage in FA and corresponding B-scan of OCT image. Also, many researches have been conducted to analyze appearance of MAs in OCT B-scans [5,46].

 

Fig. 15. MAs in B-scans of diabetic retinopathy. (a) MA in FA image. (b) corresponding B-scan of the yellow line in (a). (c) zoomed complete ring sign in image (b). (d) MA in an FA image. (e) corresponding B-scan of the yellow line in (d). (f) zoomed incomplete ring sign in image (e). (g) MA in an FA image.(h) corresponding B-scan of the yellow line in (g). (i) MA with no sign.

Download Full Size | PDF

 

Fig. 16. leakage in B-scan of diabetic retinopathy (a) FA image with leakage. Red arrow and circle show part of leakage. (b) corresponding B-scan of the yellow line in (a). yellow arrows indicate the blob shaped (cystoid) area corresponds to red circle in (a). (c) zoomed image of the blob shaped (cystoid) area in (b).

Download Full Size | PDF