Bo ' Az Klartag
email: boaz.klartag@weizmann.ac.il

Joseph Lehec
email: lehec@ceremade.dauphine.fr

Poisson processes and a log-concave Bernstein theorem

We discuss interplays between log-concave functions and log-concave sequences. We prove a Bernstein-type theorem, which characterizes the Laplace transform of logconcave measures on the half-line in terms of log-concavity of the alternating Taylor coefficients. We establish concavity inequalities for sequences inspired by the Prékopa-Leindler and the Walkup theorems. One of our main tools is a stochastic variational formula for the Poisson average.

Introduction

Let ϕ : [0, ∞) → R be a continuous function that is C ∞ -smooth on (0, ∞). Its alternating Taylor coefficients are

a t (n) = (-1) n ϕ (n) (t) n! (n ≥ 0, t > 0). (1)
A function whose alternating Taylor coefficients are non-negative is called an absolutely monotone function. Bernstein's theorem asserts that the alternating Taylor coefficients are non-negative if and only if there exists a finite, non-negative Borel measure µ on [0, ∞) with

ϕ(t) = ∞ 0 e -tx dµ(x) (t ≥ 0). (2)
In other words, ϕ is the Laplace transform of the measure µ. See Widder [START_REF] Widder | The Laplace Transform[END_REF] for proofs of Bernstein's theorem. We say that the alternating Taylor coefficients are log-concave if the sequence (a t (n)) n≥0 is a log-concave sequence for any t > 0. This means that this sequence consists of non-negative numbers and for any m, n ≥ 0 and λ ∈ (0, 1) such that λn+(1-λ)m is an integer, a t (λn + (1 -λ)m) ≥ a t (n) λ a t (m) 1-λ .

Equivalently, a t (n) 2 ≥ a t (n -1)a t (n + 1) for every n ≥ 1 and the set of integers n for which a t (n) > 0 is an interval of integers.

A measure µ on [0, ∞) is log-concave if it is either a delta measure at a certain point, or else an absolutely-continuous measure whose density f : [0, ∞) → R is a log-concave function. Recall that a function f : K → R for some convex set K is log-concave if f is non-negative and f (λx + (1 -λ)y) ≥ f (x) λ f (y) 1-λ for all x, y ∈ K, 0 < λ < 1.

Theorem 1.1 ("Log-concave Bernstein theorem"). Let ϕ : [0, ∞) → R be a continuous function that is C ∞ -smooth on (0, ∞). Then the alternating Taylor coefficients of ϕ are log-concave if and only if ϕ takes the form (2) for a certain finite, log-concave measure µ.

In Section 2 we prove Theorem 1.1 by using an inversion formula for the Laplace transform as well as the Berwald-Borell inequality [START_REF] Berwald | Verallgemeinerung eines Mittelwertsatzes von J. Favard für positive konkave Funktionen[END_REF][START_REF] Borell | Complements of Lyapunov's inequality[END_REF]. The latter inequality states that if one divides the Mellin transform of a log-concave measure on [0, ∞) by the Gamma function, then a log-concave function is obtained. It directly implies the "if" part of Theorem 1.1. Our theorem admits the following corollary: Corollary 1.2. Let µ be a finite, non-negative Borel measure on [0, ∞) and let ϕ be given by [START_REF] Borell | Complements of Lyapunov's inequality[END_REF]. Then µ is log-concave if and only if the function ϕ (n-1) (t)

-1/n is convex in t ∈ (0, ∞)

for every n ≥ 1.

In fact, in Theorem 1.1 it suffices to verify that the sequence (a t (n)) n≥0 is log-concave for a sufficiently large t, as follows from the following: Proposition 1.3. Let ϕ : (0, ∞) → R be real-analytic, and define a t (n) via [START_REF] Berwald | Verallgemeinerung eines Mittelwertsatzes von J. Favard für positive konkave Funktionen[END_REF]. Assume that 0 < r < s and that the sequence (a s (n)) n≥0 is log-concave. Then also the sequence (a r (n)) n≥0 is also log-concave.

Proposition 1.3 is proven in Section 3, alongside concavity inequalities related to logconcave sequences in spirit of the Walkup theorem [START_REF] Walkup | Pólya sequences, binomial convolution and the union of random sets[END_REF]. While searching for a Prékopa-Leindler type inequality for sequences, we have found the following:

Theorem 1.4. Let f, g, h, k : Z → R satisfy f (x) + g(y) ≤ h x + y 2 + k x + y 2 , ∀x, y ∈ Z,
where x is the lower integer part of x ∈ R and x is the upper integer part. Then

x∈Z e f (x) x∈Z e g(x) ≤ x∈Z e h(x) x∈Z e k(x) .
Our proof of Theorem 1.4 which is presented in Section 5 involves probabilistic techniques, and it would be interesting to find a direct proof. However, we believe that the probabilistic method is not without importance in itself, and perhaps it is deeper than other components of this paper. The argument is based on a stochastic variational formula for the expectation of a given function with respect to the Poisson distribution. It is analogous to Borell's formula from [START_REF] Borell | Diffusion equations and geometric inequalities[END_REF] which is concerned with the Gaussian distribution. The stochastic variational formula is discussed in Section 4.

The Berwald-Borell inequality (or Theorem 1.1) imply that when µ is a finite, log-concave measure on [0, ∞) and k ≤ ≤ m ≤ n are non-negative integers with k

+ n = + m, a t ( )a t (m) -a t (k)a t (n) ≥ 0, (3)
where a t (n) = ∞ 0 (x n /n!)e -tx dµ(x) is defined via (1) and ( 2). The following theorem shows that the left-hand side of (3) is not only non-negative, but it is in fact an absolutely-monotone function of t: Theorem 1.5. Let µ be a finite, log-concave measure on [0, ∞). Then for any non-negative integers k ≤ ≤ m ≤ n with k + n = + m there exists a finite, non-negative measure

ν = ν k, ,m,n on [0, ∞), such that for any t > 0, ∞ 0 e -tx dν(x) = a t ( )a t (m) -a t (k)a t (n), (4)
where as usual a t (n) = ∞ 0 (x n /n!)e -tx dµ(x) is defined via (1) and ( 2). Theorem 1.5 is proven in Section 3. Let us apply this theorem in a few examples. In the case where µ is an exponential measure, whose density is e -αt on [0, ∞), the measures ν from Theorem 1.5 vanish completely. In the case where µ is proportional to a Gamma distribution, also the measures ν are proportional to Gamma distributions. When µ is the uniform measure on the interval [START_REF] Berwald | Verallgemeinerung eines Mittelwertsatzes von J. Favard für positive konkave Funktionen[END_REF][START_REF] Borell | Complements of Lyapunov's inequality[END_REF], the density of the measure ν = ν 0,1,1,2 from Theorem 1.5 is depicted in Figure 1. This log-concave density equals the convex function (t -1)(t -2)/2 in the interval [START_REF] Borell | Complements of Lyapunov's inequality[END_REF][START_REF] Borell | Convex set functions in d-space[END_REF], and it equals (t -2)(4 -t) in [START_REF] Borell | Convex set functions in d-space[END_REF][START_REF] Borell | Diffusion equations and geometric inequalities[END_REF]. We suggest to refer to the measure ν from Corollary 1.5 as the Berwald-Borell transform of µ with parameters (k, , m, n). All Berwald-Borell transforms of log-concave measures that we have encountered so far were log-concave by themselves. It is a curious problem to characterize the family of measures ν which could arise as the Berwald-Borell transform of a log-concave measure on [0, ∞). Such a characterization could lead to new constraints on the moments of log-concave measures on [0, ∞) beyond the constraints given by the Berwald-Borell inequality.

Proof of the log-concave Bernstein theorem

The proof of Theorem 1.1 combines ideas of Berwald from the 1940s with the earlier Post inversion formula for the Laplace transform. The "if" direction of Theorem 1.1 follows from: Lemma 2.1. Let µ be a finite, log-concave measure on [0, ∞). Assume that ϕ is given by [START_REF] Borell | Complements of Lyapunov's inequality[END_REF]. Then the alternating Taylor coefficients of ϕ are log-concave.

Proof. In the case where µ = c δ x 0 we have that ϕ(t) = c e -tx 0 and hence

a t (n) = c e -tx 0 x n 0 n! .
Since a t (n) > 0 for every n and a t (n + 1)/a t (n) = x 0 /(n + 1) is non-increasing, this is indeed a log-concave sequence. In the case where µ has a log-concave density f , we denote f t (x) = e -tx f (x) and observe that

ϕ (k) (t) = ∞ 0 (-x) k e -tx f (x)dx = ∞ 0 (-x) k f t (x)dx (k ≥ 0, t > 0). (5)
The function f t is log-concave, and hence we may apply the Berwald-Borell inequality [START_REF] Berwald | Verallgemeinerung eines Mittelwertsatzes von J. Favard für positive konkave Funktionen[END_REF][START_REF] Borell | Complements of Lyapunov's inequality[END_REF], see also Theorem 2.2.5 in [START_REF] Brazitikos | Geometry of isotropic convex bodies[END_REF] or Theorem 5 in [START_REF] Nayar | Khinchine type inequalities with optimal constants via ultra log-concavity[END_REF] for different proofs. This inequality states that the sequence

k → ∞ 0 x k f t (x)dx k! (k ≥ 0) (6)
is log-concave, completing the proof.

We now turn to the proof of the "only if" direction of Theorem 1.1, which relies on the Post inversion formula for the Laplace transform, see Feller [START_REF] Feller | An introduction to probability theory and its applications[END_REF]Section VII.6] or Widder [START_REF] Widder | The Laplace Transform[END_REF]Section VII.1]. Suppose that ϕ is continuous on [0, ∞) and C ∞ -smooth on (0, ∞), and that the alternating Taylor coefficients a t (n) are log-concave. In particular, the alternating Taylor coefficients are non-negative. We use Bernstein's theorem to conclude that there exists a finite, non-negative Borel measure µ on [0, ∞) such that (2) holds true. All that remains is to prove the following: Proposition 2.2. The measure µ is log-concave.

The proof of Proposition 2.2 requires some preparation. First, it follows from (1) and ( 2) that for any R, t > 0,

Rt n=0 t n a t (n) = Rt n=0 t n n! ∞ 0 x n e -tx dµ(x) = ∞ 0 P (N tx ≤ Rt) dµ(x), (7)
where N s is a Poisson random variable with parameter s, i.e.,

P(N s = n) = e -s s n n! , for n = 0, 1, 2, . . .
The random variable N s has expectation s and standard-deviation √ s. From the Markov-Chebyshev inequality, for any α > 0,

P(N s ≤ αs) s→∞ -→    1 α > 1 1/2 α = 1 0 α < 1 (8)
where the case α = 1 is slightly more difficult, it is left as an exercise in Feller [START_REF] Feller | An introduction to probability theory and its applications[END_REF]ter VII], that may be solved via the central limit theorem for the Poisson distribution. The left-hand side of ( 8) is always between zero and one. Therefore we may use the bounded convergence theorem, and conclude from (7) that for any R > 0,

lim t→∞ Rt n=0 t n a t (n) = µ([0, R)) + 1 2 • µ({R}). (9)
For t > 0 define g t : [0, ∞) → R via

g t (x) = t n+1 • a t (n) x = n/t for some integer n ≥ 0 t x+1 • a t (n) 1-λ • a t (n + 1) λ x = (n + λ)/t for λ ∈ (0, 1), n ≥ 0.
Write µ t for the measure on [0, ∞) whose density is g t . We think about µ t as an approximation for the discrete measure on [0, ∞) that has an atom at n/t of weight t n a t (n) for any n ≥ 0.

Lemma 2.3. For any t > 0 the measure µ t is log-concave on [0, ∞). Moreover, if µ({0}) = 0 then for any R > 0,

µ t ([0, R)) t→∞ -→ µ([0, R)) + 1 2 • µ({R}).
Proof. We may assume that µ ≡ 0, as otherwise the lemma is trivial. Since µ({0}) = 0, for any t > 0 and n ≥ 0,

a t (n) = ∞ 0 x n n! e -tx dµ(x) > 0.
The density g t is locally-Lipschitz, and for any integer n ≥ 0 and x ∈ (n/t, (n + 1)/t),

(log g t ) (x) = log t + t log a t (n + 1) a t (n) .
The sequence

(a t (n)) n≥0 is log-concave, hence a t (n + 1)/a t (n) is non-increasing in n.
We conclude that (log g t ) (x), which exists for almost any x > 0, is a non-increasing function of x ∈ [0, ∞). This shows that µ t is a log-concave measure. In particular, the density g t is unimodular, meaning that for some x 0 ≥ 0, the function g t is non-decreasing in (0, x 0 ) and non-increasing in (x 0 , ∞). We claim that for any R, t > 0 we have the Euler-Maclaurin type bound:

R 0 g t (x)dx - Rt n=0 1 t • g t n t ≤ 3 t • sup x>0 g t (x). (10)
Indeed, the sum in ( 10) is a Riemann sum corresponding to the integral of g t on the interval I = [0, tR + 1 /t]. This Riemann sum corresponds to a partition of I into segments of length 1/t, and by unimodularity, this Riemann sum can deviate from the actual integral by at most 2/t • sup x>0 g t (x). Since the symmetric difference between I and [0, R] is an interval of length at most 1/t, the relation [START_REF] Picard | Formules de dualité sur lespace de Poisson[END_REF] follows. According to [START_REF] Picard | Formules de dualité sur lespace de Poisson[END_REF], for any R, t > 0,

µ t ([0, R)) - Rt n=0 t n a t (n) ≤ 3 t • sup x>0 g t (x) = 3 • sup n≥0 t n a t (n). (11)
Next we use our assumption that µ({0}) = 0 and also the fact that sup n e -s s n /n! tends to zero as s → ∞, as may be verified routinely. This shows that for t > 0,

sup n≥0 t n a t (n) = sup n≥0 ∞ 0 (tx) n n! e -tx dµ(x) ≤ ∞ 0 sup n≥0 (tx) n n! e -tx dµ(x) t→∞ -→ 0,
where we used the dominated convergence theorem in the last passage. The lemma now follows from ( 9) and [START_REF] Walkup | Pólya sequences, binomial convolution and the union of random sets[END_REF].

The following lemma is due to Borell, and its short proof is contained in [3, Lemma 3.3] and the last paragraph of the proof of [3, Theorem 2.1]. For A, B ⊆ R and λ ∈ R we write

A + B = {x + y ; x ∈ A, y ∈ B} and λA = {λx ; x ∈ A}.
Lemma 2.4. Let µ be a finite Borel measure on R such that for any intervals I, J ⊆ R and

0 < λ < 1, µ (λI + (1 -λ)J) ≥ µ(I) λ µ(J) 1-λ . ( 12
)
Then µ is log-concave (i.e., either µ = cδ x 0 for some c ≥ 0, x 0 ∈ R or else µ has a logconcave density).

Proof of Proposition 2.2. We may assume that µ((0, ∞)) > 0 as otherwise µ = cδ 0 and the conclusion trivially holds. Therefore a t (n) > 0 for all t and n. By log-concavity of the sequence of alternating Taylor coefficients,

∞ 0 xe -tx dµ(x) 2 ∞ 0 (x 2 /2)e -tx dµ(x) = a t (1) 2 a t (2) ≥ a t (0) = ∞ 0 e -tx dµ(x) t→∞ -→ µ({0}). (13)
Write ν t for the measure on (0, ∞) whose density with respect to µ equals x → e -tx . Then by the Cauchy-Schwartz inequality,

∞ 0 xe -tx dµ(x) 2 ∞ 0 x 2 e -tx dµ(x) = ∞ 0 xdν t (x) 2 ∞ 0 x 2 dν t (x) ≤ ν t ((0, ∞)) t→∞ -→ 0. (14)
From ( 13) and ( 14) we see that µ({0}) = 0, which is required for the application of the second part of Lemma 2.3. Let I, J ⊆ [0, ∞) be intervals and 0 < λ < 1. Thanks to Lemma 2.4, it suffices to prove the inequality [START_REF] Widder | The Laplace Transform[END_REF]. Since µ is a finite measure, it suffices to prove that µ(U ) ≥ µ(I 1 ) λ µ(J 1 ) 1-λ (15

)
where U is an arbitrary open interval containing λI + (1 -λ)J, where I 1 is a compact interval contained in the interior of I and J 1 is a compact interval contained in the interior of J. However, by Lemma 2.3,

µ(U ) ≥ lim sup t→∞ µ t (λI + (1 -λ)J) , µ(I 1 ) ≤ lim inf t→∞ µ t (I), µ(J 1 ) ≤ lim inf t→∞ µ t (J).
(16) Since µ t is log-concave, the Prékopa-Leindler inequality (see, e.g., [6, Theorem 1.2.3]) implies that for all t > 0,

µ t (λI + (1 -λ)J) ≥ µ t (I) λ µ t (J) 1-λ . ( 17
)
Now (15) follows from ( 16) and ( 17), and the proposition is proven.

The proof of Theorem 1.1 is complete.

Proof of Corollary 1.2. If µ is of the form cδ 0 , then the corollary is trivial. Otherwise, the alternating Taylor coefficients a t (n) = (-1) n ϕ (n) (t)/n! are positive for every t > 0 and n ≥ 0. The measure µ is log-concave if and only if the sequence of alternating Taylor coefficients is log-concave for any t > 0, which happens if and only if

ϕ (n) (t) n! 2 - ϕ (n-1) (t) (n -1)! • ϕ (n+1) (t) (n + 1)! ≥ 0 (n ≥ 1, t > 0). (18)
Denote by b n (t) the expression on the left-hand side of (18) multiplied by (n!) 2 . Then,

d 2 dt 2 ϕ (n-1) (t) -1/n = d 2 dt 2 (-1) n-1 • ϕ (n-1) (t) -1/n = n + 1 n 2 ϕ (n-1) (t) -(2n+1)/n • b n (t).
Hence b n (t) ≥ 0 for all t > 0 if and only if the function ϕ (n-1) (t) -1/n is convex in (0, ∞).

The log-concavity measurements are absolutely monotone

In this section we prove Proposition 1.3 and Theorem 1.5. Let k ≤ ≤ m ≤ n be nonnegative integers with k + n = + m and let ϕ be a continuous function on

[0, ∞) that is C ∞ -smooth in (0, ∞). Define c k, ,m,n (t) = a t ( )a t (m) -a t (k)a t (n) = (-1) +m ϕ ( ) (t) ! ϕ (m) (t) m! - ϕ (k) (t) k! ϕ (n) (t) n! ,
where a t (n) = (-1) n ϕ (n) (t)/n! as before. Clearly, the sequence (a t (n)) n≥0 is log-concave if and only if c k, ,m,n (t) ≥ 0 for all non-negative integers k ≤ ≤ m ≤ n with k + n = + m. We call the functions c k, ,m,n : (0, ∞) → R the log-concavity measurements of ϕ.

Lemma 3.1. The derivative of each log-concavity measurement is a linear combination with non-positive coefficients of a finite number of log-concavity measurements.

Proof. Differentiating (1) we obtain

d dt a t (n) = -(n + 1)a t (n + 1) (t > 0, n ≥ 0). Abbreviate b j = a t (j). Then, -c k, ,m,n (t) = ( + 1)b +1 b m + (m + 1)b b m+1 -(k + 1)b k+1 b n -(n + 1)b k b n+1 . ( 19
)
Assume first that < m. In this case we may rewrite the right-hand side of (19) as

(k + 1) [b +1 b m -b k+1 b n ] + ( -k) [b +1 b m -b k b n+1 ] + (m + 1) [b b m+1 -b k b n+1 ] .
Therefore, in the case < m, we expressed -c k, ,m,n (t) as a linear combination with nonnegative coefficients of three log-concavity measurements. From now on, we consider the case = m. If k = , then necessarily n = m and the log-concavity measurement c k, ,m,n (t) vanishes. If k < , then necessarily m < n and we rewrite the right-hand side of (19) as

(k + 1) [b b m+1 -b k+1 b n ] + (n + 1) [b b m+1 -b k b n+1 ] .
Consequently, in the case = m, we may express -c k, ,m,n (t) as a linear combination with non-negative coefficients of two log-concavity measurements. The proof is complete.

Corollary 3.2. Let t > 0 be such that (a t (n)) n≥0 is a log-concave sequence. Assume that k ≤ ≤ m ≤ n are non-negative integers with k + n = + m. Abbreviate f (t) = c k, ,m,n (t).
Then for all j ≥ 0, (-1) j f (j) (t) ≥ 0.

Proof. Any log-concavity measurement is non-negative at any t > 0. It follows from Lemma 3.1 that (-1) j f (j) (t) is a finite linear combination with non-negative coefficients of certain log-concavity measurements. Therefore (-1) j f (j) (t) ≥ 0.

Proof of Proposition 1.3. Write A ⊆ (0, ∞) for the set of all t > 0 for which

(a t (n)) n≥0 is a log-concave sequence. Since ϕ is C ∞ -smooth, the set A is closed in (0, ∞). From our assumption, r ∈ A. Define t 0 = inf {t > 0 ; [t, r] ⊆ A} .
Then t 0 ≤ r. Our goal is to prove that t 0 = 0. Assume by contradiction that t 0 > 0. Since A is a closed set, necessarily t 0 ∈ A. Since ϕ is real-analytic, the Taylor series of ϕ converges to ϕ in (t 0 -ε, t 0 + ε) for a certain ε > 0. Assume that k ≤ ≤ m ≤ n are non-negative integers with k + n = + m. Then also the Taylor series of f (t) = c k, ,m,n (t) converges to f in the same interval (t 0 -ε, t 0 + ε). From Corollary 3.2 we thus deduce that for all t ∈ (t 0 -ε, t 0 ], c k, ,m,n (t) ≥ 0.

Consequently, (t 0 -ε, t 0 ] ⊆ A, in contradiction to the definition of t 0 .

We proceed with yet another proof of Proposition 1.3, which is more in spirit of the Walkup theorem which we now recall: [START_REF] Nayar | Khinchine type inequalities with optimal constants via ultra log-concavity[END_REF][START_REF] Walkup | Pólya sequences, binomial convolution and the union of random sets[END_REF]). If (a n ) n≥0 and (b n ) n≥0 are log-concave sequences then the sequence (c n ) n≥0 given by

Theorem 3.3 (Walkup theorem
c n = n k=0 n k a k b n-k , (n ≥ 0)
is also log-concave.

By Taylor's theorem whenever 0 < s < t,

(t -s) k a s (k) = ∞ n=k n k (t -s) n a t (n),
assuming that ϕ is real-analytic and that the Taylor series of ϕ at t converges in (s -ε, t) for some ε > 0. Therefore Proposition 1.3 is equivalent to the following Walkup-type result:

Proposition 3.4. If (a k ) k≥0 is a log-concave sequence then the sequence (c k ) k≥0 defined by c k = ∞ n=k n k a n , (k ≥ 0)
is log-concave as well.

We do not know of a formal derivation of Theorem 3.3 from Proposition 3.4 or vice versa, yet we provide a direct proof of Proposition 3.4 which shares some similarities with the proof of Walkup's theorem and the Borell-Berwald inequality given in [START_REF] Nayar | Khinchine type inequalities with optimal constants via ultra log-concavity[END_REF][START_REF] Walkup | Pólya sequences, binomial convolution and the union of random sets[END_REF]. For integers a ≤ b we write a, b = {n ∈ Z ; a ≤ n ≤ b}. We begin the direct proof of Proposition 3.4 with the following: Lemma 3.5. Let (a n ) n≥0 be a log-concave sequence. Then for every non-negative integers k and l we have

n≥0 n k l -n k a n a l-n ≥ n≥0 n k -1 l -n k + 1 a n a l-n , (20)
where here we set n k = 0 in case where k > n or k < 0 or n < 0. Proof. Inequality (20) holds trivially if 2k > . We may thus assume that 2k ≤ . Let U be a random subset of cardinality 2k + 1 of {1, . . . , l + 1} chosen uniformly. Let X 1 , . . . , X 2k+1 be the elements of U in increasing order. Observe that the law of X k+1 is given by

P(X k+1 = n + 1) = n k l-n k l+1 2k+1 , (n ≥ 0). Therefore n n k l -n k a n a l-n = l + 1 2k + 1 E[f (X k+1 )]
where f is the function given by

f (n) = a n-1 a l+1-n , ∀n,
where we set a k = 0 for k < 0. In a similar way

n n k -1 l -n k + 1 a n a l-n = l + 1 2k + 1 E[f (X k )].
Hence the desired inequality boils down to

E[f (X k )] ≤ E[f (X k+1 )].
By Fubini it suffices to prove that

P(f (X k ) > t) ≤ P(f (X k+1 ) > t) ∀t ≥ 0. The function f satisfies f (l + 2 -n) = f (n) for all n.
The crucial observation is that because of the log-concavity of the sequence (a n ) the farther n from the midpoint (l + 2)/2 the smaller f (n). Therefore the level set {f > t} is either empty or else an interval of the form n, l + 2 -n for some integer n ≤ (l + 2)/2. Hence it is enough to prove that for every such n,

P(X k ∈ n, l + 2 -n ) ≤ P(X k+1 ∈ n, l + 2 -n ). (21)
Intuitively, since X k+1 is the middle element of U , it is more likely to be close to the center of the interval 1, l + 1 than any other element. More precisely, since X k ≤ X k+1 ,

P(X k ∈ n, l + 2 -n ) -P(X k+1 ∈ n, l + 2 -n ) = P(X k ≤ l + 2 -n; X k+1 > l + 2 -n) -P(X k < n; X k+1 ≥ n) = l+2-n k n-1 k+1 l+1 2k+1 - n-1 k l+2-n k+1 l+1 2k+1
.

In order to complete the proof we need to show that this expression is non-positive, assuming that k ≤ /2 and n ≤ ( + 2)/2. Note that

l+2-n k n-1 k+1 n-1 k l+2-n k+1 = (n -1 -k)!( -n -k + 1)! ( + 2 -n -k)!(n -k -2)! = (n -1) -k ( + 2 -n) -k . ( 22
)
We need to show that the expression in ( 22) is at most one. The denominator in ( 22) is positive, as

+ 2 -n -k = 1 + [( + 2)/2 -n] + ( /2 -k) ≥ 1.
The numerator in ( 22) is smaller than the denominator, as n -1 ≤ + 2 -n. Hence the expression in ( 22) is at most one, completing the proof of the lemma.

Direct proof of Proposition 3.4: Let k ≥ 0 be an integer. We need to prove that

n,m n k m k a n a m ≥ n,m n k -1 m k + 1 a n a m .
By grouping the terms according to the value of n + m we see that it suffices to prove that for any , k ≥ 0,

n n k l -n k a n a l-n ≥ n n k -1 l -n k + 1 a n a l-n .
This is, however, exactly the statement of the previous lemma.

When µ is a finite, log-concave measure on [0, ∞), it is well-known (e.g., [START_REF] Brazitikos | Geometry of isotropic convex bodies[END_REF]) that µ([t, ∞)) ≤ Ae -Bt for all t > 0, where A, B > 0 depend only on µ. It follows that the Laplace transform ϕ defined in ( 2) is analytic for t ∈ C with Re(t) > -C for some C > 0 depending on µ.

Proof of Theorem 1.5. By Theorem 1.1, the alternating Taylor coefficients sequence (a t (n)) n≥0 is log-concave for any t > 0. From Corollary 3.2 we thus learn that

f (t) = c k, ,m,n (t) = a t ( )a t (m) -a t (k)a t (n)
satisfies (-1) j f (j) (t) ≥ 0 for any t > 0 and j ≥ 0. The function f is real-analytic in a neighborhood of [0, ∞) and in particular it is continuous in [0, ∞). The function f is thus absolutely-monotone, and according to Bernstein theorem, there exists a finite, non-negative measure ν for which (4) holds true.

We may rewrite conclusion (4) of Theorem 1.5 as follows: For any t > 0,

∞ 0 x ! e -tx dµ(x) ∞ 0 x m m! e -tx dµ(x)- ∞ 0 x k k! e -tx dµ(x) ∞ 0 x n n! e -tx dµ(x) = ∞ 0 e -tx dν(x).
Let us now consider the Fourier transform

F µ (t) = ∞ 0 e -itx dµ(x) (t ∈ R).
By analytic continuation, Theorem 1.5 immediately implies the following: Proposition 3.6. Let µ be a finite, log-concave measure on [0, ∞). Then for any non-negative integers k ≤ ≤ m ≤ n with k + n = + m there exists a finite, non-negative measure ν = ν k, ,m,n on [0, ∞), such that for any t > 0,

F ( ) µ (t) ! F (m) µ (t) m! - F (k) µ (t) k! F (n) µ (t) n! = (-i) +m F ν (t).
Corollary 3.7. Let µ, k, , m, n, ν be as in Theorem 1.5. Write P j (µ) for the measure whose density with respect to µ is x → x j /j!. Write E t (µ) for the measure whose density with respect to µ is x → exp(-tx). Then, (i) We have ν = P (µ) * P m (µ) -P k (µ) * P n (µ) where * stands for convolution.

(ii) For any t > 0, the measure E t (ν) is the Berwald-Borell transform of E t (µ) with the same parameters (k, , m, n). The same holds for any t ∈ R for which E t (µ) is a finite measure.

Proof. We have that F (j) µ /j! = (-i) j • F P j (µ) . Proposition 3.6 thus shows that

F P (µ) F Pm(µ) -F P k (µ) F Pn(µ) = F ν . ( 23
)
The Fourier transform maps products to convolutions. Conclusion (i) therefore follows from (23). Conclusion (ii) follows immediately from the definitions.

Borell type formula for the Poisson measure

In [START_REF] Borell | Diffusion equations and geometric inequalities[END_REF], Borell gave a new proof of the Prékopa-Leindler inequality based on the following stochastic variational formula. Let γ n be the standard Gaussian measure on R n . Given a standard n-dimensional Brownian motion B and a bounded function f : R n → R we have log

R n e f dγ n = sup u E f B 1 + 1 0 u s ds - 1 2 1 0 |u s | 2 ds , ( 24
)
where the supremum is taken over all bounded stochastic processes u which are adapted to the Brownian filtration, i.e. u t is measurable with respect to the σ-field generated by {B s ; s ≤ t} for all t ∈ [0, 1].

In this section we give a discrete version of Borell's formula in which the Gaussian measure and the Brownian motion are replaced by the Poisson distribution and the Poisson process, respectively. In the following section we shall apply our formula in order to deduce a discrete version of the Prékopa-Leindler inequality. We start with some background on counting processes with stochastic intensities. Let T > 0 be a fixed number, denote R + = [0, ∞), and let (Ω, F, P) be a probability space on which our random variables will be defined.

Throughout this section, we let N be a Poisson point process on [0, T ] × R + ⊆ R 2 with intensity measure equal to the Lebesgue measure L. In particular N (F ) is a Poisson random variable with parameter L(F ) for any Borel set F ⊆ [0, T ] × R + . For a Borel subset E ⊆ [0, T ] × R + we write F E for the σ-field generated by the random variables

{N (F ); F is a Borel set, F ⊆ E} .
For t ∈ [0, T ] we set F t = F [0,t]×R + . This defines a filtration of Ω. Our next goal is to define, given a suitable stochastic process (λ t ) 0≤t≤T , a counting process (X λ t ) 0≤t≤T which satisfies

X λ t = N ({(s, u) ∈ [0, T ] × R + ; s ≤ t, u ≤ λ s }). (25)
In other words X λ t is the number of atoms of N which lie below the curve {(s, λ s ) : s ∈ [0, t]}. Let us now explain the technical assumptions regarding the stochastic intensity process (λ t ) t≤T . Recall that a process is called predictable if, as a function of t ∈ [0, T ] and ω ∈ Ω, it is measurable with respect to the σ-field P generated by the sets

{ (s, t] × A ; s ≤ t ≤ T, A ∈ F s }.
This is slightly more restrictive than being adapted, i.e., when λ t is measurable with respect to F t . We have the following standard fact: if a process is left-continuous and adapted, then it is predictable. From now on we will always assume that the stochastic intensity process λ is non-negative, bounded and predictable. The counting process X λ defined via (25) is clearly adapted and non-decreasing. Note that given M > 0, with probability 1 the process N has only finitely many atoms in the box [0, T ] × [0, M ] and no two of those lie on the same vertical line {t} × [0, M ]. Thus, with probability 1 the process X λ has finitely many jumps, all of size 1, and it is right-continuous. We sometimes refer to the jumps of X λ as atoms.

Lemma 4.1. For every non-negative predictable process (H t ) 0≤t≤T we have

E T 0 H t X λ (dt) = E T 0 H t λ t dt , (26)
where the integral on the left-hand side is a Riemann-Stieljes integral, i.e., here it is a sum of the values of H t at the atoms of X λ .

The proof of the technical Lemma 4.1 is deferred to the appendix. Equation ( 26) is frequently taken as the definition of a counting process with stochastic intensity λ. The process

X λ t = X λ t - t 0 λ s ds (0 ≤ t ≤ T )
is called the compensated process. By Lemma 4.1 it has the property that for every bounded, predictable process (H t ), the process

t 0 H s X λ (ds)
0≤t≤T is a martingale. We are now in a position to state the analogue of Borell's formula for the Poisson measure. In the following theorem π T denotes the Poisson measure of parameter T , i.e.,

π T (n) = T n n! e -T for n ∈ N = {0, 1, 2, . . .}
where we abbreviate π T (n) = π T ({n}).

Theorem 4.2. Let f : N → R be bounded and let T > 0. Then we have

log N e f dπ T = sup λ E f (X λ T ) - T 0 (λ t log λ t -λ t + 1) dt , (27)
where the supremum is taken over all bounded, non-negative, predictable processes (λ t ) 0≤t≤T , and (X λ t ) 0≤t≤T is the associated counting process, defined by (25). Moreover the supremum is actually a maximum.

Proof. Let (P t ) t≥0 be the Poisson semigroup: For every g : N → R

P t g(x) = n∈N g(x + n) π t (n).
We shall show that for every predictable bounded process (λ t ) we have

log P T (e f )(0) ≥ E f (X λ T ) - T 0 (λ t log λ t -λ t + 1) dt , (28)
with equality if λ is chosen appropriately. Let us start with the inequality. Note that for every g : N → R and t ≥ 0,

∂ t P t g = ∂ x P t g
where ∂ x g(x) = g(x + 1) -g(x) denotes the discrete gradient. Letting

F (t, x) = log P T -t (e f )(x)
we obtain ∂ t F = -e ∂xF + 1.

Let λ be a predictable, non-negative, bounded process and let

M t = F (t, X λ t ) - t 0
(λ s log λ s -λ s + 1) ds.

Almost surely, the process (M t ) 0≤t≤T is a piecewise absolutely-continuous function in t.

Hence the distributional derivative of the function t → M t is almost-surely the sum of an integrable function on [0, T ] and finitely many atoms. Namely, for any fixed t ∈ [0, T ],

M t = M 0 + t 0 ∂ x F (s, X λ s-) X λ (ds) - t 0 e ∂xF (s,X λ s ) + λ s log λ s -λ s ds (29)
where X λ s-denotes the left limit at s of X λ . Note that since X λ has only finitely many atoms we can replace X λ s by X λ s-in the second term of the right hand side. Setting α t = ∂ x F (t, X λ t-), we may rewrite (29) as follows:

M t -M 0 = t 0 α s X λ (ds) - t 0 (e αs + λ s log λ s -λ s -α s λ s ) ds. (30)
Recall that X λ t = X λ t -t 0 λ s ds is the compensated measure. Note that F (t, x) is continuous in t and that (X λ t-) is left continuous in t. Thus (α t ) is left continuous. Since it is also adapted, it is predictable. Moreover both (α t ) and (λ t ) are bounded. Consequently the first summand on the right-hand side of (30) is a martingale. Moreover, since

e x + y log y -y -xy ≥ 0 ∀x ∈ R, y ∈ R + (31)
the second integral on the right-hand side of (30) is non-negative. Therefore (M t ) 0≤t≤T is a supermartingale. In particular M 0 ≥ E[M T ], which is the desired inequality (28).

There is equality in (31) if e x = y. Hence if λ is such that

λ t = e ∂xF (t,X λ t-) , (32)
for almost every t and with probability one, then M is a martingale and we have equality in (28). Note that the function e ∂xF (t,x) is continuous in t and bounded. We prove in Lemma 4.3 below that under these conditions, a solution to (32) does indeed exist, which finishes the proof of the theorem.

Remarks.

1. It is also possible to prove Theorem 4.2 by using the Girsanov change of measure formula for counting processes. The argument presented here has the advantage of being self-contained.

2. Theorem 4.2 can probably be generalized in several ways. Firstly, up to some annoying technical details, the argument should work just the same for a function f that depends on the whole trajectory of the process rather than just the terminal point. On the left hand side, the Poisson distribution should then be replaced by the law of the Poisson process of intensity 1 on [0, T ]. In the Gaussian case, this pathspace version of the formula is known as the Boué-Dupuis formula, see [START_REF] Boué | A variational representation for certain functionals of Brownian motion[END_REF]. Then one can also replace the interval [0, T ] equipped with the Lebesgue measure by a more general measure space, leading to a Borell type formula for Poisson point processes. In a sense this was already carried out by X. Zhang in [START_REF] Zhang | Clark-Ocone formula and variational representation for Poisson functionals[END_REF], but the result in [START_REF] Zhang | Clark-Ocone formula and variational representation for Poisson functionals[END_REF] does not seem to recover our theorem.

A dual version version of

Borell's formula involving relative entropy was proved by the second-named author in [START_REF] Lehec | Representation formula for the entropy and functional inequalities[END_REF]. This can be done in the Poisson case too. The formula then reads: if µ is a probability measure on N whose density with respect to π T is bounded away from 0 and +∞, then the relative entropy of µ with respect to π T satisfies

H(µ | π T ) = inf λ E T 0 (λ t log λ t -λ t + 1) dt ,
where the infimum runs over all non-negative, predictable processes λ such that X λ T has law µ. This follows from the representation formula (27) and the Gibbs variational principle, which is the fact that the functionals ν → H(ν|π T ) and f → log e f dπ T are Legendre-Fenchel conjugates with respect to the pairing (f, ν) → f dν.

We now state and prove the technical lemma used in the proof of Theorem 4.2.

Lemma 4.3. Let G : [0, T ] × N → R + and assume that G is continuous in the first variable and bounded. Then there exists a predictable, bounded, non-negative, left-continuous process (λ t ) 0≤t≤T satisfying λ t = G(t, X λ t-), for almost every t ≤ T and with probability 1.

Proof. Consider the map H : (λ t ) 0≤t≤T → (G(t, X λ t-)) 0≤t≤T from the set of adapted, left-continuous, non-negative, bounded processes to itself. We will show that H has a fixed point. Let λ and µ be two processes in the domain of H. Since G is bounded, there is a constant C such that

E[|G(t, X λ t-) -G(t, X µ t-)|] ≤ C P(X λ t-= X µ t-). The probability that the integer-valued random variable X λ t-differs from X µ t-is dominated by E[|X λ t --X µ t -|]
, which in turn is the average number of atoms of N between the graphs of λ and µ on [0, t). Since λ and µ are predictable, it follows from Lemma 4.1 that

E|X λ t --X µ t -| = E t 0 |λ s -µ s | ds . Therefore E |G(t, X λ t-) -G(t, X µ t-)| ≤ C E t 0 |λ s -µ s | ds .
This easily implies that H is Lipschitz with constant 1/2 for the distance

d(λ, µ) = T 0 e -2Ct E[|λ t -µ t |] dt.
Thus, being a contraction, the map H has a fixed point, which is the desired result. Proof. By approximation we may assume that all four functions are bounded. Let α and β be two non-negative, bounded, predictable processes. We observe that (X α + X β )/2 coincides with the process X λ where

A discrete Prékopa-Leindler inequality

λ = (α ∧ β) χ + (α ∨ β) (1 -χ)
where χ t is the indicator function of the event that X

α t-+ X β t-is even. Similarly (X + Y )/2 = X µ where µ = (α ∧ β) (1 -χ) + (α ∨ β) χ.
Note that for every t ∈ [0, T ] either µ t = α t and λ t = β t or the other way around. In particular, for every function ϕ : [0, ∞) → R, ϕ(α t ) + ϕ(β t ) = ϕ(λ t ) + ϕ(µ t ), ∀t ∈ [0, T ]

Using the hypothesis made on f, g, h, k we get that for a continuous function ϕ, f (X α T ) + g(X β T ) - Taking the supremum in α and β and using the representation formula for f and g yields the result.

Rescaling appropriately, we obtain as a corollary a Prékopa-Leindler type inequality for the counting measure on Z. where E c denotes the complement of E. We claim that µ + and µ -coincide on I. This is in fact the statement of Theorem 1 in [START_REF] Picard | Formules de dualité sur lespace de Poisson[END_REF], we recall here the short proof for completeness. For E × A ∈ J , the random variable N (E) is independent of the set A, hence Since the class J is stable by finite intersections and generates the σ-field I the claim follows from the monotone class theorem.

Next recall the definition of the predictable σ-field P and observe that B(R + ) ⊗ P ⊆ I.

As a result, since (H t ) and (λ t ) are predictable, as a function of (t, u, ω),

H

Figure 1 :

 1 Figure 1: The density of ν 0,1,1,2 where µ is uniform on the interval [1, 2].

FollowingProposition 5 . 1 .+ y 2 + k x + y 2 ,

 512 Borell we derive in this section a Prékopa-Leindler type inequality from the representation formula (27). Recall that if x is a real number we denote its lower integer part by x and its upper integer part by x . For a, b ∈ R we denote a ∧ b = min{a, b} and a ∨ b = max{a, b}. Recall that π T denotes the Poisson distribution of parameter T . Let f, g, h, k : N → R satisfy f (x) + g(y) ≤ h x ∀x, y ∈ N. Then, N e f dπ T N e g dπ T ≤ N e h dπ T N e k dπ T .

  t ) dt.Choosing ϕ(x) = x log x + x -1, taking expectation, and using the representation formula in Theorem 4.2 for h and k we getE f (X α T ) -T 0 ϕ(α t ) dt + E g(X β T ) -T 0 ϕ(β t ) dt ≤ log N e h dπ T + log N e k dπ T .

Proof of Theorem 1 . 4 .

 14 Let Y n be a random variable having the Poisson law of parameter n and let X n = Y n -n. Applying the previous proposition to the functions f, g, h, k (translated by -n) we get E e f (Xn) E e g(Xn) ≤ E e h(Xn) E e k(Xn) .(33)On the other hand, for any fixed k ∈ Z, letting n tend to +∞ and using Stirling formula we getP(X n = k) = n n+k (n + k)! e -n = 1 √ 2πn (1 + o(1)).Hence by the dominated convergence theorem√ 2πn • E e f (Xn) n→∞ -→ x∈Z e f (x) ,and similarly for g, h, k. Hence multiplying (33) by n and letting n tend to +∞ yields the result.6 Appendix: Proof of Lemma 4.1We write B([0, T ]) for the Borel σ-field of [0, T ]. Let µ + and µ -be the measures on[0, T ] × R + × Ω equipped with the σ-field B([0, T ]) ⊗ B(R + ) ⊗ F defined by µ + (dt, du, dω) = N (ω)(dt, du) P(dω) µ -(dt, du, dω) = L(dt, du) P(dω),where L is the Lebesgue measure on the strip [0, T ] × R + while N (ω) is the measure on this strip given by the Poisson process N . Let I be the σ-field on [0, T ] × R + × Ω generated by the classJ = {E × A; E ∈ B([0, T ] × R + ), A ∈ F E c }

µ

  + (E × A) = E[N (E)1 A ] = E[N (E)]P(A) = L(E)P(A) = µ -(E × A).

  t 1 {u≤λt} is measurable with respect to I. We may therefore integrate H t 1 {u≤λt} with respect to µ + or µ -and obtain the same outcome. In other words,E [0,T ]×R + H t 1 {u≤λt} N (dt, du) = E [0,T ]×R + H t 1 {u≤λt} dtdu .H t 1 {u≤λt} N (dt, du)

	From (34) we obtain that
	T
	E	H t X λ (dt) = E
	0	[0,T ]×R +
		= E	H t 1 {u≤λt} dtdu
		[0,T ]×R +
		T
		= E	H t λ t dt ,
		0
	completing the proof of lemma 4.1.
			(34)

Acknowledgements. The first-named author was supported in part by a grant from the European Research Council (ERC).


Shamsuddin N Ladha

Piyush Yadav

Shailesh Deshpande

FEATURE SET CONSOLIDATION FOR OBJECT REPRESENTATION BY PARTS

Keywords: Object representation by parts, part pruning, SIFT, classification, local image descriptor

Image data is growing by leaps and bounds. Machine learning based applications that run on image datasets increasingly use local image feature descriptors. In a sense, we can now visualize images as objects with features as parts. Typically there are thousands of local features per image, resulting in an explosion of feature set size for already humungous image datasets. In this paper we present a feature set consolidation strategy based on two aspects: pruning of non-discriminatory features across different object types and association of matching features for the same object type. We showcase the effectiveness of our consolidation strategy by performing classification on a building dataset. Our method reduces storage space footprint (~5%) and classification runtime (~4%), and increases classification accuracy (~2%).

INTRODUCTION

Traditionally, for vision based machine learning problems, each object instance in the dataset has been cohesively represented as a single-possibly very high dimensionalvector. This approach could be attributed to the observation that certain machine learning problems are better handled in high dimensional spaces (e.g. Support Vector Machine (SVM) [START_REF] Boser | A Training Algorithm for Optimal Margin Classifiers[END_REF]). However, with the availability of discriminating local features (such as SIFT [START_REF] Lowe | Distinctive Image Features from Scale-Invariant Keypoints[END_REF]) it has been possible to conjure objects as a composition of parts (or components). For example, given an image of an animal (object), we refer as parts all the extracted local features (by a feature extraction algorithm). Here we have assumed that each image has only one object. This notion of representation by parts is different than the cognitive notion of physical parts of the animal such as head, body, etc. While locally discriminative features may provide a richer description of objects, they pose several new challenges for object representation and learning,  A fixed high dimensional vector representation of object remains no longer practical as the number and order of local features generally vary across objects.  The number of parts that represent a single object are typically in the order of hundreds (e.g. SIFT features) thus increasing the space and time complexity to store and process these features.  Both intra-class and inter-class variability shoots up significantly due to richer object representation.  The dimensionality of individual features is generally low. This could be problematic for learning algorithms where high dimensionality is desired (e.g. SVM). To address these challenges several approaches are found in the literature. These approaches can be broadly classified into three types, i.e., approaches that reduce a) the number of features; b) the dimensions of feature; c) both the number and dimensions of feature. Common approaches such as object representation by top 𝑁 features or feature clustering aim at reducing the number of features. In [START_REF] Ledwich | Reduced sift features for image retrieval and indoor localisation[END_REF], the authors reduced the number of SIFT features used for indoor scene representation based on the observation that a majority of the detected key points do not match between images that share common camera view point. Authors of [START_REF] Montazer | Content based image retrieval system using clustered scale invariant feature transforms[END_REF] represented objects by clustering SIFT features of an object with k-means to solve content based image retrieval (CBIR) problem. Methods such as PCA-SIFT [START_REF] Yan | PCA-SIFT: a more distinctive representation for local image descriptors[END_REF] and feature quantization [START_REF] Hare | Efficient Clustering and Quantisation of SIFT Features: Exploiting Characteristics of the SIFT Descriptor and Interest Region Detectors Under Image Inversion[END_REF] reduce feature dimensionality. For example, [START_REF] Foo | Pruning SIFT for Scalable Nearduplicate Image Matching[END_REF] used PCA-SIFT and applied locality sensitive hashing (LSH) for fast object retrieval. [START_REF] Hare | Efficient Clustering and Quantisation of SIFT Features: Exploiting Characteristics of the SIFT Descriptor and Interest Region Detectors Under Image Inversion[END_REF] introduced a method to quantize sift features based on inverted intensity images and then performed clustering to create a codebook for image matching. This method reduces not only feature dimensionality but also the size of feature set. A key aspect that has not been addressed, in these methods, is elimination of features that are non-discriminative across different objects. Although SIFT features are discriminative for a given image, there is a possibility that some of these features may lose their discriminatory powers in presence of similar features from altogether different objects. In other words, the distance between some of the SIFT features belonging to different objects could be small enough to cause incorrect feature matching across objects. In this paper we intend to bridge the gap in the existing methods by pruning non-discriminative SIFT features across objects. In addition, we improvise feature clustering by integrating concepts of associative memory from [START_REF] Kanerva | Sparse Distributed Memory[END_REF]. Our approach that is a combination of non-discriminatory feature elimination and associative feature clustering is termed as feature set "consolidation". In the following section, we describe in detail both aspects of consolidation.

FEATURE SET CONSOLIDATION

Consider a simplified schematic (see Figure 1), that demonstrates object representation by parts. The figure shows two object instances of two different object classes. Each object instance is represented by corresponding three parts. Consolidation has two aspects, i.e., elimination and associative clustering. To explore the impact of consolidation in greater detail we will study it in the context of object classification problem. For the schematic problem, accurate classification of objects would require learning a complex decision boundary (in this case non-linear). If we eliminate the non-discriminatory parts that lie near the decision boundary then the boundary simplifies (becomes linear in this case). In any case, elimination of non-discriminatory parts will not increase the complexity of the decision boundary because any nondiscriminatory part could be labeled as internal (not impacting decision boundary) or boundary part. In general, while building a classifier it is important to strike a good balance between the errors due to bias and variance in order to have good generalization accuracy. It is evident that pruning of non-discriminative parts helps reduce not only bias-by reducing model complexity-but also variance-by reducing unwarranted variance in the dataset. Pruning, as demonstrated in the schematic (Figure 1), helps reduce errors due to bias and variance. While pruning of parts may sound similar to pruning of object instances [START_REF] Wilson | Instance Pruning Techniques[END_REF], there are some key differences. In general, pruning of parts-unlike object pruning-does not lead to loss of representation of entire object except when all the parts that compose an object are pruned. Inter-class object parts are likely to be closer in hyper-space than interclass objects due low dimensionality of parts and large number of object parts. Large number of object parts allows us to not only be aggressive while pruning but also statistically determine the class of an object based on the underlying class distribution over parts. Thus an object will be classified accurately as long as the majority of its parts are correctly classified. Our approach can be used to complement existing approaches for feature reduction mentioned earlier.

Next, we look at the associative clustering aspect of feature set consolidation. Typical clustering methods substitute parts from same object class that belong to a single cluster with a single representative part. The strength of the cluster that is indicative of the association of cluster members is not utilized. In addition to typical clustering, we store the strength of each cluster. At the time of classification of an object, the cluster strength contributes to the class distribution that will be used to determine the class of the query object. Higher the cluster strength larger will be the contribution to the distribution and vice-versa. This concept is inspired from sparse distributed memory (SDM) [START_REF] Kanerva | Sparse Distributed Memory[END_REF]. However, it is significantly different from SDM. While SDM permits clustering of instances with same or different class labels, our approach only permits clustering of same label instances. Thus clusters generated by our method are pure, i.e., all the cluster members belong to the same object class. The two aspects of consolidation, i.e., pruning and associative clustering are complementary. Together they contribute to reduction is storage space, classification runtime and help strike a balance between model bias and dataset variance.

EXPERIMENTS AND RESULTS

To demonstrate the concept of pruning and its effectiveness we have implemented a parts based image classification system. Classification experiments consist of following three steps, 1. Store the training data set features in a database in an off-line mode. 2. At runtime match the test data set features to those stored in the database. This results in assignment of a class label to each feature in the test set. 3. Derive a class label for each query object from its feature set labels (from Step 2). Our experiments were conducted on a standard Intel core i5 processor with Windows 7 Professional operating system and 4GB RAM. The system was developed in Python with OpenCV port [START_REF]OpenCV-Python Tutorials[END_REF]. We have used the Zurich building dataset [START_REF] Zurich | Zurich Building Image Database[END_REF] for this purpose. The dataset with a total of 1005 images consists of five views of 201 different buildings in Zurich. Each view image has dimensions 640 x 480. Local features (SIFT) were extracted for each image where each SIFT feature interest point, in the context of this paper, represented a part of the image. The number of SIFT points extracted per image ranged from about 1000 to 2000. The dataset was split into training and test sets. Two sets of classification experiments were conducted, one without and another with consolidation of features. For the purpose of consolidation, the SIFT points were matched using brute force matcher. Matched points with distance less than a specific threshold were only considered as appropriate matches. This step of storing training data, without or with consolidation, was done in an off-line mode. At runtime, during classification of a query image, its local features were extracted and matched to the stored features resulting in a distribution over class labels. The query image is labeled with the mode of the distribution. A range of experiments were conducted on the dataset to assess the impact of consolidation. Experiments with varying number of SIFT features, different size training datasets, and matching threshold were performed. Table 1, Table 2, and Table 3 respectively list the number of key points, classification time, and F-score values for these experiments. Figure 1 and Figure 2 respectively show the plot of difference, without and with consolidation, in the number of SIFT key points (left vertical axis) and classification runtime (right vertical axis), as a function of the maximum number of SIFT features extracted per image for three different sizes of test data. We see that as the number of extracted SIFT points per image increase, consolidation leads to increase in savings of storage space and reduction in classification runtime. On an average our method reduces storage space requirements by 5% and shortens classification runtime by 4%. As the size of training set reduces, reduction in storage space decreases (the dashed lines). These results should be assessed in conjunction with the results of Figure 3 that show the classification performance of the two approaches (without and with consolidation) as the number of SIFT points extracted per image varies. For each experiment, assessment of the classification results was done by computing F-score. Table 3 lists the actual Fscores obtained for various experimental setups. Within a given experimental setup (e.g. with 20% dataset held out), we observe that the F-score steadily increases for both without and with pruning stored features. The F-scores for classification with consolidation is better on an average by about 2-3% particularly for experiments where the number of SIFT features are limited (highlighted entries in Table 3). For other experiments the F-scores with consolidation are similar to those without consolidation. This trend is observed for all the experimental setups (as seen in Figure 3) indicating no adverse effect of our pruning strategy on classification outcomes. A sample result of unrelated matching key points, one point from building façade and another from zebra crossings, which were pruned from two different building images, is shown in Figure 4.

SUMMARY

In this paper, we have presented a feature set consolidation strategy suitable for scenarios where objects are represented by parts. During the process of consolidation not only the non-discriminatory parts (belonging to different objects), that often complicate machine learning models (e.g. the classifier decision boundary), were eliminated but also the matching parts (across different instances of the same object class) were strengthened. Our experiments have demonstrated that such a strategy is efficient as it leads to significant savings in terms of both storage space and classification runtime without compromising classification performance. With ever exploding world of image data and greatly improved local descriptors, such consolidation strategies will be helpful in realizing ever increasing memory and time intensive applications on limited computing power devices.

Figure 1 :

 1 Figure 1: Object representation by parts

Figure 2 :

 2 Figure 2: Experiments showcasing reduction in the number of SIFT points stored and classification runtime due to pruning.

Figure 3 :

 3 Figure 3: Graph with F-scores for different experiments with and without pruning.

Figure 4 :

 4 Figure 4: Sample matching points (in blue) from two different images (separated by the yellow marker) that were pruned.

Table 1 : Number of features stored for different experimental setups.

 1

	Dataset	Experimental	Experimental	Experimental
	held out	Setup 1	Setup 2	Setup 3
	%	(20%)	(40%)	(60%)
	SIFT	Total	Pruned	Total	Pruned	Total	Pruned
	features	features	features	features	features	features	features
	per image	stored	stored	stored	stored	stored	stored
	10	8040	7892	6030	5946	4020	3974
	30	24120 23488 18090 17722 12060 11875
	60	48240 46515 36180 35183 24120 23632
	100	80400 76837 60300 58255 40200 39249
	200	160800 151861 120600 115565 80400 77981
	400	321499 309461 241099 233006 160699 155940

Table 2 : Classification runtime (in seconds) for different experimental setups.

 2

	Dataset	Experimental	Experimental	Experimental
	held out	Setup 1	Setup 2	Setup 3
	%	(20%)	(40%)	(60%)
	SIFT	Time	Time	Time	Time	Time	Time
	features	(No	(With	(No	(With	(No	(With
	per image	pruning)	pruning)	pruning)	pruning)	pruning)	pruning)
	10	74.16	67.16	135	134.32 201.34 201.64
	30	72.78	73.57 142.38 143.35 210.01 209.87
	60	91.93	91.08 171.84 170.73 236.74 235.68
	100	133.84 129.38 237.14 234.89 297.72 295.54
	200	352.59 341.84 573.09 504.33 573.69 548.1
	400	1007 985.38 1542.26 1450.8 1618.37 1574

Table 3 : Classification F-scores for different experimental setups.

 3

	Dataset	Experimental	Experimental	Experimental
	held out	Setup 1	Setup 2	Setup 3
	%	(20%)	(40%)	(60%)
	SIFT	F-score	F-score	F-score	F-score	F-score	F-score
	features	(No	(With	(No	(With	(No	(With
	per image	pruning)	pruning)	pruning)	pruning)	pruning)	pruning)
	10	0.11	0.10	0.11	0.10	0.04	0.05
	30	0.30	0.32	0.30	0.34	0.19	0.20
	60	0.51	0.55	0.53	0.54	0.36	0.40
	100	0.68	0.69	0.75	0.74	0.54	0.57
	200	0.84	0.83	0.84	0.86	0.73	0.72
	400	0.92	0.90	0.92	0.90	0.81	0.81


Georges-Pierre Bonneau

Stefanie Hahmann

e p r i n t Polyhedral modeling

Keywords: triangular meshes, visual continuity, arbitrary topology, visualization

Polyhedral meshes are used for visualization, computer graphics or geometric modeling purposes and result from many applications like iso-surface extraction, surface reconstruction or CAD/CAM. The present paper introduces a method for constructing smooth surfaces from a triangulated polyhedral mesh of arbitrary topology. It presents a new algorithm which generalizes and improves the triangle 4-split method [7] in the crucial point of boundary curve network construction. This network is then filled-in by a visual smooth surface from which an explicit closed form parametrization is given. Furthermore, the method becomes now completely local and can interpolate normal vector input at the mesh vertices.

INTRODUCTION

Polyhedral meshes consisting of a collection of vertices, edges and triangular faces which describe an oriented 2-manifold in IR 3 are dealt with in many applications. They result f.ex. from iso-surface extraction algorithms like marching cubes. Surface reconstruction methods used for virtual reality purposes of CAD/CAM applications produce a fine triangulated polyhedral mesh approximating an object from which a 3D scanner has sampled millions of unorganized data points. In medical imaging where CT or MRI scanner produce a set of slices with a grid of gray level values a 3D volume is reconstructed. But contour line approximation on each slide followed by contour reconstruction also produces a polyhedral mesh. They also occur in CAD/CAM applications where the classical tensor product surface models fail. Topological complex situations like inner corner blends need an extraneous treatment. Surfaces defined on polyhedral meshes can however represent surfaces of arbitrary topological type.

In particular in CAD/CAM applications or surface reconstruction but also for computer graphics or visualization purposes, it is often desired to get a visual smooth representation of the object which is roughly approximated by a polyhedral mesh. On one hand surface subdivision algorithms converge to a smooth surface starting from a polyhedral mesh, but they generally don't provide an explicit closed form expression of the resulting surface. Therefore, computations on such surfaces like intersections or evaluations of intrinsic quantities like curvature can become tricky. On the other hand, and this is the subject the present paper deals with, parametric polynomial surfaces fulfill all these requirements. The polyhedral mesh serves as a control mesh. Recall that, since every polyhedra can be triangulated, it is sufficient to consider as input a triangular mesh. Triangular surface patches which join with geometric continuity in order to produce an overall smooth surface while interpolating the mesh vertices can then be computed.

The main difficulty in developing a triangular interpolation scheme for polyhedral meshes of arbitrary topological type consists of ensuring the G 1 continuity (visual smoothness, tangent plane continuity). Three kinds of schemes exist for polynomial representations: the Clough-Tocher-like domain splitting methods [START_REF] Farin | A construction for visual C 1 continuity of polynomial surface patches[END_REF], [START_REF] Piper | Visually smooth interpolation with triangular Bézier patches[END_REF] the G 2 continuous constructions [START_REF] Peters | Smooth interpolation of a mesh of curves[END_REF]; or the twist compatible boundary curve schemes [START_REF] Loop | A G 1 triangular spline surface of arbitrary topological type[END_REF], [START_REF] Hahmann | Triangular G 1 interpolation by 4-splitting domain triangles[END_REF]. Further non-polynomial methods also exist [START_REF] Gregory | N-sided surface patches[END_REF], [START_REF] Hagen | Geometric surface patches without twist constraints[END_REF], [START_REF] Nielson | A transfinite, visually continuous, triangular interpolant[END_REF], [START_REF] Bajaj | Smoothing polyhedra using implicit algebraic splines[END_REF]. In addition to an explicit closed form parametrization, the scheme should be local. This means that a local modification of the polyhedral control mesh should only locally affect the surface. This allows an interactive modeling of smooth surfaces by interactive modeling of polyhedral meshes.

The present paper deals with such a polyhedral modeling method. It is based on a recently developed method [START_REF] Hahmann | Triangular G 1 interpolation by 4-splitting domain triangles[END_REF], called 4-split method (chapter 2), and presents a new algorithm (chapter 3) which generalizes and improves the 4-split scheme in such a way that the number of degrees of freedom for the construction of the boundary curves can significantly be increased. A more flexible and more intuitive control of the surface shape is now possible. Further features of the scheme are localness and normal vector interpolation at the mesh vertices.

THE 4-SPLIT METHOD

The triangular 4-split method interpolates a polyhedral mesh of arbitrary topology by a visual smooth (G 1 continuous) piecewise quintic Bézier surface. It is local, affine invariant and provides an explicit closed form parametrization. Several degrees of freedom are available for shape control.

The following subsections briefly present first the method's basic idea and why this is benefit in comparism to other existing schemes. Then the algorithm itself is outlined without going into too much detail, but enough for reproducing this method.

Let the polyhedral control mesh M be a set of vertices, edges and triangular faces that describe an oriented 2-manifold in IR 3 .

The surface S which interpolates the vertices of M is composed of triangular surfaces M i which are in one-to-one correspondence to the mesh facets and which meet with tangent plane (G 1 ) continuity along the common boundary curves. The M i 's are piecewise A u t h o r P r e p r i n t polynomial images of the unit domain triangle, we call them macropatches.

Four split of the mesh triangles.

The basic idea of the present method is to subdivide each domain triangle into four subtriangles by joining together the edge midpoints. Each macro-patch therefore consists of four Bézier patches. A convenient parametrization of the macro-patches concentrates on the mesh vertices, see fig. 1. -it allows to introduce more degrees of freedom for each macro-patch, -the boundary curves and cross-boundary tangents are piecewise polynomial functions (2 pieces),

u i-1 M i . u i+1 u i M i M i-1 p 0 1/2 1
-the surface scheme is local.

It will be seen later that the requirements on the boundary curves concern their first and second derivatives at the endpoints. In order to keep the surface scheme local the derivatives should be independent from one vertex to another. Instead of degree 5 curves, piecewise curves need only to be of degree 3. Similar considerations hold for the cross-boundary tangents and finally the macro-patches are composed of four quintic Bézier triangles, which is the lowest possible degree for this kind of methods.

The polyhedral 4-split mesh interpolation method consists of three main steps:

-generating first and second derivative informations at the mesh vertices, -computing first order cross derivatives along the curves, -computing inner Bézier points for each macro-patch.

Generating first and second derivative informations at the mesh vertices.

The first step of the algorithm consists of constructing a network of piecewise cubic C 1 boundary curves of the macro-patches interpolating the mesh vertices. It turns out that the G 1 conditions for joining a network of patches together imposes quite restrictive conditions on the set of boundary curves incident to a common mesh vertex.

Let p be a mesh vertex of order n, i.e. p is common to n mesh facets (or edges) and therefore has n neighbour points p i . The n first and second derivatives of the n boundary curves in p are not free. The boundary curve network is constructed from these derivative informations by computing for each vertex the curve pieces incident to that vertex. Since each piece of curve is a cubic Bézier parameterized over 0 1 2 ], this means that the control points b0, b1, b2 must be computed. They are directly related to the first and second derivatives at the curve's endpoints by x(0) = b0 x 0 (0) = 6(b1 ; b0) x 00 (0) = 24(b2 ; 2b1 + b0):

The complementary curve pieces are computed when the neighbour mesh vertices are treated. Finally, the curve network is closed by computing for each boundary curve the middle control point b3 which is common to both curve pieces and which should be calculated such that the two pieces join with C 1 continuity.

Let us consider a mesh vertex p and its n curve pieces as illus- Then the first and second derivative informations of the n boundary curves have to be calculated in order to determine b i 1 and b i 2 .

trated in fig. 2. b =p 0 i b 1 i b i b 3 i p p i-1 p i+1 i 2 tangent plane
These derivatives have to satisfy a set of equations (see ( 4) and [START_REF] Hahmann | Triangular G 1 interpolation by 4-splitting domain triangles[END_REF] of chapter 3) and are therefore not free. The following solution has been found [START_REF] Hahmann | Triangular G 1 interpolation by 4-splitting domain triangles[END_REF], which consists of

b i 0 = p b i 1 = b i 0 + n P n j=1 cos ; 2 (j;i) n p j b i 2 = 1 ; 2 3 2 ]b i 0 + 1 n P n j=1 cos ; 2 (j;i) n p j + 2 P n j=1 p j ( 1=6 if j = i ; 1 i + 1 1=3 if j = i 0 otherwise (1)
where i = 1 : : : n is taken modulo n. A few degrees of freedom f 1 2 g are nevertheless available for each vertex. The control points b i 1 form an affine transformation of a regular n-gon. They all lie in the same plane, which is the tangent plane at p. controls the magnitude of the first derivatives. The second derivative control points b i 2 lie in a plane containing the 3 points b i 0 , b i 1 and di = 2 6 (2p + p i;1 + 2 p i + p i+1 ). More details can be found in [START_REF] Hahmann | Triangular G 1 interpolation by 4-splitting domain triangles[END_REF].

Computing first order cross derivatives along the edges.

The second step of the algorithm is also related to the G 1 continuity conditions. It consists of computing the network of tangent planes along the boundary curves. The boundary curves determine already the tangent planes of the macro-patches at the mesh vertices. At each curve point ui 2 0 1] the tangent plane is spanned by the curve derivative @M i @u i (ui 0) and one of the cross-boundary derivatives @M i @u i+1 (ui 0).

A u t h o r P r e p r i n t The cross-boundary derivatives are defined as @Mi @ui+1 (ui 0) = i(ui) @Mi @ui (ui 0) + i(ui)Vi(ui) (2) where i and i are scalar functions and Vi is a vector valued function. These functions should -be of minimal degree, and -interpolate the first and second derivatives at the curve endpoints.

p i-1 p i+1 p i p M i (u ,0) u i M i (u ,0) u i+1 i i
It turns out that Vi is a piecewise continuous degree two polynomial, which can be written in Bézier form. As in the case of the boundary curves, only the pieces incident to the mesh vertex p are computed:

Vi(ui) = P 2 j=0 v i j B 2 j (ui) ui 2 0 1 2 ] with v i 0 = 6 n P n j=1 sin 2 (j;i) n p j v i 1 = 6 n P n j=1 1 0 i ( 1 ; 8 0 + 4 0 ) t a n ( n ) ; 1 i sin ; 2 (j;i) n + 4 0 i 2 0 n 1 if j = i + 1 ;1 if j = i ; 1 v i 2 is free subject to Vi( 1 2 + ) = Vi( 1 2 ; ) (3)
where 0 = cos 2 n , 1 = 1 + 2 sin 2 n , 0 = s i n 2 n , 1 i = sin 2 n i ; sin 2 n , and ni is the order of p i i is a linear function, therefore the cross-boundary tangents are therefore piecewise cubic. However, the surface must be of degree 5 in order to ensure C 1 continuity between the 4 inner Bézier patches. The ribbons of cross-boundary tangents are now determined, and the first inner row of control points of the macro-patches can be computed from ( 1) and (3).

Computing inner Bézier points.

Each macro-patch is composed of four quintic triangular Bézier patches. The boundary curves of the macro-patch are the twice degree elevated curves of section 2.2. The cross-boundary tangents of sect. 2.3 determine the first inner row of control points after one degree elevation. The remaining 15 inner control points, which are highlighted in fig. 4-left, are now computed by joining the four inner patches with C 1 continuity. The necessary and sufficient C 1continuity conditions between two internal Bézier patches inside one macro-patch are shown in fig. 4-right: all pairs of adjacent triangles must form a parallelogram, otherwise C 1 continuity would not be possible. In [START_REF] Hahmann | Triangular G 1 interpolation by 4-splitting domain triangles[END_REF] it was shown, that the first and last pairs of adjacent triangles in fig. 4 thus free for shape control of the macro-patch. Minimization of an energy functional produces a well shaped smooth surface.

GENERATING DERIVATIVE INFORMA-TION AT THE MESH VERTICES

The polyhedral mesh interpolation method which we now want to describe in detail has the same skeleton as the 4-split method sketched out in chapter 2. It preserves all the attractive features but generalizes the method in an important point, that is the construction of the boundary curve network. A more flexible approach will be presented offering more degrees of freedom which are necessary in order to obtain a well shaped curve network. It is well known that visual smoothness, which stands for G 1 continuity, generally doesn't always guarantee a "nice shape" [START_REF] Mann | A survey of parametric scattered data fitting using triangular interpolants[END_REF]. "Bad" shape features, like wiggles, bumps or self-intersections, should be avoided and the shape inherent to the polyhedral input mesh should be reproduced smoothly. It is therefore obvious that the construction of the boundary curve network is the crucial point of such an interpolation method; it is predominant for the final result.

The generalized method we propose for the construction of a G 1 compatible curve network interpolating a polyhedral mesh is due to a detailed study of the G 1 and the twist conditions which are imposed at the mesh vertices. This leads to a new and quite different approach as what has been done until now in related works [START_REF] Hahmann | Triangular G 1 interpolation by 4-splitting domain triangles[END_REF], [START_REF] Loop | A G 1 triangular spline surface of arbitrary topological type[END_REF]. The benefit is threefold: more freedom and intuitive control for the boundary curve network is available; one gets a more local interpolant; normal vector interpolation becomes possible.

G 1 Continuity.

Let us start by remembering that two patches X and Y sharing a common boundary curve meet with first order geometric continuity, denoted G 1 , if they share a common tangent plane at all points of their common boundary, i.e. if there exist three scalar functions such that (s)Ys(s 0) = (s)Xr(0 s ) + (s)Yt(s 0) : When joining polynomial patches together to a network of patches with G 1 continuity special attention has to be paid to the patch corners. The idea is to construct the triangular surface patches in oneto-one correspondence to the polyhedral mesh faces. Let us consider the n patches parameterized as in fig. 1, each triangular patch A u t h o r P r e p r i n t

Y (s,0) t X (0,s) r Y (s,0)
M i = M i (ui u i+1)
is the piecewise polynomial image of the unit triangle. i = 1 : : : n is taken modulo n.

In order to obtain at the end an explicit expression for each patch some simplifying and symmetric assumptions on the G 1 conditions have to be made as follows:

i(0) = i;1(0) =: 0 0 i (0) = 0 i;1 (0) =: 1 i(ui) i(ui)
1 2 i = 1 : : : n :

The G 1 condition between a pair of surface patches M i and M i;1 now looks as follows:

i(ui) @M i @ui (ui 0) = 1 2 @M i;1 @ui;1 (0 u i) + 1 2 @M i @ui+1 (ui 0): (4)

At the patch corner p, i.e. for ui = 0 , the condition (4) can be written in matrix form as follows:

P r 1 = 0

:

Note, that r 1 is a (n 3) vector containing the first partial derivatives of the patch boundary curves at the vertex ui = 0 . The function i can be determined from the fact that det P = 0 should be satisfied and that the derivative of i at ui = 0 should not depend on the value [START_REF] Bajaj | Smoothing polyhedra using implicit algebraic splines[END_REF]. i should also be of lowest possible degree which leads to the following choice:

i(ui) = ( cos 2 n (1 ; 2ui) + ui for ui 2 0 1 2 ]
(1 ; ui) + ( 1 ; cos 2 n i )(2ui ; 1) for ui 2 1

2 1] (6)

Twist Compatibility.

When joining n polynomial patches around a common vertex with G 1 continuity an additional condition has to be satisfied at the common vertex. It is called twist compatibility condition and given by T t = 1 r 1 + 0 r 2 @ 2 M 1 @u 1 @u 2 (0 0) . . . @ 2 M n @un@u 1 (0 0) :

It is obtained by differentiating (C) with respect to ui and evaluation at ui = 0 . The condition [START_REF] Hahmann | Triangular G 1 interpolation by 4-splitting domain triangles[END_REF] must hold because polynomial patches lie in the continuity class C 2 , implying that the twist terms (mixed partial derivatives) are identical. Different approaches exist in order to solve the twist compatibility problem [START_REF] Piper | Visually smooth interpolation with triangular Bézier patches[END_REF], [START_REF] Nielson | A transfinite, visually continuous, triangular interpolant[END_REF], [START_REF] Peters | Smooth interpolation of a mesh of curves[END_REF], [START_REF] Loop | A G 1 triangular spline surface of arbitrary topological type[END_REF], [START_REF] Hahmann | Triangular G 1 interpolation by 4-splitting domain triangles[END_REF]. We solve the system [START_REF] Hahmann | Triangular G 1 interpolation by 4-splitting domain triangles[END_REF] which is singular if n is even, by constructing a twist compatible boundary curve network i.e. boundary curves are determined such that the vectors r 1 and r 2 lie in the image space of T . It will turn out in the following subsections that the kernel of P and the image space of T is explicitly known [START_REF] Davis | Circulant Matrices[END_REF]. This will enable us to characterize the vectors r 1 and r 2 with a maximum number of degrees of freedom while simultaneously satisfying equations (4), ( 5), [START_REF] Hahmann | Triangular G 1 interpolation by 4-splitting domain triangles[END_REF].

Once the first and second order derivatives of the boundary curves at the vertices, r 1 , r 2 , are known, the whole curve network is also known. The construction of the cross-boundary tangents and the inner patches can then be continued as in the 4-split method, see chapter 2 in order to obtain a smooth polyhedral mesh modeling method.

First order derivatives.

Let us now explain how the boundary curves can be determined. They will consist of piecewise (two pieces) cubic curves in order to keep the scheme local by separating first and second derivatives from one vertex to its neighbouring vertices. The final triangular patches will be given explicitly in Bézier form, let us therefore use the Bernstein-Bézier notations [START_REF] Farin | Curves and Surfaces for Computer Aided Geometric Design[END_REF] for the boundary curves as well.

Let x(t) = For the boundary curve pieces incident to the vertex p = b i 0 we adopt the notations as illustrated in fig. 2. The joining curve pieces are determined from the neighbouring vertices p i (locality of the scheme). b i 0 = p is already known (interpolation of the mesh vertices), the vector r 1 is therefore determined by the control points b i 1 and r 2 is determined by the control points b i 2 , i = 1 : : : n . From ( 5) and [START_REF] Hahmann | Triangular G 1 interpolation by 4-splitting domain triangles[END_REF] The vector r 1 can therefore be written as a linear combination of k a and k

b " r 1 # n 3 = " k a # n 1 ( a )1 3 + " k b # n 1 ( b )1 3 (9)
where the vectors a and b are the degrees of freedom and can be chosen arbitrarily.

Geometrical interpretation:

Equations ( 8) and ( 9) can be interpretated geometrically in order to indicate the correct way for choosing a and b. The components r 1 i are the tangents of the boundary curves at the common vertex p which are all lying in the same plane, which is the tangent plane T p at p (see G 1 -condition (4)). fa bg is the basis of the tangent plane and k a i k b i are the coordinates of r 1 i in this basis, see fig. 6. These coordinates, k a i = c o s 2 i n , k b i = sin 2 i n imply that the first derivative control points b i 1 = b i 0 + 1 6 r 1 i of the boundary curves form an affine transformation of a plane regular n-gon. This affine transformation is explicitly given by the choice of the tangent plane basis a and b.

For practical use of these theoretical considerations it is now important to choose a and b such that the n tangents of the boundary curves in p fit the best. For this purpose we first determine "optimal tangents", denoted r 1 opt , by a simple heuristic rule and then approach them in a least squares sense. The heuristic rule [START_REF] Piper | Visually smooth interpolation with triangular Bézier patches[END_REF] consists of choosing the optimal tangent vectors r 1 iopt as lying in the plane spanned by the edge pp i and a mean normal vector of the mesh facets at p, see fig. 7.

These are the tangents which we want to reach by using a linear least squares approximation

Second order derivatives.

A further requirement on the boundary curves is to have second order derivatives at the end points lying in the image space of the matrix T . We can therefore choose any vector r 2 , such that r 2 = T r 2 [START_REF] Nielson | A transfinite, visually continuous, triangular interpolant[END_REF] where r 2 is an arbitrary (n 3) vector. T is a (n n) cyclic matrix. If n (order of the vertex p) is odd, then T is invertible, and thus n second derivatives for the n boundary curves r 2 can be chosen arbitrarily. And if n is even, T has only rank n ; 1, therefore some chosen n second derivatives, denoted r 2 opt can only be approximated, for example through a linear least squares approximation kr 2 ; r 2 opt k 2 ! min :

Since the matrix T has nearly full rank, this least squares minimization approaches quite good the desired second derivatives. As an example, see fig. 12, where the input mesh consists only of vertices of even order.

Note, that the total amount of n resp. n ; 1 vector valued degrees of freedom are available at each vertex. Once these second order derivatives are fixed, i.e. the control points b i 2 are fixed, the whole boundary curve network is entirely fixed. Remember, that each boundary curve is a piecewise cubic curve, which is required to be C 1 continuous, therefore the control points b i 3 and b k 3 of a boundary curve between the vertex p and its neighbour vertex p i have to be identical and b i 3 = b k 3 = 1 2 (b i 2 + b k 2 ).

Target middle points & target derivatives.

It's not quite evident how to choose the "optimal" second order derivatives at a vertex p for the boundary curves. On one hand, the scheme should keep a local one, i.e. both curve pieces of a boundary curve have to be constructed independent from each other. On the other hand, fixing b i 2 and b k 2 (k is the index of p in the neighbourhood of p i ) determines finally the shape of the boundary curve between p and p i .

In order to cope with both problems we propose the following solution. For each boundary curve we choose a middle control point, called target point and a derivative at that point, see fig. 8, which we want to approach or interpolate by the curve. This is equivalent to choosing b i 2 and b k 2 . It keeps the scheme local, but allows at the same time a control of the whole curve. Furthermore, these control handles offer an intuitive way for controlling the curve network. Optimal values can be found either by using a heuristic rule, or by minimizing some energy norm.

Note, that for odd order vertices, these optimal target handles can be interpolated exactly, while otherwise they have to approximated in a least squares sense. Figure 9 illustrates how important it is to have a maximum number of degrees of freedom for the second order derivatives. On the right some optimal target points and target derivatives have been approached. They are shown in red with respect to a cubic parametrization. The smooth surface is also shown together with the control polygon of the boundary curves (blue). They are shown as degree 5 curves (two times degree elevated). The sensitivity of the interpolation scheme is shown on the left figures, where some "bad" target points and tangent have been chosen.

Twist information.

Once the vectors r 1 and r 2 are determined the system (T) can be solved for the twist vectors t. An explicit expression can be given

as follows t = 1 ( k a a + k b b) + 0 r 2 : (11)
This expression is used in the next step of the algorithm, which is the construction of the cross-boundary tangent ribbons along the boundary curve network. Here it is important to know the twists explicitly, because otherwise it would not be possible to provide an explicit expression of the final surface. For this purpose and the last setp of filling-in the macro-patches C 1 -continuously one proceeds analogously to the 4-split method, see sections 2.3 and 2.4.

Remark 1: normal interpolation

Since we are able to exactly choose the tangent plane basis vectors, all one has to do in order to interpolate a given normal vector at the vertex is to choose these two tangent plane vectors orthogonal to that given normal.

Remark 2: Locality of the scheme The scheme is local in the sense that moving one vertex will only affect the macro-patches around that vertex, and more over will leave the tangent ribbons of the opposite edges to that vertex unchanged.

RESULTS

The first example stresses the ability to interpolate meshes with arbitrary topology. The input mesh is available at ftp:// ftp.cs.washington.edu/pub/graphics/meshes. Fig. 10-left shows the input mesh. Fig. 10-right shows the resulting G 1 surface together with the boundary curves. The order of the vertices ranges from 4 to 8 in this example. Fig. 11-left shows the control polygons of each of the triangular Bézier patches. Different colors have been used in order to distinguish the four Bézier patches per macro-patch. The curve control polygons are also highlighted, while the interpolated points are shown in red. Fig. 11-right shows the surface together with the boundary curve control polygons. It should be pointed out that we used a very crude treatment for the surface border and some additional work should be done on how to handle open meshes. Fig. 12, 13 show at the top an input mesh, in the middle the resulting surface together with the interpolated points in red, and the bottom the control polygons. For the torus input mesh in Fig. 14 a heuristic rule in order to find optimal second derivatives r 2 opt (see sect. 3.4) has been used. Even though all vertices have an even order, which means less degrees of freedom in the choice of the second derivatives, the torus is nearly reproduced.

A u t h o r P r e p r i n t

Figure 1 :

 1 Figure 1: Parameterization of n patches around a common vertex. This triangle 4-split has several advantages:

Figure 2 :

 2 Figure 2: Boundary curves incident to a mesh vertex. The control points b i 0 are fixed by the interpolation condition.

Figure 3 :

 3 Figure 3: Cross-boundary tangents.

Figure 4 :

 4 Figure 4: C 1 -conditions between two adjacent quintic Bézier patches.

Figure 5 : G 1

 51 Figure 5: G 1 continuity between two adjacent patches.

  ), t 2 a b] be a cubic Bézier curve.Its derivatives are given by x r (t) n;r i (t). At the end point t = a one obtains for the first and second derivativesx 0 (a) = 3 b;a (b1 ; b0)x 00 (a) = 6 (b;a) 2 (b2 ; 2b1 + b0):

  follows that the vector r 1 should lie in the kernel of P and in the image space of T . P and T are cyclic matrices [DAV79]. The kernel of P is spanned by two vectors

  k b lie also in the image space of T , because two vectors k a and k b exist such that T k a = k a and T k b = k b where

Figure 6 :

 6 Figure 6: Geometrical interpretation of the vectors a and b.

kr 1 ( 1 Figure 7 :

 117 Figure 7: "Optimal" first derivatives.

Figure 8 :

 8 Figure 8: Target point and target vector.

Figure 9 :

 9 Figure 9: Torus mesh, left with bad target points, right with correct target points.

Figure 10 :

 10 Figure 10: Polyhedral mesh -smooth surface with boundary curve network.

Figure 11 :

 11 Figure 11: control polygon -smooth surface with boundary curves.

Figure 12 :

 12 Figure 12:Figure 13: Top to bottom: polyhedral mesh, smooth surface, control polygon.

Figure 14 :

 14 Figure 14: Polyhedral mesh -smooth surface -smooth surface with boundary curves.