Thibault Bourgeron
email: bourgeron@ann.jussieu.fr

Carlos Conca
email: cconca@dim.uchile.cl

AND Rodrigo Lecaros

DETERMINING THE DISTRIBUTION OF CALCIUM CHANNELS IN THE OLFACTORY SYSTEM *

Keywords: Distribution of Calcium Channels, Olfactory System, Inverse Problem, Fredholm Integral Equation, Mellin Transform AMS subject classifications. 35Q92, 45Q05, 65R20

.

 for an approximated kernel. We also introduce a better approximated kernel for which we prove continuity and observability results, and for which we perform numerical computations using experimental data. Lastly, we prove that no observability inequality, in certain weighted L 2 spaces, can be found for the original model. Throughout this paper we present numerical simulations illustrating the theoretical results.

1. Introduction. Identification of detailed features of neuronal systems is a major issue in the biosciences for the coming years. In this respect, inverse problem methods have already proved to be efficient. As a contribution to this field, making use of a model developed in [START_REF] French | Numerical approximation of solutions of a nonlinear inverse problem arising in olfaction experimentation[END_REF], this article proposes an improved method to determine the spatial distribution of CNG ion channels along the length of a cilium. We also apply our method to experimental data published in [START_REF] French | Integral equation models for the inverse problem of biological ion channel distributions[END_REF].

Cilia are long thin cylindrical structures that extend from an olfactory receptor into the nasal mucus. When an odorant molecule binds to an olfactory receptor in a cilium membrane, it successively activates two enzymes, which results in an increase in the concentration of cyclic adenosine monophosphate (cAMP) concentration within the olfactory receptor neuron. Some of the cAMP binds to cyclic nucleotide-gated (CNG) ion channels, causing them to open. This allows a depolarizing influx of Na + ions to flow into the cell, which causes the neuron to depolarize. Although the single-channel properties have been well described, the distribution of these channels along the cilia remains still widely unknown and may well turn out to be crucial in determining the kinetics of the neuronal response.

Experimental procedures to isolate a single (grass frog) cilium have been developed, [START_REF] Kleene | Origin of the chloride current in olfactory transduction[END_REF][START_REF] Kleene | Transmembrane currents in frog olfactory cilia[END_REF][START_REF] Kleene | Calcium-activated chloride conductance in frog olfactory cilia[END_REF][START_REF] Kleene | An electrophysiological survey of frog olfactory cilia[END_REF][START_REF] Flannery | Clustering of cyclic-nucleotide-gated channels in olfactory cilia[END_REF]. One olfactory cilium is drawn into a pipette which is then moved to a pseudointracellular bath which contains no cAMP. The pipette containing the cilium is then transferred to a bath containing cAMP. Contact with the bath initiates the diffusion of cAMP into the cilium. The current through the cilium is recorded.

A very natural issue is whether it is possible to determine the CNG channel distribution along the length of a cilium from current experimental data. In [START_REF] French | Numerical approximation of solutions of a nonlinear inverse problem arising in olfaction experimentation[END_REF] the authors propose a mathematical model for the dynamics of cAMP concentration in this experiment, consisting of two nonlinear differential equations and a constrained Fredholm integral equation of the first kind. Numerical methods to compute the channel distribution were proposed in [START_REF] French | Numerical approximation of solutions of a nonlinear inverse problem arising in olfaction experimentation[END_REF], [START_REF] French | Perturbation approximation of solutions of a nonlinear inverse problem arising in olfaction experimentation[END_REF]. Some computations indicated that this mathematical problem is highly ill-conditioned.

To tackle the problem mathematically, [START_REF] French | Integral equation models for the inverse problem of biological ion channel distributions[END_REF] introduced some simplifications in the model of the experiment. The simplified inverse problem, which interests us in this work, consists in determining the function ρ, representing the distribution of the CNG channels, from the measurement of the electrical current I 0 [ρ], given by

I 0 [ρ](t) = L 0 ρ(x)H 0 (c(t, x)) dx, (1.1)
for t 0, where the Hill function H 0 is defined by

H 0 (x) = x n x n + K n 1/2
.

(1.2)

In (1.1), c(t, x) denotes the concentration of cAMP which diffuses along the cilium. The Hill function appearing here is widely used in biochemistry to model the binding of a ligand to a macromolecule, here the binding of cAMP molecules to the cilium.

After [START_REF] French | Integral equation models for the inverse problem of biological ion channel distributions[END_REF][START_REF] Conca | Determination of the calcium channel distribution in the olfactory system[END_REF] we will assume that c(t, x) is known, analytically, and given by:

c(t, x) = c 0 erfc x 2 √ Dt . (1.3)
Here, erfc = 1 -erf and erf is the Gauss error function.

In (1.2), K 1/2 > 0 is a constant corresponding to the half-bulk concentration, and c 0 the concentration of cAMP in the solution with which the pipette comes into contact.

All the parameters have been experimentally estimated, see [START_REF] French | Integral equation models for the inverse problem of biological ion channel distributions[END_REF][START_REF] Badamdorj | Identification of Cl(Ca) channel distributions in olfactory cilia[END_REF] and the references therein, to be: n = 1.7, K 1/2 = 0.17 µM. For the numerical simulations throughout this paper we will also use c 0 = 1 µM, D = 1 m 2 s -1 . This value of c 0 has been chosen so as to show the interesting behaviours in the figures. The value of D > 0 does not play an important role in the analysis we propose, cf. (2.6).

In [START_REF] French | Integral equation models for the inverse problem of biological ion channel distributions[END_REF] an explicit formula for ρ is given if H is taken to be a single Heaviside function (H(x) = a1 x>α ). Paper [START_REF] Conca | Determination of the calcium channel distribution in the olfactory system[END_REF] presents identifiability, stability and reconstruction results when H is a finite sum of Heaviside functions.

The aim of this article is to give a complete solution of this inverse problem. We make use of the Mellin transform to invert (1.1), which allows us to complete the results obtained in [START_REF] Conca | Determination of the calcium channel distribution in the olfactory system[END_REF], in a simple way. We then introduce a better approximation of H 0 and we prove continuity and observability results for the corresponding operator I. For this operator we also perform numerical computations using experimental data. Lastly we prove that for the original model (1.1), (1.2), (1.3) presented above, no observability inequality can be found, in certain weighted L 2 spaces.

As noted in [START_REF] Conca | Determination of the calcium channel distribution in the olfactory system[END_REF], (1.1) is a Fredholm integral equation of the first kind. The associated inverse problem is in general ill-posed. For instance, if the kernel of an integral operator is sufficiently smooth, then this operator is compact from L 2 to L 2 , thus, even if it is injective, its inverse cannot be continuous. The mathematical results contained in this article can be understood in this way. Depending on the regularity of the kernel of the operator (1.1) (i.e. depending on the choice for the approximation of the Hill function H 0 ), the integral operator (1.1) has a more or less strong regularizing effect. Therefore the inverse problem is more or less ill-posed, see [START_REF] Engl | Regularization of inverse problems[END_REF]. Precisely in Propositions 3.4, 3.5, the kernel is at most continuous and the degree of ill-posedness of the inverse problem is 1. Whereas in Proposition 5.7 the kernel belongs to some space of smooth functions (denoted S[0, ∞), see Definition 5.4) and its degree of ill-posedness is infinite.

The paper is organized as follows. Section 2 introduces the main idea for using the Mellin transform to invert the integral operator I. Section 3 states continuity and observability inequalities for two different approximations of the Hill function and a non observability result for the genuine Hill function. Section 4 presents some numerical simulations performed with experimental data. Sections 5 and A (in the Appendix) give the proofs of these results. Finally the Mellin transform is briefly recalled Section B, in the Appendix.

2. Notations and use of the Mellin transform. The aim of this section is to present the strategy used to find the function ρ > 0, which is the distribution of the CNG channels, from the measurement of the current I 0 [ρ]. Before studying this problem we introduce two approximations of I 0 [ρ], obtained by approximating the Hill function H 0 .

The first approximation of the Hill function H 0 is linear. This case has already been studied in [START_REF] Conca | Determination of the calcium channel distribution in the olfactory system[END_REF]. We focus on it here both to explain our method and because all the computations can be performed explicitly. Definition 2.1. For (α k ) a finite increasing sequence of elements of (0, c 0 ) and a k positive numbers, let:

H 1 (x) = m k=1 a k 1 x α k .
(2.1)

As t tends to ∞ the factor x/ √ Dt appearing in (1.3) tends to 0. Therefore we consider the following approximation of H 0 . Definition 2.2. For positive real numbers a < c 0 , K, let:

H (a) 2 (x) = x n x n + K n 1/2 1 x a + K1 a<x c0 .
(2.

2)

The function

H (a) 2
is represented in Figure 2.1. For the sake of simplicity we will often drop the superscript in the notations H

I l [ρ](t) = L 0 ρ(x)H l (c(t, x)) dx, (2.3)
where c(t, x) is still defined by (1.3).

Appendix B contains some details regarding the Mellin transform and the definition of multiplicative convolution. As the Mellin transform of a multiplicative convolution is the pointwise product of Mellin transforms, it can be used to invert the operator I l [ρ]. Before going further, let us introduce some convenient notations.

Notations. For l = 0, 1, 2, let:

G l (x) = H l (c 0 erfc(2 -1 D -1/2 x -1 ))
and J l (x) = H l (c 0 erfc(x)).

(2.4)

The functions G l and J l and their Mellin transforms are linked by the relations:

G l (x) = J l (2 -1 D -1/2 x -1
) and MG l (s) = 2 s D s/2 MJ l (-s).

(2.5)

The relation (2.3) defining I l can then be inverted:

ρ(x) = 2 -q-1 D -q/2 M -1 q M I l [ρ] s 2 MJ l (-s) . (2.6)
The main idea of this article is that the quantities I l [ρ] can be written as multiplicative convolutions, and the following formal identity holds. For l = 0, 1, 2:

I l [ρ](t 2 ) = ∞ 0 xρ(x)1 x∈[0,L] G l t x dx x = xρ(x)1 x∈[0,L] * G l , (2.7)
where * is the multiplicative convolution.

Therefore finding continuity or observability inequalities for the operator I l [ρ] is reduced to bounding MJ l (s), from above or from below respectively, on the vertical line on which the inverse Mellin transform is taken. The next section presents the bounds obtained in each case l = 0, 1, 2, and the continuity/observability results they imply.

3. Main results. Before stating our main results, let us introduce some weighted L p spaces for which the Mellin transform is well-suited.

Notations. For real numbers p 1, q and for a function f : R → R we denote f L p q = R |f (x)| p |x| q dx 1/p . We denote L p q the Banach space L p ([0, ∞), x q dx) endowed with the norm

L p q .
For a function f : q +iR → C we denote f L p (q+iR) = R |f (q + it)| p dt 1/p . We also denote L p (q + iR) the Banach space L p (q + iR, dx) with the norm L p (q+iR) . As stated in the previous section, continuity/observability results for the operators I l (l = 0, 1, 2), are consequences of the behaviours of MJ l on some vertical lines. The following lemma makes this link precise. Lemma 3.1. For l = 0, 1, 2, the operators I l are defined by (1.2), (1.3), (1.1), and (2.1), (2.2), (2.3), and the functions J l are defined by (2.1), (2.2), (2.4)

. Let k ∈ N, q ∈ R, C > 0 such that k-1 j=0 (s + 2j) MJ l (s) C (resp. C) on - q + 1 2 + iR.
Then:

(I l [ρ]) (k) L 2 2k-2+ q+1 2 C2 -k+q/2 D (q-1)/4 ρ L 2 q resp. C2 -k+q/2 D (q-1)/4 ρ L 2 q .
As a consequence of Lemma 3.1 the strategy to prove continuity/observability results for the operators I l , in a weighted L 2 space, is to find an integer power k ∈ N such that s k MJ l (s) is bounded from above/below on a vertical line, whose abscissa depends on the weight.

For the approximated operators I 1 and I 2 the integer k appearing in the lemma can be taken to be 1, whereas for the non approximated operator I 0 , MJ 0 (s) decays faster than any s k . The converse of Lemma 3.1, which is true under an extra assumption on MJ l (cf. Proposition 3.6 and its proof, Section 5.5), implies that no observability inequality of the previous type holds. These results are proved in Section 5. We now state the results for each case.

β k = erfc -1 α k c 0 , q 1 = ln 1 a1 m k=2 a k ln β 1 -ln β 2 > 0, (3.1)
where the numbers a k , α k are defined by (2.1).

Lemma 3.3. Let the β k 's and q 1 be defined by (2.1), (3.1) and, let J 1 be defined by (2.1), (2.4). For s in the right half-plane, the Mellin transform of J 1 is given by :

MJ 1 (s) = 1 s m k=1 a k β s k .
-

∈ 1 + iR, a = 0.1573, K = c n 0 K n 1/2 +c n 0 . ∀q > q 1 ∃c > 0 ∀s ∈ q + iR : |sMJ 1 (s)| c
This lemma can be rephrased. On vertical lines of the right half-plane the function s k MJ 1 (s) is bounded from above for k ∈ [0, 1] and is bounded from below for k 1, if the vertical lines are far enough to the right in the complex plane. This allows us to establish Proposition 3.4.

Proposition 3.4. Let I 1 [ρ] be defined by (1.3), (2.1), (2.3). • Let q < -1. If ρ is in L 2 q
, there exist constants C, C > 0 such that:

I 1 [ρ] L 2 q-3 2 C ρ L 2 q and (I 1 [ρ]) L 2 q+1 2 C ρ L 2 q .
• Let q 1 be defined by (2.1), (3.1) and, let q < -2q 1 -1.

If (I 1 [ρ]) is in L 2 q+1 2
, there exists a constant C > 0 such that:

(I 1 [ρ]) L 2 q+1 2 C ρ L 2 q .
The proofs of these statements can be found in Section 5.2. This proposition is the same as the one obtained in [START_REF] Conca | Determination of the calcium channel distribution in the olfactory system[END_REF], except that the weight q appearing in that paper is positive.

Second, let us consider the non linear approximation

I 2 [ρ] of the operator I 0 [ρ].
In this case the computations are no longer explicit. Nevertheless an argument in the spirit of the stationary phase method proves that MJ (a) 2 has the same behaviour as MJ 1 (s), if the real number a has a small enough value. This is illustrated in Figure 3.1 and this allows to establish Proposition 3.5. Proposition 3.5. Let a > 0 and I (a) 2 [ρ] be defined by (2.2), (1.3), (2.3). Let q < -1.

• If ρ is in L 2 q , there exists a constant C > 0 such that: • If ρ is in L 2 q , and a has a small enough value, then there exists a constant C > 0 such that:

I (a) 2 [ρ] L 2 q-3 2 C ρ L 2 q . -10 -8 -6 -4 - 2
I (a) 2 [ρ] L 2 q+1 2 C ρ L 2 q . • If I (a) 2 [ρ] is in L 2 q+1 2
, and a has a small enough value, then there exists a constant C > 0 such that:

I (a) 2 [ρ] L 2 q+1 2 C ρ L 2 q .
The proofs of these statements can be found in Section 5.4.

The smallness condition on a is needed for the observability inequality to be true. Indeed for

J (c0) 2 = J 0 (because H (c0) 2
= H 0 ) no inequality of the last type holds. This fact is a consequence of Proposition 3.6.

Finally let us consider the case of the Hill function defined by (1.2). Explicitly computing the derivatives of J 0 shows that the Mellin transform MJ 0 has a fast decay on vertical lines. Consequently on vertical lines of the right half-plane, and for any k ∈ N, the function s k MJ 0 (s) is bounded from above and is not bounded from below. This fact is illustrated in Figure 3.2. As a consequence, no observability inequality of the previous type holds. Proposition 3.6. Let q < -1.

• If ρ is in L 2 q , there exists C > 0 such that:

I 0 [ρ] L 2 q-3 2 C ρ L 2 q .
• Furthermore, let us assume that µ s ∈ q+1 2 + iR | Mρ(s) = 0 = 0 where µ is the Lebesgue measure on q+1 2 + iR. For n > 0, there is no k ∈ N, C > 0 such that: For instance, the condition µ s ∈ q+1 2 + iR | Mρ(s) = 0 = 0 holds if Mρ is non-zero and holomorphic on a neighborhood of q+1 2 + iR in C and this is satisfied if ρ is in a non-zero L 1 w function for w in a neighborhood of q-1 2 , see Section B. The proofs of these statements can be found in Section 5.5.

(I 0 [ρ]) (k) L 2 2k+ q-3 2 C ρ L 2 q .

Numerical simulations.

In this section we present some numerical simulations based on experimental data. The aim is to use formula (2.6) to compute the distribution ρ of the CNG channels from an experimentally measured current I exp .

We consider that the experimentally observed current is given by the following formula (see [START_REF] French | Integral equation models for the inverse problem of biological ion channel distributions[END_REF] and the references therein).

I exp (t) =    0 if t t D I max 1 + K I t-t D n I -1 if t > t D (4.1)
The parameters are given by: I max = 150 pA, n I = 2.2, K I = 100 ms, t D = 30 ms. The function I exp is shown in Figure 4.1.

Note that the derivative of I exp is in L 2 q , for some q < 0 so Proposition 3.5 (or Proposition 3.4) can be applied directly. If the data were noisy, a standard regularization method for an inverse problem with a finite degree of ill-posedness (see [START_REF] Engl | Regularization of inverse problems[END_REF] for instance) could be applied to the data before applying Proposition 3.5. In this extend our numerical method consisting of Proposition 3.5 can be applied to any noisy experimental data. For the approximated operator I For the operator I 0 [ρ] the distribution ρ cannot be found from the current I exp , even after a regularization, as shown in Figure 4.3. This is coherent with Proposition 3.6, which states that no observability inequality can be found for I 0 in some weighted L 2 q spaces.

Technical aspects

As seen in Section 5 the functions MJ 0 , MJ 1 , MJ 2 are holomorphic on the right half-plane. As t D > 0, n I > 0, I exp is in L 1 k for every k < -1, then, M I exp is holomorphic on the left half-plane. Therefore the quotient (-s) is meromorphic on the left half-plane. The inverse Mellin transform of this function has been computed on different vertical lines, whose abscissas are q < 0. Fig. 4.2. For s ∈ q + iR, for some values of q < 0: Fig. 4.3. For s ∈ q + iR, for some values of q < 0: -(left)

M Iexp( s 2 ) MJ l
-(left) M Iexp( s 2 ) MJ (a) 2 (-s) , -(right) ρexp = M -1 q M Iexp( s 2 )

MJ

M Iexp( s 2 ) MJ 0 (-s) , -(right) M -1 q M Iexp( s 2 ) MJ 0 (-s) .
As the quotient

M Iexp( s 2 ) MJ (a) 2 (-s)
does not vanish in the left half-plane, the distribution ρ does not depend q < 0. More details regarding this fact can be found in [START_REF] Bourgeron | Estimating the division rate of the growthfragmentation equation with a self-similar kernel[END_REF]. Figure 4.2 illustrates Proposition 3.5. For the numerical simulations, the parameter a is taken to be the biggest for which the quotient defined above does not vanish in the strip {x + it | x ∈ [-4, -2], t ∈ R} (that is on the vertical -4 + iR thanks to ), that is a = 0.1573. If the negative part of ρ exp is removed then its profile is similar to the one exhibited in [START_REF] French | Numerical approximation of solutions of a nonlinear inverse problem arising in olfaction experimentation[END_REF], [START_REF] French | Perturbation approximation of solutions of a nonlinear inverse problem arising in olfaction experimentation[END_REF], [START_REF] French | Integral equation models for the inverse problem of biological ion channel distributions[END_REF]. 5. Proofs. Notation. (Pochhammer symbol) For a real number x and a nonnegative integer n we write:

(x) n = x • • • (x -n + 1) = n-1 j=0 (x -j).
(5.1)

For instance ∀q ∈ R and k ∈ N: (x q ) (k) = (q) k x q-k .

1/2 M I[ρ](s/2) = M(xρ(x))(s) MG l (s) = Mρ(s + 1) 2 s D s/2 MJ l (-s).
This leads to:

M I[ρ](s) = 2 2 2s D s MJ l (-2s) Mρ(2s + 1) (s -k) k M I[ρ](s -k) = 2 2 2(s-k) D s-k (s -k) k MJ l (-2(s -k)) Mρ(2(s -k) + 1).
As the Mellin transform is an isometry (up to the factor (2π

) -1/2 ) of L 2 2k-1+ q-1 2 onto L 2 k + q-1
4 + iR , for s in k + q-1 4 + iR the previous relation implies:

(I[ρ]) (k) L 2 2k-2+ q+1 2 = (2π) -1/2 (-1) k (s -k) k M I[ρ](s -k) L 2 (k+ q-1 4 +iR) = (2π) -1/2 2 2 (q-1)/2 D (q-1)/4 (s -k) k MJ l (-2(s -k))Mρ(2(s -k) + 1) L 2 (k+ q-1 4 +iR)
= (2π) -1/2 2 2 (q-1)/2 D (q-1)/4 (s) k MJ l (-2s)Mρ(2s + 1) L 2 ( q-1 4 +iR) .

The inequality

|(s) k MJ l (-2s)| C (resp. C) on q -1 4 + iR is equivalent to: -s 2 k MJ l (s) C (resp. C) on - q -1 2 + iR.

And:

-s

2 k = k-1 j=0 - s 2 -j = (-1) k 2 -k k-1 j=0 (s + 2j).
Therefore, as M is an isometry of L 2 q onto L 2 q+1 2 + iR , the assumption implies:

(I[ρ]) (k) L 2 2k-2+ q+1 2
(resp. ) C2 -k (2π) -1/2 2 2 (q-1)/2 D (q-1)/4 Mρ(2s

+ 1) L 2 ( q-1 4 +iR) = C2 -k (2π) -1/2 2 2 (q-1)/2 D (q-1)/4 2 -1/2 Mρ(s + 1) L 2 ( q-1 2 +iR) = C2 -k (2π) -1/2 2 q/2 D (q-1)/4 Mρ(s) L 2 ( q+1 2 +iR)
= C2 -k+q/2 D (q-1)/4 ρ L 2 q . As erfc is a decreasing function, for c, a > 0 we have: 1 c erfc(t)-a 0 = 1 erfc -1 ( a c )-t 0 . The explicit formula for MJ 1 is then a consequence of: M(1 0 x a )(s) = 1 s a s for Re s > 0; and of the linearity of M .

Approximation of the

As a k > 0 and 0 < α

1 < • • • < α n the β k satisfy: 0 < β n < • • • < β 1 < +∞. As q > 0: m k=1 a k β q+it k a 1 β q 1 - m k=2 a k β q k a 1 β q 1 -β q 2 m k=2 a k .
This last term is bounded from below by a positive constant for q > q 1 . Proof. [Proof of Proposition 3.4]

• Let q > 0, s ∈ q + iR, k ∈ [0, 1].
The explicit formula of Lemma 3.3 implies:

s k MJ 1 (s) |s| k-1 m k=1 a k β q k C |s| k-1 C.
Now let q < -1, that is -q+1 2 > 0. Applying Lemma 3.1 for k = 0, 1 on the vertical line -q+1 2 + iR leads to the continuity inequalities. • Let q > q 1 > 0, s ∈ q + iR, k 1. The explicit formula of Lemma 3.3 leads to the bounds:

s k MJ 1 (s) |s| k-1 a 1 β q 1 -β q 2 m k=2 a k C > 0. Now let q < -1 -2q 1 , that is -q+1 2 > q 1 . Lemma 3.

applied for k = on the vertical line -q+1

2 + iR concludes the proof.

5.3. Lemmas. Lemma 5.1. Let A and B be two elements of [0, ∞], k ∈ N an integer and f a function such that f (j) is in L 1 j (A, B) for every j = 0, . . . , k. For every real number t, we have:

B A f (x)x it dx = k-1 j=0 (-1) j Q j (t) x j+1 f (j) (x)x it B A +(-1) k Q k-1 (t) B A x k f (k) (x)x it dx, where Q j (t) = j l=0 1 + l + it -1 .
Proof. The proof is by induction on k ∈ N. For k = 0 there is nothing to prove. We assume that the formula is true for an integer k ∈ N.

As (k + 1 + it)Q k = Q k-1
it remains to be proved that:

(k + 1 + it) B A x k f (k) (x)x it dx = x k+1 f (k) x it B A - B A x k+1 f (k+1) (x)x it dx.
As d dx x it = it x x it , the previous relation follows by integration by parts:

it B A x k f (k) (x)x it dx = B A x k+1 f (k) (x)(x it ) dx = x k+1 f (k) (x)x it B A -(k + 1) B A x k f (k) (x)x it dx - B A x k+1 f (k+1) (x)x it dx. Corollary 5.2. Let f : [A, B] → R with A, B ∈ [0, ∞] be a piecewise C 1 function. If f is non-negative, f is non-positive, f ∈ L 1 (A, B), f ∈ L 1 1 (A, B
) and for all t ∈ R: [xf (x)x it ] B A = 0, then:

t B A f (x)x it dx B A f (x) dx,
where t = (1 + t 2 ) 1/2 . Proof. From Lemma 5.1 with k = 1 one obtains:

∀t ∈ R (1 + it) B A f (x)x it dx = - B A xf (x)x it dx.
As A, B 0 and f 0, using the formula for t = 0:

t B A f (x)x it dx B A |xf (x)| dx = B A f (x) dx.
Lemma 5.3. Let n, K > 0 and f = erfc n erfc n +K . Let q ∈ R. There exists x q > 0 such that the function g q : x ∈ [x q , ∞) → f (x)x q-1 is decreasing. Let q = inf E q where E q = {c 0 | ∀x c : g q (x) < 0}. The function q → q is increasing and q ∼ (q/(2n)) 1/2 as q → +∞.

Proof. Let n, K > 0. Let q ∈ R. As f > 0, the inequality g q (x) 0 is equivalent to

f (x) f (x) - q x .
(5.2)

Let us compute f f . To do so, let u = erfc n , so that: f = u u+K . We have:

f f = u u K u + K = n erfc erfc K u + K . (5.3)
The derivative of erfc is given by: erfc (x) = -2π -1/2 e -x 2 . And, as x tends to +∞: erfc(x) ∼ π -1/2 e -x 2 1

x . Thus, as x tends to +∞:

f (x) f (x) ∼ n erfc (x) erfc(x) ∼ -2nx. (5.4)
This asymptotics proves that the inequality (5.2) is satisfied for large enough values of x. As a consequence for every q in R the set E q is not empty, which justifies the definition of q. Note that the definition of q implies: g q ( q) = 0, that is: f ( q) f ( q) = -q q , using (5.2).

Let q 1 q 2 be two real numbers. In order to show that q 2 q 1 , it is enough to prove that g q1 ( q 2 ) 0. This holds because:

g q1 ( q 2 ) = q 2
q1-1 (f ( q 2 ) q 2 +f ( q 2 )q 1 ) q 2 q1-1 (f ( q 2 ) q 2 +f ( q 2 )q 2 ) = q 2 q1-q2 g q2 ( q 2 ) = 0.

To find the asymptotics on q, let us recall a lower bound on erfc(x) for x 0:

x + (x 2 + 2) 1/2 -1 1 2 π 1/2 exp(x 2 ) erfc(x).
As the function u = erfc n takes its values in (0, 1], we have nK 1+K nK u+K n. Consequently, thanks to (5.3):

-n x + (x 2 + 2) 1/2 f (x) f (x) .
(5.5)

Let q > 0 and set

x q = q (2n) 1/2 (n+q) 1/2 . The inequality -q x -n x + (x 2 + 2) 1/2 is equivalent to x x + (x 2 + 2) 1/2 q n .
A simple computation shows that this inequality is satisfied for x = x q (and becomes and equality). Thanks to (5.5), we conclude that x q satisfies:

f (xq) f (xq) -q
xq , which leads to q x q , by definition of q and by (5.2). This last inequality implies that q tends to +∞ as q tends to +∞. Finally we get the asymptotics for q, using (5.4):

-2n q ∼ f ( q) f ( q) = - q q .
5.4. Non linear approximation of the Hill function. Proof. [Proof of Proposition 3.5] From the estimate for erfc at ∞, given in the proof of Lemma 5.3, the function J 2 (s) , from above or from below, on the vertical lines q + iR, for q > 0.

The function J (a) 2 is written:

J (a) 2 (x) = H 2 (c 0 erfc(x)) = c n 0 erfc(x) n c n 0 erfc(x) n + K n 1/2 1 x α + K1 0<x<α ,
with: α = erfc -1 (a/c 0 ). Then the Mellin transform of J (a) 2 at s = q + it is given by:

MJ (a) 2 (s) = K α 0 x s-1 dx + ∞ α f (x)x s-1 dx = K α s s + ∞ α f (x)x q-1 x it dx, where f (x) = erfc(x) n erfc(x) n +c -n 0 K n 1/2 .
For any a 0, q > 0 and s ∈ q + iR we have:

MJ (a) 2 (s) K α q q + ∞ α f (x)x q-1 dx,
which is finite. Let q > 0. According to Lemma 5.3 the function x → f (x)x q-1 is decreasing for

x x 0 . Let a < c 0 erfc(x 0 ). Then α = erfc -1 (a/c 0 ) x 0 . Let g(x) = f (x)x q-1 1 x α .
For every t ∈ R : f (x)x it ∞ x0 = 0 because f vanishes for x α and x 0 α, and g(x) ∼ +∞ π -n/2 x -n+q-1 e -nx 2 . Then Corollary 5.2 can be applied to the function g, with A = α, B = +∞, for s ∈ q + iR, to give:

sMJ (a) 2 (s) K |α s | + |s| t t ∞ α f (x) x s-1 dx Kα q + max(1, q) ∞ α f (x)x q-1 dx < ∞, because: |s| t ∈ [q, 1] ∪ [1, q],
either q 1 or q 1. For small values of a, the first part dominates the second one. The same calculation as above leads to:

sMJ (a) 2 (s) Kα q -max(1, q) ∞ α f (x)x q-1 dx < ∞,
The last expression is equivalent to Kα q as α tends to ∞, therefore, it is positive for large values of α.

The Hill function

. Definition 5.4. Let S[0, ∞) be the space of functions f in C ∞ ([0, ∞), C) which statisfy: ∀j ∈ N, k ∈ N lim x→∞ f (j) (x)x k = 0.
If f is a function of the Schwartz space S(R) then f

1 x 0 is in S[0, ∞).
Remark. The converse is also true. Let f be in S[0, ∞). By Borel's lemma there exists a function g ∈ C ∞ (R) such that for every k ∈ N: f (k) (0) = g (k) (0). Then the function ψ(g1 x<0 + f 1 x 0 ), where ψ is a smooth function vanishing for x < -1 and equal to 1 for x 0, is in S(R) and its restriction to [0, ∞) is f . Lemma 5.5. Let f ∈ S[0, ∞). Its Mellin transform Mf is holomorphic on the right half-plane, and:

∀q > 0 ∀k ∈ N ∃C 0 ∀t ∈ R : |MJ 0 (q + it)| C t -k , where t = (1 + t 2 ) 1/2 . Proof. [Proof of Lemma 5.5] Let f ∈ S[0, ∞), q > 0. By the definition of S[0, ∞) for every l in N and k > -1, the function x → x k f (l) (x) is in L 1 . Proposition B.
3 of the Appendix implies that Mf is holomorphic on the right half-plane. Lemma 5.1 with g(x) = f (x)x q-1 is written:

Mf (q+it) = ∞ 0 f (x)x q-1 x it dx = k-1 j=0 (-1) j Q j (t) x j+1 g (j) x it ∞ 0 +(-1) k Q k-1 (t) ∞ 0 x k g (k) (x)x it dx, where Q j (t) = j l=0 1 + l + it -1 .
To prove the result it is enough to show that the terms between brackets vanish and that the last integral is finite.

Let l, k ∈ N. By the Leibniz rule:

x l g (k) (x) = k j=0 k j f (k-j) (x)(x q-1 ) (j) x l = k j=0 k j (q -1) (j) f (k-j) (x) x q+l-1-j .
For l = k +1 and for x = 0 this expression vanishes because f (k-j) (0) is finite and: q+k-j q > 0. As x tends to ∞ the expression tends to 0 as: f (k-j) (x) x q+k-j → 0.

For l = k this expression shows that the integral ∞ 0 x k g (k) (x) dx is finite because for every j ∈ {0, . . . , k}: x → x q-1+j f (j) (x) is in L 1 because q -1 + j q -1 > -1. Then:

|Mf (q + it)| C |Q k-1 (t)| ∼ t→∞ C t -k .
Lemma 5.6. Let J 0 be the function defined by (1.2), (2.4). For n > 0, the function J 0 is in S[0, ∞).

The proof of Lemma 5.6 is given in Appendix A. Proposition 3.6 can be strenghened in Proposition 5.7 below (in the sense that Lemma 5.6 and Proposition 5.7 imply Proposition 3.6).

Proposition 5.7. Let q < -1, l in {0, 1, 2} such that J l , defined by (2.1), (2.2), (2.4), is in S[0, ∞).

• If ρ is in L 2 q , there exists C > 0 such that:

I l [ρ] L 2 q-3 2 C ρ L 2 q .
• Furthermore assume that µ s ∈ q+1 2 + iR | Mρ(s) = 0 = 0 where µ is the Lebesgue measure on q+1 2 + iR. There is no k ∈ N, C > 0 such that:

(I l [ρ]) (k) L 2 2k+ q-3 2 C ρ L 2 q .
Proof. [Proof of Proposition 5.7] It is enough to show that the two inequalities:

(I[ρ]) (k) L 2 2k+ q-3 2 C ρ L 2 q resp. C ρ L 2 q (5.6) and k-1 j=0 (s + 2j) MJ l (s) C (resp. C) on - q -1 2 + iR (5.7)
are equivalent.

This equivalence having been proved, Lemma 5.5 implies that, for any k ∈ N, the inequality is true and the inequality is not satisfied.

The implication (⇐) is Lemma 3.1. For the other implication (⇒), let us rewrite equations (5.6), (5.7) in equivalent forms (cf. proof of Lemma 3.1). Respectively, we have:

(s -k) k MJ l (-2(s -k)) Mρ(2(s -k) + 1) L 2 ( q-1 4 +k+iR) C Mρ(2(s -k) + 1) L 2 ( q-1 4 +k+iR) (5.8) resp. C Mρ(2(s -k) + 1) L 2 ( q-1 4 +k+iR) and |(s -k) k MJ l (-2(s -k))| C (resp. C) on q -1 4 + k + iR. (5.9)
If (5.8) is satisfied then (5.9) is satisfied on the set:

s ∈ q -1 4 + k + iR | |Mρ(2(s -k) + 1)| > 0 ,
that is (5.7) is satisfied on the set:

s ∈ - q -1 2 + iR | |Mρ(-s + 1)| > 0 .
As q < 1, and as

J l is in S[0, ∞), the function J l is in L 1 -q+1 2
. By Theorem B.2 (in the Appendix), MJ l is continuous on -q-1 2 + iR. The left-hand side of (5.7) is then a continuous function of s. As µ s ∈ q+1 2 + iR | Mρ(s) = 0 = 0, the inequality (5.7) holds for every s in -q-1 2 + iR, which concludes the proof. 6. Conclusion. We focussed on the problem of finding the spatial distribution of CNG ion channels from the experimental current data, following [START_REF] Conca | Determination of the calcium channel distribution in the olfactory system[END_REF][START_REF] French | Perturbation approximation of solutions of a nonlinear inverse problem arising in olfaction experimentation[END_REF][START_REF] French | Integral equation models for the inverse problem of biological ion channel distributions[END_REF][START_REF] French | Numerical approximation of solutions of a nonlinear inverse problem arising in olfaction experimentation[END_REF][START_REF] Badamdorj | Identification of Cl(Ca) channel distributions in olfactory cilia[END_REF]. The self-similar structure of the integral inverse problem (1.1), (1.2), (1.3) allowed us to use the Mellin transform, and to obtain a thorough comprehension of it. It permited us to reduce its study to estimating some Mellin transform, on vertical lines of the complex plane. To do so, explicit computations were carried out using techniques inspired by the stationary phase method.

As a result the inverse problem studied has been shown not to be controlable in some weighted L 2 spaces. The kernel of the integral operator is smooth and the associated inverse problem has an infinite degree of ill-posedness. This conclusion could probably be linked to the fact that for the original problem, introduced in [START_REF] French | Numerical approximation of solutions of a nonlinear inverse problem arising in olfaction experimentation[END_REF], certain numerical computations indicate that it is ill-conditioned. We also introduced a better approximation than the one already studied, for which we performed numerical simulations and provided estimates. In this case the kernel of the integral operator is at most continuous and the associated inverse problem has a finite degree of ill-posedness (which is 1). The profiles obtained from the experimental current data consolidate the ones obtained in [START_REF] French | Numerical approximation of solutions of a nonlinear inverse problem arising in olfaction experimentation[END_REF], [START_REF] French | Perturbation approximation of solutions of a nonlinear inverse problem arising in olfaction experimentation[END_REF], [START_REF] French | Integral equation models for the inverse problem of biological ion channel distributions[END_REF].

To go further, one could study the inverse problem (1.1), (1.2) where c is defined as the solution of the linear heat equation

       ∂ t c -D ∂ xx c = 0, t > 0, x ∈ (0, L), c(t, 0) = c 0 , t > 0, ∂ x (t, L) = 0, t > 0, c(0, x) = 0, x ∈ (0, L).
If the Hill function (1.2) is changed into its Taylor polynomial extension around c 0 , and if this polynomial has a degree m 8, then [START_REF] Conca | Determination of the calcium channel distribution in the olfactory system[END_REF] proves that the resulting integral operator I is identifiable in L 2 . To go further, a first step, could be to study this problem without approximating the Hill function. Studying the complete inverse problem for the original model could be a last step. It would seem that these two problems cannot be solved using the Mellin transform technique alone, as no self similar structure is involved.

Appendix A. Proof of Lemma 5.6. We recall the Faà di Bruno formula, first proved by L. F. A. Arbogast see [START_REF] Arbogast | Du calcul des dérivations[END_REF].

(f • g) (k) = (m1,...,m k )∈N k k j=1 jmj =k k! k j=1 m j ! j! mj f ( k j=1 mj ) • g k j=1 g (j) mj .
This formula allows us to make explicit computations, to prove that the functions appearing in Lemmas 5.6 are in S[0, ∞).

Lemma A.1. The function erfc is of the class C ∞ on [0, ∞). For every k ∈ N * :

erfc (k) (x) = P k (x) exp(-x 2 ),
where P k is a polynomial which has a degree of k -1 whose leading coefficient is π -1/2 (-2) k . In particular, as x tends to ∞ :

erfc (k) (x) ∼ π -1/2 (-2) k x k-1 e -x 2 .
Proof. The proof is by induction on k ∈ N * . By the definition of erfc: erfc (x) = -2π -1/2 exp(-x 2 ) and with simple calculations we have:

P k+1 = -2XP k + P k .
Lemma A.2. The function H defined by (1.2) is of the class C ∞ on [0, ∞). For n 0 it satisfies for every k ∈ N * :

H (k) (x) = P k (x) (x n + K n 1/2 ) k+1 , where P k (x) = (-1) k-1 (n -k -1) k K n 1/2 x (n-1)k + • • • + (n) k K nk 1/2 x n-k . In particular: H (k) (x) ∼ 0 K nk 1/2 (n) k x n-k .
Proof. The proof is by induction on k ∈ N * . Simple calculations lead to: P 1 (x) = nK n 1/2 x n-1 and:

P k+1 (x) = x n P k (x) -(k + 1)nx n-1 P k (x) + K n 1/2 P k (x).
It follows that the leading term in P k (x) as x tends to +∞ is of the form a k x (n-1)k where a k satisfies:

a 1 = nK n 1/2 , a k+1 = a k ((n -1)k -n(k + 1)) = -a k (n + k).
As n 0, the leading term in P k (x) as x tends to 0 is of the form b k x n-k where b k satisfies:

b 1 = nK n 1/2 , b k+1 = b k (n -k)K n 1/2 .
Lemma A.3. For n > 0, the function J 0 defined by (1.2), (1.3), (2.4) is in S[0, ∞).

Proof. The function is of the class C ∞ (0, +∞) because erfc is of this class and n > 0. For every integer k ∈ N, H (k) (c 0 ) and erfc (k) (0) are finite, thus for every k ∈ N, J (k) 0 (0) is finite. Applying the Faà di Bruno formula one gets, for x 0:

J (k) 0 (x) = (m1,...,m k )∈N k k j=1 jmj =k k! k j=1 m j ! j! mj H ( k j=1 mj ) (c 0 erfc(x)) k j=1 c 0 erfc (j) (x) mj .
As x tends to +∞, from the previous lemmas, we have, for S = k j=1 m j :

H (S) (c 0 erfc(x)) ∼ +∞ K -n 1/2 (n) S (c 0 erfc(x)) n-S ∼ +∞ K -n 1/2 (n) S (c 0 π -1/2 x -1 exp(-x 2 )) n-S .
Then, for

k j=1 jm j = k and k j=1 m j = S: H ( k j=1 mj ) (c 0 erfc(x)) k j=1 c 0 erfc (j) (x) mj ∼ +∞ K nS 1/2 (n) S (c 0 π -1/2 x -1 exp(-x 2 )) n-S k j=1 c mj 0 (π -1/2 (-2) j x j-1 exp(-x 2 )) mj = (-2) k K nS 1/2 (n) S c n 0 π -n/2 x k-n exp(-nx 2 )
Let us denote:

C(n, k) = (-2) k k!c n 0 π -n/2 (m1,...,m k )∈N k k j=1 jmj =k (n) k j=1 mj K n k j=1 mj 1/2   k j=1 m j ! j! mj   -1 .
If C(n, k) = 0, the previous calculations lead to:

J (k) 0 (x) ∼ C(n, k)x k-n exp(-nx 2 ).
If C(n, k) = 0: J (k) 0 (x) = o x k-n exp(-nx 2 ) . In both cases, for every j ∈ N: J The construction of the Mellin transform on iR can be done in the general context of the Fourier transform on a locally compact abelian group, see [START_REF] Rudin | Fourier Analysis on Groups[END_REF]. Here we consider the multiplicative abelian group G = (0, ∞) (with unit 1), equipped with the topology inherited from R and with the Haar measure dx

x (that is the unique measure on G, up to a positive multiplicative constant, which is translation-invariant). It is easy to show that the dual group Γ of all the characters of G with the Gelfand topology is isomorphic to iR with the topology inherited from C via iR → Γ it → (x → x -it ). In this general context the L 1 and L 2 theories can be built, with the same results, as the Fourier transform on the topological group (R n , +, dx). The extension of the Mellin transform to a vertical strip q + iR of the complex plane C is simply obtained by defining Mf (q + it) = Mg(it) with g(x) = x q f (x).

Definition B.1. Let f be in L 1 q . The Mellin transform of f is a complex valued function defined on the vertical strip q * + iR, where q * = q + 1; writing s -1 = q + it for a real number t it is defined by:

Mf (s) = ∞ 0 x s f (x) dx x .
Theorem B.2 (Riemann-Lebesgue). The Mellin transform is a linear continuous map of L 1 q into C 0 (q * + iR) ⊂ L ∞ (q * + iR), where q * = q + 1, its norm is 1.

Proposition B.3. If f is in L 1 q for every real number q in (a, b) then its Mellin transform Mf (s) = • If H 0 is the Hill function H 0 (x) =

x n x n +K n q and if Mf L 1 (q * +iR) is finite, where q * = q + 1, then almost everywhere x > 0 we have:

f = M -1 q * (Mf ),
where M -1 q ϕ(x) = 1 2π R ϕ(q + it)x -(q+it) dt.

Definition B.6. For two functions f, g we define the multiplicative convolution f * g as:

(f * g)(x) = ∞ 0 f (y)g
x y dy y .

Proposition B.7.

M(f * g)(s) = Mf (s)Mg(s), whenever this expression is well defined.

Proposition B.8. For a function f in L 1 q-1 ∩ L 2 2q-1 we have:

f L 2 2q-1 = (2π) -1/2 Mf L 2 (q+iR) .
As the subspace L 1 q-1 ∩ L 2 2q-1 is dense in L 2 2q-1 this identity allows us to extend the Mellin transform to L 2 2q-1 .

Theorem B.9 (Plancherel transform). According to the previous formula the Mellin transform can be extended, in a unique manner, to an isometry (up to the multiplicative constant (2π) -1/2 ) of L 2 q onto L 2 (q * + iR), where q * = q+1 2 .

  0

First, let

  us consider the linear approximation I 1 [ρ] of the operator I 0 [ρ]. In this case the Mellin transform MJ 1 can be computed explicitly. Lemma 3.3 shows explicit bounds. Definition 3.2. Define:

  Fig. 3.1. MJ (a) 2 (s) (left) and MJ (a) 2 (s)s k , for some values of k (right), for s ∈ 1 + iR, a = 0.1573, K =

Fig. 3 . 2 .

 32 Fig. 3.2. |MJ 0 (s)| (left) and MJ 0 (s)s k , for some values of k (right), for s ∈ 1 + iR.

Fig. 4 . 1 .

 41 Fig. 4.1. Functions Iexp defined by (4.1) and t → |M Iexp(-q + it)| for some values of q.

  ] the distribution ρ can be found from the current I exp , as shown in Figure4.2.

5. 1 .

 1 Proof of Lemma 3.1. Let us apply the Mellin transform to (2.7) and use Proposition B.4 (of the Appendix):

  Hill function by a step function. Proof. [Proof of Lemma 3.3] The function J 1 is in L 1 k for every k > -1 thus MJ 1 is holomorphic on the right half-plane, see Proposition B.3 of the Appendix.

1 k

 1 for every k > -1. Thus MJ (a) 2 is holomorphic on the right half-plane, see Proposition B.3 of the Appendix. Using Lemma 3.1 on the vertical line -q+1 2 +iR with -q+1 2 > 0, as for the proof of Proposition 3.4, it amounts to bounding sMJ (a)

  x j ), which concludes the proof.Appendix B. The Mellin transform.

∞ 0 f 1 √ π 1 s Γ s+1 2 .

 0112 (x)x s-1 dx is holomorphic in the strip S = {s ∈ C | a + 1 < Re s < b + 1}. Proposition B.4. function Mellin transform f (at), a > 0 a -s Mf (s) f (t a ), a = 0 |a| -Mf (a -1 s) f (k) (t) (-1) k (s -k) k Mf (s -k) Example. • For erfc(x) = 2 √ π +∞ xexp(-t 2 ) dt one obtains for Re s > 0: M erfc(s) = 1

1 / 2

 12 then for s ∈ nZ: MH 0 (s) = 5 (inversion theorem). If f is in L 1

This work has begun during a visit by the first author to the CMM in Santiago, Chile. This was made possible thanks to ECOS no. C11E07. The authors received partial support from Ecos-Conicyt Grant C13E05. The author was partially supported by PFBasal-01, PFBasal-03 projects, Fondecyt Grant N • 1140773 and Basal-CMM project.

§ Center for Mathematical Modelling (CMM), UMI 2807 CNRS-Chile (rlecaros@dim.uchile.cl). The author was partially supported by PFBasal-03 and Basal-CMM project. 1


A Gademer

F Mainfroy

L Beaudoin
email: beaudoin@esiea.fr@www.esiea-recherche.eu

L Avanthey

V Germain

C Chéron

S Monat

J P Rudant

FAUCON NOIR UAV PROJECT DEVELOPMENT OF A SET OF TOOLS FOR MANAGING, VISUALIZING AND MOSAICING CENTIMETRIC UAV IMAGES

Keywords: UAV, high resolution imagery, flight-plan, intuitive database navigation, image mosaicing

Eager for greater flexibility, and lesser costs, we are often tempted to acquire high-resolution data with home-made acquisition systems. However, this approach cannot benefit from usual tools to tackle problems such as image navigation within a large flow of very high resolution pictures, data completion and image mosaicing of small patches of the overflown area. This article presents the tools suite our laboratory has developed. In our case, we have developed homemade quad-rotor for low altitude imagery, equipped with up to three compact cameras on a stabilized platform. This configuration represents 2800 images at a ground resolution of 2 to 5 centimeters for each 20 minute flight. In order to take advantage of this flow of data, we needed a relevant set of tools for optimal acquisition and exploitation of the images.

INTRODUCTION

With the spreading of Geographic Information Systems and free visualization tools such as Google Earth TM , more and more users tend to use satellite and aerial photography in their work to save time and reduce financial cost of their field study.

Nowadays, one can see growing need for very high resolution imagery in various fields such as crisis management (risk map, unreachable place exploration, damage assessment), biodiversity study (land-cover map, extended ground truth), radar images interpretation (roughness estimate) or even recreative cartography (virtual reality panorama).

From this point of view, UAVs (Unmanned Aerial Vehicles) are very interesting tools for low altitude imagery. Quickly deployable, they can fly when you need them, as often as you need without extracost. On the downside, the area covered are limited to a few hundred meters, trajectory is more sensitive to wind and thus less rectilinear and at this resolution the flow of data can grow huge rapidly. Therefore it is vital to develop a relevant set of tools to optimize the data acquisition and exploitation from the UAVs.

We will first present the Libellule Violette software which is an important part of our ground control station. This software applies optimization algorithm to the flight path to reduce the flight time and This work was supported by ONERA and DGA in the Challenge Minidrones contest: http://minidrones.fr/ thus increase the coverage of the UAV. It also has a real-time completion detection algorithm which allows the pilot to know during the flight if there is no missing data and allows him to correct it before the landing, saving the cost of a backup flight.

Then we will present the Drone Eye software, our data exploitation library, which allows the end-users to navigate instinctively in the image database, just after the UAV landing, by an accurate indexation of the geolocalized images.

We will end this article with JIM, a powerful home-made opensource image mosaicer, that we have used to make planar mosaics of the scene.

UAS

Unmanned Aerial System acronym cover all the elements orbiting around the UAV : the quad-rotor itself and its control ground station, the payload and even the exploitation station, that allow the visualization of the data. It is the all-inclusive product the end-user will receive.

UAV

The UAVs history began in the early XX th century where the military discovered that they needed an unmanned aircraft for two main missions : to be able to self-destruct on precise location or to be used as mobile target for the training of jet pilots. With the improvement of technology, UAVs are entrusted with more and more missions both military (intelligence, weapon, cargo, etc.) and civil (security, rescue, cartography, etc.). They come in all shapes and sizes, from automatized real scale plane to micro coaxial rotors UAV [START_REF] Martnez De Dios | Multiple Heterogeneous Unmanned Aerial Vehicles[END_REF][START_REF] Valavanis | Advances in Unmanned Aerial Vehicles[END_REF].

In our case we have made the choice of a small scale quad-rotor. Indeed, this fixed structure with four fixed blade rotors has a lot of advantages, as it is simple from the mechanical point of view, sturdy and very maneuverable. Its degree of freedom comes from the differential speed of the pair of rotors, thus needing quite efficient control algorithms as it is intrinsically unstable. (cf. Fig 1)

The last technical hindrances for the democratization of civil UAV have been solved recently : higher energetic density with the Lithium-Polymer battery, higher efficiency with brushless motors and a global miniaturization of all the electronic components (internal sensors, computing chip, detection sensors, etc.) with the microelectromechanical systems (MEMS) and the surface mounted components (SMC). All of these reduce drastically the weight of the system and thus allow the building of small UAV, fitted for operational missions, at a reasonable cost.

Payload

The payload corresponds to the weight of material the UAV can lift besides itself ; 1 kilogram, in our case, without limiting our maneuverability. Therefore the payload is also the material dedicated to the UAV's mission. Depending on our mission, it can be a photo camera, a video camera, a thermal camera, a passive or active sensor in all the frequency spectrums. It can be a single sensor or in group of sensors within reasonable size and weight. [START_REF] Gademer | Réalisation d'un drone pour des missions de télédétection THR environnementales à la demande[END_REF] shows the interest of simultaneous multi-acquisition of high resolution images for environmental applications and for quick cartography in crisis situation. For situation assessment, operational users need information as much as possible ; they want recent, reliable images of their interest regions from several points of view to help the evaluation.

On the one hand, simultaneous acquisition reduces the number of flights needed to cover an area and thus increases the range of the UAV, it also permits to take a scene from different angles (cf. Fig 2), opening new data exploitation possibilities (intuitive navigation in the image database, stereoscopy, etc.). On the other hand the weight and size constraints become drastic and need more efforts on sensor integration.

As it can be seen in (cf. Fig 2), we use on-the-shelf photo camera in our system. Camera have the advantage of an autonomous storage capacity, all the processing of the image being done on board. They also have a big CCD captor, optics of professional quality, a mechanical shutter and an included stabilization system. However, these cameras are bigger and heavier than industrial CCD cameras and are made for human manipulation and thus are not easy to integrate in on-board systems.

Fortunately, after a slimming process and the development of an electronic card interface we are now able to control perfectly our three cameras and acquire simultaneous images from the sky.

Geolocalization

Relatively simple information is sufficient for the UAV stabilization, but the geolocalization of the data produced by the payload needs a more precise measurement of the trajectory and of the orientation of the UAV.

The problem is not simple as UAV have more complex trajectories than satellites or planes, they are more sensible to environment and their positions and heights can rapidly change due to a gust of wind. The main sensor on the UAV is the Inertial Measurement Unit (IMU) ; with accelerometers, gyroscopes and magnetometers it assesses the orientation of the UAV and thus allows to stabilize it. The orientation is also used to determine in which direction the payload's Fig. 2. The simultaneous multi-acquisition payload cameras look. The position and altitude is given by the Global Positioning System (GPS) which can be accurate enough on flying systems. Another important information is the height above the ground, which is obtained by the data fusion of ultrasound, barometric and accelerometer systems.

To make things a little more complex, all these sensors measure physical property from which it derives the state of the system and thus present noisy and biased results. To make sure to have accurate knowledge of the true trajectory we must use complex filtering and data fusion [START_REF] Chatfield | Fundamentals of high accuracy inertial navigation[END_REF].

All the tools we will present in the next section depend on the quality of this measurment and that is why one never stops working on improving it.

DEDICATED TOOLS

Acquiring images from UAV can be a bit challenging (limited flight time, huge flow of data produced, user craving for usable data, etc.) and therefore we need a relevant set of tools to help us during this process.

From our experience, obtain a fully operational UAS implies to place importance on the UAV and on the ground control station which allows the user to interface with the system, give it orders, check its status in real-time. The first two softwares that we will present are important part of the ground station : the mission programmer and the mission completion surveyor.

Libellule Violette: Mission programmer

The mission programmer is a tool which allows us to plan the UAVs missions. Among other things it offers to register, before take off, all the pieces of information we need in order to perform the mission in good conditions.

A flight-plan is composed of a sequence of geographic coordinates, which the UAV will follow during the mission. We thought that it would be easier to set them knowing approximately the scene we wanted to cover. Therefor we add a functionality to load any geolocalized image as a background map. The image can be loaded from the hard-drive or directly from Internet giving its centers geographical coordinates, a radius and a resolution. You can have the names of the streets and other items on this map or not according to your wishes. The map is directly downloaded and cached for later use.

Once loaded, the map is shown in background and you can move on it or change the zoom factor at any moment. On this image we can add manually some vector drawings, if we already possess some field information. The possibility to add areas and points of interest on the map can prove very useful during the planning and even during the flight as we will show in paragraph 3.2. For example it is important to note any uneven relief or obstacle. It can also be used to mark any change seen between the map and the reality as seen on the field. Every object can be customized with an identifier and a color so it is easy to track it down instantaneously on the map. (cf. Fig 3)

Once aware of the characteristics of the area, the user can plan the flight by adding, deleting and moving way-points. The program will calculate in real-time the optimized path to cover all the waypoints in the shortest time. A traveling salesman algorithm has been integrated for optimal calculation. In the future we hope to add other algorithm based on multiple criteria (distance, priority, sensor ground-print, etc.) All the raster and vector informations can be saved in XML format to be reloaded later and to be use by the other programs, as the mission completion surveyor.

Libellule Violette: Mission completion Surveyor

The ground station gives in real-time the progress of the UAV on the flight-plan map to check his status and the correlation between the flight and the planned mission. We wanted to do more.

As cartographic UAV users, the most important thing for us is to avoid incomplete data. If we have two or more pictures of the same scene we can choose the best one but if we lack one, we are not able to give a correct assessment and that is not desirable. Most of the systems (satellite, aerial or UAVs) take a security margin : all the images cover a great percentage of the others, just in case. We wanted to optimize our flight and reduce this margin but how can we assure that we will not forget any important image ?

It is the job of the mission completion surveyor. This program, active during the flight use the data from the flight-plan programmer to show the progress of the UAV. But it takes advantage of all the information about the images we received in near real time from our UAV.

Firstly, as we have real-time retransmission of an on-board mini video camera , we add the functionality to edit during the flight all the vector information of the flight-plan. Therefore, it is possible, as in the mission programmer phase to add areas and points of interest, to tag them, to modify them, and so on. So for example, we may want to indicate that an area is flooded, or that there are rock slides and mass of fallen rocks at this place, that a detail requires ones interest here, that this area is marshy and this one forested, or whatever else you want to be highlighted on your map. Secondly, we use the trajectory information paired with the time of acquisition of the pictures to calculate the approximate ground print of the pictures. Similar to the "fog-of-war" you can see in video games, all the places that have been seen by the UAV's payload are highlighted. All remaining areas are darkened, so the ground station user can see rapidly if there is any lack of data or not and possibly edit the flight plan as a result.

Finally, we made the data structure of these maps simple to be able to share with all the operational users on the ground, and that again in near real time. This map will allow them to get an overall view of all the data gathered and then they can act accordingly. Now we can assure before the landing of the UAV that we have all the images we want, we will need to transfer and process them to be integrated in our image database. It is the role of our next tool, named Drone Eye.

Drone Eye: Generating pertinent image database

As seen in section 2.2, our small quad-rotor UAV can carry up to three compact cameras on a stabilized platform. From the operational point of view, that represents a total of 2800 images at a resolution of 2 to 5 centimeters for each 20 minute flight, which would be completely useless without powerful tools to index and visualize those pictures. We have developed an innovative solution to transform these raw images into highly valuable scenes.

The first step of our solution is to geolocalize the acquired data, using our embedded and GPS and Inertial System. We are thus able to compute the footprint of every picture and infer topological relationships between them. To ease the process, we use the WGS84 datum for GPS compliance (longitude, latitude, altitude coordinates), and a local tangent plane (north, east, up coordinates) which makes local navigation in our images database fast and easy.

The geographic metadata are then stored as separated XML files. Our first idea was to store these metadata into EXIF, but this data format does not provide a field for all the information we need, for example the attitude of the camera, which is obviously needed to reproject the pictures on the ground. Storing metadata in separated XML files, and not in the images files, has another benefit: we are dealing with very small text files and do not have to manipulate large binary images files when looking for a specific image. XML format also ensures temporal durability of those critical links between pictures and their geographic coordinates. These data are then stored in an XML Native Database (XND). The choice of using an XND was driven by the need of a fast access to metadata, in order to make navigation and visualisation seamless. This database scheme fully brings the benefits or XML trees, and only pertinent data is requested using the XQuery language. It is necessary make the distinction between XML-enabled databases and XML Native databases. The former are conventional databases able to parse and request XML text files, as XND are designed specifically for this task. XND are much faster though when dealing with XML files than conventional databases.

Drone Eye: Visualization & Navigation

As seen in section 3.3, we have designed a way to access our images through geographic requests.

We developed a software solution to access our pictures in a natural way, through the use of Google Earth TM .

Users navigate through the database with the help of their keyboard, using only a few keys to access desired pictures. Valid requests are "next image north from this one", "same location but with a 45 point of view", "same location but viewed from the south". Figure 4 shows Drone Eye usage : an image from the database is currently displayed, and the user can navigate through images around the area using the keyboard, as entered keys are intercepted by Drone Eye, which displays requested pictures.

Technically speaking, every adjacent image to the current one is requested, using XQuery language, from the XND. When the user accesses another image, it can be displayed very fast as we don't have to search the entire database.

Drone Eye is the glue between the XND, Google Earth TM and the user, which hides complex geographic projections and matrices mathematics, and lets the end-user concentrate on his needs: "show me the next image on the left", with a very intuitive navigation system.

JIM: Open-source automatic mosaicker

As we acquire multiple images from a very low altitude, we have multiples patches, each covering a very narrow area.

This approach gives centimetric resolution, but denies the user the possibility of getting the global picture.

We have developed our open-source solution to automatically generate mosaics from our images based on [START_REF] Capel | Image Mosaicing and Super-Resolution[END_REF].

The first and critical step of our algorithm is to find repeatable interest points on acquired scenes. Different point detectors are suitable for this kind of application, but [START_REF] Bauckhage | Evaluation of keypoint detectors[END_REF] recommend using Harris corner detector [START_REF] Harris | A combined corner and edge detector[END_REF] for its high repeatability. To be more precise, we use an evolution of this detector described by [START_REF] Shi | Good features to track[END_REF]. We use Harris's [START_REF] Harris | A combined corner and edge detector[END_REF] corner detector to ensure free software compliance. We then use Zero-mean Normalized Sum of Square Differences (ZNSSD) to match interest points and derive the most likely homography from these matches using RANdom SAmple Consensus (RANSAC) [START_REF] Fischler | Random sample consensus: a paradigm for model fitting with applications to image analysis and automated cartography[END_REF].

More matches are then found using guided matching based on the homography found previously, which makes the software able to refine the transformation between images.

If overlapping between two (or more) non consecutive images is detected, Jim is able to recompute only the necessary transformations using the N-views graph described by Capel, taking advantage of Dijkstra algorithm [START_REF] Dijkstra | A note on two problems in connexion with graphs[END_REF] to minimize error propagation.

Transformed images are then mapped on a plane, using the previously computed homographies.

But image mosaicing would only be stitching without radiometric rearrangement. Jim is here helped by Enblend, another very efficient and open-source tool, which detects radiometric seams and moves separations between images in high frequencies areas, to make the transitions between images almost invisible.

Jim can be downloaded under GPL v3.0 on SourceForge.net.

CONCLUSION

To conclude we will say that the development of a set of tools to optimize acquisition and exploitation of data from UAV is very important. It is the missing link between the UAV makers and the operational users around the world. We showed how promising UAV would be in the near future for the civil users and how important a precise geolocalization was. All ground station software allows the optimization of the flight to increase the range of the system. Exploitation software should give the operational users the ability to make a relevant use of this new type of data.

We will end, saying that we are currently running field test with end-users that will be the subject for a future publication.

Fig. 1 .

 1 Fig. 1. The quad-rotor maneuverability is obtained by the differential speed of the rotors.

Fig. 3 .

 3 Fig. 3. The background map, with the flight-plan and the vector information.

Fig. 4 .

 4 Fig. 4. Images can be requested with the keyboard and are immediatly displayed in Google Earth TM .

Fig. 5 .

 5 Fig. 5. A mosaic realized by JIM and Enblend using 6 12Mpx images.