Dinh-Hoan Trinh

Christian Daul
email: christian.daul@univ-lorraine.fr

On Illumination-Invariant Variational Optical Flow for Weakly Textured Scenes

This paper deals with variational optical flow approaches for motion estimation under varying illumination conditions in weakly textured scenes. It proposes a systematic and complete study on descriptor-based data-terms that lead to a robust variational optical flow model. Unlike the literature which most often only experimentally shows that a descriptor is illumination invariant, this contribution gives a theoretical proof of this invariance. First, a local illumination change model is proposed and used to mathematically check whether a descriptor is invariant or not with respect to illumination variations between images. Then, this contribution proposes two general mathematical formulations which can be used to design a wide variety of new illumination-invariant descriptors. To illustrate the interest of the proposed approach, two novel illumination-invariant descriptors are constructed using the proposed general formulations. Moreover, the performance of the descriptors was evaluated on numerous datasets with known ground truth optical flow, while the robustness of the variational optical flow approach was highlighted using complex medical image sequences without ground truth. These experimental results have shown that data-terms based on the proposed descriptors led to accurate and constant optical flow under varying illumination conditions.

variational methods [START_REF] Xu | Motion detail preserving optical flow estimation[END_REF][START_REF] Drulea | Motion estimation using the correlation transform[END_REF][START_REF] Mohamed | Illumination-robust optical flow using a local directional pattern[END_REF][START_REF] Ali | Illumination invariant optical flow using neighborhood descriptors[END_REF].

As an example, Fig. 1 illustrates the OF result obtained by the CPM-Flow patch matching method [START_REF] Hu | Efficient coarse-to-fine patchmatch for large displacement optical flow[END_REF] and the proposed variational method for a pair of gastroscopic images with very weak textures. Figs. 1(a)-1(b) are two consecutive video-frames of the inner stomach wall obtained with the endoscope axis which is relatively perpendicular to the epithelium surface and with a camera translation between the two acquisitions (i.e. with a camera movement without in-and out-plane rotations). In such a situation the flow field vectors should all be almost parallel and their magnitude should only slightly change according to the distance between the camera and the scene points. It is clearly visible in Fig. 1(c) that the vector field obtained with the CPM-Flow method [START_REF] Hu | Efficient coarse-to-fine patchmatch for large displacement optical flow[END_REF] is not uniform. Indeed, the three red image rectangles highlight zones in which there is a large disparity in the vector magnitudes and orientations reveal rather large errors in the OF determination. On the contrary, as shown in Fig. 1(d), with the variational method proposed in this paper, the vector field is quite uniform, both in orientation and magnitude 1 .

As illustrated by the previous example, variational methods remain valuable approaches for difficult scenes with few textures and with changing illumination conditions. Moreover, the variational methods can also be used as a refinement step in other approaches such as the FlowNet [START_REF] Dosovitskiy | Flownet: Learning optical flow with convolutional networks[END_REF]. Contrary to recent publications which highlighted the interest of CNN-based or patch matching methods, this contribution focuses on the contributing role of variational methods for weakly textured scenes under strong illumination changes. The aim of this study is not to compare the performance of different OF approaches (learning based methods, patch matching approaches or variational methods) since they have all optimal performances for different scene types, and their applicability depends on the scene type (e.g., in medical endoscopy it is not possible to obtain ground 1 The difference in uniformity of the OF fields of the IIOF-NLDP method and the CPM-Flow method is also shown in the 120 images of the gastroscopic video-sequence (pyloric antrum region of the stomach) provided as a supplementary material of this paper. OF results of the CPM-Flow method [START_REF] Hu | Efficient coarse-to-fine patchmatch for large displacement optical flow[END_REF] and the proposed method on a pair of gastroscopic images.

truth flow fields for the learning step of CNN-based methods).

The goal of this paper is rather to deeply investigate patch-based illumination invariant descriptors leading to robust data-terms in a variational OF model. From the theoretical point of view, the major contribution of this paper lies in the general mathematical formulation of illumination-invariant descriptors that facilitate the design of robust data-terms.

Previous works

The importance of variational OF estimation has been demonstrated since the pioneering work of Horn and Schunck [START_REF] Horn | Determining optical flow[END_REF]. Since this publication, numerous improvements were presented in the field of variational OF. Interested readers may refer to [START_REF] Baker | A database and evaluation methodology for optical flow[END_REF][START_REF] Sun | A quantitative analysis of current practices in optical flow estimation and the principles behind them[END_REF][START_REF] Fortun | Optical flow modeling and computation: A survey[END_REF] for a comprehensive overview of OF methods and the general principles behind them.

Variational models for OF estimation can be generally formulated as follows. Given a source image I s and a target image I t , the dense flow field u = (u x , u y ) between I s and I t is computed by minimizing

E(u) = E reg (u) + λE data (I s , I t , u), (1)
where E reg is a regularization term that assumes smoothness of solution u, E data stands for the data-term that measures the similarity of the pixels in I s and I t while parameter λ > 0 controls the relative importance of the data and regularization terms.

The data-terms E data (I s , I t , u) proposed in the literature correspond to an energy with one or several components E i (I s , I t , u):

E data (I s , I t , u) = K i=1 λ i E i (I s , I t , u). (2)
The scalars λ i > 0 act as weights controlling the relative importance of the K different data-term parts. Each component E i (I s , I t , u) in ( 2) is based on a certain constancy assumption.

The classical variational OF methods [START_REF] Horn | Determining optical flow[END_REF][START_REF] Nagel | An investigation of smoothness constraints for the estimation of displacement vector fields from image sequences[END_REF][START_REF] Zach | A duality based approach for realtime TV-L 1 optical flow[END_REF][START_REF] Werlberger | Anisotropic Huber-L1 optical flow[END_REF] define the data-term using the well-known brightness constancy assumption (BCA):

I t (x + u x ) = I s (x). (3)
More specifically,

E data = E BCA = x∈Ω Ψ(I s (x) -I t (x + u x )), (4)
where Ψ(•) is a penalty function and Ω stands for the image domain. It is clear that the BCA is most often not an appropriate assumption since, in many real scenes, the illumination of homologous pixels changes between two images.

Brox et al. [START_REF] Brox | High accuracy optical flow estimation based on a theory for warping[END_REF] have shown that the robustness against illumination changes can be improved by combining the BCA with the gradient constancy assumption (GCA, energy E GCA ). The two-component data-term becomes:

E data = λ 1 E BCA + λ 2 E GCA , (5)
with

E GCA = x∈Ω Ψ( ∇I t (x + u x ) -∇I s (x) ). (6)
This popular BCA-GCA data-term was implemented in several optical flow methods [START_REF] Yang | Dense, accurate optical flow estimation with piecewise parametric model[END_REF][START_REF] Zhang | Robust non-local TV-L1 optical flow estimation with occlusion detection[END_REF]. However, Xu et al. [START_REF] Xu | Motion detail preserving optical flow estimation[END_REF] showed that a simultaneous use of BCA and GCA in all pixels is not optimal in terms of robustness. In fact, an efficient simultaneous use of BCA and GCA is difficult to reach because determining the weights (which give the appropriate relative importance of the BCA and GCA terms) is not easy and is strongly scene dependent. Moreover, the simultaneous use of both constancy assumptions leads also to a computational complexity in solving the optimization problem. Consequently, Xu et al. [START_REF] Xu | Motion detail preserving optical flow estimation[END_REF] proposed to use a binary mapping method that locally selects (according the pixel values in space Ω) E data based either on BCA or on GCA. Although GCA is less sensitive than BCA to illumination variations, it is only invariant under additive illumination changes. As pointed out in [START_REF] Hafner | Why is the census transform good for robust optic flow computation?[END_REF] and [START_REF] Kumar | A decoupled approach to illumination-robust optical flow estimation[END_REF], GCA is not able to compensate for the multiplicative illumination changes.

In particular, in scenarios with large illumination changes, GCA is usually not fulfilled since multiplicative illumination changes are significant in such scene conditions.

Other constancy assumptions based on high order derivatives such as the Laplacian and the Hessian have been investigated

in [START_REF] Papenberg | Highly accurate optic flow computation with theoretically justified warping[END_REF]. The data-term proposed in [START_REF] Papenberg | Highly accurate optic flow computation with theoretically justified warping[END_REF] is obtained by a linear combination (as in ( 2)) of BCA, GCA, the Laplacian constancy assumption, and the Hessian constancy assumption. However, the use of higher order derivatives significantly increases both the computational complexity and the sensitivity with respect to noise (data-term based on GCA is less sensitive to noise than that in [START_REF] Papenberg | Highly accurate optic flow computation with theoretically justified warping[END_REF] since it contains lower order derivatives). More generally, a major drawback of a linear combination of constancy assumptions lies in the difficulty to find the weights λ i in [START_REF] Ilg | Flownet2.0: Evolution of optical flow estimation with deep networks[END_REF] ensuring a robust OF determination for different scenes.

Wedel et al. [START_REF] Wedel | An improved algorithm for TV-L1 optical flow[END_REF] proposed an alternative to higher-order constancy assumptions that preprocesses the images using a structure-texture decomposition. The input image for OF estimation is a linear combination of the structure and texture components, with an emphasis on the texture components. Although this method was successfully used in [START_REF] Sun | A quantitative analysis of current practices in optical flow estimation and the principles behind them[END_REF], its computation complexity is high and the weight adjustment (balance between texture and structure) is not trivial for different scenes.

In other approaches [START_REF] Golland | Montion from color[END_REF][START_REF] De Weijer | Robust optical flow from photometric invariants[END_REF][START_REF] Mileva | Illumination-robust variational optical flow with photometric invariants[END_REF][START_REF] Zickler | Color sub-spaces as photometric invariants[END_REF][START_REF] Zimmer | Optical flow in harmony[END_REF][START_REF] Park | Illumination robust optical flow estimation by illumination-chromaticity decoupling[END_REF], the data-term was computed with photometric invariants to improve the robustness of the OF against illumination changes in color image sequences. The overall idea behind photometric invariants is detailed in Mileva et al. [START_REF] Mileva | Illumination-robust variational optical flow with photometric invariants[END_REF]. These authors proposed a dataterm combining several photometric invariants obtained from normalized RGB channels, log-derivatives, from the HSI color space, or from the rφθ color space. However, this model is only applicable to color images. Zimmer et al. [START_REF] Zimmer | Optical flow in harmony[END_REF] proposed to use constraint normalization, and a HSV color representation with higher order constancy assumptions.

Kumar et al. [START_REF] Kumar | A decoupled approach to illumination-robust optical flow estimation[END_REF] assume that the image intensity I is proportional to the product of illumination L with reflection R, i.e. I ∝ L × R. The authors proposed to decouple the intensities of images I s and I t in a illumination component L and reflection component R. Then, the OF field between I s and I t is estimated from the two images I s = β log L s + log R s and I t = β log L t +log R t , with β ∈ (0, 1) since reflectance is invariant to illumination changes. It was demonstrated that the decoupled method combined with the BCA and GCA assumptions achieves improved accuracy under high illumination variation.

However, the performance of this method strongly depends on the accurate estimation of the illumination component L. No theoretical basis that proves that this method is invariant to illumination changes was provided.

Haussecker and Fleet [START_REF] Haussecker | Computing optical flow with physical model of brightness variation[END_REF] introduced a framework that tackle the illumination issue by using explicit models of the underlying physical process that causes illumination changes. In this framework, the OF field and the parameters of the illumination model have to be simultaneously estimated. Unlike [START_REF] Haussecker | Computing optical flow with physical model of brightness variation[END_REF],

Negahdaripour [START_REF] Negahdaripour | Revised definition of optical flow: Integration of radiometric and geometric cues for dynamic scene analysis[END_REF] proposed the generalized dynamic image model (GDIM) that describes the illumination changes as a combination of a multiplicative and an additive component regardless of the underlying physical events. Accordingly, the author replaced the classical BCA in (3) by a more general constancy constraint:

I t (x + u x ) = a x I s (x) + b x , (7)
where the brightness of corresponding pixels in two consecutive images is related via the motion vector u x , as well as by the radiometric parameters a x and b x . It is clear that the GDIM model is suited to many scenes due its general nature. Kim et al. [START_REF] Kim | Robust motion estimation under varying illumination[END_REF] used this model to construct a robust variational method for OF computation. However, the optimization of the OF with this method is quite complex for a limited effectiveness.

Recently, descriptor-based variational OF methods became very popular because of their effectiveness in tackling the problem of illumination changes. The underlying idea of these methods is to construct a robust data-term by defining illumination-invariant descriptors at each pixel. Some popular descriptor-based approaches are the census transform-based methods [START_REF] Zabih | Non-parametric local transforms for computing visual correspondence[END_REF][START_REF] Stein | Efficient computation of optical flow using the census transform[END_REF][START_REF] Müller | Illuminationrobust dense optical flow using census signatures[END_REF][START_REF] Ranftl | Pushing the limits of stereo using variational stereo estimation[END_REF], the rank transform-based method [START_REF] Demetz | The complete rank transform: A tool for accurate and morphologically invariant matching of structures[END_REF],

the correlation-based methods [START_REF] Molnár | Illumination-robust variational optical flow using cross-correlation[END_REF][START_REF] Werlberger | Motion estimation with non-local total variation regularization[END_REF][START_REF] Drulea | Motion estimation using the correlation transform[END_REF][START_REF] Fortun | Aggregation of patch-based estimations for illumination-invariant optical flow in live cell imaging[END_REF], the MLDP method of Mohamed et al. [START_REF] Mohamed | Illumination-robust optical flow using a local directional pattern[END_REF], and the NND method of Ali et al. [START_REF] Ali | Illumination invariant optical flow using neighborhood descriptors[END_REF].

Motivation of the paper

As demonstrated in [START_REF] Ali | Illumination invariant optical flow using neighborhood descriptors[END_REF], the descriptor-based methods are able to preserve accuracy of the OF under changing illumination conditions. However, up to now, most of the contributions describing descriptor-based methods only introduce a formulation of the descriptors and use experiments to highlight their performances, without giving a mathematical justification of their invariance according to illumination change models.

Moreover, although many of the descriptors are based on similar principles, there is no unified theoretical basis for defining these descriptors.

In preliminary works we proposed a new illuminationinvariant descriptor [START_REF] Trinh | Illumination-invariant optical flow: Application to endoscopic image mosaicing[END_REF] and a general form [START_REF] Trinh | A general form of illuminationinvariant descriptors in variational optical flow estimation[END_REF] from which descriptors can be derived. However, no general and detailed mathematical discussion was given about the illumination invariance of descriptors and their design. Moreover, no thorough assessment and comparison of existing and new descriptors were made on reference benchmarks (KITTI, MPI-Sintel, etc.). This paper proposes further developments in the definition of illumination-invariant data terms in variational OF model. The major contributions of this paper are as follows:

1. We use a local illumination change model to describe the illumination changes between homologous images regions. Complex scene illumination changes can be considered when the corresponding neighborhoods in images I s and I t are small enough.

2. We propose a unified theoretical basis for defining illumination-invariant descriptors in the data-term of variational approaches. This theoretical basis is used to gives answers to following questions : What explains the robustness against illumination changes of some well-known models as the census transform, the correlation transform or the complete rank transform? What are their underlying mathematical concepts?

3. We introduce two generalized formulations facilitating the design of illumination-invariant descriptors (a completely novel formulation and a more detailed description of the formulation in [START_REF] Trinh | A general form of illuminationinvariant descriptors in variational optical flow estimation[END_REF]).

4. Based on the generalized formulations, we propose two illumination-invariant descriptors (a novel descriptor and an improved mathematical justification of the descriptor in [START_REF] Trinh | A general form of illuminationinvariant descriptors in variational optical flow estimation[END_REF]).

5. We also present a comparison and evaluation of illumination-invariant descriptors using the same variational model.

Paper organization

The remainder of this paper is organized as follows. Section 2 presents the global and general framework for the descriptor-based variational OF methods. Section 3 proposes an illumination change model, as well as a criterion to determine whether a descriptor is invariant or not to illumination changes.

This section also recalls some existing descriptors which are robust to varying illumination, and mathematically proves their illumination invariance. Two generalized forms of illuminationinvariant descriptors are introduced in Section 4. The general forms given in Section 4 are then used to propose two novel descriptors in Section 5. In Section 6, the performance of the new descriptors is compared to that of reference descriptors of the literature. Section 7 shows the interest of the proposed OF approach in the frame of complex scene mosaicing. The last section gives a general conclusion for this contribution.

Descriptor-based Variational Optical Flow Model

In the general energy given in [START_REF] Dosovitskiy | Flownet: Learning optical flow with convolutional networks[END_REF], descriptor-based dataterms can be formulated as follows:

E data = x∈Ω Ψ(D(P I s (x)) -D(P I t (x + u x ))), (8)
where I s and I t : Ω → R are source and target images respec-

tively, Ω = {x = (x, y) | 1 ≤ x ≤ N, 1 ≤ y ≤ M} ⊆ N 2 , Ψ(•)
is a penalty function and D is a descriptor that locally characterizes image region similarities in images I s and I t . D(P I (x))

denotes the feature descriptor of the pixels in the neighborhood of pixel x in image I, this descriptor being computed with the pixel values of patch P I (x) centered on x.

While the data-term in (8) will take various forms, the regularization term in (1) is kept constant in this contribution. We use the non-local total variation as in [START_REF] Werlberger | Motion estimation with non-local total variation regularization[END_REF][START_REF] Sun | Secrets of optical flow estimation and their principles[END_REF][START_REF] Krähenbühl | Efficient nonlocal regularization for optical flow[END_REF][START_REF] Drulea | Motion estimation using the correlation transform[END_REF][START_REF] Zhang | Robust non-local TV-L1 optical flow estimation with occlusion detection[END_REF] to define the regularization term:

E reg (u) = x∈Ω x ∈N x w x x u x -u x 1 , (9)
where N x is the set of neighbor pixels centered on pixel x and the weights w x x depend on the similarity between pixels x and x . Similarly to [START_REF] Drulea | Motion estimation using the correlation transform[END_REF], the weights w x

x are defined as follows:

w x x = exp       - x -x 2 2σ 2 1 - L(x) -L(x ) 2 2σ 2 2       , (10)
where σ 1 and σ 2 are parameters controlling the similarity measure, and L(x) is the color vector in the CIE Lab color space.

According to [START_REF] Drulea | Motion estimation using the correlation transform[END_REF], by defining a weight matrix W in space

R |Ω|×|N| , W =                   w x ∈N (1,1) (1,1) • • • w x ∈N (1,1) (1,1)
. . . . . . . . .

w x ∈N (N,M) (N,M) • • • w x ∈N (N,M) (N,M)                   , (11)
and a linear operator

K : R 2×|Ω| → R 2×|Ω|×|N| Ku =                   u x ∈N (1,1) -u (1,1) • • • u x ∈N (1,1) -u (1,1) . . . . . . . . . u x ∈N (N,M) -u (N,M) • • • u x ∈N (N,M) -u (N,M)                   ,
where |N| is the size of the neighborhoods N x , the regularization term (9) can be rewritten as

E reg = F(Ku), (12)
with F : R 2×|Ω|×|N| → R is a function defined by

F(z) = W • z 1 . (13)
Therefore, flow field u is computed by solving following optimization problem:

min u F(Ku) + λE data (I s , I t , u). (14)
This optimization can be effectively implemented using the projected-proximal-point algorithm [START_REF] Chambolle | A first-order primal-dual algorithm for convex problems with applications to imaging[END_REF]. The details of the algorithm are beyond the scope of this work. Interested readers can refer to [START_REF] Chambolle | A first-order primal-dual algorithm for convex problems with applications to imaging[END_REF][START_REF] Drulea | Motion estimation using the correlation transform[END_REF] for more details on the mathematical concepts.

The following section details the underlying principles of the descriptor-based data-terms that are robust to illumination changes.

Descriptor-based Illumination-Invariant Data-Terms

The descriptor-based data-term in ( 8) is robust against illumination changes when local descriptor D is illuminationinvariant, i.e., if pixel x in I s and pixel (x + u x ) in I t are two corresponding pixels then D(P I s (x)) = D(P I t (x + u x )).

Local descriptors D are computed with small patches P I (x 0 ) centered at pixel x 0 in image I and with a size of (2k+1)×(2k+1) pixels, where k is a positive integer. P I (x 0 ) can be represented by a vector taking the pixel intensity values as components:

P I (x 0 ) = [I(x 0 ), I(x 1 ), . . . , I(x n )] T ∈ R n+1 , (15)
where n = (2k + 1) 2 -1 is the number of neighbor pixels of x 0 in the patch. The general formulation of descriptors D for any pixel x 0 in image I can be written as:

D(P I (x 0 )) =                           f 0 (P I (x 0 )) f 1 (P I (x 0 )) . . . f m (P I (x 0 ))                           ∈ R m , (16)
where f i : R n+1 → R, i = 0, . . . , m are real functions.

The main issue arising when designing an illumination invariant data-term lies in the definition of the functions f i (•) leading to an accurate OF field.

Illumination change model and illumination-invariance criterion for descriptors

In numerous scenes (e.g. in outdoor or medical scenes) the illumination changes between consecutive images of a videosequence are complex and cannot be represented by a global model. For instance, gradient constancy over complete images (or large image parts) is most often not a realistic illumination change assumption since the illumination variation is barely constant for all homologous pixels of two images. Illumination changes between images should rather be locally modeled to be as general as possible in terms of scene types.

Negahdaripour [START_REF] Negahdaripour | Revised definition of optical flow: Integration of radiometric and geometric cues for dynamic scene analysis[END_REF] simulated the illumination changes by an affine transformation (see [START_REF] Revaud | Epicflow: Edgepreserving interpolation of correspondences for optical flow[END_REF]). This pixel-wise model is able to capture complex variations. However, as shown in [START_REF] Kim | Robust motion estimation under varying illumination[END_REF], integrating this model into a variational calculus algorithm leads to a complex optimization problem and high computation times since, besides a flow vector, parameters a x and b x have to be estimated for each pixel.

An interesting fact is that two consecutive video frames have usually the property of local stationarity such that all pixels in a small image region may share the same parameters in the illumination variation model. Therefore, the assumption is made that the illumination changes of all corresponding pixels in small homologous neighborhoods in two consecutive images can be accurately represented by a model with constant parameters.

Thus, instead of using the linear model in [START_REF] Revaud | Epicflow: Edgepreserving interpolation of correspondences for optical flow[END_REF] at pixel-level, a patch-based model is used to describe locally the illumination changes between homologous neighborhoods:

P I t (x + u x ) = a x P I s (x) + b x , (17)
where a x ∈ R >0 , b ∈ R, P I s (x) and P I t (x + u x ) are two corresponding patches in I s and I t , and u x is the displacement vector at pixel x of image I s . Constant values of parameters a x and b x in small image regions allows for modelling complex illuminations changes between images.

Unlike in the GDIM model in [START_REF] Revaud | Epicflow: Edgepreserving interpolation of correspondences for optical flow[END_REF], parameter a x in the proposed model is assumed to be greater than 0, since the intensity values of pixels are non-negative. Rather than computing the values of parameters a x and b x , the model in [START_REF] Horn | Determining optical flow[END_REF] is only used to design illumination independent descriptors.

The criterion allowing to confirm or to contradict the illumination invariance property of descriptor D is first given by exploiting [START_REF] Horn | Determining optical flow[END_REF].

Definition 1. Descriptor D is illumination-invariant if it is in-variant with respect to the illumination variation model given in [START_REF] Horn | Determining optical flow[END_REF]. In other words, D has to satisfy following constraint:

D(P I (x 0 )) = D(a x 0 P I (x 0 ) + b x 0 ), ( 18
)
for all a x 0 ∈ R >0 , b x 0 ∈ R.

Equation [START_REF] Sun | A quantitative analysis of current practices in optical flow estimation and the principles behind them[END_REF] gives the necessary condition for descriptors to be illumination-invariant.

The literature presents descriptor-based variational OF methods that have only been experimentally proven to be robust to illumination changes [START_REF] Zabih | Non-parametric local transforms for computing visual correspondence[END_REF][START_REF] Stein | Efficient computation of optical flow using the census transform[END_REF][START_REF] Müller | Illuminationrobust dense optical flow using census signatures[END_REF][START_REF] Ranftl | Pushing the limits of stereo using variational stereo estimation[END_REF][START_REF] Demetz | The complete rank transform: A tool for accurate and morphologically invariant matching of structures[END_REF][START_REF] Molnár | Illumination-robust variational optical flow using cross-correlation[END_REF][START_REF] Werlberger | Motion estimation with non-local total variation regularization[END_REF][START_REF] Drulea | Motion estimation using the correlation transform[END_REF][START_REF] Mohamed | Illumination-robust optical flow using a local directional pattern[END_REF][START_REF] Ali | Illumination invariant optical flow using neighborhood descriptors[END_REF].

The data-terms of these methods use patch-based descriptors under the form given in [START_REF] Hu | Efficient coarse-to-fine patchmatch for large displacement optical flow[END_REF]. In sub-sections 3.2 and 3.3, we rewrite some well-known descriptors and explain how it can be mathematically checked whether these descriptors are invariant or not with respect to the illumination change model proposed in [START_REF] Horn | Determining optical flow[END_REF].

Existing Local descriptors

Census transform

The Census transform [START_REF] Zabih | Non-parametric local transforms for computing visual correspondence[END_REF] has been used as an effective solution for computing OF under outdoor lighting conditions [START_REF] Stein | Efficient computation of optical flow using the census transform[END_REF][START_REF] Müller | Illuminationrobust dense optical flow using census signatures[END_REF][START_REF] Ranftl | Pushing the limits of stereo using variational stereo estimation[END_REF]. Writing the Census transform under the form given in [START_REF] Ali | Illumination invariant optical flow using neighborhood descriptors[END_REF] leads to:

D Census (P I (x 0 )) =                           sgn(I(x 0 ) -I(x 1 )) sgn(I(x 0 ) -I(x 2 )) . . . sgn(I(x 0 ) -I(x n ))                           (19)
with

sgn(υ) =          1, υ > 0 0, otherwise. (20)

Complete rank transform (CRT)

The CRT descriptor [START_REF] Demetz | The complete rank transform: A tool for accurate and morphologically invariant matching of structures[END_REF] is computed as follows:

D CRT (P I (x 0 )) =                           n j=0 sgn(I(x 0 ) -I(x j )) n j=0 sgn(I(x 1 ) -I(x j )) . . . n j=0 sgn(I(x n ) -I(x j ))                           . ( 21
)
The i-th component of D CRT gives the number of pixels in P I (x 0 ) having a smaller intensity than pixel x i .

Local Directional Pattern (LDP)

The LPD descriptor [START_REF] Kabir | A local directional pattern variance (LDPv) based face descriptor for human facial expression recognition[END_REF] corresponding to pixel x 0 is a binary code computed for patch P I (x 0 ) of size 3 × 3 pixels and based on the eight Kirsch edge kernels M i , i ∈ {1, 2, . . . , 8} given in Fig. 2. The LDP descriptor is mathematically defined as:

D LDP (P I (x 0 )) =                           sgn(|M 1 ⊗ P I (x 0 )| -m k ) sgn(|M 2 ⊗ P I (x 0 )| -m k ) . . . sgn(|M 8 ⊗ P I (x 0 )| -m k )                           (22)
where operator ⊗ gives the sum of the element-wise product of two matrices, and m k is the k-th largest element in array {|M 1 ⊗ P I (x 0 )|, |M 2 ⊗ P I (x 0 )|, . . . , |M 8 ⊗ P I (x 0 )|}. As shown in Fig. 2, the edge response in direction i is given by M i ⊗ P I (x 0 ) and i = 1 corresponds to the East direction. The edge magnitudes along the i-th direction (defined by pixel pair x 0 and x i ) is used to determine the binary code for pixel x 0 .

Modified Local Directional Pattern (MLDP)

Mohamed et al. [START_REF] Mohamed | Illumination-robust optical flow using a local directional pattern[END_REF] used the descriptor of Section 3.2.3 to obtain the modified local directional pattern (MLDP). Similarly to the LDP descriptor, the MLDP vector corresponds to an 8-bit information generated from the signs of the eight Kirsch edge responses. This descriptor is defined by:

D MLDP (P I (x 0 )) =                           sgn(M 1 ⊗ P I (x 0 )) sgn(M 2 ⊗ P I (x 0 )) . . . sgn(M 8 ⊗ P I (x 0 ))                           , (23)
where sgn(•) is a binary function as defined in [START_REF] Nagel | An investigation of smoothness constraints for the estimation of displacement vector fields from image sequences[END_REF].

Correlation transform

For the correlation transform [START_REF] Drulea | Motion estimation using the correlation transform[END_REF], the descriptor is defined with the mean µ P and variance σ 2 P of the n + 1 pixels intensities of the patch:

D Corr (P I (x 0 )) =                           I(x 0 )-µ P σ P I(x 1 )-µ P σ P . . . I(x n )-µ P σ P                           , (24)
with

µ P = n i=0 I(x i )
n+1 , and

σ 2 P = n i=0 (I(x i )-µ P ) 2 n+1
.

Normalized Neighborhood Descriptor (NND)

In order to compute their descriptor, Ali et al. [START_REF] Ali | Illumination invariant optical flow using neighborhood descriptors[END_REF] used a patch P I (x 0 ) with a minimal size of 5 × 5 pixels (i.e., k ≥ 2).

At central pixel x 0 , a sub-patch P 0 ⊂ P I (x 0 ) of size (2k + 1)

× (2k + 1) is defined (k > 0).
Here, x 0 is also the central pixel of sub-patch P 0 , and the size of P I (x 0 ) and P 0 is constrained by

k = 2k . Likewise, there exists a set of m = (2k + 1) 2 -1 sub- patches {P j | P j ⊂ P I (x 0 )} m j=1
where P j has the same size as P 0 . The pixels within sub-patch P j are denoted by x 1 j , x 2 j , . . . , x m+1 j .

Using theses patches, the NND descriptor is defined as:

D NND (P I (x 0 )) =                             exp -P 1 -P 0 2 h 2 x 0 exp -P 2 -P 0 2 h 2 x 0 . . . exp -P m -P 0 2 h 2 x 0                             , (25)
where

P j -P 0 2 = m+1 i=1 (I(x i j ) -I(x i 0 )) 2 , (26)
and

h 2 x 0 = 3 k=0 P 2k+1 -P 0 2 /4. (27)
In [START_REF] Kumar | A decoupled approach to illumination-robust optical flow estimation[END_REF], P 1 , P 3 , P 5 , and P 7 are 4 sub-patches centered at 4

neighbor pixels of x 0 , namely x 1 , x 3 , x 5 and x 7 , respectively (refer to Fig. 2 for the relative pixel positions).

Descriptor Illumination-Invariance Assessment

Until now, the effectiveness in estimating OF under changing illumination conditions has only been shown experimentally for the descriptors listed in Section 3.2. To the best of our knowledge, none of the corresponding contributions have discussed the mathematical concepts explaining why such descriptors can cope with illumination changes. This sub-section mathematically justifies the appropriateness of these descriptors. Proof. Denote D i as the i-th component of descriptor D and P I (x 0 ) = a x 0 P I (x 0 ) + b x 0 . One needs to prove that for all a x 0 ∈ R >0 and b x 0 ∈ R,

D i (P I (x 0 )) = D i ( P I (x 0 )). (28)
-For D Census , the i-th component can be written as:

D i Census (P I (x 0 )) = sgn(I(x 0 ) -I(x i )). ( 29
)
Due to a x 0 > 0,

D i Census ( P I (x 0 )) = sgn(a x 0 (I(x 0 ) -I(x i ))) = D i Census (P I (x 0 )). (30)
-The i-th D CRT component is:

D i CRT (P I (x 0 )) = n j=0 sgn(I(x i ) -I(x j )). (31)
Because

a x 0 > 0, sgn[(a x 0 I(x i ) + b x 0 ) -(a x 0 I(x j ) + b x 0 )] = = sgn[a x 0 (I(x i ) -I(x j ))] = sgn(I(x i ) -I(x j )). (32)
Therefore, D i CRT ( P I (x 0 )) = D i CRT (P I (x 0 )). -Similarly, each component of D LDP , can be written as:

D i LDP (P I (x 0 )) = sgn(|M i ⊗ P I (x 0 )| -m s k ) (33)
where m s k is the k-th largest element in array

A s = {|M 1 ⊗ P I (x 0 )|, . . . , |M 8 ⊗ P I (x 0 )|}. (34)
Thus,

D i LDP ( P I (x 0 )) = sgn |M i ⊗ P I (x 0 )| -m t k , (35)
where m t k is the k-th largest element in array

A t = {|M 1 ⊗ P I (x 0 )|, . . . , |M 8 ⊗ P I (x 0 )|}. ( 36
)
As it can be seen in Fig. 2, the sum of the elements in the Kirsch matrix M i equals zero for all i = 1, . . . , 8. Thus, it is easy to conclude that

M i ⊗ a x 0 P I (x 0 ) + b x 0 = a x 0 (M i ⊗ P I (x 0 )) . (37)
Therefore, ∀i, j ∈ {1, . . . , 8}, ∀a x 0 > 0,

M i ⊗ P I (x 0 ) ≥ M j ⊗ P I (x 0 ) ⇔ |M i ⊗ P I (x 0 )| ≥ M j ⊗ P I (x 0 ) . ( 38
)
It means that the element order in array A t is the same as the element order in the corresponding array A s . More precisely, if m s k = |M j ⊗ P I (x 0 )| is the k-th largest element in array A s , then its corresponding m t k = |M j ⊗ P I (x 0 )| is also the k-th largest element in array A t . Therefore,

D i LDP ( P I (x 0 )) = D i LDP (P I (x 0 )). (39)
-For D MLDP , the i-th component can be written as:

D i MLDP (P I (x 0 )) = sgn(M i ⊗ P I (x 0 )). (40)
From ( 37) it follows:

sgn(M i ⊗ a x 0 P I (x 0 ) + b x 0 ) = sgn(a x 0 (M i ⊗ P I (x 0 ))).
Therefore, with a x 0 > 0

D i MLDP ( P I (x 0 )) = D i MLDP (P I (x 0 )). (41)
-Component i of descriptor D Corr becomes:

D i Corr (P I (x 0 )) = (I(x i ) -µ P ) /σ P . (42)
One has

D i Corr ( P I (x 0 )) = a x 0 I(x i ) + b x 0 -µ P /σ P , (43)
a x 0 I(x i ) + b x 0 -µ P = a x 0 I(x i ) + b x 0 -a x 0 µ P -b x 0 = a x 0 (I(x i ) -µ P ), (44)
and

σ P = 1 n + 1 n i=0 (a x 0 I(x i ) + b x 0 -µ P ) 2 = 1 n + 1 n i=0 a 2 x 0 (I(x i ) -µ P ) 2 = a x 0 σ P . (45)
By combining [START_REF] Demetz | The complete rank transform: A tool for accurate and morphologically invariant matching of structures[END_REF], [START_REF] Molnár | Illumination-robust variational optical flow using cross-correlation[END_REF], and (45) we deduce

D i Corr ( P I (x 0 )) = D i Corr (P I (x 0 )). (46)
-Finally, the i-th component of descriptor D NND is:

D i NND (P I (x 0 )) = exp -P i -P 0 2 /h 2 x 0 . ( 47
)
Let us consider the NND descriptor on P I (x 0 ). From ( 26) and ( 27) one gets

D i NND ( P I (x 0 )) = exp -P i -P 0 2 / h2 x 0 , (48)
where

P i -P 0 2 = m+1 i=1 (a x 0 I(x i j ) + b x 0 -a x 0 I(x i 0 ) -b x 0 ) 2 = m+1 i=1 a 2 x 0 (I(x i j ) -I(x i 0 )) 2 = a 2 x 0 P i -P 0 2 , (49)
and

h2 x 0 = 3 k=0 P 2k+1 -P 0 2 4 = a 2 x 0 3 k=0 P 2k+1 -P 0 2 4 = a 2 x 0 h 2 x 0 . (50)
By combining [START_REF] Trinh | A general form of illuminationinvariant descriptors in variational optical flow estimation[END_REF], [START_REF] Sun | Secrets of optical flow estimation and their principles[END_REF], and (50), one obtains

D i NND ( P I (x 0 )) = exp -a 2 x 0 P i -P 0 2 a 2 x 0 h 2 x 0 = D i NND (P I (x 0 )). (51)
In this section it was theoretically proven why descriptors

D

Generalized Formulations for Illumination-Invariant

Descriptors

Consider two vectors v s = [v s 0 , v s 1 , . . . , v s n ] T and v t = [v t 0 , v t 1 , . . . , v t n ] T in R n+1 such that v t i = a x v s i + b x , ∀i = 0, 1, . . . , n (52)
with a x ∈ R >0 and b x ∈ R. These vectors correspond to the descriptors of patches P I s and P I t of the source (I s ) and target (I t ) images, respectively. The relationship between the two vectors given in ( 52) is based on the illumination change model given in [START_REF] Horn | Determining optical flow[END_REF].

To propose general forms of illumination-invariant descriptors, let us start with the following lemma.

Lemma 1. Suppose v s and v t are two vectors in R n+1 satisfying [START_REF] Kabir | A local directional pattern variance (LDPv) based face descriptor for human facial expression recognition[END_REF].

If {α 0 , α 1 , . . . , α n } is a sequence of real numbers such that              α 2 0 + α 2 1 + • • • + α 2 n 0 α 0 + α 1 + • • • + α n = 0, ( 53
)
then n i=0 α i v t i = a x n i=0 α i v s i , and (54)
sgn n i=0 α i v t i = sgn n i=0 α i v s i . (55)
Proof. We have,

n i=0 α i v t i = n i=0 α i (a x v s i + b x ) v t i =a x v s i +b x = a x n i=0 α i v s i + b x n i=0 α i =0 . ( 56
)
From n i=0 α i = 0 and a x > 0, one can deduce that

       n i=0 α i v t i = a x n i=0 α i v s i sgn n i=0 α i v t i = sgn n i=0 α i v s i . (57)
Using Lemma 1, we can propose two general forms of illumination-invariant descriptors as follows.

Descriptor Form based on Sign-Invariance

Using Lemma 1, a generalized form of illumination-invariant descriptors can be defined as:

D(P I (x 0 )) =                          
sgn n j=0 α 0, j I(x j ) sgn n j=0 α 1, j I(x j ) . . .

sgn n j=0 α m, j I(x j )                           , ( 58
)
where α i = {α i,0 , α i,1 , . . . , α i,n }, i = 0, 1, . . . , m are the sequences of real numbers satisfying [START_REF] Ali | Anisotropic motion estimation on edge preserving riesz wavelets for robust video mosaicing[END_REF].

With [START_REF] Lurie | 3D reconstruction of cystoscopy videos for comprehensive bladder records[END_REF], it is possible to generate new descriptors by choosing appropriate values for α i . For instance, when

α i = {α i,0 , α i,1 , . . . , α i,n } with α i,0 = 1, α i,i = -1, α i, j = 0 with j {0,
D i (P I (x 0 )) = D i (a x 0 P I (x 0 ) + b x 0 ), ∀i. (59)
To this end, the i-th component in D is defined by:

D i (P I (x 0 )) = Ψ g 1,i (P I (x 0 )) g 2,i (P I (x 0 )) (60)
where Ψ(•) is a non-constant function and g 1,i , g 2,i : R n+1 → R are functions such that

g 1,i (a x 0 P I (x 0 ) + b x 0 ) = h(a x 0 )g 1,i (P I (x 0 )) (61)
g 2,i (a

x 0 P I (x 0 ) + b x 0 ) = h(a x 0 )g 2,i (P I (x 0 )), (62)
with function h(.) taking a x 0 as unique variable. Thus, one has for each component i of a descriptor:

D i (a x 0 P I (x 0 ) + b x 0 ) = Ψ g 1,i (P I (x 0 )) g 2,i (P I (x 0 )) = D i (P I (x 0 )),
and, consequently, D(a x 0 P I (x 0 ) + b x 0 ) = D(P I (x 0 )).

As an example, the following function can be used both as g 1,i and g 2,i in ( 61) and (62) :

g(P I (x 0 )) = γ         L i=1 n j=0 α i, j I(x j ) τ         η , (63)
where the sequences {α i, j } n j=0 , i = 0, 1, . . . , L satisfy condition [START_REF] Ali | Anisotropic motion estimation on edge preserving riesz wavelets for robust video mosaicing[END_REF], and γ, τ and η are fixed positive numbers. In this case, g(a x 0 P I (x 0 ) + b x 0 ) = (a x 0 ) τη g(P I (x 0 )).

Descriptors D Corr and D NND are two specific cases of the generalized form [START_REF] Ben-Hamadou | Flexible calibration of structured-light systems projecting point patterns[END_REF]. For example, for descriptor D Corr , one can see that

D i Corr (P I (x 0 )) = I(x i ) -µ P σ P = Ψ g 1,i (P I (x 0 )) g 2,i (P I (x 0 )) (64)
in which functions Ψ, g 1,i and g 2,i are given by:

                     Ψ(x) = x g 1,i (P I (x 0 )) = I(x i ) -µ P g 2,i (P I (x 0 )) = σ P = 1 √ n+1 n i=0 (I(x i ) -µ P ) 2 1 2 . ( 65
)
We have, I(x i ) -µ P = n j=0 α i, j I(x i ) with the sequence {α i, j } n j=0 is given by

α i, j =              -1 n+1 , if j i n n+1 , otherwise.
(66) It is remarkable that sequence {α i, j } n j=0 satisfies condition ( 53) for all i = 0, 1, . . . , n. Therefore, functions g 1,i and g 2,i in (64) can be rewritten as:

g 1,i (P I (x 0 )) = n j=0 α i, j I(x i ) (67) g 2,i (P I (x 0 )) = 1 √ n + 1 n i=0 n j=0 α i, j I(x i ) 2 1 2 (68)
where α i, j is defined by (66). Referring to (63), it can be seen that in (67), the parameters γ, L, τ, and η are set to 1, while in (68), γ = 1 √ n+1 , L = n, τ = 2, and η = 1 2 .

From the Generalized Formulations to New Descriptors

As shown before, the existing descriptors given in Section 3. As visible in [START_REF] Horn | Determining optical flow[END_REF], a linear model is used to represent the illumination changes. This model leads to high accuracy if it is applied on small enough neighborhoods. Thus, to define the descriptor for pixel x 0 in image I, we consider a 3 × 3 patch P I (x 0 ) = [I(x 0 ), I(x 1 ), . . . , I(x 8 )] T centered at x 0 (see Fig. 3).

First Proposed Descriptor

The first descriptor (D 1 ) based on the elimination of the illumination change model parameters a x and b x is designed as follows:

D 1 (P I (x 0 )) = A x 0 A x 0 2 , ( 69
)
where

A x 0 = [M 1 ⊗ P I (x 0 ), . . . , M 8 ⊗ P I (x 0 )] T is a vector in R 8
and M 1 , M 2 , . . . , M 8 are eight Robinson compass kernels (see Fig. 3).

Note that the sum of the elements in the Robinson compass kernels is 0. Thus, a Robinson kernel corresponds to a sequence of coefficients {α j } 9 j=1 satisfying condition [START_REF] Ali | Anisotropic motion estimation on edge preserving riesz wavelets for robust video mosaicing[END_REF]. The descriptor given in (69) relates to local images structures since each kernel gives a response in one of eight directions around pixel x 0 .

Referring to [START_REF] Ben-Hamadou | Flexible calibration of structured-light systems projecting point patterns[END_REF], one can see that the i-th component in D 1 is given by D

i 1 (P I (x 0 )) = Ψ g 1,i (P I (x 0 )) g 2,i (P I (x 0 )) in which                      Ψ(x) = x g 1,i (P I (x 0 )) = M i ⊗ P I (x 0 ) = 9 j=0 α i, j I(x i ) g 2,i (P I (x 0 )) = 8 i=0 9 j=0 α i, j I(x j ) 2 1 2 . ( 70
)
Functions g 1,i , g 2,i are defined using (63), and {α i, j } 9 j=1 is the sequence of the elements of matrix M i satisfying [START_REF] Ali | Anisotropic motion estimation on edge preserving riesz wavelets for robust video mosaicing[END_REF].

Descriptor D 1 simulates a star-shaped structure which gives normalized grey-level variations in eight directions, along horizontal, vertical and diagonal line segments originating all from the patch center (pixel x 0 ). It captures a local 2D intensity variation information which relates to the shape and sharpness of the textures overlapped by the patch. For homologous points in two images, the descriptor response is similar when the local grey-level distributions are similar. The similarity measurement with this descriptor is robust since intensity differences attenuate the effect of additive intensity changes (elimination of b x 0 ), whereas the normalization with A x 0 limits the effect of multiplicative intensity changes (elimination of a x 0 ).

Second Proposed Descriptor

In this example, the functions in the general form (60) are defined as:

                     Ψ(x) = exp(x)
g 1,i (P I (x 0 )) = I(x i )min(P I (x 0 )) g 2,i (P I (x 0 )) = max(P I (x 0 )min(P I (x 0 )).

(71) This leads to descriptor D 2 :

D 2 (P I (x 0 )) =                          
exp I(x 0 )-min(P I (x 0 )) max(P I (x 0 )-min(P I (x 0 )) exp I(x 1 )-min(P I (x 0 )) max(P I (x 0 )-min(P I (x 0 )) . . . exp I(x n )-min(P I (x 0 )) max(P I (x 0 )-min(P I (x 0 ))

                          . ( 72
)
The idea behind descriptor D 2 is similar to that which motivated the design of D 1 . The main difference is that the structure corresponding to D 2 is not necessarily "star-shaped" because the structure center from which the line segments originate is not x 0 but the pixel with the smallest intensity min(P I (x 0 )) in patch P I (x 0 ). With this minimum, the intensity changes relating to textures are maximized. Effects of additives and multiplicative illumination terms are attenuated by subtracting min(P I (x 0 )) from the patch intensities and by normalization with (max(P I (x 0 )) -min(P I (x 0 ))).

The next section compares the performance of all descriptors discussed in this paper (the descriptors of the literature and the two proposed ones).

Descriptor Performance Comparison

As mentioned in Section 1, the aim of this work is not to propose a competitive OF method. The main objective of this contribution is to propose a theoretical study on the illuminationinvariant descriptors used in variational OF models. For this reason, this paper focuses on the evaluation and comparison of the performance of descriptors by estimating the OF in scenes with and without strong illumination changes.

The proposed descriptors are compared to the descriptors presented in Section 3.2: Census [START_REF] Zabih | Non-parametric local transforms for computing visual correspondence[END_REF], CRT [START_REF] Demetz | The complete rank transform: A tool for accurate and morphologically invariant matching of structures[END_REF], LDP [START_REF] Kabir | A local directional pattern variance (LDPv) based face descriptor for human facial expression recognition[END_REF],

MLDP [START_REF] Mohamed | Illumination-robust optical flow using a local directional pattern[END_REF], Corr [START_REF] Drulea | Motion estimation using the correlation transform[END_REF], and NND [START_REF] Ali | Illumination invariant optical flow using neighborhood descriptors[END_REF]. For a fair evaluation, all descriptors are placed in the same variational OF model as presented in Section 2. The configuration of the variational model is presented in subsection 6.1.

The performance of the descriptors were evaluated on five datasets with known ground truth OF and including images corresponding to five scene conditions:

1. Weak illumination changes combined with small displacements (Section 6.2).

2. Weak illumination changes combined with large displacements (Section 6.3).

3. Strong illumination changes associated with small displacements (Section 6.4).

4. Strong illumination changes associated with large displacements (Section 6.5).

All possible illumination change and displacement combinations (weak or strong illumination changes combined

with small, or large displacements, Section 6.6).

Moreover, this section also proposes a subjective comparison using real endoscopic images. In addition, the effects of some parameters on the performance of each descriptor are reported in this section.

Experimental configuration

In this performance study, all descriptors are placed in exactly the same optimization scheme such that the only factors influencing the differences in the OF field results are the descriptors themselves. Thus, to evaluate the descriptors, the same regularization term (see [START_REF] Baker | A database and evaluation methodology for optical flow[END_REF]) and the same dataterm (see [START_REF] Hu | Efficient coarse-to-fine patchmatch for large displacement optical flow[END_REF]) with a quadratic penalty function

Ψ(v) = v 2
were used as energy in the minimization process (the quadratic penalty function facilitates the optimization problem solving).

Similarly to [START_REF] Drulea | Motion estimation using the correlation transform[END_REF], the minimization in ( 14) is performed with the projected-proximal-point algorithm [START_REF] Chambolle | A first-order primal-dual algorithm for convex problems with applications to imaging[END_REF] and the classical coarse-to-fine warping strategy is used to cope with large displacements. This pyramidal approach uses a classical bilinear interpolation for building the images at different levels, and for In the data-term, the default size of the descriptor patch P I (x 0 ) is 3 × 3 pixels, except for the patch in the NND descriptor where this number is 5 × 5 (a size lower than 5 × 5 is not possible for this descriptor). In the regularization-term, the size of neighborhood N x in ( 9) is systematically set to 5 × 5.

The remaining parameters consisting of σ 1 and σ 2 in (10), λ in [START_REF] Drulea | Motion estimation using the correlation transform[END_REF], and the pyramid scale factor Py s (parameter Py s corresponds classically to the ratio of both the image width and height when passing from level n to level n + 1) are specifically adjusted for each descriptor. The optimal parameter quadruplet (σ 1 , σ 2 , Py s , λ) is obtained for each descriptor by computing the OF for all combinations of following parameter values: σ 1 and σ 2 ∈ {1, 3, 5, 7}, Py s ∈ {0.5, 0.6, 0.7, 0.8, 0.9}, and λ ∈ (0, 140]. The effects of parameters Py s and λ will be studied because they are the most crucial for ensuring robust and accurate flow field estimation.

Weak illumination changes associated with small displacements

The Middlebury training benchmark [START_REF] Baker | A database and evaluation methodology for optical flow[END_REF] The KITTI 2012 [START_REF] Geiger | Are we ready for autonomous driving? the KITTI vision benchmark suite[END_REF] and KITTI 2015 [START_REF] Menze | Object scene flow for autonomous vehicles[END_REF] 1 for all descriptors.

Table 4 shows a comparison of the proposed descriptors D 1 and D 2 with the existing descriptors. It is noticeable that the Corr descriptor is ranked at the first position according to all the quality criteria. Descriptors D 1 and D 2 take the remaining two positions in the top three. The 4th and 5th places were obtained by NND and MLDP. However, there is no significant difference in the quality indexes between the descriptors of the top five.

Globally, the average errors of descriptors Census, CRT and LDP are significantly higher than those of the five first places. Table 4. Results for the KITTI 2012 [START_REF] Geiger | Are we ready for autonomous driving? the KITTI vision benchmark suite[END_REF] and the KITTI 2015 [START_REF] Menze | Object scene flow for autonomous vehicles[END_REF] training datasets. The reported error measures are the mean of average end-point errors (AEE), and the percentage of erroneous pixels over a threshold of 3 pixels (BP3). Noc and Occ stand for non-occluded and occluded, respectively. The illumination changes between the source and target images are controlled using:

Descriptors

I out (x, y) =                      0 if M(x, y).I in (x, y) + C ≤ 0 255 if M(x, y).I in (x, y) + C ≥ 255 [M(x, y).I in (x, y) + C] otherwise (73)
where It is noticeable that, compared to the optimal values obtained for the weak illumination changes given in Table 1, the opti- mal values of the pyramid scale parameter Py s corresponding to descriptors CRT, MLDP and Corr increased from 0.5 to 0.6.

I
The impact of the parameter pair (Py s , λ) on the performance of each descriptor will be discussed in detail in Section 6.6.1. The experimental results are reported in

The results of these experiments with strong illumination

Descriptor comparison on a more general dataset

In subsections 6.2 to 6.5, the descriptor performances were evaluated on images corresponding to four distinct situations obtained by using images with either small or large displacements combined with either weak or strong illumination changes. The aim of this subsection is to assess the descriptor performances on a dataset including simultaneously all the four situations tested separately in previous subsections.

This "more general" dataset consists of nine image pairs involving two benchmarks, namely Middlebury and MPI Sintel.

The first three image pairs are those given in Fig. 4 important parameters on the performance of the descriptors.

Influence of parameters

This subsection investigates the impact of two crucial OF scheme parameters on the descriptor performances. These parameters are the pyramid scale parameter Py s and the tradeoff parameter λ in [START_REF] Drulea | Motion estimation using the correlation transform[END_REF]. It is usually difficult to adjust these two parameters under different displacements between image and for different illumination changes. Therefore, a descriptor should not only be accurate, but (i) its accuracy should as less as possible be affected by scene variations and (ii) the crucial OF scheme parameters should be easily adjustable.

The experiments in Sections 6.2, 6.3, 6.4 and 6.5 have shown that setting the values of σ 1 and σ 2 in (10) to 3 and 5 respectively is appropriate for the four different scene conditions.

These values were also chosen in these experiments. The results in this section were obtained by testing in a systematic way all combinations of Py s and λ. The pyramid scale parameter Py s took values 0.5, 0.6, 0.7, 0.8 and 0.9 successively, while the trade-off parameter λ in ( 14) varied from 0.05 to 140. However, to show more clearly the optimal value range of λ, we only present here, for each descriptor, the experimental results with λ-values around its optimal range. For each descriptor, the AEE and AAE curves according to λ are given for every values of Py s . Thus, for the Census descriptor λ belongs to The OF of every image pair in the "general" dataset was computed for each descriptor and for all (Py s , λ) pairs. The two mean values of the nine AEE and AAE criteria obtained for each (Py s , λ) pair act as descriptor accuracy measure. Fig. 8 shows, for each descriptor and for each value of Py s , the behavior of the mean AEE and AAE errors according to the λ values.

As visible in Fig. 8 for all descriptors, the mean AEE and AAE curves with Py s = 0.5 (light blue plots) are always above the curves corresponding to Py s = 0.6, 0.7, 0.8 and 0.9. These • LDP : [START_REF] Horn | Determining optical flow[END_REF][START_REF] Zhang | Robust non-local TV-L1 optical flow estimation with occlusion detection[END_REF] (interval width : 8, see Fig. 8(c)),

• MLDP : [START_REF] Sun | Pwc-net: Cnns for optical flow using pyramid, warping, and cost volume[END_REF][START_REF] Baker | A database and evaluation methodology for optical flow[END_REF] (interval width : 6, see Fig. 8(d)),

• Corr : [START_REF] Sun | Pwc-net: Cnns for optical flow using pyramid, warping, and cost volume[END_REF][START_REF] Menze | Object scene flow for autonomous vehicles[END_REF] (interval width : 9, see Fig. As visible in Fig. 9, even if the differences between the eight descriptors are very small, the values of the objective criteria By considering globally the results of section 6.6, one can conclude that the eight descriptors can be separated into four groups as follows:

• CRT and LDP are among the group of four descriptors with the least accuracy when the parameters are optimally tuned and their OF is rather less λ independent (scene condition independent).

• Census and NDD are the other two descriptors of the group with the least accuracy when the parameters are optimally tuned, but their OF is more λ independent.

• MLDP and Corr are among the group of four descriptors with the best accuracy when the parameters are optimally tuned. However their OF is rather less λ independent.

• D 1 and D 2 are both in the descriptor group with the best accuracy when the parameters are optimally tuned and in the group of descriptors with a λ independent OF.

Descriptor comparison in very complex scenes

The aim of this section is to test the descriptors on a scene type with extrem conditions (very inhomogeneous scene illumination, strong illumination changes between images, lack of textures and large camera displacements). A test on such scenes is of particular interest since dense OF remains a useful technique in many applications which need to precisely register images under complicated conditions. For instance, the mosaicing of image sequences acquired in gastroscopy [START_REF] Ali | Anisotropic motion estimation on edge preserving riesz wavelets for robust video mosaicing[END_REF] or in white light [START_REF] Miranda-Luna | Mosaicing of bladder endoscopic image sequences: Distortion calibration and registration algorithm[END_REF] or fluorescence [START_REF] Hernández-Mier | 2D panoramas from cystoscopic image sequences and potential application to fluorescence imaging[END_REF] cystoscopy involves such complex conditions.

As visible in Figs. 10( get images in Fig. 10). Objective quality criteria such as AAE or AEE cannot be estimated since no ground truth is available for these images. Thus, to present OF results, we use the classical color code representation [START_REF] Baker | A database and evaluation methodology for optical flow[END_REF] for which the hue and the saturation correspond to the vector orientation and modulus respectively. In regions with constant colors the OF is also constant.

Computation time

This subsection evaluates the run-times requested by each descriptor to determine the descriptor vectors at all pixels in an image (descriptor time in millisecond, ms), as well as the total time for computing the OF using that descriptor (OF time in second, s). ally not corrected to show the image superimposition. However, these discontinuities can be corrected with an approach as described in [START_REF] Weibel | Contrast-enhancing seam detection and blending using graph cuts[END_REF].

Fig. 13 gives another mosaicing result for textureless images. This panorama was built with 161 images and represents scheme parameters. This illustrates the variety of scenes with can targeted by the proposed variational OF method.

In the cystoscopic data of Fig. 14, the endoscopic bladder images include textures whose aspect and contrast vary inside and between the sequences. White light cystoscopy (see Fig. 14(a))

is the standard examination for cancerous bladder lesion diagnosis and patient follow-up. Wide field of view mosaics have several advantages in cystoscopy: the organ inspection is facilitated since the endoscope and potential lesions can be located with respect to anatomical landmarks (this is most often not the case when simply observing video-sequences), two mosaics built at a some week or month interval facilitate the lesion evolution assessment (urologists are unable to assess the evolution by comparing two video-sequences), and the mosaics can be used by different specialists as a medium of concertation (i.e., urologists, surgeons, radio-physicians and oncologists can use a mosaic as a basis of discussion). In the white light modality, although more or less contrasted textures are present, the images are affected by strong illumination changes between different viewpoints and the endoscopic displacement can be large and saccadic, as in gastroenterology. Fig. 14(a) shows that the proposed variational OF scheme can be used to significantly increase the field of view in white light cystoscopy (the dark and bright regions in this mosaic are due to illumination changes). Some cystoscopic set-ups also enable to switch between the white light modality and the fluorescence modality.

In the latter (see Fig. The components of the vectors of both groups are all related to a description of the texture overlapped by a patch. However, the thresholding of real values in the binary descriptors generally leads to a loss of valuable information, especially when the textures are weakly contrasted. For this reason, the real-value descriptor performance is, in average, higher for images with few or poorly contrasted textures. This observation must be taken into account when developing a descriptor.

The proposed method is quite appropriate for scenes without textures or with very weakly contrasted textures. However, this method is not limited to such scenes since it can also be applied to textured scenes. For textured scenes, the proposed OF method can be improved in terms of computation speed and accuracy by combining it with matching-based approaches.

As a perspective for applications in endoscopy, the results of this paper will be used to optimize the descriptor design for the registration of endoscopic images with the aim to build wide-field of view 2D mosaics of the internal bladder wall (cystoscopy) and the internal stomach epithelium (gastroscopy).

The OF will be used to determine the non-linear geometrical transformations between homologous image points with the aim to superimpose the images of the video-sequences.

A natural and further extension of this work is to build 3D mosaics (internal bladder or stomach wall surfaces superimposed by the image textures/colors). To do so, the homologous image points of image pairs will be associated to structure from motion (SfM) or simultaneous localization and mapping (SLAM) techniques. The feasibility of such approaches was shown in [START_REF] Soper | Surface mosaics of the bladder reconstructed from endoscopic video for automated surveillance[END_REF][START_REF] Lurie | 3D reconstruction of cystoscopy videos for comprehensive bladder records[END_REF] for the bladder. However, these approaches use feature extraction methods (as SIFT) in the point correspondence step. Such correspondences can often not be robustly determined due to a lack of pronounced textures in the cystoscopic and gastroscopic images [START_REF] Ali | Fast mosaicing of cystoscopic images from dense correspondence: Combined SURF and TV-L1 optical flow method[END_REF]. The proposed OF method not only determines homologous points in a more robust way, even in the stomach images, but provides also a more dense correspondence field as that given by a SIFT approach (for few and weakly contrasted textures a SIFT method provides few matches). It is noticeable that SfM or SLAM techniques are only based on images. This avoids the use of active vision methods [START_REF] Ben-Hamadou | Flexible calibration of structured-light systems projecting point patterns[END_REF] which require the modification of standard endoscopes [START_REF] Ben-Hamadou | Construction of extended 3D field of views of the internal bladder wall surface: a proof of concept[END_REF] or the use of non conventional devices [START_REF] Soper | Surface mosaics of the bladder reconstructed from endoscopic video for automated surveillance[END_REF].

  Fig. 1. OF results of the CPM-Flow method [8] and the proposed method

Fig. 2 .

 2 Fig. 2. Kirsch edge kernels for the LDP and MLDP descriptor computation. Kernel M i , i ∈ {1, 2, . . . , 8} is used to determine the gradient component along the i-th direction defined by the line passing through points x 0 and x i .

Theorem 1 .

 1 Descriptors D Census , D CRT , D LDP , D MLDP , D Corr , and D NND are illumination-invariant.

  Census , D CRT , D LDP , D MLDP , D Corr and D NND are effective forOF estimation under changing illumination conditions. However, experimental results in[START_REF] Drulea | Motion estimation using the correlation transform[END_REF][START_REF] Mohamed | Illumination-robust optical flow using a local directional pattern[END_REF][START_REF] Ali | Illumination invariant optical flow using neighborhood descriptors[END_REF] demonstrated that there is a significant difference in the accuracy of OF fields estimated by these descriptors. Thus, the design of illumination invariant descriptors leading to accurate flow fields and involving a limited computational complexity remains an open challenge. The next section provides two generalized types of illumination-invariant descriptors that will facilitate the design of new descriptors that are robust to illumination changes.

Fig. 3 .

 3 Fig. 3. Robinson kernels used to define the first descriptor. The patch corresponds to a 3 × 3 neighborhood around pixel x 0 .

2 are

 2 particular cases of the two generalized formulations given in Section 4. Using these two formulations, one can design new illumination-invariant descriptors. Note that the generalized formulations are the starting points for the construction of illumination invariant descriptors. The main issue now lies in the choice of appropriate functions and parameters in[START_REF] Lurie | 3D reconstruction of cystoscopy videos for comprehensive bladder records[END_REF] for the sign invariance descriptors and in[START_REF] Ben-Hamadou | Flexible calibration of structured-light systems projecting point patterns[END_REF] for the illumination model parameter elimination approach. This section introduces, as examples, two novel descriptors derived from the second general form (60)-(63). Similarly, other descriptors could be derived from the first general form[START_REF] Lurie | 3D reconstruction of cystoscopy videos for comprehensive bladder records[END_REF]. The assessment of the performance of the proposed descriptors is given in Section 6.

  datasets are composed of real-world images taken from a driving platform and consist of 194 and 200 image pairs, respectively. All sequences are with known ground truth and, even if the illumination changes are rather weak in average, these changes are varying along the image sequences since they are arising in real scenes. However, even for the large displacements between images, the illumination changes remain moderate. The AEE and BP3 (the percentage of bad-pixels which have an end-point error above 3 pixels) are the two quality criteria used for an objective comparison. Four sequences in the KITTI 2012 training set, namely sequences 11, 15, 44, and 74, were chosen to determine the optimal (default) values of parameters σ 1 , σ 2 , Py s and λ for each descriptor. The adjusted parameter combinations are those leading to the lowest mean BP3 value in the non-occluded areas. The default values for the parameters are shown in the KITTI column of Table

6. 4 .

 4 Strong illumination changes associated with small displacements In order to compare descriptor performance in this case, we report experimental results on image pairs with simulated illumination changes. The RubberWhale image pair of the Middlebury training dataset was used in these experiments. As visible in Fig. 4(d) and Fig. 4(e), this image pair exhibits few illumination changes and small displacements. Illumination variations were simulated and applied on this pair to obtain images with strong illumination changes. The ground truth OF of these image pairs being available, the AEE and AAE criteria are again used to quantify the descriptor performance.

• Case 2 :Fig. 4 (• Case 3 : 2 ,

 2432 Fig. 4(g), where the illumination of the target image is changed by the multiplier M illustrated in Fig. 4(a), and by additive factor C = 20.

Fig. 4 .

 4 Fig. 4. Experiments with simulated illumination changes. Images (a), (b) and (c) represent the multiplicative masks used to generate the synthetic illuminations in (g), (h) and (i), respectively. Sub-figures (d) to (i) present the source-target image pairs used to evaluate the performance of the descriptors and their stability against illumination changes in terms of OF accuracy. (d) and (e) : image pair simulating weak illumination changes. (f) and (g) : strong illumination change in the target image center. (h) and (i) strong vertical illumination intensity gradient between the source and target images. (j) Intensity values of the two 3 × 3 homologous regions marked by squares in Figs. 4(h) and 4(i). The central grid gives the intensity ratios of homologous pixels of these regions.

Fig. 5 . 2 Fig. 6 . 6 . 5 .Fig. 7 )

 526657 Fig. 5. Performances and stability in terms of AEE and AAE values obtained for different descriptors with the three different illumination change simulations represented in Fig. 4.

  (modified Middleburry images) and including either weak or strong illumination changes combined with small displacements. The other six image pairs are extracted from two final pass sequences of the MPI Sintel training dataset [10]. More specifically, four consecutive images from frame-0039 to frame-0042 of the Bamboo-2 sequence (see Figs. 7(a)-7(d)) represent three image pairs, and four other consecutive images from frame-0001 to frame-0004 of the Shaman-3 sequence (see Figs. 7(e)-7(h)) correspond to the last 3 image pairs. These MPI-Sintel data were chosen since they include both small and strong illumination changes associated to large displacements. The AEE and AAE metrics are used to objectively evaluate the accuracy of the OF field. Besides the OF accuracy tests, this custom-made dataset is also used to evaluate the effects of some

Fig. 7 .

 7 Fig. 7. Four consecutive images in the final pass of the Bamboo-2 and Shaman-3 sequences of the MPI Sintel training dataset [10]. The first row includes the Bamboo images and the second row shows the Shaman images.

[ 1 : 1 :

 11 40], where [1:1:40] means that λ takes the values ranging from 1 to 40 and are changing with a step of 1. Similarly, λ ∈ [0.05:0.05:5] for CRT, λ ∈ [1:1:40] for LDP, λ ∈ [0.05:0.05:1] ∪ [2:1:100] for MLDP and Corr, λ ∈ [1:1:140] for NND, λ ∈ [1:1:100] for D 1 , and λ ∈ [1:1:40] for D 2 .

2 Fig. 8 .

 28 Fig. 8. OF accuracy according the values of parameters Py s and λ. The mean AAE and AEE values are computed for the dataset described at the beginning of Section 6.6 and plotted according λ for different values of parameter Py s .

  Fig.8(f), D 1 in Fig.8(g) and D 2 in Fig.8(h)) that the plateaus are preserved when Py s increases from 0.6 to 0.9. One can also notice that an increase of the Py s parameter value also increases the number of levels in the pyramid scheme. For this reason, a solution to minimize the computation time while preserving the CPME in the AEE and AAE plots is to set Py s to 0.6.For a given descriptor, the optimal value of λ should be in

6 . 6 . 2 .

 662 8(e)),• NND : [75, 120] (interval width : 45, see Fig.8(f)),• D 1 :[START_REF] Nagel | An investigation of smoothness constraints for the estimation of displacement vector fields from image sequences[END_REF][START_REF] Ben-Hamadou | Flexible calibration of structured-light systems projecting point patterns[END_REF] (interval width : 40, see Fig. 8(g)), and • D 2 : [12, 31] (interval width : 19, see Fig. 8(h)). It is noticeable that, for Py s = 0.6, descriptors Census, NND, D 1 , and D 2 lead to constantly high accuracy for large ranges of λ-values. In practice it means that the optimal λ value is easily adjustable and that these descriptors are the less sensitive to illumination changes. Descriptor LDP has also high accuracy in a large range of λ-values when Py s = 0.9 (Dark blue curve in Fig. 8(c)). Descriptors CRT, MLDP, and Corr have small λ intervals with high accuracy and seem less appropriate when the illumination and displacement magnitudes change from one scene type to another. Descriptor performance with optimal parameter settings The mean AAE and AEE values over the nine image pairs of the "more general" dataset (as described at the beginning of Section 6.6) are used to compare the performance of the eight descriptors. These values are computed using the parameter set optimized in Subsection 6.6.1 (i.e., σ 1 = 3, σ 1 = 5, Py s = 0.6 for all 8 descriptors, while λ = 20, 0.8, 25, 8, 12, 100, 50, and 15 for descriptors Census, CRT, LDP, MLDP, Corr, NND, D 1 , and D 2 , respectively).

Fig. 9 .

 9 Fig. 9. Objective comparison between the eight descriptors on the synthetic dataset.

2 Fig. 10 .

 210 Fig. 10. Experimental result on a gastroscopic image pair. (a) and (b) are two endoscopic images of the pyloric antrum region (without textures and with illumination changes) used to compute OF fields. These images were provided by Pr. Dominique Lamarque (Ambroise Paré Hospital, Boulogne-Billancourt, France). Figures from (c) to (j) give the OF results for each descriptor using the classical color code.

  a)-10(b) and Figs. 11-12, endoscopic images often include regions without textures, inhomogeneous scene illumination, varying lighting conditions according the viewpoint and/or specular reflections. In order to demonstrate the potential and robustness of the proposed descriptors in such medical scenes, this section gives OF results for an endoscopic image pair of the pyloric antrum region (see the source and tar-

  The homologous point displacements between the two images in Figs. 10(a)-10(b) are caused both by inner stomach movements and the endoscopic camera movements. In practice, due to the varying depth of the surface, the modulus of the flow vectors should be larger for tissue which is close to the camera (e.g. pixels in the down left image corner in Fig.10(a)), and smaller for tissue which is farther from the camera (e.g. pixels in the image center). Although there is a difference in magnitude and direction between the motion vectors of different pixels, the motion fields are smooth for such scenes. Moreover, the motion vectors inside the black "hole" (duodenum) and the flow vectors at the circular border of this "hole" (pyloric sphincter limits) have to be equal since the movement is constant in this region. As noticeable in Fig.10, the OF results of Census, NND, D 1 and D 2 globally correspond to realistic flow fields since i) their OF is smooth, and ii) the flow field is constant onto and around the duodenum black disc (in Figs. 10(g)-10(j) the absence of the disc shows that flow field is effectively constant in this region). These results are globally less coherent for the CRT, LDP, MLDP and Corr descriptors. In particular, the light and dark blue discs in Figs. 10(c)-10(f) indicate a flow discontinuity which does not correspond to a real displacement continuity.

Fig. 11 (

 11 Fig. 11(f) shows a panoramic image computed with 21 images, while Fig. 12(e) shows the mosaicing result computed for a sequence of 45 images. The precise mosaics (without structure discontinuities) of Fig. 11 and Fig. 12 confirm the potential and robustness of the proposed decriptor 2 . In these mosaicing examples, the illumination discontinuities were intention-

Fig. 11 .Fig. 12 .

 1112 Fig. 11. Mosaic built with 21 images of a gastroscopic sequence (5 images of the sequence are also given in the first row). The OF between the images was directly used to place the pixels of the 21 frames into the coordinate system of the first image which acts as mosaic reference.

Fig. 13 .

 13 Fig. 13. Mosaic of a landscape of the red planet. This mosaic was built with 161 images almost weak textures.

  14(b)), the natural aspect of the bladder colours and textures is lost, but the detection of cancerous lesions located under the epithelium becomes possible. In urology, white light and fluorescence endoscopy are two complementary image modalities. The mark delineated by a circle in the mosaic of Fig. 14(b) is due to a transurethral resection of a lesion. The mosaic is used to check whether the whole tumoral tissue was removed or not during the surgical intervention. It is worth noticing that the same descriptor and OF scheme parameters (that used for the gastroscopic data and given at section beginning) were again used to build both the white light and fluorescence mosaics.

  This contribution shows how illumination invariant dataterms can be constructed in the frame of OF algorithms. The illumination changes between images are modelled with affine transformations which locally link the intensities in small images regions (in this paper the affine transformation parameters are constant in a 3 × 3 neighborhood). These local affine transformations enable to deal with complicated (i.e. strong and/or inhomogeneous) illumination changes between images. It was also shown how it can be verified whether a descriptor used to construct a data-term is illumination invariant or not with respect to the proposed model. Until now, this invariance was most often experimentally highlighted (i.e. without mathematical proof).The major contribution of this paper lies in the development of two general formulations of illumination invariant descriptors, one formulation being based on a descriptor sign invariance, while in the other formulation the idea is to eliminate the parameters of the affine transformation. These two formulations help to understand the appropriate way to construct illumination invariant descriptors and act as a starting point for constructing accurate and robust data-terms for complex scenes.

Fig. 14 .

 14 Fig. 14. Two mosaics of the inner bladder wall epithelium. (a) Wide field of view map constructed with 200 images acquired in the white light modality. In this reference modality in urology, a clinician (urologist or surgeon) can, due to the natural colors and textures, interpret the scene content in order to (mentally) localize the endoscope with respect to anatomical landmarks. This modality allows for the diagnosis of multifocal cancerous lesions which pierce the epithelium. (b) Panoramic mosaic constructed with 95 images acquired in the fluorescence modality. This complementary modality facilitates the detection of tumoral lesions potentially located under the epithelium and which are not visible in the white light modality.

Table 1 .

 1 Default values of parameters σ 1 , σ 2 , Py s and λ obtained with the Middlebury and KITTI training datasets.

			Middlebury		KITTI
	Desc.	σ 1	σ 2	Py s λ	σ 1	σ 2	Py s λ
	D 1	3	5	0.8 50	1	5	0.9 50
	D 2	3	5	0.7 15	1	5	0.9 30
	Census 3	5	0.8 20	1	5	0.9 30
	CRT	5	7	0.5 0.8 1	5	0.9 1.1
	LDP	5	7	0.8 17	1	5	0.7 25
	MLDP 3	5	0.5 9	1	5	0.8 7
	Corr	3	5	0.5 12	1	5	0.9 10
	NND 3	5	0.7 100 1	5	0.9 75

up-sampling the flow field to the finer level. At each pyramid level, 5 warps and 40 iterations per warp are used to optimize the energy. As recommended in

[START_REF] Sun | A quantitative analysis of current practices in optical flow estimation and the principles behind them[END_REF]

, a median filter of size 5×5 is applied to the intermediate flow results after every warping iteration.

Table 2 .

 2 AEE and AAE values obtained for the Middlebury training database. Parameters σ 1 , σ 2 , λ and Py s were set on the default values given in Table1.

		Dimetrodon	Grove2	Grove3	Hydrangea RubberWhale	Urban2	Urban3	Venus
	Des.
		AEE AAE AEE AAE AEE AAE AEE AAE AEE AAE	AEE AAE AEE AAE AEE AAE
	D 1	0.11 2.09 0.13 1.80 0.46 4.79 0.17 2.07 0.08 2.68	0.35 3.30 0.48 3.55 0.25 3.88
	D 2	0.20 4.37 0.13 1.76 0.48 5.08 0.18 2.16 0.08 2.68	0.35 3.48 0.71 5.06 0.25 3.75
	Census 0.26 5.51 0.15 2.01 0.52 5.34 0.19 2.21 0.09 2.73	0.37 3.50 0.68 4.79 0.91 7.48
	CRT	0.22 4.82 0.17 2.24 0.56 5.83 0.19 2.17 0.23 2.99	0.44 3.98 0.71 5.02 0.30 4.06
	LDP	0.14 2.83 0.26 3.32 1.03 9.31 0.22 2.49 0.11 3.36	0.64 4.94 2.8	11.32 7.25 18.48
	MLDP 0.13 2.35 0.13 1.84 0.48 5.02 0.17 2.06 0.08 2.71	0.33 3.18 0.57 4.18 0.26 3.68
	Corr	0.23 4.86 0.20 2.53 0.49 5.21 0.17 2.01 0.08 2.60	0.34 3.65 0.75 5.67 0.28 4.26
	NND	0.22 4.62 0.19 2.77 0.63 6.50 0.18 2.29 0.09 3.06	0.52 3.93 0.59 4.28 0.34 5.3

Table 3 .

 3 Descriptor ranking according to the accuracy (the average errors) on the eight image pairs of Table2.

	Rank	Average AEE	Average AAE
	1	D 1	0.255 D 1	3.02
	2	MLDP 0.268 MLDP 3.12
	3	D 2	0.297 D 2	3.54
	4	Corr	0.320 Corr	3.85
	5	NND	0.345 CRT	3.89
	6	CRT	0.352 NND	4.09
	7	Census 0.396 Census 4.19
	8	LDP	1.556 LDP	7.00
	6.3. Weak illumination changes associated with large dispal-
	cements

Table 5 .

 5 By consid-

ering only this table, one can conclude that four descriptors (D 1 , D 2 , Corr and MLDP) outperform the remaining four descriptors which are almost systematically among the least accurate for all images. One can notice that Corr notably exhibits high accuracy in this table.

Table 5 .

 5 Results for the Bamboo 2 and Shaman 3 sequences of the MPI Sintel[START_REF] Butler | A naturalistic open source movie for optical flow evaluation[END_REF] training dataset. These sequences include illumination changes.

	Clean

  The evaluation was performed on the Middle-bury Urban2 image pair with a size of 640 × 480 pixels and using MATLAB R2016a running under Windows 10, 64bits on a HP Desktop Core i7-4910MQ (2.90 GHz, 16 GB RAM). Descriptor vectors are completely computed with the Matlab code, while the core of the algorithm optimizing the OF (see[START_REF] Drulea | Motion estimation using the correlation transform[END_REF]) is available as a .mex file (it is first written with the C++ language, and then compiled into a .mex code to be called in the Matlab environment). To compute the OF time, parameters σ 1 , σ 2 , Py s , and λ are set based on the results in Section 6.6.1. Specifically, the three parameters σ 1 , σ 2 , Py s are set to 3, 5 and 0.6, respectively. The λ parameters of descriptors D 1 , D 2 , Census, Subsections 6.2 to 6.5). However, when considering the criterion relating to OF accuracy independence towards the λ parameter, the two proposed descriptors can be considered as being more robust than MLDP and Corr. Moreover, in the experiment on gastroscopic images shown in Fig.10, the pro-

		posed descriptors gave subjectively more coherent results than
		Corr and MLDP.
		The ability to deal both with textureless scenes and with
		scenes with various textures is illustrated and confirmed in the
		next section on different real-data scenes which are not repre-
		sented in the standard benchmarks (Middlebury, KITTI 2012,
	CRT, LDP, MLDP, Corr, and NND are set to 50, 15, 20, 0.8, 25,	KITTI 2015 and MPI-Sintel).
	8, 12, 100, respectively.
	The descriptors are ordered according to increasing compu-	7. Applications to Image Mosaicing in Gastroscopy and
	tation time:	Other Scenes
		A major challenge in an image mosaicing application is
	MLDP (91.68 s), Corr (92.14 s), NND (94.22 s), LDP
	(94.98 s), and CRT (96.72 s).
	6.9. Global discussion
	The experimental results have shown that when a descriptor
	fulfills condition (18) it is indeed illumination invariant. With
	the variational OF model presented in this work, all the eight
	descriptors can reach their best OF accuracy with the pyramid
	scale parameter Py s ≥ 0.6.
	Considering globally Sections 6.2, 6.3, 6.4, 6.5 and 6.6, de-
	scriptors Census, CRT and LDP led to the less accurate OF ac-
	curacy in comparison to the other descriptors. The NND de-
	scriptor accuracy is nearly equivalent to that of the proposed
	descriptors (D 1 and D 2 ) on the KITTI training datasets (see Ta-
	ble 4), whereas its ranking was not too high on the Middlebury

-Descriptor time : MLDP (23.5 ms), Census (42 ms), D 1 (44.5 ms), Corr (57 ms), D 2 (60 ms), LDP (81 ms), NND (104 ms), and CRT (119 ms).

-OF time : D 2 (87.22 s), D 1 (87.67 s), Census (89.12 s), training set (see Table

3

). Moreover, the computation time of NND is higher (104 ms) than that of the proposed descriptors.

Descriptors MLDP, Corr, D 1 and D 2 can be considered as being globally the most accurate for particular illumination change conditions (to register accurately and robustly a long sequence of images with poor quality. This subsection shows the mosaicing results for two image sequences (two patients) of the pyloric antrum region (see Figs.

11-12

). The video-sequences were provided by Pr. Dominique Lamarque (Ambroise Paré Hospital, Boulogne-Billancourt, France). The image registration is based on the flow field obtained with descriptor D 1 (i.e. the IIOF-NLDP scheme) using the parameter values as in Table

1:

(σ 1 , σ 2 , Py s , λ) = (3, 5, 0.8, 70).

Two videos illustrating the image mosaicing process from which Figs. 11 and 12 have been built are provided as supplementary material.

Acknowledgments

This work is funded by the Agence nationale de la recherche (ANR) in the framework of the EMMIE (Endoscopie Multi-Modale pour les lésions Inflammatoires de l'Estomac) project (project number: ANR-15-CE17-0015).