François Dubois
email: francois.dubois@cnam.fr

Zeno Toffano
email: zeno.toffano@centralesupelec.fr

Eigenlogic: a Quantum View for Multiple-Valued and Fuzzy Systems

Keywords: ; ed. July 7, 2016.] finite elements, quantum gates, Boolean functions. AMS classification: 03B52, 81Q99, 94-04

We propose a matrix model for two-and many-valued logic using families of observables in Hilbert space, the eigenvalues give the truth values of logical propositions where the atomic input proposition cases are represented by the respective eigenvectors. For binary logic using the truth values {0, 1} logical observables are pairwise commuting projectors. For the truth values {+1, -1} the operator system is formally equivalent to that of a composite spin 1 /2 system, the logical observables being isometries belonging to the Pauli group. Also in this approach fuzzy logic arises naturally when considering non-eigenvectors. The fuzzy membership function is obtained by the quantum mean value of the logical projector observable and turns out to be a probability measure in agreement with recent quantum cognition models. The analogy of many-valued logic with quantum angular momentum is then established. Logical observables for three-value logic are formulated as functions of the Lz observable of the orbital angular momentum = 1. The representative 3-valued 2-argument logical observables for the Min and Max connectives are explicitly obtained.

Introduction

Quantum logic developed by Birkhoff and von Neumann in their seminal article in 1936 [START_REF] Birkhoff | The Logic of Quantum Mechanics[END_REF] considers logical propositions as subspaces of a quantum state Hilbert space. As will be shown hereafter and also underlined in [START_REF] Ying | Foundations of Quantum Programming[END_REF], these subspaces can be viewed as eigenspaces of projectors, the projectors corresponding to logical propositions. A true proposition is then associated to the eigenvalue +1. The representation of logical propositions in a vector space could be of interest in modern semantic theories such as distributional semantics, for example using the "Hyperspace Analogue to Language" algorithm as was done in [START_REF] Barros | Contextual Query Using Bell Tests[END_REF], or in connectionist models of cognition [START_REF] Busemeyer | Quantum models of cognition and decision[END_REF]. In this work we show that a proposition in a logical system can be represented by an observable in Hilbert space. When interpreted in the context of quantum mechanics this model uses finite dimensional projectors and angular momentum observables. Conversely, a quantum system when considered in its eigenspace is formally equivalent to a logical propositional system. The view here, which comes under the name of "Eigenlogic" (for the original motivation and more detailed discussion see [START_REF] Toffano | Eigenlogic in the spirit of George Boole[END_REF]), considers that the eigenvalues of the logical observables are the truth values of a proposition and the associated eigenvectors correspond to the different input atomic propositional cases. When considering vectors outside of the eigensystem this view leads to a "fuzzy" measure of the degree of truth of a logical proposition.

In our model for binary valued logic, using numbers {0, 1}, the logical observables are pairwise commuting projectors. The model is extended to the other binary system using numbers {+1, -1}, differences reside in the symmetry of the corresponding logical observables. In the latter case the observables are equivalent to quantum spin 1 /2 observables, no more idempotent projectors but isometric self-inverse reflection observables squaring to 1. These are equivalent to the recently proposed "quantum Boolean functions" [START_REF] Montanaro | Quantum Boolean functions[END_REF] developed in the context of the research topic "Fourier analysis of Boolean functions" having many applications in theoretical computer science, information theory and also in social decision and voting theory. We then propose an algebraic generalization, based on the finite-elements method, that can be applied to whatever m-value n-argument logical system.

The paper is organized as follows: we start with Boolean two-valued {0, 1} logic and we demonstrate important expressions for the projector observables in the 2-argument case indicating also the general method for n-arguments. The case for binary values {+1, -1} is then presented. Then we consider the case for fuzzy logical propositions and give the method for calculating fuzzy membership functions by using the Born rule and show that these functions can be identified with probabilities. The last section is devoted to the many-valued systems (m > 2) the case of 3-valued 2-argument logic is discussed with some examples of applications.

2 Two-valued Eigenlogic

Projector two-valued logic

We will consider a two-dimensional rank-1 projector Π acting on a single set. What are the expected outcomes when applying this projector? If, for example, vector |a > corresponds to an element of the set, the following matrix equation will be verified: Π • |a >= 1 • |a >. The value 1 being the eigenvalue of the projector associated with the eigenvector |a >. Interpretable results [START_REF] Toffano | Eigenlogic in the spirit of George Boole[END_REF] considered in a two-value {0, 1} logical system will correspond to the possible eigenvalues 0 and 1, where 0 is the result for elements not belonging to the set. So in this way a question concerning the proposition of belonging or not to a particular set, will have as an answer one of the two eigenvalues. The "true" value 1 will correspond to the eigenvector |a >, now named |1 >, and the "false" value 0 will correspond to the complementary eigenvector |a >, named |0 >. When these properties are expressed in matrix form: vectors |1 > and |0 > become 2 dimensional orthonormal column vectors and the projection operators 2 × 2 square matrices. This gives:

|1 > = 0 1 , |0 > = 1 0 .
The choice of the position of the value 1 in the column vectors is arbitrary, here it follows the quantum information convention for a "qubit-1" [START_REF] Nielsen | Quantum Computation and Quantum Information[END_REF]. As usual in Quantum Mechanics we can find the set of projectors that completely represent the quantum system, in particular by lifting the eventual degeneracy of the eigenvalues. Here eigenvalues are always equal to 0 or 1 and the question about the multiplicity of eigenvalues is natural. In this contribution we focus on different projective structures that completely define the logical system. In the very simple case where 0 and 1 are both not degenerate eigenvalues, the projectors relative to the eigenvector basis take the form:

Π 1 = Π = 0 0 0 1 , Π 0 = I -Π = 1 0 0 0 . ( 1
)
We systematically consider all the possible structures of such projectors. When representing logic with n atomic propositions using projectors various possibilities are intrinsically present in a unique structure with 2 2 n different projectors.

Once the eigenbasis is chosen the remaining structure is intrinsic.

For example the two projectors shown in equation ( 1) are complementary and idempotent. One can give a general expression of a one-argument "logical observable" as an expansion over the commuting projectors Π 0 and Π 1 spanning the vector space:

F = f (0) Π 0 + f (1) Π 1 = f (0) 0 0 f (1) (2)
the coefficients f (0) and f (1) in the expansion are the truth values of the corresponding {0, 1} Boolean logical connective. Equation (2) represents the spectral decomposition of the operator and because the eigenvalues are real the logical operator is Hermitian and can thus be considered as a quantum observable. In this way, in Eigenlogic, the truth values of the logical proposition are the eigenvalues of the logical observable. One can then construct the 4 logical observables corresponding to the 4 one-argument Boolean connectives: A = Π 1 is the "logical projector" and A = I -Π 1 = Π 0 its complement. The "True" operator corresponds here to the identity operator I . The "False" observable corresponds to the null operator. These four observables form a complete family of commuting projectors. The extension to more arguments is obtained by using the Kronecker product ⊗ in the same way as for the composition of quantum systems (for technical details on this operation see for example [START_REF] Nielsen | Quantum Computation and Quantum Information[END_REF]).

In the case of n = 2 arguments we will have an expansion over 4 commuting orthogonal rank-1 projectors . Some properties of the Kronecker product on projectors have to be specified: (i) The Kronecker product of two projectors is also a projector; (ii) If projectors are rank-1 projectors (a single eigenvalue is equal to 1, all the others are 0) then their Kronecker product is also a rank-1 projector. Using these properties, the 4 commuting orthogonal rank -1 projectors Π 00 , Π 01 , Π 10 , and Π 11 , spanning the 4 dimensional vector space are calculated in a straightforward way:

Π 00 = (I -Π) ⊗ (I -Π) , Π 01 = (I -Π) ⊗ Π , Π 10 = Π ⊗ (I -Π) , Π 11 = Π ⊗ Π .
So one can write the logical observable for n = 2 arguments:

F = f (0, 0) Π 00 + f (0, 1) Π 01 + f (1, 0) Π 10 + f (1, 1) Π 11 . (3)
In an explicit way:

F =     f (0, 0) 0 0 0 0 f (0, 1) 0 0 0 0 f (1, 0) 0 0 0 0 f (1, 1)     .
Equation ( 3) represents a spectral decomposition with the eigenvalues being the truth values, in this case we will have a family of 16 possible different observables. All these observables are pairwise commuting projectors and in general their product (matrix product) is not equal to zero. This last point is essential in the model, because not only mutually exclusive projectors are representative for a logical system, the complete family of projectors must be used. For example the observables for conjunction, AND, and disjunction, OR, which have in common the truth value, (1, 1), for the input combination (True ≡ 1 , True ≡ 1), have their matrix product different from zero. This method can be extended to whatever number of arguments n using the "seed" projector Π, its complement (I -Π) and by applying the Kronecker product. So given the number of input arguments n and knowing the truth table of the logical connective one directly obtains the corresponding binary Eigenlogic observable. Now let's develop the case for n = 2 arguments: one can express the connectives corresponding to a "logical projector" according to the composition rule, thus obtaining two commuting projector observables:

A = Π ⊗ I , B = I ⊗ Π , A • B = Π ⊗ Π (4)
the conjunction, AND, observable becomes simply the product of these two logical projectors A • B. The disjunction, OR, and exclusive disjunction, XOR, observables are shown on Table 1, where the algebraic expansions for Boolean connectives explicitly derived in [START_REF] Toffano | Eigenlogic in the spirit of George Boole[END_REF] are used. Negation (complementation) is obtained by subtracting from the identity operator for projective logical observables and by multiplying by -1 for isometric logical observables (see hereafter). Useful transformations are obtained by De Morgan's theorem (for general theorems in logic see for example Knuth [START_REF] Knuth | The Art of Computer Programming[END_REF]), for the negative conjunction, NAND one has the identity A ∧ B = A ∨ B in the same way one can obtain NOR with the identity A ∨ B = A ∧ B. Implication observables are also shown on Table 1.

connective for truth table {0, 1} projective {+1, -1} isometric Boolean {F, T} : logical logical A, B {0, 1} ; {+1, -1} observable observable False F F F F F 0 +I NOR ; A ∨ B F F F T I -A -B + A • B 1 2 (+I -U -V -U • V ) A B F F T F B -A • B 1 2 (+I -U + V +U • V ) A F F T T I -A -U A B F T F F A -A • B 1 2 (+I + U -V + U • V ) B F T F T I -B -V XOR ; A ⊕ B F T T F A + B -2A • B U • V = Z ⊗ Z NAND ; A ∧ B F T T T I -A • B 1 2 (-I -U -V + U • V ) AND ; A ∧ B T F F F A • B = Π ⊗ Π 1 2 (+I + U + V -U • V ) A ≡ B T F F T I -A -B + 2A • B -U • V B T F T F B = I ⊗ Π V = I ⊗ Z A ⇒ B T F T T I -A + A • B 1 2 (-I -U + V -U • V ) A T T F F A = Π ⊗ I U = Z ⊗ I A ⇐ B T T F T I -B + A • B 1 2 (-I + U -V -U • V ) OR ; A ∨ B T T T F A + B -A • B 1 2 (-I + U + V + U • V ) True T T T T T I -I
Table 1. The sixteen two-argument two-valued logical connectives and the respective Eigenlogic observables for eigenvalues {0, 1} and {+1, -1}.

Isometric reversible two-valued logical observables

There is a linear bijection (isomorphism) from the projector logical observables F towards reversible observables G:

G = I -2F .
The two families of observables commute and have the same system of eigenvectors. Practically to obtain G from F one just has to substitute the eigenvalue 0 with +1 and 1 with -1. The observables G are "isometries": unitary reflection operators. From projector Π in equation ( 4) one obtains the observable Z:

Z = I -2Π = +1 0 0 -1 = σ z
which is actually one of the Pauli matrices σ z and corresponds in quantum mechanics, to the z component of a spin 1 /2 observable S z = ( /2) σ z where is the reduced Planck's constant. In the field of quantum information this operator is also named the "Pauli-Z" gate or "phase-π" gate [START_REF] Nielsen | Quantum Computation and Quantum Information[END_REF]. Here, U = Z designates the "logical projector" connective and U = -Z its complement (negation), nota bene in this case the connective "logical projector" is not a projection operator, in order to avoid ambiguity it is often named [START_REF] Montanaro | Quantum Boolean functions[END_REF] "dictator".

For n = 2 arguments one can then write directly the expression for a logical isometric observable by using its spectral decomposition. The logical "dictators" U and V become:

U = Z ⊗ I , V = I ⊗ Z , U • V = Z ⊗ Z .
The exclusive disjunction XOR observable is here simply given by the product of the dictators: U • V . Negation is obtained by multiplying by the number -1.

From table 1 one sees that there are more complicated relations, for example the conjunction, AN D, observable is:

1 2 (I + U + V -U • V ) =     +1 0 0 0 0 +1 0 0 0 0 +1 0 0 0 0 -1     = C Z .
Those familiar with the domain of quantum information can easily recognize the reversible logical gate "control-Z" or simply named C Z [START_REF] Nielsen | Quantum Computation and Quantum Information[END_REF].

3 From deterministic logic to fuzzy logic Fuzzy logic deals with truth values that may be any number between 0 and 1, here the truth of a proposition may range between completely true and completely false. It is generally considered that probability theory and fuzzy logic are related to different forms of uncertainty, the first is concerned with how probable it is that a variable belongs to a given set and the second one uses the concept of fuzzy set membership, intended as the degree of membership. This was the first motivation of fuzzy logic [START_REF] Zadeh | Fuzzy sets[END_REF]. But this distinction when considering the quantum probabilistic Born rule is not so strict from a formal point of view. We will start the discussion by giving the interpretation of a vector state in Eigenlogic.

In the preceding sections we considered operations on the eigenspace of a logical observable family. For example for n = 2 arguments a complete family of 16 commuting logical observables represents all possible logical connectives and becomes "interpretable" [START_REF] Toffano | Eigenlogic in the spirit of George Boole[END_REF] when applied to one of the four possible canonical eigenvectors of the family. These vectors, corresponding to all the possible atomic input propositional cases, are represented by the vectors |00 >, |01 >, |10 > and |11 > forming a complete orthonormal basis. When applying a logical observable on one of these vectors the resulting eigenvalue will correspond to the truth value for the considered input.

Now what happens when the state-vector is not one of the eigenvectors of the logical system? In quantum mechanics, where vectors operate in Hilbert space, one can always express a state-vector as a decomposition on a complete orthonormal basis. In particular we can express it over the canonical eigenbasis of the logical observable family. For two-arguments this vector can be written as:

|Ψ > = C 00 |00 > + C 01 |01 > + C 10 |10 > + C 11 |11 > .
We can interpret this in the following way: when only one of the coefficients is non-zero (in this case its absolute value must take the value 1) then we are back in the preceding situation of a determinate input atomic propositional case.

But when more than one coefficient is non-zero we are in a "mixed" or "fuzzy" propositional case. Such a state could also possibly be interpreted as a quantum superposition of atomic propositional cases.

We can then calculate the "mean value" of a logical observable. In particular the logical projector observables F will give a "fuzzy measure" of the logical proposition in the form of the "fuzzy membership function" µ. Let's show this on some examples: in the case of one argument one can express an arbitrary 2-dimensional quantum state as: |ϕ >= sin α |0 > +e iβ cos α |1 > where the "angles" α and β are real numbers. The quantum mean value of the "logical projector" observable A = Π can then be calculated using the Born rule:

µ(a) = < ϕ|Π|ϕ > = cos α e -iβ < 1|1 > < 1| cos αe iβ |1 > = cos 2 α ;
in the same way one can calculate the complement

µ(a) = < ϕ|I -Π|ϕ > = sin 2 α = 1 -µ(a) .
This verifies one of the requirements of fuzzy logic for the complement (negation) of a fuzzy set.

According to standard notations for spin 1 /2 quantum states, or qubits, on the Bloch sphere [START_REF] Nielsen | Quantum Computation and Quantum Information[END_REF] we use the transformation α = (π-θ) /2 and β = ϕ. A quantum compound state can be built by taking the tensor product of two elementary states: |ψ >= |ϕ p > ⊗|ϕ q >, where |ϕ p >= cos θp 2 |0 > +e iϕp sin θp 2 |1 > (for |ϕ q > we have a similar expression). Now sin 2 θp 2 = p and sin 2 θq 2 = q represent the probabilities of being in the "True" state |1 > for spins 1 /2 oriented along two different axes θ p and θ q .

One can calculate the fuzzy membership function of the corresponding "logical projector" for the two-argument case using equation [START_REF] Busemeyer | Quantum models of cognition and decision[END_REF].

µ(a) =< ψ|Π ⊗ I|ψ >= p(1 -q) + p • q = p , µ(b) =< ψ|I ⊗ Π|ψ >= q .
This shows that the mean values correspond to the respective probabilities. Now let's "measure" for example the conjunction and the disjunction, using the observables in table 1, this gives:

µ(a∧b) = < ψ|Π ⊗ Π|ψ > = p • q = µ(a) • µ(b) , µ(a∨b) = p + q -p • q = µ(a) + µ(b) -µ(a) • µ(b) .
Similar results for conjunction and disjunction have been outlined recently, also using projector operators, when considering concept combinations [START_REF] Aerts | Quantum structure of negation and conjunction in human thought[END_REF] for quantum-like experiments in the domain of quantum cognition.

What happens when the state-vector cannot be put in the form of a tensor product, that is when it corresponds to an entangled state? The problem is outside the scope of this paper but an interesting result can be shown: the mean value of whatever logical observable of the type F on an arbitrary quantum state |Ψ > will always verify the inequality:

< Ψ |F |Ψ > = Tr (ρ Ψ • F ) ≤ 1 , with ρ Ψ ≡ |Ψ > < Ψ | ,
and can thus be interpreted as a probability measure.

From two-valued to multi-valued logic

Multi-valued logic requires a different algebraic structure than an ordinary binaryvalued one. Many properties of binary logic do not support set of values that do not have cardinality 2 n . Multi-valued logic is often used for the development of logical systems that are more expressive than Boolean systems for reasoning [START_REF] Miller | Multiple Valued Logic: Concepts and Representations[END_REF]. Particularly three and four valued systems, have been of interest with applications to digital circuits and computer science. The total number of possible logical connectives for an m-valued n-argument system is the combinatorial number m m n , so in particular for a binary 2-valued 2-argument system, as shown above, the number of connectives will be 2 2 2 = 16, the complete list indicated on table 1. For a binary three-argument system, the number increases to 2 2 3 = 256. For a 3-valued 1-argument system the number of connectives will be 3 3 1 = 27 and for a 3-valued 2-argument system: 3

3 2 = 19683.
So it is clear that by increasing the values from two to three the possibilities of new connectives becomes intractable for a complete description of a logical system, but some special connectives play important roles and will be illustrated hereafter. We will proceed by showing the general algebraic method.

Interpolation with finite elements

The finite element method (see for example [START_REF] Zienkiewicz | The Finite Element Method in Engineering Science[END_REF]) allows one to interpolate a function, id est to make explicit the values f (x) from the given values of specific numbers, the (so-called) degrees of freedom. Let's consider the following simple example: given the values f (+1), f (0) and f (-1) of a function f at the particular points x = +1, 0, -1, and using the appropriate Dirac linear forms, we can write: < δ +1 , f >= f (+1) , < δ 0 , f >= f (0) and < δ -1 , f >= f (-1), where Σ ≡ { δ +1 , δ 0 , δ -1 } is called the set of degrees of freedom. This linear structure shows that it is natural to consider a three-dimensional space. The so-called "basis function" ϕ i associated to the set of degrees of freedom Σ and to the polynomial space solves this problem. Using the degrees of freedom, and second-degree polynomials, we obtain the three basis functions

ϕ +1 (x) = 1 2 x (x + 1) , ϕ 0 (x) = 1 -x 2 , ϕ -1 (x) = 1 2 x (x -1) . (5)
So in general, an arbitrary function f can be written:

f (x) = i=+1,0,-1 f (i) ϕ i (x) , i=+1,0,-1 ϕ i (x) ≡ 1 (6)
where the completeness of the basis functions is verified by their sum being 1.

Formalization of three-valued Eigenlogic

We use an operator system which is equivalent to the one of orbital angular momentum = 1. In general angular momentum is characterized by two quantum numbers: j the angular momentum number and m j the magnetic momentum number. Both these numbers must be integer or half integer. The rules are: j ≥ 0, and attached to this value we have the condition: -j ≤ m j ≤ j. The value j = 0 is possible and gives a single value m j = 0 the next is j = s = 1 /2 giving two values m s = ± 1 /2 corresponding to the two-valued spin system. The value j = 1 gives three possible values m j = {+1, 0, -1} and so on. We consider for j = = 1 the z-component orbital angular momentum observable [START_REF] Schiff | Quantum Mechanics[END_REF]

L z = Λ =   +1 0 0 0 0 0 0 0 -1   . (7)
In the above matrix the three eigenvalues {+1, 0, -1} will be considered as the logical values. A convention for these values, extending binary logic, is the following: False : F ≡ +1 , Neutral : N ≡ 0 , True : T ≡ -1 .

We can now express the three-value logical observables as spectral decompositions over the rank-1 projectors spanning the vector space: Π +1 , Π 0 and Π -1 . These operators correspond to the pure state density matrices of the three eigenstates | + 1 > , |0 > and | -1 > of L z . The three projectors can be expressed as a function of the dimensionless observable Λ, using directly the expressions given above in [START_REF] Toffano | Eigenlogic in the spirit of George Boole[END_REF] where the basis functions ϕ i become the projectors and the symbol x the observable Λ given in [START_REF] Nielsen | Quantum Computation and Quantum Information[END_REF]:

Π +1 = 1 2 Λ (Λ + I) Π 0 = I -Λ 2 Π -1 = 1 2 Λ (Λ -I) (8)
Then every one-argument "local projector" F (Λ) can be obtained using the relation (6).

Three-valued, two-argument examples: Min, Max

When considering a 2-argument 3-valued system we find the expansion by using the Kronecker product in the same way as for the binary system in equation ( 3):

F = i, j =+1, 0, -1 f ij Π i ⊗ Π j , f ij ∈ {+1, 0, -1} . ( 9
)
these observables are now 9 × 9 matrices. We can define the two argument "dictators", U and V , simply by the rule of composition, this leads to:

U = Λ ⊗ I V = I ⊗ Λ U • V = Λ ⊗ Λ . (10)
In trivalent logic (see e.g. [START_REF] Miller | Multiple Valued Logic: Concepts and Representations[END_REF]) popular connectives are Min and Max, defined in the maps on table 2.

Here the connectives Min and Max are symmetric, they are equivalent for a complete inversion of signs on inputs and outputs. Using the relations ( 8), ( 9) and [START_REF] Aerts | Quantum structure of negation and conjunction in human thought[END_REF] in conjunction with reduction rules we obtain the following observables:

     Min(U , V ) = 1 2 U + V + U 2 + V 2 -U • V -U 2 • V 2 Max(U , V ) = 1 2 U + V -U 2 -V 2 + U • V + U 2 • V 2 (11) Min U \ \ V F N T F ≡ +1 +1 +1 +1 N ≡ 0 +1 0 0 T ≡ -1 +1 0 -1 Max U \ \ V F N T F ≡ +1 +1 0 -1 N ≡ 0 0 0 -1 T ≡ -1 -1 -1 -1
Table 2. The Min and Max maps for a three-valued two-argument logic.

The proof of the relations ( 11) is a direct consequence of relations ( 9) and ( 5).

We have on one hand:

Min (U, V ) = ϕ 1 (U ) ⊗ ϕ 1 (V ) + ϕ 1 (U ) ⊗ ϕ 0 (V ) + ϕ 1 (U ) ⊗ ϕ -1 (V ) + ϕ 0 (U ) ⊗ ϕ 1 (V ) + ϕ -1 (U ) ⊗ ϕ 1 (V ) -ϕ -1 (U ) ⊗ ϕ -1 (V ) = ϕ 1 (U ) + ϕ 1 (V ) -ϕ 1 (U ) ⊗ ϕ 1 (V ) -ϕ -1 (U ) ⊗ ϕ -1 (V ) due to (6) = 1 2 U (U + I) + 1 2 V (V + I) -1 4 U (U + I) V (V + I) -1 4 U (U -I) V (V -I) = 1 2 U 2 + U + V 2 + V -U 2 V 2 -U V
and the first relation of ( 11) is proven. On the other hand, we have Max (U,

V ) = ϕ 1 (U ) ⊗ ϕ 1 (V ) -ϕ 1 (U ) ⊗ ϕ -1 (V ) -ϕ 0 (U ) ⊗ ϕ -1 (V ) -ϕ -1 (U ) ⊗ ϕ -1 (V ) -ϕ -1 (U ) ⊗ ϕ 1 (V ) -ϕ -1 (U ) ⊗ ϕ 0 (V ) = ϕ 1 (U ) ⊗ ϕ 1 (V ) -ϕ -1 (U ) -ϕ -1 (V ) + ϕ -1 (U ) ⊗ ϕ -1 (V ) due to (6) = 1 4 U (U + I) V (V + I) -1 2 U (U -I) -1 2 V (V -I) + 1 4 U (U -I) V (V -I) = 1 2 U 2 V 2 + U V -U 2 -V 2 + U + V
and the second relation of ( 11) is proven. The proof presented above exploits the properties of the Kronecker product and reduction rules due to the completeness of the finite projection space. Reduction of logical expressions is an important topic in logic. In binary logic it is formalized by using Karnaugh maps which represent canonical SOP (Sum Of Products) disjunctive normal forms [START_REF] Knuth | The Art of Computer Programming[END_REF]. Binary logic is "included" in ternary logic, we want to verify this by eliminating the "neutral" state, N ≡ 0, and considering only the two logical values {+1, -1}. In this case we have: U 2 = V 2 = I and so [START_REF] Miller | Multiple Valued Logic: Concepts and Representations[END_REF] reduces to:

Min (U , V ) = 1 2 I + U + V -U • V , Max (U , V ) = 1 2 -I + U + V + U • V considering that
for binary logic the Min connective becomes the conjunction, AND, and the Max connective the disjunction, OR, we find the previous results given on table 1 for binary {+1, -1} observables.

Discussion and Conclusion

We have presented an operational formalism named "Eigenlogic" using observables in Hilbert space. The original feature being that the eigenvalues of a logical observable represent the truth values of the corresponding logical connective, the associated eigenvectors corresponding to one of the fixed combination of the inputs (atomic propositions). This approach differs from other geometric formalizations of logic (for references and discussion see [START_REF] Toffano | Eigenlogic in the spirit of George Boole[END_REF]). Here the outcome of a "measurement" or "observation" on a logical observable will give the truth value of the associated logical proposition, and becomes "interpretable" when applied to the eigenspace leading to a natural analogy with the measurement postulate in quantum mechanics. One of the referees proposed the following diagram to summarize the point of view presented in this contribution: eigenvectors in Hilbert space -→ atomic propositional cases projectors -→ logical connectives eigenvalues -→ truth values.

At first sight this method could be viewed as "classical" because exactly the same results are obtained in Eigenlogic as in ordinary propositional logic. This is in itself an important result demonstrating a new method in logic based on linear algebra, the method being also developed in multivalued logic. But when considering vector states, id est input propositions, that are not eigenvectors, the measurement outcomes are governed by the quantum Born rule, and interpretable results are then given by the mean values. This fact led us to apply the method to Fuzzy logic.

Another important point is the general algebraic method, based on classical interpolation framework suggested by the finite-element method. Our method can be employed for whatever m-valued n-argument logical system and in each case the corresponding logical observables can be defined. Some observables can be formally compared with angular momentum observables in quantum mechanics.

Because of the exponential increase of complexity, an analytical formulation is only tractable for a low number of logical values and arguments. We treated the two-argument binary case completely and the three-valued case using the logical observables Min and Max. An algorithmic approach for logical connectives with a large number of arguments could be interesting to develop using Eigenlogic observables in high-dimensional vector spaces. But because the space grows in dimension very quickly, it may not be particularly useful for practical implementation without logical reduction. It would be interesting to develop specific algebraic reduction methods for logical observables inspired from actual research in the field. For a good synthesis of the state of the art, see e.g. [START_REF] Yanushkevich | Introduction to Logic Design[END_REF].

Eigenlogic could create a new perspective in the field of quantum computation because several of the observables turn out to be well-known quantum gates.

Here we represent them as diagonal matrices, id est in their eigenbasis, other "normal" forms being easily recovered by unitary transformations. It would be interesting to operate quantum gates in our framework. Many-valued logic is being investigated in quantum computation for example with ternary-logic quan-tum gates using "qutrits". Our formulation of multivalued logical observables could be used for the design of new quantum gates. Dynamical evolution of the logical system could be included in the model by identifying the appropriate Hamiltonian operators. Standard procedures for expressing interaction Hamiltonians as a function of angular momentum observables could be used [START_REF] Schiff | Quantum Mechanics[END_REF].

More generally we think that this view of logic could add some insight on more fundamental issues. Boolean functions are nowadays considered as a "toolbox" for resolving many problems in theoretical computer science, information theory and even fundamental mathematics. In the same way Eigenlogic can be considered as a new "toolbox" and could be of interest for the "Quantum Interaction" community where quantum-like approaches in human and social sciences need to be founded on a logical basis.

Acknowledgments

The authors thank both referees for their precise and constructive remarks and suggestions. Some of them have been included in the present version of this contribution.


R Hegedus

A Lucat

J Redon

R Pacanowski

Workshop on Material Appearance Modeling (2016) H. Rushmeier and R. Klein (Editors) Isotropic BRDF Measurements with Quantified Uncertainties

Keywords: (according to ACM CCS): I, 3, 3 [Computer Graphics]: Picture/Image Generation-BRDF Measurements Analysis 1, Motivation and Related Work

Image-based BRDF measurements on spherical material samples present a great opportunity to shorten significantly the acquisition time with respect to more traditional, non-multiplexed measurement methods for isotropic BRDFs. However, it has never been analyzed deeply, what measurement accuracy can be achieved in such a setup; what are the main contributing uncertainty factors and how do they relate to calibration procedures. In this paper we present a new set of isotropic BRDF measurements with their radiometric and geometric uncertainties acquired within such an imaging setup. We discuss the most prominent optical phenomena that affect measurement accuracy and pave the way for more thorough uncertainty analysis in forthcoming image-based BRDF measurements. Our newly acquired data with their quantified uncertainties will be helpful for comparing the quality and accuracy of the different experimental setups and for designing other such image-based BRDF measurement devices.

Since the seminal work of Marschner [Mar98] or Matusik et al. [MPBM03] on BRDF acquisition, a lot of devices and measurement setups have been developed to acquire the appearance of materials.

For the last two decades, most researchers from the Vision or CG communities have focused on acquiring acquire more complex appearance functions, such as spatiallyvarying BRDF (e.g., [McA02, LBAD * 06]) to Bidirectional Texture function (e.g., [START_REF] Dana | Reflectance and texture of real-world surfaces[END_REF]) or tried to reduce setups complexity or cost (e.g., [RWS * 11]). Unfortunately, most setups or available measurements data are hard to compare because they rely on qualitative validation instead of a quantitative one. For example, many papers (e.g., [NDM05]) that try to fit BRDF analytical models discard measurements at arbitrary grazing angles because they «seem» unreliable.

We believe that new available data with their quantified uncertainties will be helpful for the research community, especially in order to compare the quality and accuracy of the different experimental setups and even to help other research teams to develop their devices. Furthermore, uncertainty values can be useful to guide approximation techniques (e.g., Levenberg-Marquardt, SQP, Weighted Least-Square),which are necessary to obtain parameters of analytical models. This paper presents a new set of isotropic BRDF measurements with their radiometric and geometric uncertainties acquired from our imaging device (cf. Section 2). The uncertainties are derived from a model (cf. Section 3) that take into account different source of errors (e.g., focal length, mechanical, ...). In section 4 we illustrate through various plots the different uncertainties obtained with our setup and also illustrate the acquired BRDFs. Finally, we discuss the limitations of our current uncertainty model and offer possible improvements for it.

BRDF Setup

Our setup shares similarities from the one introduced by Matusik et al.. As shown in Figure 1 in our setup we measure the BRDF of materials on a spherical sample by acquiring multiple images using a Jenoptik scientific camera (ProgRes C14 Plus). The camera is firmly mounted on an arm that is fixed to a rotation stage (from Zaber T-RS60). The material sample (i.e., a target sphere of approx. 6cm) is placed on a support also fixed to the rotating arm. Contrary to previous setups, our light source is a fixed collimated beam that lits the entire target sphere. In another words the camera and the target sphere are always facing each other, regardless of the rotation step. However, they are moving wrt. to the light source thanks to the rotating arm, which makes possible to c 2016 The Author(s)

Lens Camera

Material Sample

Collimation Setup (3 lenses + aperture)

Inco min g colli mat ed light

Liquid light guide

Optical Table

Rotation Stage Light Source

Figure 1: Schematic view of our acquisition setup.

have the three degrees of freedom required to measure an isotropic BRDF.

The goal of our setup is to measure a spatially uniform isotropic BRDF. As defined by Nicodemus [NRH * 77], the BRDF f r is a positive quantity expressed in sr -1 (inverse steradian) representing how a material reflects angularly the light. For a given point p S on a surface S f r is defined as the ratio of the differential reflected radiance L r and the differential irradiance dE:

f r (θ i , θ o , ∆ φ ) = L r (θ o , φ o ) dE(θ i , φ i )
where all angles (θ i , θ o , ∆ φ = |φ o -φ i |) are expressed in the local (tangent) frame of the surface orientated with a normal vector n s . When using a broadband white light source with an RGB-sensor, the measured BRDF can be seen as an angular function returning an RGB HDR value:

f r (θ i , θ o , ∆ φ ) → (ρ r , ρ g , ρ b ) . (1)
The geometric calibration step is responsible to establishing the angular values (θ i , θ o , ∆ φ ), whereas the radiometric calibration is necessary for the BRDF values themselves.

Geometric Calibration As in Matusik et al. setup, the fact to have a fixed camera with respect to the target sphere presents the advantage to calibrate (geometrically) the camera only once. Calibrating the camera is done using the Caltech library [Bou] with the associated thin-lens model. This model is also able to take into account distortion, but the calibration process did not compute meaningful values for our camera and therefore we simply discarded these parameters. The camera calibration allows us to recover the perpixel view direction as well as the focal length.

The second calibration step is to recover the target sphere position in the camera frame, assuming its radius is known (e.g., provided by the manufacturer). This can be done easily by noticing that a perspective projection of a sphere provides a perfect elliptical contour. Fitting this ellipse directly gives the sphere center, with a very good accuracy. space), we use a mirror sphere as target and record multiple images of the specular highlight visible on the sphere. By assuming a perfect mirror, we compute the reflected ray impacting from the light spot for various camera-light configurations. This process gives a 3D locus of the light direction, which can be accurately described by a rotation axis and a latitude.

We deduce from the camera calibration and the radius of the target sphere its position wrt. to the camera. Then, we ray-trace a virtual sphere in order to estimate the normal vector in each pixel. Finally, combining the light direction, the view vectors and the per-pixel normal vectors we get an estimate, per pixel, of the three angles (θ i , θ o , ∆φ) of an isotropic BRDF measurement.

Radiometric Calibration

We analyzed the radiometric response of our camera, which is an important factor to consider in order that accurate radiance values be measured.

Spectralon White Teflon

Copper Carbon-Steel Since our setup accommodates a stabilized light source and the camera has an electronic shutter that allows for precise exposure time control, we have chosen a simple, yet efficient way of acquiring the response curve. The light was directed towards the camera lens, in front of which a thin diffuser sheet was fixed, thus the same image of the diffused source light was recorded multiple times with different exposure times. Then the median pixel value in each of R,G,B channels throughout the whole image was selected. In this way median value vs. exposure time directly yields the response curve, since exposure time is directly proportional to incoming total radiance. These data showed that radiometric response is totally linear as long as pixel value, encoded with 14 bits, is less than 15800 (values over 15800 were simply discarded). The other common correction related to color measurements is to make RGB values independent from the sensor spectral sensor sensitivity and the spectral light composition. This correction procedure can be done using a reference object, in our case a spherical Spectralon that is designed precisely to have a uniform spectral reflectance. Imaging this object allows to find the correct multiplicative factors to apply on each color channel in order to recover correct relative (inter-color-channel) values.

Material Samples As shown in Figure 3, we acquired four materials ranging from a very diffuse (Spectralon) to an almost perfect mirror sphere made of carbon-steel. We used the carbon-steel ball as reference for the radius, since the manufacturer provides its radius (6cm) with a 50 µm uncertainty. For the other material samples, we compute the radius using the carbon steel sphere as a reference. If we replace it by another sphere of unknown radius, we can recover its position and radius by the ellipse fitting method previously described, knowing the depth in the camera frame instead of the radius. We assumed bilateral symmetry (i.e., symmetry of the BRDF wrt. to the incident plane) of an isotropic BRDF and therefore our setup acquired 221 images (1360x1024 pixels) by moving the rotation arm 0.75 degrees in each step. This resolution gives about 60 million samples for each measured material. More precise measurements can be achieved by reducing the rotation step of the stage (up to. 0.0002 degree). However, it is important to keep in mind that the angular resolution of the BRDF measurement also depends on the camera sensor resolution, as it ranges from 0.1 • resolution for a normal viewing direction to 2 • at 50 • and reaches 12 • for grazing angles. θ i is ranging from 0 to 90 • , θ o from 0 to 85 • , ∆φ from 0 to 180 • . Finally, all material measurements will be released on the ALTA [BCP * 15] project website.

Uncertainty Model

In this paper we are interested in establishing the uncertainties σ associated with the BRDF values (σ ρ r , σ ρ b , σ ρ b ) as well as with the geometrical positions of the measurements (σ θ o , σ θ i , σ ∆φ ). Figure 4 shows how the uncertainties of a BRDF measurement are related to the calibration parameters of the different parts of the setup. Some uncertainties such as the radius of the sphere or the pixel size are directly given by manufacturers. Radii of the spheres were given up to 50µm, whereas the pixel size is known up to 0.005µm. For the mechanical uncertainty of the rotation stage we verified that the repeatability and uncertainty numbers given by the manufacturer were correct. In order to increase the precision of the measurement, the stage is always set to rotate only by a discreet number of steps. All other uncertainties are derived from either radiometric or geometric calibration procedure. More formally, the uncertainties σ of a given parameter is a function of the following parameters:

σ θ i ,∆φ = U(F, s ccd , r s , c s , α, A, φ stage ) (2) σ θ o = U(F, s ccd , r s , c s ) (3) σ ρ r , ρ g , ρ b = U(E r,g,b (p ccd ), ω o , F, θ i ) (4)
Instead of deriving an analytical formula we choose to numerically compute the partial derivatives of U for each parameter. For our setup we obtain the following uncertainties: the focal length of the camera uncertainty is σ F = 0.3360mm and the sphere position σ c s = (0.0051225, 0.0026681, 4.879) mm.

When considering a uniform transmitting lens, a uniform illumination light-field, and neglecting convolution effects (diffraction, aberrations, and light divergence), one can show that the BRDF is governed to the irradiance on the sensor by the following equation:

f r (θ o , θ i , ∆φ) = K E r,g,b (p ccd ) F 2 cos θ i (cos θ o ) 4 (5)
where K is scaling factor that can be determined by using a known target such as our Spectralon. Therefore, the radiometric related uncertainties σ ρ r , ρ g , ρ b depend on the HDR values, for which the main source of error is due to photon c 2016 The Author(s)

Figure 4: Relationship and connection of the different uncertainties. The different colors illustrate the source of errors that are taken into account in our uncertainty model. The full analytical formula that relates every setup parameter to the BRDF value is in reality more complex and, particularly, non-invertible. For this reason, we have to use an algorithm that gives an approximate result of the measurement, in our case, the zero-order approximation (equation 5). Much complex approximations would lead to more accurate BRDF values, but we will always have to consider a reconstruction bias. shot noise as well as other sources such as the square of the focal length, the cosine factor cos θ i vignetting effect (or fall-off) (cos θ o ) 4 . For our setup the fall-off factor is very low since its maximum value is 0.9957 at the silhouette of the sphere (compared to 1.0 at the center). Figures 5 presents the uncertainties values related to the geometrical calibration. As one could expect, geometrical calibrated values become less and less reliable (up to 4 order of magnitude) when approaching grazing angle. The full analytical formula that relates every setup parameter to the BRDF value is in reality more complex and, particularly, non-invertible. For this reason, we have to use an algorithm that gives an approximate result of the measurement, in our case (equation 5), the zero-order approximation.

Measurements Results

Figure 6 presents the BRDF measurements in the incident plane with uncertainties for the different material samples.

Figure 7 shows the measured values of BRDF for Teflon material in 2D polar plot. Overall, the error on the BRDF values varies from 3% at normal incidence to 5% at 50 • and reaches 30% at 85 •

Discussion

The presented measurements with their uncertainties are primarily valid for our precise setup. Because of the complex interrelations among uncertainties of various parameters pertaining to such measurements, some general conclusions can be hard to find, nevertheless some hints can be deduced from our study. Below we discuss three important optical effects Spectralon Teflon Copper Mirror Diffraction from the Camera Aperture An important optical factor that needs proper attention in a camera-based BRDF acquisition setup, and is often overlooked, is diffraction originating from the aperture of the camera. The more stopped-down is the aperture the more diffraction will arise. Since we use a well stopped down aperture (f/8) in the camera lens and with good reason the effect cannot be neglected.

Although opening the aperture would practically eliminate problems caused by diffraction, it is not a good idea, because then lens aberrations would quickly become dominant and depth of field (which is needed for spherical samples) would also reduce, causing inaccuracies of an other nature that are actually more difficult to handle. The main effects of diffraction on the accuracy of obtained BRDF values are twofold (although in reality they arise from the same, single physical effect causing a convolution in the image). Firstly, there is a general blurring effect that can be modeled by taking into account the first order of a circular Airy disk pattern as convolution kernel. The second effect, is the so-called starburst effect on highlights, the cause of which is that light gets diffracted through the aperture (hexagonal blades in our setup). In our experience the general blurring effect in our setup is negligible compared to other convolution-type effects (depth of field confusion, imperfect collimation of the light source, and surface roughness), and that is for this reason that we set the aperture to a large value (f/8). However, the star-burst pattern around highlights is very much apparent and practically unavoidable. This forces us to drastically reduce the time of exposure range of the HDR imaging, and consequently the dynamic range of the BRDF measurement.

In reality, the sharper the highlights on the material surface the stronger this effect. Thus, in the worst-case scenario of BRDF with a very specular component, we can only reach a 20× factor between the shortest and the longest time of exposure. The more diffuse the material, the higher this factor, so that this restriction is eased. In this work, we have not addressed this problem yet and we simply stop recording data when the effect arises.

Diffraction White Light Speckle Another last effect of the diffraction is due to the surface roughness, it is generally known as the white light speckle [McK76]. It arises from the fact that on a small enough surface the light becomes coherent and produces grating-like light pattern when bouncing off a rough surface. If the surface is smooth or rough enough, this effect disappears but in between it becomes prominent.

For incident light angles up to 60 • , the main effect is changing locally color of the reflected light, and for more grazing angles, it also introduces a shift in the reflection angles (similarly to the effect of a grating). It then falsifies the material spatial uniformity and introduces a structured color noise. A major consequence is that, in reality, this coherence disappears with the divergence of the incoming light (as is the case in usual daily life when observing a scene), so using BRDF data from such a measurement introduces a perception bias, which can be inconvenient for a rendering purpose.

Out-of focus effect of common lenses One common issue with lenses is their lack of large depth of field, which introduces an error when measuring the BRDF with an imaging system. Indeed, one pixel of the CCD sensor recovers an irradiance value that is the convolution of the illumination with a portion (and not only a single point p S ) of the sphere and the BRDF. Deriving the uncertainty is a complex task that we did not take into account in this paper. However, we simulated precisely how large the different areas on the sphere are, across which the convolution applies. As shown in Figure 8 the area size ranges from 0.001mm 2 (part of the sphere which are in focus) to 2.1mm 2 (grazing angles). Note that we set the focus distance at 45 • latitude as a compromise.

Future Work

As future work we would like to improve our uncertainty model by taking into account the potential (spatial and directional) non-uniformity of the light source, or the lens and also remove on-the-fly pixels where the aperture diffraction appears. Furthermore, we would also like to average temporally our measurements to avoid white light speckle effects.We would also like to take into account correlation between the different uncertainties by computing numerically the co-variance matrix of the different factors. Furthermore, we would also like to study how the different BRDF models approximate these data with respect to the uncertainties provided. Finally, we will also provide additional measurements for other materials to improve our qualified database.

Acknowledgments

R. Hegedus is grateful to the Alexander von Humboldt Foundation, and acknowledges their support for his fellowship for experienced researchers.

Finally, in orderFigure 2 :

 2 Figure 2: Mathematical notations used in this paper. Please refer to Figure 4 to see how these quantities relate to each other in terms of uncertainties.

Figure 3 :

 3 Figure 3: The four materials measured for this paper.

Figure

  Figure Uncertainty values for the BRDF angles i , θ o , ∆φ) when φ stage = 7 • .

Figure 6 :Figure 7 :

 67 Figure 6: Incident plane BRDF measurements with uncertainty bars for the four materials.

Figure 8 :

 8 Figure8: Illustration of the out-of-focus effect when imaging a sphere. The picture shows for each pixel of the acquired image how much the optical rays from the lens camera cover the target sphere. The larger this spread the larger the error on the BRDF measurement. The spread is minimal only for the part of the target sphere which is at perfect focus. For our setup, the sphere radius is 6cm and the focus is made at 4.2cm from the sphere center (corresponding to an in-focus plane at 45 • latitude). The black area starts at θ o = 68 • and finishes at θ o = 85 • where the maximum area is about 2.1mm 2 .

c 2016 The Author(s)

Physical Elements of the Setup Calibration in the camera frame Measured and Calibrated Values


Xin Wang

Nicolas Thome

Matthieu Cord

GAZE LATENT SUPPORT VECTOR MACHINE FOR IMAGE CLASSIFICATION

Keywords: weakly supervised learning, latent SVM, eye-tracking, gaze features, CCCP

This paper deals with image categorization from weak supervision, e.g. global image labels. We propose to improve the region selection performed in latent variable models such as Latent Support Vector Machine (LSVM) by leveraging human eye movement features collected from an eye-tracker device. We introduce a new model, Gaze Latent Support Vector Machine (G-LSVM), whose region selection during training is biased toward regions with a large gaze density ratio. On this purpose, the training objective is enriched with a gaze loss, from which we derive a convex upper bound, leading to a Concave-Convex Procedure (CCCP) optimization scheme. Experiments show that G-LSVM significantly outperforms LSVM in both object detection and action recognition problems on PASCAL VOC 2012. We also show that our G-LSVM is even slightly better than a model trained from bounding box annotations, while gaze labels are much cheaper to collect.

INTRODUCTION

In the era of statistical machine learning, an overwhelming amount of images is available. For example, Ima-geNet [START_REF] Deng | ImageNet: A Large-Scale Hierarchical Image Database[END_REF] contains more than 10 million labeled images. The number of images also grows rapidly as multimedia communications and mobile devices become essential in our decades. Ever since 2013, 350 million new photos are posted on Facebook everyday 1 . We may never be able to label all images in a reasonable time by human labor.

The success of deep learning in computer vision largely improved the image classification accuracy on ImageNet [START_REF] Krizhevsky | Imagenet classification with deep convolutional neural networks[END_REF]. Also, the deep models trained on large scale datasets were found to adapt well to other specific vision tasks [START_REF] Wang | Recipe recognition with large multimodal food dataset[END_REF][START_REF] Chevalier | LR-CNN for Fine-Grained Classification With Varying Resolution[END_REF]. However, the problem is exacerbated when detailed annotations are used, e.g. bounding-boxes or pixel-wise labels, which are much onerous and expensive to collect than image-level annotations [START_REF] Blaschko | Tutorial: Visual learning with weak supervision[END_REF]. Additionally, using such accurate annotations during training can significantly boost categorization performances [START_REF] Durand | Image classification using object detectors[END_REF].

One option to get the best of both worlds is to develop 1 http://www.businessinsider.com/facebook-350-million-photos-eachday-2013-9?IR=T weakly supervised learning (WSL) frameworks. Basically, WSL consists in designing accurate models able to predict detailed annotations, e.g. region localization, while being trained from coarse labels, e.g. global image labels. In machine learning, the Latent Support Vector Machine (LSVM) [START_REF] Andrews | Support Vector Machines for Multiple-Instance Learning[END_REF] introduces a theoretically sound formalism for WSL. Several attempts have been devoted to applying LSVM for weakly supervised object detection and scene recognition [START_REF] Felzenszwalb | Object Detection with Discriminatively Trained Part Based Model[END_REF][START_REF] Pandey | Scene recognition and weakly supervised object localization with deformable partbased models[END_REF][START_REF] Russakovsky | Object-centric spatial pooling for image classification[END_REF][START_REF] Bilen | Object classification with latent window parameters[END_REF][START_REF] Juneja | Blocks that shout: Distinctive parts for scene classification[END_REF][START_REF] Sun | learning discriminative part detectors for image classification and cosegmentation[END_REF]. One challenge with LSVM is due to the introduction of latent variables, which makes the resulting optimization problem non-convex. To improve training performances, different solutions have recently been explored to apply the curriculum learning idea, i.e. how to find easy samples to incrementally train the model [START_REF] Kumar | Self-paced learning for latent variable models[END_REF][START_REF] Russakovsky | Object-centric spatial pooling for image classification[END_REF][START_REF] Bilen | Object classification with latent window parameters[END_REF]. When using sliding window approaches, the size of the latent space becomes huge. To overcome this issue, incremental exploration strategies have been proposed in [START_REF] Durand | Incremental learning of latent structural svm for weakly supervised image classification[END_REF][START_REF] Russakovsky | Object-centric spatial pooling for image classification[END_REF][START_REF] Bilen | Object classification with latent window parameters[END_REF]. Finally, recent works focus on enriching the prediction function, by using several (top) instance scores instead of using a single max [START_REF] Li | Multiple instance learning for soft bags via top instances[END_REF], or by incorporating negative evidence [START_REF] Durand | MANTRA: Minimum Maximum Latent Structural SVM for Image Classification and Ranking[END_REF][START_REF] Durand | WEL-DON: Weakly Supervised Learning of Deep Convolutional Neural Networks[END_REF].

In this paper, we propose to incorporate gaze features extracted from an eye-tracker device, to improve the training of LSVM models. Gaze features are appealing since they can be generated by humans at almost zero-cost when performing a recognition task. Collecting gazes is more user-friendly and less time-consuming than collecting traditional annotations: it takes about 1 second to collect gazes for one image [START_REF] Dim | Training object class detectors from eye tracking data[END_REF], comparing to 26s for drawing a bounding-box [START_REF] Su | Crowdsourcing Annotations for Visual Object Detection[END_REF] and 15-60 min for labeling the segmentation mask for an image [START_REF] Kohli | Robust Higher Order Potentials for Enforcing Label Consistency[END_REF]. For different purposes, people design different collection protocols to acquire gazes [START_REF] Lopez | One gaze is worth ten thousand (key-) words[END_REF][START_REF] Dim | Training object class detectors from eye tracking data[END_REF][START_REF] Mathe | Action from still image dataset and inverse optimal control to learn task specific visual scanpaths[END_REF][START_REF] Karthikeyan | From Where and How to What We See[END_REF]. The collection protocols can be grouped into two main groups: task-driven and free-viewing. Task-driven means the annotators are given a specific semantic to look at, e.g. dog or a group of actions. Free-viewing means the annotators view the image freely. In this paper, we use two task-driven datasets [START_REF] Dim | Training object class detectors from eye tracking data[END_REF][START_REF] Mathe | Action from still image dataset and inverse optimal control to learn task specific visual scanpaths[END_REF]. Recently, attempts have been devoted to incorporating gazes as weak supervision signals [START_REF] Fathi | Learning to Recognize Daily Actions using Gaze[END_REF][START_REF] Dim | Training object class detectors from eye tracking data[END_REF][START_REF] Samaras | Action Classification in Still Images Using Human Eye Movements[END_REF] for improving the performance of classification or segmentation systems. In [START_REF] Dim | Training object class detectors from eye tracking data[END_REF], objects detectors are trained from gaze features instead of accurate bounding boxes, showing promising results.

This paper introduces G(aze)-LSVM, a new weakly supervised learning model for image classification. G-LSVM generalizes latent SVM (LSVM) by exploiting human gazes for localizing objects. Fig. 1 illustrates the rationale of our model, where latent variables correspond to all possible regions in the image. The goal of LSVM is to select semantically meaningful regions, e.g. those containing the target object class (region Z in Fig. 1a). To improve the quality of the region selection, G-LSVM supports regions with high density of gazes (region Z in Fig. 1b) with respect to the region with the highest density of gazes (region Z i in Fig. 1b), by assuming that gaze features are related to regions relevant for the recognition task. Unlike [START_REF] Shcherbatyi | GazeDPM: Early Integration of Gaze Information in Deformable Part Models[END_REF][START_REF] Fathi | Learning to Recognize Daily Actions using Gaze[END_REF], G-LSVM only exploits gazes during training phase, and uses the pure visual information at test time without gazes. While [START_REF] Dim | Training object class detectors from eye tracking data[END_REF][START_REF] Sigal | Action is in the eye of the beholder: Eye-gaze driven model for spatio-temporal action localization[END_REF] focuses respectively on object detection and action recognition in the video, our targeted goal is to improve classification performances in the still images by using improved predicted latent regions.

In section 2, we formally define G-LSVM and its training procedure. Experiments conducted in section 3 show that G-LSVM significantly outperforms LSVM on both object dataset and action dataset of PASCAL VOC 2012.

GAZE-LSVM (G-LSVM) MODEL

We consider the problem of learning with weak supervision in a binary classification context. The interpretation is in the section 1.

Our prediction function f : X → Y takes as input an image x, and outputs a binary y ∈ {+1, -1}. Each image x is associated to latent variables z ∈ Z(x), which corresponds to sub-regions, as illustrated in Figure 1. For each region z in image x, we extract a feature vector Φ(x, z) ∈ R d . Our model is linear with respect to Φ, i.e. each region z is assigned the score w, Φ(x, z) . The problem is weakly supervised since the region-specific labels are unknown during training. Our prediction takes the maximum score over latent variables:

f w (x) = max z∈Z(x) w, Φ(x, z) (1)
Note that our prediction function in Eq. ( 1) is the same as in LSVM [START_REF] Felzenszwalb | Object Detection with Discriminatively Trained Part Based Model[END_REF], such that the label prediction at test time does not need gaze information.

G-LSVM Training

Our training scheme, however, penalizes the selection of latent regions based on gaze information. The general expres-sion of G-LSVM training objective is as follows:

L G (w) = 1 2 w 2 + n i=1 ∆ c (ŷ i , y i ) + γ • ∆ g (ẑ i , z i ) (2)
where y i is the true label of image x i , z i is the region with the maximum number of gazes, ŷi = sgn (f w (x i )) is the label predicted by our model, ẑi = arg max z∈Z(xi) w, Φ(x i , z) is the selected region, and 1 2 w 2 is the standard max margin regularization term. For each training example, Eq. ( 2) includes a classification loss ∆ c , and a gaze loss ∆ g , with a trade-off parameter γ. A standard classification metric is the 0/1 loss, which is, however, difficult to optimize. As in SVM, we use the hinge loss as upper-bound, so that ∆ c (ŷ i , y i ) = max(0, 1 -

y i f w (x i )).
The novelty in our training scheme is the introduction of a gaze loss. Its preliminary definition is δ g :

δ g (ẑ i , z i ) = 1 - g(x i , ẑi ) g(x i , z i ) (3)
where g(x i , z) is the number of gazes in the region z for image x i . Fig. 1b illustrates the proposed gaze loss, with blue circles representing gaze annotations, and z i is shown in red. In this example the region z contains 18 gazes out of 20 for z i , so that the gaze loss is 0.1, leading to a small penalization.

δ g ( ẑi , z i ) in Eq. ( 3) is difficult to optimize, because the dependency on w is complex and non-smooth. To overcome this issue, we derive a convex upper-bound ∆ g , inspired from margin-rescaling [START_REF] Joachims | Cutting-plane training of structural svms[END_REF]:

∆g(ẑi, zi)= max z∈Z(x i ) [δg(z, zi) + w, Φ(xi, z) ] -w, Φ(xi, zi) (4)
Our training objective in Eq. ( 2) is thus biased by the gaze loss ∆ g , so that G-LSVM learns different w parameters compared to LSVM.

Optimization

To minimize our training objective function, we first show that Eq. ( 2) can be rewritten as a difference of convex functions, i.e. u(w) -v(w), where v(w

) = C np ip=1 f w (x ip ), u(w) = 1 2 w 2 + C np ip=1 [γ∆ g (ẑ i , z i ) + max(1, f w (x ip ))] +C nn in=1 max(0, 1+f w (x in ))
, where n p (n n ) is the number of positive (negative) examples. Note that in the previous decomposition, the non-convex classification loss of every positive example is first decomposed into a difference of two convex functions: max(0,

1 -f w (x)) = max(1, f w (x)) -f w (x).
We then optimize u(w) -v(w) by CCCP (algo.1). The CCCP algorithm is guaranteed to decrease the objective function at every iteration and to converge to a local minimum or saddle point [START_REF] Yuille | The Concave-Convex Procedure (CCCP)[END_REF]. In algo 1 the line 3 involves linearizing the concave part -v(w). We calculate the supergradient v t of -v(w) at the point w t , where v t = -np ip=1 Φ(x i , ẑi ). At line 4, the problem becomes convex, we can use any convex optimization tool for solving this problem, e.g. SGD.

Algorithm 1: Concave-Convex Procedure

Output: w * 1 Set t = 0, stopping criterion and initialize w by w 0 2 repeat 3 Find hyperplane v t to linearize -v(w):

-v(w) ≤ -v(w t ) + (w -w t ) • v t , 4 Solve w t+1 = argmin w u(w) + w • v t , 5 Set t = t+1, 6 until [u(w t ) -v(w t )] -[u(w t-1 ) -v(w t-1 ))] < ;

EXPERIMENTAL RESULTS

Two datasets, PASCAL VOC 2012 action and object, are used for evaluation. For both datasets, [START_REF] Mathe | Action from still image dataset and inverse optimal control to learn task specific visual scanpaths[END_REF][START_REF] Dim | Training object class detectors from eye tracking data[END_REF] collected gaze annotations in task-driven manners.

Statistical consistency of gaze information

Before evaluating G-LSVM, we first provide a detailed analysis of the gaze data consistency. We compute statistics for the proportion of gazes falling into or outside the bounding boxes and compare it to the proportion of image pixels (Fig. 2). Statistically, for action dataset, 68.8% of the gazes fall into the ground-truth bounding-box, while the score of pixels is only 30.6%. Similarly, the scores of object dataset is 77.3% vs 36.9%. This preliminary study provides a quantitative validation that human gazes are highly related to object localization, and convey relevant features for classification.

Experimental results

Setup: In this paper, latent regions correspond to square image regions extracted with a multi-scale sliding window strategy. Region sizes vary from 90% to 30% of the whole image area, with a stride of 10%. Each region is described by the state-of-the-art deep features extracted from the pre-trained imagenet-vgg-m-2048 deep model 2 , which are subsequently L2-normalized.

2 http://www.vlfeat.org/matconvnet/pretrained/ Performance comparison: We compare G-LSVM to the baseline LSVM, with C = 10 4 for both models, and γ = 0.2 for G-LSVM. G-LSVM and LSVM are trained independently for each scale. We perform scale combination using a simple object bank representation [START_REF] Li | Object bank: A high-level image representation for scene classification & semantic feature sparsification[END_REF], leading to a 8-dimensional vector for each image. We also report performances of a SVM classifier trained on deep features computed on the whole image, denoted as wSVM.

The results are gathered in Table . 1, using 5 random folds on the train+val sets [START_REF] Mark Everingham | The Pascal Visual Object Classes Challenge a Retrospective[END_REF], and evaluating performances with the standard mAP metric. We show that G-LSVM outperforms LSVM by a margin of 2.1% for action (resp. 0.4% for object). Paired T-tests reveal that the improvement is statistically significant for a risk of less than 0.5% for action (resp. 2% for object). Both methods largely outperform wSVM, which clearly validate that training WSL models is able to capture local information.

G-LSVM

LSVM wSVM action 70.5 ± 0.8 68.4 ± 1.0 60.8 ± 1.2 object 92.4 ± 1.0 92.0 ± 1.1 88.2 ± 1.2 Table 1: mAP(%) of combination multi-scale model. Fig. 3 shows the performance evolution for LSVM and G-LSVM when varying the region scale s. We observe that the improvement of G-LSVM is more pronounced at small scales. This is expected: for large scales, all regions are informative, whereas at smaller scales, the model has to focus on relevant localized features. Note that s =100% corresponds to wSVM, for which the mAP for action and object is 60.8% and 88.2% . G-LSVM thus outperforms SVM at all scales of action dataset, as well for scales in [50, 90] on object dataset. . 2 gives per-class performances at the smallest scale 30%. G-LSVM outperforms LSVM by a margin of 3.1% and 0.8% for respectively action dataset and object dataset. Paired T-tests show that G-LSVM is significant than LSVM for a risk less than 0.5% and 1% for action and object datasets. For action dataset, the performance gain of G-LSVM is especially large for the categories phoning, reading, walking. We note that these actions are usually associated with tiny objects, i.e. cellphone, book and small person (e.g Fig. 5b). For object dataset, G-LSVM performs well at cow and motorbike and improves over LSVM for most categories.

Further analysis: The impact of the parameter γ in Eq. ( 2) is shown in Fig. 4 for scale 50%. We can see that we show that G-LSVM selects areas with more gaze features than LSVM. On the second row, we present results for test images, for which gaze features are not available. Interest-ingly, we can see that G-LSVM extracts regions which are more semantic than LSVM for the classification task. We validate this idea by measuring the detection performances of G-LSVM vs LSVM by computing the Intersection over Union (IoU) metric between the predicted region and the ground-truth bounding boxes. The results in Table . 3 at every scale show that G-LSVM always outperforms LSVM. Table 3: IoU (%) between predicted region and ground-truth bounding boxes.

Finally, we perform the last experiment using bounding box annotations during training, leading to a model denoted as G-LSVM*. We replace the gaze loss in Eq. ( 3) by a ground-truth loss computed as 1 -IoU (z, z gt ), where z gt is the ground-truth region in the dataset. The experiment reveal that G-LSVM is even slightly better than G-LSVM* (↑ 0.4 % (0.2%) mAP for the action (object) dataset). This shows that gaze features contain as relevant information as bounding box annotations, while being much cheaper to collect.

CONCLUSION

In this paper, we introduce G-LSVM, a new latent variable model which leverages human gaze features during training. We derive a concave-convex upper bound of the non-convex problem and solve it by the CCCP. When gaze annotations are scarce, an appealing feature of G-LSVM is that the model only uses the gazes for training, whereas only visual information is used for prediction. Experimental results show that G-LSVM significantly outperforms LSVM on classification and localization tasks, and that the model achieves similar performances as a model trained with expensive bounding box annotations. In the future, we plan to research gazes as timeseries for mining more information.

Fig. 1 :

 1 Fig. 1: Gazes bias the selection of latent regions for LSVM. The interpretation is in the section 1.

Fig. 2 :

 2 Fig. 2: Proportions of gazes and pixel numbers in (outside) the ground-truth bounding boxes.

Fig. 3 :

 3 Fig.3: mAP(%) at different scales. Table.2 gives per-class performances at the smallest scale 30%. G-LSVM outperforms LSVM by a margin of 3.1% and 0.8% for respectively action dataset and object dataset. Paired T-tests show that G-LSVM is significant than LSVM for a risk less than 0.5% and 1% for action and object datasets. For action dataset, the performance gain of G-LSVM is especially large for the categories phoning, reading, walking. We note that these actions are usually associated with tiny objects, i.e. cellphone, book and small person (e.g Fig.5b). For object dataset, G-LSVM performs well at cow and motorbike and improves over LSVM for most categories.Further analysis: The impact of the parameter γ in Eq. (2) is shown in Fig.4for scale 50%. We can see that

  performances of LSVM, corresponding to γ = 0, can be improved for most values in γ ∈ ]0, 1.0]. It is worth noticing that the performances in Fig.4are shown on average for all classes. We can further substantially boost the performances by cross-validating γ. For example, on the action dataset, a class-wise cross validation (γ ∈ [0, 1; 0.1]) at scale 50% leads to nearly 1% improvement compared to γ = 0.2.

Fig. 4 :

 4 Fig.4: For scale = 50%, the effect of parameter γ.

Fig. 5 :

 5 Fig. 5: Localization results. (a)(b): training results, (c)(d): test results. red: LSVM, blue: G-LSVM, yellow: groundtruth bounding-box. cyan: gazes.

  21.4 25.8 27.6 28.3 29.0 29.3 28.1 LSVM 14.5 20.4 24.3 26.7 27.9 28.9 28.

Table 2 :

 2 AP(%) at scale 30%

	Action Dataset	mAP	jump	phone instru' read	bike	horse	run	photo	comp'	walk
	G-LSVM	61.29	69.20	50.51 79.49 50.57 78.86 83.88 53.62 38.64	72.08	36.03
	LSVM-Standard 58.17	68.93	41.95 79.21 39.11 79.26 84.20 55.11 36.74	73.69	23.52
	Object Dataset	mAP aeroplane	cow	dog	cat	motor boat	horse	sofa	din'table bike
	G-LSVM	85.39	96.76	76.78 91.71 90.77 88.15 88.08 82.82 71.14	82.08	85.59
	LSVM	84.59	96.72	71.97 91.27 90.03 86.30 87.84 84.05 71.19	81.83	84.75


Introduction

One of the goals of GPR measurements is to determine the geometries of fine structures by imaging the shallow subsurface. In general, the GPR profiles are recorded on nearly flat surfaces and in this case, if highly dipping reflections and/or diffractions are present in data, a standard migration is needed in order to precisely determine the geometries of shallow structures [START_REF] Zeng | 3D imaging of a reservoir analogue in point bar deposits in the Ferron Sandstone, Utah, using ground-penetrating radar[END_REF]. For variable topographic relief, a standard processing procedure includes the application of static shifts [START_REF] Annan | Topographic Corrections of GPR data[END_REF] followed by a classical migration commonly performed with a flat datum plane. Nevertheless, this processing technique does not give good results for high topographic variation. In addition, the inadequacy of conventional elevation static corrections in accounting for gentle to rugged topographic relief was shown to be a particular problem [START_REF] Lehmann | Topographic migration of georadar data: Implications for acquisition and processing[END_REF]. To obtain reliable images from GPR data acquired on the areas showing irregular topography might require a special processing which accounts for the topography. [START_REF] Lehmann | Topographic migration of georadar data: Implications for acquisition and processing[END_REF] adapted a topographic migration for GPR data based on the Kirchhoff algorithm proposed by [START_REF] Wiggins | Kirchhoff integral extrapolation and migration of nonplanar data[END_REF] for the seismic data collected in mountainous areas. According to these authors, the topographic migration should be considered when the surface gradient exceeds 10%. This migration method has been successfully used, in 3D, by [START_REF] Mcclymont | Visualization of active faults using geometric attributes of 3D GPR data: An example from the Alpine Fault Zone[END_REF] for the GPR data acquired on active fault areas showing rugged topography.

In this paper, we first present an overview of the Kirchhoff topographic migration and demonstrate the diffraction equation used in this method as presented by [START_REF] Lehmann | True Topographic 2-D Migration of Georadar Data[END_REF]. To show the efficiency of the method we use synthetic data from a single diffraction point model, and compare the results of migration with flat datum and topography, respectively. Next, we present different examples of real GPR data recorded on areas with local and high topographic variations. The data are from the Chadian desert and Mongolia. Finally, we compare the result of GPR profiles processed with static shift and topographic migration, and discuss the supremacy of the later one.

The Kirchhoff topographic migration

We assume a zero offset survey with transmitting and receiving antennas which move on a flat horizontal surface at z = 0. The goal of the migration is to find the geological model (in the x-z plane) from the zero offset GPR profile (in the x-t plane). For a resistive medium, we can use the Kirchhoff method which gives the wave field at the location of the diffractor (x d , z d ) from the zero-offset wave field measured at the surface z = 0 [START_REF] Schneider | Integral formulation for migration in two and three dimensions[END_REF][START_REF] Leparoux | Adaptation of prestack migration to multi-offset groundpenetrating radar (GPR) data[END_REF]. Practically, the Kirchhoff migration will calculate the diffraction hyperbola (migration template) for each point of the GPR profile, and by adding the amplitudes along the template will place it at the top of the template, in the migrated profile [START_REF] Claerbout | Imaging the Earth's Interior[END_REF][START_REF] Yilmaz | Seismic Data Analysis[END_REF]. For the standard Kirchhoff migration, at a location x on the surface z = 0 (i.e. the antennas move on the flat datum plane), the two way travel time t(x) along the gray line path in Figure 1a is given by: ( )

2 2 2 0 4 ) ( V x x t x t d - + = , ( 1
)
where t 0 = 2z d /V is the two way travel time above the diffraction point (black point), x d is the horizontal position of the diffraction point, z d is the depth of diffraction point from the flat datum surface (see Figure 1a) and V is the EM wave velocity.

When GPR measurements are performed over a surface with topography, the migration template is not any more a diffraction hyperbola, instead it will be a distorted diffraction curve. This is shown by the thick line in Figure 1a where the topography is chosen to be a circle whose centre is on the diffraction point (x d , z d ), and on both sides of this circle the topography is flat. The distance between the diffraction point and the antennas (in zero offset), moving on the surface along the circle, is constant. Therefore the migration template, shown in Figure 1b, is flat on the top, and on both sides it is presented by two flanks of a diffraction hyperbola. In this case the imaging result of the classical Paris, France, 3-5 September 2012

Kirchhoff migration with a flat datum plane will be spurious (see Figure 1c). For that reason we need absolutely to take into account the topography of the GPR acquisition surface. Thus, we have to choose for migration template the diffraction curve shown in Figure 1b instead of using the diffraction hyperbola given by equation ( 1). This will allow the template to follow exactly the real travel path of the GPR data. Indeed, for the same x location (Figure 1a) the z position of the antennas (moving on the rugged surface) has been changed and the two way travel time t(x) is now calculated along the thick line path in Figure 1a to obtain:

( ) 2 2 2 4 ) ( V x x t x t d top - + = , with ( ) V z x z t t top 0 0 ) ( 2 - - = (2)
Where z(x) is the topography of the acquisition surface, and now t(x) is given by:

( ) ( ) ( ) V z x z t V z x z x x t x t d 0 0 2 2 0 2 2 0 ) ( 4 ) ( 4 ) ( - - - + - + = (3)
This equation is the same as the one shown without any demonstration by [START_REF] Lehmann | True Topographic 2-D Migration of Georadar Data[END_REF] in the case of 2D migration. In Figure 1d we present the result of the topographic migration appropriately weighted by an amplitude factor proportional to cos(α) = t top /t(x) which also depends on the topography z(x). The data are well imaged and, as expected, focused to a bright single point located at its real depth of 1.5 m.

Figure 1 a) The model of the diffraction point, with topography (thick line); b) Zero offset GPR data corresponding to a survey over this area; c) classical Kirchhoff migration with a flat surface at z = 0 (dashed line); d) The result of the topographic migration, the thick line corresponds to the real topography.

Near Surface Geoscience 2012 -18 th European Meeting of Environmental and Engineering Geophysics Paris, France, 3-5 September 2012

Real GPR data examples

The GPR profile, shown in Figure 2a, has been obtained over an Aeolian dry dune in the Chadian desert by using a 450 MHz shielded antenna. The acquisition mode was a constant offset of 0.25 m, the antennas were moved by steps of 0.125 m with a stack of 64 to improve the signal to noise ratio.

The following standard processing sequence was used: constant shift to adjust the time zero followed by NMO corrections; DC filter to remove the low frequency; flat reflections filter to remove continuous flat reflections noise; a band-pass filter and finally a time varying gain function was applied. The GPR profile of Figure 2a shows complex geometries, with imbricate reflections corresponding to different phases of deposits. The undulating reflection indicated by four white arrows in this figure represents the base of the dune, which in fact is flat and consists of pebbles (>2.0 mm diameter). In order to apply the topographic static shift and/or migration we need to know the velocity of GPR waves. In Figure 2a we also observe a nice diffraction hyperbola of 10 m wide situated just under the Near Surface Geoscience 2012 -18 th European Meeting of Environmental and Engineering Geophysics Paris, France, 3-5 September 2012 base of the dune (see black circle). After analyzing this diffraction, with different velocities, we have found that it can be fitted very well with a constant velocity model of 0.18 m/ns, which is in good agreement with values for dry sand. Figure 2b shows the same GPR profile as in Figure 2a, but the topographic corrections and the depth conversion of the time axis are performed by using a velocity of 0.18 m/ns. The topography shows a local variation (at profile coordinate 38 m) and its global variation of about 5 m (5%) is comparable to the investigation depth. The reflection from the base of the dune is roughly flattened while the diffraction hyperbola (at profile coordinate 25m) is still present and the dipping reflectors (see black arrow) are crossing down-dip leftwards the base of the dune instead of starting from the top of the base and going up-dip rightwards. The measured dips for the two reflectors indicated by the white arrows (from left to right) are 23.5° and 16°, respectively.

Figure 2c presents the topographic migration with the same velocity (0.18 m/ns) and a specific migration template width of 13 m, due to the number of traces per meter (8). The base of the dune is flattened and the diffraction at 25 m is now correctly focused on a single point inside the black circle, which justifies our choice of 0.18 m/ns for the GPR velocity. In addition, the dipping reflector shown by the black arrow has undergone a vertical and horizontal shift of 1.1 and 3.8 m. It starts at the base of the dune and goes up-dip rightwards as expected. The measured dips on the migrated section of the same reflectors shown by white arrows are higher than the ones found on the section with topographic corrections. Their values are now 26.5° and 19.5° on the migrated section, which corresponds to an error varying from 12% to 18%. Although the low value of topographic gradient (5%), the error of dips is very large due to the high EM velocity medium. This highlights the importance of the topographic migration even in the case where the surface gradient does not exceed 5%.

Conclusion

In the presence of relief variations of the same order as the investigation depth of GPR data, a topographic migration is necessary to correctly locate the dipping reflectors and correctly focus the diffractions. The topographic migration, presented in this paper, is based on the Kirchhoff algorithm similar to the method shown by [START_REF] Lehmann | Topographic migration of georadar data: Implications for acquisition and processing[END_REF]. We demonstrate the equation of migration template, as a function of the topography, along which the amplitudes are added together to give a single point on the migrated section. By comparing processed section obtained from GPR data measured over media of high EM velocity (dry sand) having a topographic gradient of 5%, we show that reflectors obtained by standard processing (static shift corrections) have dip angles that deviate from the angles in topographically migrated profile by 12% to 18%. Thus, for high velocity media, even in the case where the surface gradient does not exceed 5%, the application of the topographic migration is necessary and efficient.

Figure 2

 2 Figure 2 GPR profile acquired over a Chadian dry dune with a 450 MHz antennae. a) After the standard processing described in the text. b) The same as in a) but the topographic corrections and depth conversion of time axis are performed by using a velocity of 0.18m/ns. c) After Kirchhoff topographic migration with the same velocity of 0.18m/ns.


J Guillemoteau

P222 Influence of Grain Size Distribution and Shape on GPR Waves -Study of Aeolian Dunes

Ground penetrating radar (GPR), a geophysical method based on electromagnetic (EM) wave propagation, can provide very detailed and continuous images of the internal structures of aeolian dunes. In order to model and explain the origin of observed reflections, we build a model of electric permittivity which accounts for the grain size distribution and shape. By modelling the propagation of GPR waves in frequency domain, we have shown that grain size and grain shape have an influence on GPR reflectivity and we confirmed it by adjusting the modelled data to the real GPR data acquired in arid zones.

Introduction

Ground penetrating radar (GPR), a geophysical method based on electromagnetic (EM) wave propagation, can provide very detailed and continuous images of the internal structures of aeolian dunes [START_REF] Bristow | Internal structure of Aeolian sand dunes in Abu Dhabi determined using ground penetrating radar[END_REF][START_REF] Bano | Radar reflections and water content estimation of Aeolian sand dune[END_REF]. In these studies, change of water content is believed to be the origin of reflections within the aeolian dunes. However, the present day dunes are located on arid zones and the probability to observe change of the humidity in the first few meters of the dune is very low. On the other hand, some studies have shown that the properties which might change are the grain size distribution [START_REF] Watson | Grain size variations on a longitudinal dune and a barchan dune[END_REF][START_REF] Thomas | Analysis of linear dune sediment-form relationships in the Kalahari dune desert: Earth Surface Processes and Landforms[END_REF][START_REF] Lancaster | The Namib sand sea: Dunes forms, processes and sediments: Rotherdam[END_REF][START_REF] Dong | Grain size characteristics of dune sands in the central Taklimakan Sand Sea[END_REF] and/or the grain shape [START_REF] Sen | A self-similar model for sedimentary rocks with application to dielectric constant of fused glass beads[END_REF] in the case of a sand/sandstone contact. These two structural parameters are not able to cause a big change of permittivity. However, since dry earth materials have a small relative permittivity (between 2 and 6), the formula of the GPR wave reflectivity suggests that small changes of permittivity are able to cause a strong variation of the reflection coefficient.

First, we propose to explain these reflections by using a model of dry materials in which we can vary the grain size distribution. We developed a method to predict the permittivity of a mix of two grain sizes depending on their respective volumetric fraction. Afterwards, in order to model the reflectivity of GPR waves due to the change of grain shape, we simulated a sand/sandstone contact.

Finally we computed the propagation and the reflection of GPR waves for a two layers model varying the grain size distribution and/or the grain shape. The results of the forward modelling using reasonable geological predictions are in good accordance with real data acquired on present day dunes of the Chadian desert.

Influence of the grain size and grain shape on the electric permittivity of sands

The propagation and the reflectivity of GPR waves depend on the electric permittivity of the medium.

The permittivity of a mix ( n ) of several grains sizes (of permittivity  q , quartz in our case) can be computed by the recursive formula introduced by [START_REF] Robinson | Effect of particle size distribution on the effective dielectric permittivity of saturated granular media[END_REF]:

) ( 2 ) ( 3 1 1 1 1 1 1 1                   n q n n n n n q n q n n n n n n n f f f f            (1)
n f is the volumetric fraction of the grains of size "n" which has to validate the condition:

     1 1 N n n f (2)
Where N is the total number of sizes and  is the porosity. Generally, the data from granulometry are given in percentage of weight. If A n is the percentage of weight of the size "n", the volumetric fraction of the present grain size is given by:

   N n n n A A f 1 ) 1 (  (3)
This model depends on the porosity of the medium. However, one has to keep in mind that the porosity depends also on the grain size distribution. Thus, the two parameters and n f  are joined, and we suggest to compute, first, a realistic porosity associated to the grain size distribution. It is very difficult to predict the porosity for a mix of more than two sizes. We used the empirical "fractional packing model" [START_REF] Koltermann | Fractional packing model for hydraulic conductivity derived from sediment mixtures[END_REF] which is adapted to bimodal mixtures. Therefore, in order to predict the permittivity using the right porosity, we propose to combine the formula (1) and the fractional packing model as it is shown in figure 2. The variation of the permittivity due to the change of grain shape can be predicted by the formula [START_REF] Sen | A self-similar model for sedimentary rocks with application to dielectric constant of fused glass beads[END_REF]:

                        L q q * 0 0 * (4)
The parameter L is called depolarisation factor and depends on the grain shape [START_REF] Stratton | Electromagnetic theory[END_REF]. Its value varies between 0 and 1.

Table 1 Depolarising factor for different kind of grain shape.

Modelling of GPR reflections observed inside an aeolian dry sand dune

Figure 2 GPR section acquired over an aeolian dune with an antenna of 450 MHz. The two vertical black lines refer to the traces which have been studied.

A GPR profile (Fig. 2) acquired on a dry dune in Chadian desert shows many strong reflections in its interior. Owing to the dry climate of the area (arid zone) and also to the geometry of the reflections (which are not flat), one can think that these reflections have no chance to be caused by changes of the water content inside the dunes. We modelled the first reflection of the interior of the dune considering 72 nd EAGE Conference & Exhibition incorporating SPE EUROPEC 2010 Barcelona, Spain, 14 -17 June 2010 a change of grain size distribution. We considered a two layers model in which the permittivity has been computed using small grains for the first layer and a mix of small and large grains (with volumetric fraction of 20 % and 80 %, respectively) for the second layer. The porosity of the two layers is calculated by the fractional packing model described above. This gives a porosity of 0.4 for the first layer and 0.24 for the second one. The resulting permittivity for the two sands is between 2 and 3, which is in good agreement with values given in the literature (e.g. [START_REF] Davis | Ground-penetrating radar for high resolution mapping of soil and rock stratigraphy[END_REF].

The GPR response related to this model has been computed by using a 1D modelling method in frequency domain [START_REF] Bano | Radar reflections and water content estimation of Aeolian sand dune[END_REF][START_REF] Bano | Modelling of GPR waves for lossy media obeying a complex power law of frequency for dielectric permittivity[END_REF]. These synthetic data are computed for a large number of arrival time The same experience has been performed in order to interpret the GPR data acquired on a sand/sandstone contact. Figure 4 shows GPR data (without topographic corrections) acquired on a dry aeolian dune whose base (the bent reflection) is constituted of sandstone. We computed the reflectivity from the base of the dune (sand/sandstone contact) using the formula (4). We modelled the sand by a layer composed of spherical grains (L=0.33) with a porosity equal to 0.36, while the sandstone is considered as a medium composed of oblate grains (L=0.5) with a lower porosity equal to 0.1. The results of this modelling for the two soundings, shown by two red lines in figure 4, are presented in figure 5.

Conclusions

In this study, we tried to relate physics parameters (like grain size or shape) to the reflectivity of GPR waves. By using analytical models to compute the permittivity, we have shown that grain size and grain shape have an influence on GPR reflections and we confirmed it by adjusting the modelled data to the real GPR data acquired in arid zones. This kind of modelling is restricted to the case of dry media in which changes of grain size distribution or shape induces small change on the permittivity, but enough variations on the reflectivity. In addition, the presence of the water (that causes strong change on the permittivity) might hide the effects of the granular composition of the media. However, by this study, we show that one has to integrate these effects to model the GPR data in which water content variations and change of lithology are combined. Extending the application of the proposed method to radar data acquired on controlled dry materials in laboratory is the logical next step.

Figure 1

 1 Figure 1 Schematic presentation showing how to predict the permittivity from the data of granulometry. A n is the grain size distribution.

  (propagation time) and compared to the measured signals. The right propagation time is matched by minimizing the error between the synthetic and real data. The comparison between modelled and measured data is presented in figure 3. The real signals are from traces pointed out by two black lines on figure 2.

Figure 3 Figure 4

 34 Figure 3 Comparison between modelled (in red) and observed (in black) data for a contact between two sands with different grain size distribution. The real wavelets are located at x=70 m and x=43 m on the section of figure 2.

Figure 5

 5 Figure 5 Comparison between modelled (in red) and observed (in black) data for a sand/sandstone contact. The real wavelets are located at x=27 m and x=8.5 m on the section of figure 4.

nd EAGE Conference & Exhibition incorporating SPE EUROPEC 2010 Barcelona, Spain, 14 -17 June 2010