"MAKING ROOM": A PRELIMINARY STAGE TOWARD IMPROVING PRISON AND CORRECTIONAL MANAGEMENT? REPRESENTATIONS AND MANAGEMENT OF PRISON OVERCROWDING IN FRANCE AND CANADA

Keywords:

published or not. The documents may come

followed from the French media coverage and presence in the political debate of prison overcrowding, particularly after two parliamentary inquiry committees pointed to a "discourse of shame" in 2000.

Of course, the French prison administration had for decades deplored the lack of resources responsible for the dilapidation of buildings and the chronic shortage of beds in prisons in numerous reports. At the same time, international bodies had regularly emphasized the disastrous consequences of prison overcrowding in terms of inmates' fundamental rights. Increasing the number of prisons became a preferred -if not the sole response to the evils of the prison institution even if it failed to resolve the overcrowding issue. Why build more prisons? Is this building policy related to the objective situation of prisons? Is it possible to do without new prisons or to delay building them?

In a former study, (3) we showed that other contemporary societies made opposite choices. In Canada, for example, increasing the number of prison beds was not a public policy priority until the 2010s. The first Canadian report on prison overcrowding was published in May 2014. In that report, the Auditor General of Canada aimed to "determine whether the Correctional Service of Canada (CSC)(4) had increased the capacity of its correctional facilities in a manner that met its needs and was costeffective." (5) He reported the installation of more than 2,500 beds since 2009. That year, CSC had anticipated that "changes in criminal justice legislation would result in longer sentences for many offenders, leading to an increased offender population." However, the offender population did not grow as much as expected. This at least had the virtue to bring CSC to update its forecast on its institutional capacity, something it had not done in the last thirty years despite a rapid population growth rate similar to that observed in France, particularly in the 1990s. In short, moderation prevailed on both federal and provincial levels, with the exception of Ontario in the early 2000s. Our hypothesis is that the way in which the political and administrative elites perceived, and more generally, dealt with the prison overcrowding issue partly explains the sharp contrast between the two countries' political choices.

CONSTRUCTION OF "PRISON OVERCROWDING" INTO A PUBLIC ACTION ISSUE

Starting in the 1980s, prison overcrowding was gradually turned into a genuine public action issue by French governments, becoming a favorite justification ground for the adoption of prison building legislation. There was fear that worsened detention conditions due to overpopulation might call into question the reputation of the "country of human rights." Increasing media coverage of monitoring bodies of prison conditions -the European Council, the comptroller general of custodial centers, NGOs -made things worse. Refusal to increase the capacity of correctional facilities in a situation of overcrowding would then amount to a denial of humanity: the modernization of facilities offered a way out of this shameful situation.

Without negating the reality of prison overcrowding or the fact that historically, the needs of the French prison system were not sufficiently addressed, one has to admit that the overcrowding argument gradually became an essential part of the discourse of successive French governments. However, one is stricken by the fact that the moment when prison population growth was at its highest (the 1980s and 1990s) does not coincide with the moment when the term "overcrowding" gained greater prominence in national debate (the 2000s). In 1986, promiscuity was already one of the explicit reasons given by the Ministry of Justice to justify the "13,000 Program" (6) -which shows that the authorities had taken up overcrowding to turn it into a political issue long before its dissemination in public space. At the time, the term "overcrowding" was not widely used in the Ministry's discourse.

Thus, there was a significant time lag between the reality of the situation and the visibility it acquired in political discourse and preoccupations. Overcrowding was brought to the fore in 2000, when it was largely highlighted by two inquiry committees and set off a wave of emotion among deputies. The increased number of parliamentary items dealing with it at the French National Assembly attests to that: from none before 1985, they passed to three during the eighth Parliament (1986)(1987)(1988), which marked the launch of the first important building program, and to more than a hundred during the twelfth Parliament (2002)(2003)(2004)(2005), whereas the number of missing beds had been considerably reduced. In 1987, for instance, the number of inmates exceeded the capacity of correctional facilities by 17,000, and so did it in 2000.

In Canada, however, the prison overcrowding issue was often underestimated or even put into perspective by authorities and caused less sensation than in France. Does this mean that Canadian correctional institutions did not exceed their rated capacities? Far from it. In the mid 1990s, overcrowding and its consequences were exposed in the correctional investigator's public reports, but policymakers did not take up the issue. Institutional archives documents show that provincial correctional administrations, particularly the most populous ones (Ontario, Quebec, and British Columbia) also coped with the problem. However, they did not turn overcrowding into a first-order public action issue, judging that it would be short-lived. The same was true at the federal level, where the word "overcrowding" seldom appeared in the political debate, emphasis being rather on "shared accommodation" and "double bunking" -to use CSC's euphemistic but eloquent expressions. ( 7)

The fact that prison overcrowding did not concern all the correctional facilities of the Canadian territory may explain why political attention on the issue was limited. Its existence may also have been called into question by civil service officers -who did not submit it as evidence to apply for funding -or by association representatives. The latter considered that while new cells would relieve an a priori short-lived pressure, they would also be automatically and lastingly filled, and this would result in a de facto growing offender population. Another reason is the policymakers' representation of overcrowding and ability to solve it. While French political and administrative authorities viewed it as inexorable, their Canadian counterparts considered it a transitory phenomenon, an "accident" in Canadian penal history rather than something inevitable, thinking that the correctional system had its own relief mechanisms in temporary absences (8) and parole.

OVERCROWDING VS. OVER-SENTENCING

Another issue, over-sentencing, caught the attention of Canadian policymakers to such a point that it warranted priority public action. Indeed, the Canadian justice system was widely thought to be excessively punitive. The overrepresentation of certain specific categories of the population (women, juveniles, ethnic minorities) led Canadian policymakers across the political spectrum to realize that incarceration rates were inordinately high. At the federal level, even ministers who were the most conservative about criminal justice issues were compelled publicly to lament the overrepresentation of Aboriginals in prisons -their incarceration rate is six times (7) "Shared accommodation cells" are designed for two inmates. "Double bunking" is the practice of holding two offenders in a cell designed for one. They are contrary to the intent of CSC policy.

(8) Temporary absences "provide inmates with opportunities to access the community for medical, administrative, parental responsibility, compassionate reasons, community service, family contact, personal development, and for rehabilitative purposes," http://www. csc-scc.gc.ca/lois-et-reglements/710-3-cd-eng.shtml, retrieved November 1, 2016. that of the national average -and to propose measures addressing a situation unacceptable in a society so deeply concerned with the rule of law.

It must be said that starting in the 1990s, many actors from the academic world, organizations, or independent administrative authorities, became "civic leaders" who brought to the attention of government and parliament officials the socially unjust character of problems engendered by over-incarceration rather than life conditions in overcrowded facilities. The argument on over-sentencing was twofold. On the one hand, the development of a "differentiated" statistical approach to the prison population had shown that the penal system imprisoned a disproportionate number of certain offender types, which could damage the reputation of a system that had been striving for justice and equity, particularly since the 1982 implementation of the Canadian Charter of Rights and Freedoms. On the other hand, it turned out that the incarceration rate in Canada exceeded that of other western countries, except for the US, whose penal policy acted as a foil to Canadian political and administrative officials. For these reasons, Canadian political authorities chose not to settle for informal and piecemeal measures (9) to relieve prisons and adopted substantial legislative reforms which included revamping of the legislation on juvenile delinquents.

In France, however, increasing concern about delinquency and security in public debate during the same years limited political authorities' ability to talk about reducing incarceration rates. The idea that overcrowding could result from over-sentencing went unheeded. Organizations involved in the protection of the rights of inmates were reluctant to lose focus by taking a stand on the questions of reductionism or abolitionism. As a matter of fact, political measures toward reduction only concerned remand.

While it is true that generally and historically, Canadian correctional institutions have less qualitative inadequacies than their French counterparts, the fact is that both countries fundamentally differed from each other on the perception and choice of which problem to solve -this explains why different policies were implemented. While French elites often tend to consider that the improvement of detention conditions is inseparable from the state of infrastructures, in Canada, improvement is achieved less by material modernization than by access to rights.

  These practices include placing mattresses on the ground or accommodating supernumerary inmates in specialty cells -segregation, medical, or observation cells. TEPSIS PAPER 17 / SEPTEMBER 2016 Laboratoire d'Excellence TEPSIS, 190-198 av. de France, 75013 Paris, https://tepsis.fr

(6) The 13,000 Program [programme 13 000] is the name of a real estate project instituted by the law 87-432 of June 22, 1987 on the prison public service.

Laboratoire d'Excellence TEPSIS, 190-198 av. de France, 75013 Paris, https://tepsis.fr


Romuald Josserand

GALILEO IN THE JUVENILE SECTION: A PRISON TEACHER'S PERSPECTIVE

Keywords:

School attendance at the penitentiary of Varennes le Grand is mandatory for all juveniles. In this prison environment which often interferes with education, our task is to change the binding and discontinuous nature of the time spent in class into a constructive and linear moment conducive to learning and providing each pupil with clear objectives: pursuing a training course, preparing for an exam, changing one's attitude towards knowledge and school, getting high school credits.

Classes consist of groups of three to four juveniles. Pupils have different levels, with most of them being dropout or special needs students. We try to tailor learning objectives to individual needs without neglecting group work. The structure of the course is modular: a stable group of pupils share common or partially common objectives over a short, specific time. Contents are interdisciplinary.

USING SCIENTIFIC INQUIRY TO INTERROGATE THE WORLD

In all the disciplinary fields that we cover, the pupils regularly invoke the God hypothesis to explain the world. In the wake of the tragic events of January and November 2015, we consider it fundamental to conduct educational sequences leading the pupils to observe the world without exclusively resorting to this hypothesis.

During scientific modules that are short and organized on a just-in-time basis, structure time and address the permanence of humanity's great questions, we hypothesize that introducing complex experimentation and the epistemology of science as a remedy would allow pupils to construct their own knowledge without resorting to textbooks or prerequisites that they would have to accept without knowing their justification.

Using Galileo's works and the epistemological break they constitute, we bring the pupils to conduct hands-on experiments, isolate parameters in order to interpret their results, and emit, verify, or invalidate hypotheses. In this way, they come to observe and test phenomena that radically contradict their representations and are totally counter-intuitive: they thus come to acknowledge the fact that all objects freefall at the same speed regardless of their mass. The pupils design their experiments in writing or using a diagram and systematically describe the results in individual essays. Writing engages them in the work at hand and allows them to get some distance from their manipulations and formalize their results. This is necessary for an eventual transfer.

At the outset, pupils must discover, just like Galileo did, that the period of oscillation of a simple pendulum solely depends on the length of the cord to which the object is attached, regardless of the mass of the object. But does this experiment actually work? Can one reach clear-cut conclusions through manipulation? And if the oscillations of a pendulum do not depend on mass, would the same object, if detached from the cord, free-fall at a speed independent of mass?

The pupils of the penitentiary are not ready to abandon their representations of the world, particularly if it means accepting counter-intuitive physical phenomena. Nor are they ready for concessions regarding experimental facts. Going back to the works of Galileo, who did not compromise on facts either, helps us out of this predicament. We can use his observation notes and discover that it is possible to solve the question with a thought experiment, an experiment that is as structured and dependable as real ones, but conducted in the imagination, and that can thus be realized an infinite number of times, anywhere. Why not on the Moon? Galileo did it, and reached major conclusions. All objects free-fall at the same speed regardless of their mass. The American astronauts of the Apollo 15 mission reproduced it on the Moon. "For real?" the pupils ask. We can then watch a video realized in 1971. (1) The astronauts' conclusion, 450 years later, is final: "Galileo was correct!" The pupils realize that one can ask the same questions as Galileo did five centuries ago, to devise similar experimental protocols, obtain the same results, and experience the same doubts. This stability of phenomena over time, from the Earth to the Moon, gives them a feeling of security. They experience this permanence because they have built it through their experiments.

Throughout these experiments, the epistemology of science encourages the pupils to observe and interrogate the world without having recourse to the God hypothesis. Even so, the purpose is not to organize into a hierarchy or set against each other two ways of apprehending the world, but to allow oneself to interrogate the world without resorting to God during an experiment in which everything can become an object of research. This, in fact, is what surprises the pupils most: one can and is allowed to seek to explain any "scientific phenomenon," and if one does not reach understanding, one can always rely on someone else to reach an explanation. However, if one uses the God hypothesis, one cuts the inquiry short. For this reason, science does not allow itself to invoke God to explain a phenomenon.

INTERROGATING THE BELIEF/ KNOWLEDGE DISTINCTION…

Here, we directly deal with the epistemology of science with the pupils, justifying the scientific community's methodological and systematic refusal of recourse to the God hypothesis by citing four essential characteristics of scientific knowledge. It is important not to approach and address the legitimacy of belief.

First, scientific knowledge has this specificity that it is refutable. As Karl Popper wrote, "a theory that is not refutable by any conceivable event is non-scientific."

(2) According to Popper, this criterion helps differentiate scientific knowledge from nonscientific discourse. Our task is to make the pupils understand how scientific knowledge can be intellectually structured. In order to remain valid and structure reasoning, it must be liable to be called into question and discussed at any time. This implies, of course, that it can, as a last resort, be proved wrong. This is referred to as falsifiability. Beyond the epistemology of science, this scientific character criterion enables one to apply reasoning when analyzing everyday discourse and distinguish between what pertains to common thought and what is a refutable proposition.

Next, in the modern period -which is traditionally thought to start with Galileo's works -experimental science has the particularity of being elaborated on the basis of observation and analysis of reality. In order to produce a verifiable, and consequently refutable, knowledge, scholars must elaborate valid and thorough hypotheses on the modes of study of phenomena before attempting any experimental exploration. Once these hypotheses have been discussed and sorted out, the pupils, just like scholars before them, can devise how they will set up structured experiments in order to test their relevance.

Scientific experimentation is reproducible. This is the third characteristic of the elaboration of scientific knowledge. An experimentation whose results cannot be replicated should then be excluded from scientific inquiry. Concretely, the same experiment conducted by two different groups of pupils must yield identical results.

Comparison of these results can then necessitate that groups negotiate in order to agree on what has actually been tested. This experimentation should also yield the same results as those conducted by experimentalists remote from them in space or time. This criterion establishes that knowledge constructed in such a way is endowed with spatial and temporal stability which, in the future, will permit to draw upon it to predict the behavior of objects submitted to the same phenomenon.

Here, we touch on the fourth characteristic of scientific knowledge, prediction of the behavior of an object submitted to known influences. It differentiates scientific knowledge from traditional knowledge.

The completion of processes allowing the pupils to appropriate these epistemological characteristics takes part in the development of a scientific and reasoned approach, since the stability and predictability of the covered scientific content make this approach relevant to them.

Recourse to the history of science and the use of narratives also lead our pupils to feel the necessity of interrogating written sources and differentiating between myth, primary and secondary sources. Last, the history of science allows us to make a pedagogical detour: Galileo stood against the God of the Roman Apostolic Church 450 years ago; this allows us not to address the existence of a God present in our pupils' daily culture frontally.

A SECULAR PEDAGOGY

It seems to us that this educational sequence on Galileo's discoveries engages our pupils in acquiring an inquisitive approach and free will. In February 2016, a pupil told us: "My dad used to say there are things you know and things you must learn." Scientific inquiry takes our pupils on the complex path to distinction between belief and knowledge. It allows them at least to bring that question which has no definite answer to the foreground. We do not seek the stability and permanence of the answer, but that of the inquiry, regardless of the environment. Instilling critical thinking, a scientific mind, a whole new attitude in our pupils on a long-term basis offers a possible way towards the construction of secularism and citizenship.


Cécile Barbachoux

Jérome Gariel

Gérard

G Le Denmat

Cosmological energy in a thermo-horizon and the first law

may come

Cosmological energy in a thermo-horizon and the first law

Introduction

The generalization of the Thermodynamics of black holes (BH) to cosmological horizons represents an important stake to understand different issues in cosmology such as the nature of the dark energy (DE) in relation with the problems of the cosmological constant (CC) and of the vacuum energy, the acceleration of the present universe, the coincidence problem and the early inflation.

This generalization was first introduced for de Sitter spacetime [1]. Thereafter, it was tentatively extended to quasi-de Sitter FRW spacetimes in different frameworks (see [2]- [6]). In an interesting approach, Bousso [7,8] considers the flow of energy through the horizon as a null surface. He interprets the variation of the entropy of the horizon through the variation of its surface as the response of the horizon to the flux of energy, in the same way as the "first law" of the BH.

Following this approach, several authors (see for example [9,10]) have estimated that the apparent horizon (a.h.) is the only limit surface (excluding other horizons such as the event horizon (e.h.)) having coherent thermodynamical properties to address problems such as the nature of the DE.

Our main goal is to shed some light on the contradiction between the amount of energy calculated from the first law as defined in [7,8] and the definition of the energy contained inside the horizon, independently of the choice of the thermo-horizon (t.h.).

We restrict our study to a spatially flat FRW spacetime, which is the starting point of other studies (non spatially flat spacetimes, cases with interactions,...).

After a brief review of the definition of a t.h. in a Q-space introduced in [7,8], we show that any t.h. obeys the second law (Section 2). In section 3, we present the contradiction between the amount of energy derived from the first law and the definition of the energy inside the horizon. We then show that this contradiction is resolved in a thermodynamical model for a DE [4,5] based on the e.h. (Section 4).

Definition of a thermo-horizon

In a spatially flat FRW spacetime

ds 2 = -dt 2 + a(t) 2 (dr 2 + r 2 dΩ 2 ), (1)
the dynamical evolution of the scale factor a(t) is given for a perfect fluid with energy density ρ and pressure P by

ȧ a 2 = H 2 = χ ρ 3 , (2)
ä a = - χ 6 (ρ + 3P ), (3)
where χ = 8π is the Einstein constant, with G = 1 and c = 1. The equation of state (EoS) ω of the fluid is given by P = ωρ and we introduce the parameter ε = 3 2 (1 + ω). In the following, we restrict our study to the Q-space [8], namely accelerated universes, for which 0 < ε < 1. Extending the reasoning of [8], we consider an horizon (null surface) with a given radius L. According to the first law, the flow of energy through this surface is given by

- • E = 4πL 2 ρ(1 + ω) = T • S. (4)
We assume that we can associate a temperature T and an entropy S to the dynamical horizon of radius L, given by the relations of Bekenstein-Hawking for a BH or a de Sitter horizon

T = 1 2πL
, and S = πL 2 .

(

) 5
Any horizon of radius L with a temperature and a entropy given by ( 5) and which obeys the first law ( 4) is called thermo-horizon (t.h.).

Using (5), we obtain directly T

• S =

• L and Eq. ( 4) becomes

εL 2 H 2 = • L. (6)
2

With our notations, the Eq. ( 3) is given by

• 1 H = ε, (7)
and the first law (6) rewrites

• H = • ( 1 L ). ( 8
)
After integration, this equation leads to

HL -1 = CL, (9)
where C is a constant. Eq. ( 9) establishes a general relation between the t.h.

L and the a.h. R A = 1 H which is satisfied by any thermo-horizon of radius L without restriction on ε (in particular without assuming ε = constant). With the constant C, this relation is more general than the Eq. ( 28) of [8].

If L is the a.h., then

L = 1 H = R A , implying C = 0. Conversely, only C = 0 leads to L = 1 H
. Therefore the a.h. obeys the first law (4) if and only if C = 0.

This special case only is considered by [8]. More generally, any horizon L defined by ( 9) with a temperature and an entropy given by ( 5) is a t.h. and it obeys the first law (4).

The equation ( 6) can be rewritten with the help of ( 9)

• L = ε(1 + CL) 2 . ( 10
)
Any t.h. verifies this equation. Using Eq. ( 10), L is strictly increasing in the Q-space (accelerated universe) where 0 < ε < 1. The same result can be derived for the entropy S given by (5).

3 Energy in a thermo-horizon and the first law

On one side, the total amount of energy contained inside the a.h. for a spatially flat FW space-time is (e.g. [9] before Eq. ( 20))

E = ρ 4π 3 R 3 A = R A 2 . ( 11
)
Let us remark that in [10] this relation is used for a non-spatially flat FRW spacetime, albeit no more valid in this case.

On the other side, using ( 4) and ( 5), the first law applied to the a.h. considered as a t.h. leads to -

• E = • R A , (12)
where -

• E is the total amount of energy crossing the a.h. by unit of time. According to the conservation of the energy, this amount of energy is equal to the variation of the total energy (11) per unit time,

• E = • R A 2
. This result is in contradiction with (12) except when R A is constant, which corresponds to a de Sitter spacetime where the a.h. identifies with the e.h.. This result is not restricted to the a.h. and can be extended to any t.h.. Using (7), the left hand side of (4) becomes for a t.h. of radius L -

• E = L 2 H 2 ε = -L 2 • H, (13)
while the total energy inside the horizon is

E = 1 2 H 2 L 3 . (14)
Differentiating ( 14) and equating with (13), we obtain with (8)

3 2 (HL) 2 -HL + 1 = 0, (15)
where • H = 0 has been assumed. No real root can be found for this equation.

In particular, L = 1 H is not a solution. Therefore, the above contradiction can only be removed for

• H = 0, namely for a de Sitter spacetime.

4 Thermodynamical model of the event horizon for the dark energy

The preceding results are independent of the underlying model for the DE. They depend only on the assumptions of the existence of a temperature and an entropy associated to a t.h. through the relation ( 5) and of the validity of the extension of the first law of BHs (4) to cosmological t.hs.. With these assumptions, we obtain (9) (assuming C = 0), as demonstrated by [7,8] and by [9], in sections II-A and II-B.

In section II-B of [9], the authors consider only the specific model for the DE developed in [6]. Let us emphasize that their results can be obtained independently of any model for the DE (see Section 2) because the demonstration involves only the density of the total energy ρ. Consequently, the reasoning developed in [9] cannot question the validity of the model assumed for the component DE and in particular the approach proposed in [6]. The first law ( 4) is a relation between the density of energy ρ and the entropy. It does not involve the density of energy of the DE ρ Λ and therefore cannot be used to discuss or refute its expression.

In [9], the authors emphasize the apparent discrepancy between the horizon chosen as IR cut-off in the expression of ρ Λ in the holographic model of Li [6], which is the e.h. r χρ Λ = 3c 2 r 2 , (16) and the t.h. L which must be the a.h. R A in order to satisfy the first law. They suggest that the e.h. r should be identical to the t.h. L obeying the first law. So, they implicitly assume that in any holographic model for ρ Λ , i.e. χρ Λ = 3/l 2 , we have to choose for cut-off l the same horizon as the t.h. L obeying the first law.

Let us discuss the full consequences of the assumption l = L. On one hand, by setting L = l = r in the expression of the flux of energy in the first law (4), the first law is no more satisfied (see Section 3). On the other hand, if we choose L = l = R A , we obtain for the holographic model:

χρ Λ = 3 R 2 A = 3H 2 . ( 17
)
which with (2) implies ρ = ρ Λ , excluding any other contribution to the total energy (dark matter (DM), dust, radiation,...) in contradiction with the present observations where the DM takes a non negligible part (about 1/3) of the energy of the universe. Two points of view can be followed to solve this dilemma:

i) First, we can choose to preserve the first law ( 4) and propose an holographic model of the form (17) albeit not compatible with the observations [6,12].

ii) Secondly, we can consider the holographic model of the DE (16) compatible with the present observations (acceleration and EoS today ω Λ ≃ -1 [5,6,11]) and modify the first law (4) to be compatible with the chosen model for the DE.

Because compatible with the observational features, the second alternative is more reasonable . This leads naturally to question the first law (4) which seems to fail because, as seen in Section 3, it contradicts the definition (11) of the total energy in the horizon. The first law must be modified in order to include the DE through a model linking the DE with the chosen horizon. To achieve this goal, we propose a Gibbs' equation describing the thermodynamics of the component DE instead of the first law (4). This approach was introduced in the model [4,5], where a DE component with an energy density of type (16) (with c = 1) and an EoS of vacuum-type were considered

Λ χ = ρ Λ = -P Λ = 3/χr 2 , (18)
with r the radius of the e.h.. In this model, the Gibbs' equation relative to this component DE is given at the specific level by

T Λ ds Λ = dε Λ + P Λ dv Λ , (19)
where ε Λ ≡ Λ n Λ χ is the specific energy, v Λ ≡ 1 n Λ the specific volume and n Λ the number density. With the EoS (18), the Eq. ( 19) becomes

• Λ n = χT Λ • S Λ = - • r, (20)
with

T Λ = 1 2πr
. After integration, we obtain for the specific entropy

S Λ = -πr 2 + K, ( 21
)
where K is a constant. In [4,5] the local equation of the conservation of the energy u β ∇ α T αβ = 0 is used instead of the expression of the global energy E and without assuming any peculiar expression for the entropy.

Considering the global energy of the DE component inside the horizon as in [9,10],

E Λ = ρ Λ 4π 3 r 3 = r 2
, the Gibbs'equation associated to this component becomes

dE Λ = T Λ dS Λ -P Λ dV Λ , or • E Λ = T Λ • S Λ -P Λ • V Λ . (22)
With ( 18) and ( 21), this Eq. leads to

• E Λ = • r 2
, in full agreement with the previous definition of E Λ = r 2 . The previous results show that the contradiction described in Section 3 is essentially related to the mainly questionable assumption of the validity of the expression of the static entropy S = πL 2 of the BH for the cosmological horizon [7,8]. On the contrary, our expression (21) of S Λ is not postulated but deduced from the Gibbs' equation for the model (18) of the DE [4,5], which is supported by the holographic approach [6].

By comparison of the two Eqs. ( 5) and ( 21), let us first note the difference of sign in the expression of the entropy. This difference is in accordance with the predictions because in the first case, the observer is outside the limit-surface of the BH and looses information, while in the second case (cosmological case), the observer is inside the limit-surface. This difference of sign is pointed out in [7,8] where it is attributed to the exchange of energy -dE rather than to T dS, while we interpret it as a lost (or a gain) of information from the limit-surface.

Secondly, let us remark the presence of a pressure term in the Gibbs' equation ( 22) albeit not appearing in the first law (4). This term takes into account the fact that the DE has a (negative) pressure, and that the transformation is not isochore because the surface varies. This term is responsible for the sign "minus" in the expression of the entropy Eq. ( 21). This point relates the two preceding remarks and strengthens the consistency of this approach. The existence of such a term was recently considered and discussed in another context for the first law of the BH by [13].

Conclusion

By a general demonstration independent of the underlying model for the DE, we show in this article a contradiction (except in de Sitter) between the first law introduced in [7,8] for the thermo-horizon that leads to To solve this contradiction, we propose to replace the first law by a Gibbs' equation for the DE component, which is naturally associated to the e.h. in the model [4,5], later supported by an holographic model [6] in an independent approach.

  the definition of the total energy in the horizon E = R A 2 .


Armine Karami
email: armine.karami@lip6.fr

Dimitri Galayko

Philippe Basset

teaching and research institutions in France or abroad, or from public or private research centers.

Characterization of the capacitance variation of electrostatic vibration energy harvesters biased following rectangular charge-voltage diagrams

Characterization of the capacitance variation of electrostatic vibration energy harvesters biased following rectangular charge-voltage diagrams 1. Introduction Recent research on micro-machined electrostatic vibration energy harvesters (e-VEHs) have shown that the electrical conditioning and load interfacing part of such systems is critical to the performances in terms of harvested power. In the light of this consideration, researchers have started to investigate new types of conditioning circuits [START_REF] Queiroz | IEEE 5th Latin American Symposium on Circuits and Systems[END_REF][START_REF] Karami | [END_REF][START_REF] Wei | Simple and efficient interface circuit for vibration electrostatic energy harvesters[END_REF]. These new conditioning circuits have in common that their biasing of the transducer, across one cycle of its capacitance variation, can be summarized by a rectangular charge-voltage characteristic diagram (QV diagram) [4].

However, the way this type of biasing affects the device's dynamics via the electromechanical coupling effect is still unclear: although semi-analytical approximations are possible for a limited range of input conditions [5], no measurement method allowed to measure the electromechanical coupling effect on the dynamics of a given e-VEH. In particular no method was reported for the measurement of the extremal values of capacitance variation for a device submitted to harmonic input excitation and biased following a rectangular QV diagram.

In this paper, a simple method is proposed to carry out such a measurement. This method only involves the use of characterization equipment that is typically found in MEMS e-VEHs characterization labs, i.e., generic electronics components, characterization gear, and a mechanical shaker. This method is based on the dynamics of two simple circuits implementing biasing schemes described by rectangular QV diagrams. After a brief presentation and analysis of these circuits, the method is presented and carried out as an example on a state-of-the-art micro-machined, symmetrical gap-closing geometry e-VEH. Consider the two conditioning circuits depicted in Fig. 1 and Fig. 2. The voltage source elements E 1 and E 2 model external voltage sources, eventually superposed with the effect of electret charging. The QV diagrams of the corresponding circuits are depicted in Fig. 3 and Fig. 4. The diodes are supposed to follow an ideal diode voltage-current law, with threshold voltage V T ≥ 0. Defining a variation cycle by the transducer's variable capacitance varying from a maximal value C max to a minimal value C min and then from C min to C max , it comes for the local evolution laws of each of the circuits that:

V #1 n+1 = C S V #1 n + (C max -C min )E 1 -(C max + C min )V T C min + C S , (1)
V #2 n+1 = C S V #HW 2 n + (C max -C min )E 2 -(C max + C min )V T C max + C S , (2)
where

V #1 n and V #2
n denote the value of the voltage across C S at cycle index n for the circuits depicted in Fig. 1 and Fig. 2, respectively.

Presentation of the measurement method

From the preceding equations, denoting ∆

1 = V #1 n+1 -V #1 n , ∆ 2 = V #1 n+1 -V #1 n , V 1 = V #1 n and V 2 = V #2
n , it comes:

C min = C S ∆ 2 (E 1 -V T ) -∆ 1 (E 2 -V T -∆ 2 -V 2 ) (E 1 + ∆ 1 + V 1 + V T )(E 2 -∆ 2 -V 2 -V T ) -(E 1 -V T )(E 2 + V T ) , (3)
C max = C S ∆ 2 + C min (E 2 + V T ) E 2 -V 2 -∆ 2 . ( 4
)
Figure 5. 3D schematic of the MEMS transducer that is characterized as an example of application of the presented measurement method. Let's consider an e-VEH for which the value of maximal capacitance C max is to be measured, when the device is submitted to an harmonic external excitation, and when its biasing across one cycle of variation is described by a rectangular charge-voltage diagram of extremal voltages V L and V R (see Fig. 6). Suppose the slope of the non-horizontal portions of the QV cycles of the two conditioning circuits is large compared to C max , i.e., C S C max (see Fig. 3 and Fig. 4). By choosing the values of the external voltages sources and initial voltages across C S such that:

E 1 = V L , (5)
E 2 = V R , (6) V 1 = V 2 = V R -V L , (7)
it comes that evaluating (3) and ( 4) with the values of ∆ 1 and ∆ 2 measured with these parameters yields the values of C max and C min as affected by the electromechanical coupling induced by the rectangular biasing depicted in Fig. 6.

Note that the values of C max and C min result from a dynamical process, and hence the measurement method is accurate only if the electromechanical coupling impact on the e-VEH's dynamics can be considered as a quasi-static process: otherwise, a systematic error is introduced in the measurement. This hypothesis is verified if C S is chosen such that the values of ∆ 1 and ∆ 2 do not change abruptly one cycle to the next, i.e., that C S is chosen large enough. Also, this results in the QV diagrams in Fig. 3 and Fig. 4 to become more alike, which improves the measurement accuracy. However, a too large value of C S results in ∆ 1 and ∆ 2 becoming indiscernible one from another because of a non-null noise floor due to random noise sources. Thus, prior to the measurement, a study should be done to chose the optimal value for C S given those constraints.

Example of transducer characterization using the presented method

As an illustrating example, the following describes the application of the method to the measurement of the value of C max for a practical transducer, that has been previously reported in [6]. The schematic of the transducer is depicted in Fig. 5. The measurement is carried out for the transducer biased following a rectangular QV diagram, with V L = 10V and V R = 15V, and mechanically forced by an external excitation of 1.5g amplitude, 150Hz frequency. Because of the symmetrical gap-closing geometry of the device, the value of C min is fixed throughout the e-VEH's operation, and is measured by standard capacitance measurement means, as C min = 55pF. In both circuits, the fixed capacitor is chosen such that C S = 500pF.

Because C min is fixed, only the circuit in Fig. 2 is necessary. The parameters for the measurements are E 2 = 15V, and V 2 = 5V (see ( 6) and ( 7)). Several measurements are carried out for the characterization, with approximatively equal values of V 2 , to statistically estimate the uncertainty due to type A error sources on the measurements. With the obtained results for (∆ 2 , V 2 ), applying (4) yields a value of C max = 94.7pF ± 2.8pF. Equation (3) can be used to check the results as C min is fixed, by carrying out the same measurement with the circuit in Fig. 1, choosing E 1 = 10V, and at V 1 = 5V (see ( 5) and ( 7)). If the geometry of the device is such that C min is not unalterable, this latter measurement has to be carried out to measure the value of C min as impacted by the electromechanical coupling resulting from the bias depicted in Fig. 6.

Conclusion

This work presented the first method to measure the capacitance variation of electrostatic MEMS transducers of an e-VEHs, when electrically conditioned following a rectangular charge-voltage diagram. This enables more thorough studies on the effect of the electromechanical coupling on the dynamics of e-VEH using electrical interfaces implementing these types of conditionings schemes, such as charge-pump conditioning circuits.

Figure 1 .

 1 Figure 1. First conditioning circuit used for the characterization of C max and C min .

Figure 2 .Figure 4 .

 24 Figure 2. Second conditioning circuit used for the characterization of C max and C min .

Figure 6 .

 6 Figure 6. Generic rectangular QV diagram, with extreme voltages V L and V R .

Figure 7 .

 7 Figure7. Complete time-evolution of the voltage across C S for the circuit depicted in Fig.1, with E 1 = 10V, with the transducer submitted to input acceleration of 1.5g amplitude at 150Hz frequency. The quantities V 1 and ∆ 1 are highlighted.

Figure 8 .

 8 Figure8. Complete time-evolution of the voltage across C S for the circuit depicted in Fig.2, with E 2 = 15V, with the transducer submitted to input acceleration of 1.5g amplitude at 150Hz frequency. The quantities V 2 and ∆ 2 are highlighted.


J Wolfmann

FROM NEAR-BENT TO BENT: A SPECIAL CASE

Keywords:

We consider the construction of (2t)-bent functions from two (2t -1)-near-bent functions in a special case. We deduce new families of bent-functions.

Introduction

F 2 is the finite field of order 2. A m-Boolean function (or Boolean function in m dimensions) is a map F from F m 2 to F 2 . Bent functions are the Boolean functions whose Fourier coefficients have constant magnitude and were introduced by Rothaus in [6]. It means that a m-Boolean function F is bent if all its Fourier coefficients are in {-2 m/2 , 2 m/2 }. Since Fourier coefficient are in Z then bent functions in m dimensions exist only when m is even.

Bent functions are of interest for Coding Theory, Cryptology and well-correlated binary sequences. For example, they have the maximum Hamming distance to the set of affine Boolean functions. It is easy to prove that the set B(m) of m-bent functions is invariant under the action of the product of the affine linear group of F m 2 with the group of translations of affine Boolean functions. The corresponding partition of B(m) under the action of this group defines an equivalence.

Two main infinite families of bent functions are known (see [START_REF] Wolfmann | Bent Functions and Coding Theory in Difference Sets, Sequences and their Correlation properties[END_REF]) but in general it is very difficult to decide if any bent function is equivalent to a member of these families.

Bent functions have been the topic of a lot of works but the complete classification of bent functions is still open. In order to improve the knowledge on bent functions, it is convenient to find new properties and constructions. This is the goal of this work. By definition, a m-Boolean function F is near-bent (sometimes called semi-bent [START_REF] Khoo | Stinson A new characterisation of semi-bent and bent functions on finite fields Des[END_REF]) if all its Fourier coefficients are in {-2 (m+1)/2 , 0, 2 (m+1)/2 } Near-bent functions exist only when m is odd.

It is known that the restrictions of a (2t)-bent function to any hyperplan and to the complement of this hyperplan (view as (2t-1)-Booleans functions) are near-bent. Properties of these near-bent functions are investigated in [START_REF] Wolfmann | Cyclic Code Aspects of Bent Functions in Finite Fields Theory and Applications[END_REF] in connexion with the theory of cyclic codes. In the present paper we consider the question of the construction of (2t)-bent functions from two (2t -1)-near-bent functions. A first approach of this problem appears in [START_REF] Leander | Construction of Bent Functions from Near-Bent Functions Journal of Combinatorial Theory[END_REF].

We restrict the problem to the special situation where the sum of the two (2t -1)-near-bent functions is an affine linear form of

F 2t-1 2 .
The study in this case was initiated by Leander and McGuire in [START_REF] Leander | Construction of Bent Functions from Near-Bent Functions Journal of Combinatorial Theory[END_REF] where they present several properties and constructions. One of their results is Theorem 14 of the present paper. In particular they obtain a non-weakly-normal bent function proving that this special point of view could produce new classes of bent functions. In this work we introduce new properties of the (2t -1)-near-bent functions involved in the construction of bent functions and we deduce new infinite families of bent functions.

This paper is a continuation and a generalisation of [START_REF]Wolfmann Special Bent and Near-Bent Functions in Advances in Mathematics of Communication[END_REF]. Lemma 5, Lemma 6 and Theorem 7 are already in [START_REF]Wolfmann Special Bent and Near-Bent Functions in Advances in Mathematics of Communication[END_REF] and are recalled for convenience, Theorem 14 is in [START_REF] Leander | Construction of Bent Functions from Near-Bent Functions Journal of Combinatorial Theory[END_REF] while all the other Theorem, Proposition, Corollary, Lemma and the entire subsection 3.5 are new.

Preliminaries

Classical definitions and results.

We gather here some definition and well-known results that we will use in the paper. Properties (P 1 ) to (P 5 ) are classical and follow immediately from the definitions or from straightforward calculations. The distribution given in (P 6 ) is a special cases of Proposition 4 in [START_REF] Canteaut | Decomposing Bent Functions[END_REF]. See [START_REF] Canteaut | On Cryptographic Properties of the Cosets of R(1,m[END_REF], [START_REF] Canteaut | Decomposing Bent Functions[END_REF], [START_REF] Dillon | Elementary Hadamard Difference Sets Ph[END_REF], [START_REF] Wolfmann | Bent Functions and Coding Theory in Difference Sets, Sequences and their Correlation properties[END_REF] for details.

Boolean functions.

F 2 m is the finite field of order 2 m . A m-Boolean function is a map F from F m 2 to F 2 . Its weight is the number of X in F m 2 such that F (X) = 1 and is denoted by w(F ). If e ∈ F m
2 then the Derivative of F with respect to e is the m-Boolean function D e F defined by: D e F (X) = F (X) + F (X + e)

The Fourier transform (or Walsh transform) F of F is the map from F m 2 into Z defined by:

F (v) = X∈F 2 m (-1) F (X)+<v,X>
where <, > denotes any inner product of

F 2 m over F 2 . F (v) is called the Fourier coefficient of v. Notation: T v is the Boolean function defined by T v (X) =< v, X >.
It comes immediately:

(P 1 ) : F (v) = 2 m -2w(F + T v ).
Remark: the set of F (v) when v runs through F 2 m is independent of the choice of the inner product <, >. The dual F of a (2t)-bent function F is the (2t)-Boolean function F defined by: F (v) = (-1) F (v) 2 t where F is the Fourier transform of F . It is easy to prove that F is bent and that the dual of F is F .

(P 4 ): Let F be the dual of a (2t)-bent function F . Then:

F (v) = 1 if and only if F (v) = -2 t .
2.1.3. Near-bent functions.

A (2t-1)-Boolean function is "near-bent" if all its Fourier coefficients are in {-2 t , 0, 2 t }.

In the litterature, near-bent functions are sometimes called semi-bent functions (see [START_REF] Khoo | Stinson A new characterisation of semi-bent and bent functions on finite fields Des[END_REF]).

(P 5 ) Let f be a (2t -1)-boolean function and let l be an affine linear form of F 2t-1 2 . f is a near-bent function if and only if f + l is a near-bent function.

(P 6 ): The distribution of the Fourier coefficients of a (2t -1)-near bent function f is well known (see Proposition 4 in [START_REF] Canteaut | Decomposing Bent Functions[END_REF]). f (v) = 2 t number of v:

2 2t-3 + (-1) f (0) 2 t-2 f (v) = 0 number of v: 2 2t-2 f (v) = -2 t number of v: 2 2t-3 -(-1) f (0) 2 t-2 .

2.2.

A two-variable representation.

Special description of F 2t

2 . We identify F 2t 2 with the finite field F 2 2t and F 2 2t with:

F 2 2t-1 × F 2 = {X = (u, ν) | u ∈ F 2 2t-1 , ν ∈ F 2 }.
If m = 2t-1 the inner product used to calculate the Fourier coefficients is defined by < a, x >= tr(ax) where tr is the trace function of F 2 2t-1 . Notation: For every a ∈ F 2 2t-1 the (2t -1)-Boolean function t a is defined by t a (x) = tr(ax).

If m = 2t, a special inner product adapted to the above special description of F 2t 2 will be defined in 2.2.3. 2.2.2. Special representation of (2t)-Boolean functions. Using the description of

F 2 2t as F 2 2t-1 ×F 2 then a (2t)-Boolean function F now is a map from F 2 2t-1 × F 2 to F 2 . ∀(u, ν) ∈ F 2 2t-1 × F 2 : F (u, ν) = 0 or 1
Let F be such a function. Define two (2t -1)-Boolean functions f 0 and f 1 by f 0 (u) = F (u, 0) and f 1 (u) = F (u, 1). Now let φ be the (2t)-Boolean function defined by: φ(x, y)

= (y + 1)f 0 (x) + yf 1 (x). If ν = 0 then φ(u, 0) = f 0 (u) = F (u, 0). If ν = 1 then φ(u, 1) = f 1 (u) = F (u, 1). Therefore, for all (u, ν) then F (u, ν) = φ(u, ν) whence φ = F . F (x, y) = (y + 1)f 0 (x) + yf 1 (x)
This is the two-variable representation of F . The (2t)-Boolean function F is completely defined by the two (2t -1)-Boolean functions f 0 and f 1 .

Notation: F is denoted by [f 0 , f 1 ] f 0 and f 1 are called the components of F .
From the definitions of f 0 and f 1 we have:

(P 7 ) : w(F ) = w(f 0 ) + w(f 1 )
2.2.3. Representation of (2t)-linear forms. The purpose of this part is to express (2t)-linear forms and the inner product <, > over F 2 2t used in the calculation of the Fourier coefficients in such a way which is consistent with the identification of

F 2 2t as F 2 2t-1 × F 2 .
With this identification it is easy to check that <, > defined by < (a, η), (x, ν) >= tr(ax)+ην is an inner product of F 2 2t (non-degenerate symmetric bilinear form). We use it to calculate the Fourier coefficients.

Definition 1. The inner product T v of F 2 2t such that T v (X) =< v, X > is now defined by T (a,η) (x, ν) = tr(ax) + ην.
Consequently, every affine linear form of F 2 2t is of the kind T (a,η) + ω with (a, η) ∈ F 2 2t and ω ∈ F 2 . We immediately obtain :

Proposition 2. (*) T (a,η) = [t a , t a + η]. Let F = [f 0 , f 1 ] be a (2t)-boolean function. (**) F + T (a,η) = [f 0 + t a , f 1 + t a + η].

Representation of bent functions.

The next proposition is a special version of a well known result which appears in several papers ( [START_REF] Canteaut | Decomposing Bent Functions[END_REF], [START_REF] Wolfmann | Bent Functions and Coding Theory in Difference Sets, Sequences and their Correlation properties[END_REF],...). A proof is given in [START_REF]Wolfmann Special Bent and Near-Bent Functions in Advances in Mathematics of Communication[END_REF].

Problem:

The previous proposition leads to the following problem: construction of (2t)-bent-functions by means of (2t -1)-near-bent functions f 0 and f 1 satisfying (b).

In order to study this point we have an additional information on f 0 and f 1 .

Proposition 4. If [f 0 , f 1 ] is bent then f 0 + f 1 is balanced (considered as a Boolean function on F 2 2t-1 ). Proof. D 0,1 (F )(u, ν) = F (u, ν) + F (u, ν + 1).
The two-variable representation of F is:

F (x, y) = (y + 1)f 0 (x) + yf 1 (x). Hence F (u, ν) = (ν + 1)f 0 (u) + νf 1 (u) F (u, ν + 1) = (ν + 1 + 1)f 0 (u) + (ν + 1)f 1 (u). Finally D 0,1 (F )(u, ν) = f 0 (u) + f 1 (u).
On the other hand, we know from (P 2 ) that D 0,1 (F ) is balanced which means that the number of (u, ν)

such that D 0,1 (F )(u, ν) = 1 is 2 2t-1 . Now, note that D (0,1) (F )(u, 1) = f 0 (u) + f 1 (u) = D (0,1) (F )(u, 0). Since for a given u such that f 0 (u) + f 1 (u) = 1 there are two (u, ν) such that D (0,1) (F )(u, ν) = 1 then the number of u such that f 0 (u) + f 1 (u) = 1 is 1 2 2 2t-1 = 2 2t-2
and this proves that f 0 + f 1 is balanced. Special case: Because of the above proposition, a possible approach to attack our problem is to search f 0 and f 1 such that f 0 + f 1 is balanced. For every ∈ F 2 the Boolean function tr + where tr is the trace function of F 2 2t-1 over of F 2 , is balanced and is therefore a good candidate for f 0 + f 1 . This leads us to focus in this work on the case f 1 + f 0 = tr + .

Results in this case are given by McGuire and Leander in [START_REF] Leander | Construction of Bent Functions from Near-Bent Functions Journal of Combinatorial Theory[END_REF].

Results

Previous results.

The three next results have already been published in [START_REF]Wolfmann Special Bent and Near-Bent Functions in Advances in Mathematics of Communication[END_REF] with their proofs. We need them for the rest of the paper.

Lemma 5. Let F = [f 0 , f 1 ] be a (2t)-Boolean function. Then: a) F (u, 0) = f0 (u) + f1 (u). b) F (u, 1) = f0 (u) -f1 (u).
Proof. see [START_REF]Wolfmann Special Bent and Near-Bent Functions in Advances in Mathematics of Communication[END_REF], Lemma 13. Lemma 6.

If f 0 and f 1 are the components of a bent function F . and ω ∈ F 2 then:

D 1 f 0 = ω if and only if D 1 f 1 = ω + 1.
Proof. see [START_REF]Wolfmann Special Bent and Near-Bent Functions in Advances in Mathematics of Communication[END_REF], Proposition 16.

The next theorem is a fundamental result.

Theorem 7. Let f 0 be a (2t -1)-near-bent function. If the derivative D 1 f 0 is a constant function then the (2t)-Boolean function F = [f 0 , f 0 + tr] is bent.
Proof. see [START_REF]Wolfmann Special Bent and Near-Bent Functions in Advances in Mathematics of Communication[END_REF], Theorem 1.

3.2. New results.

3.2.1.

Results on near-bent functions.

Definition 8. If f is a (2t -1)-near-bent function then Îf is the indicator of the set {x ∈ F 2 2t-1 | f (x) = 0}
where f is the Fourier transform of f . (In other words, Îf (x) = 1 if and only if f (x) = 0).

Lemma 9.

Let f be a (2t -1)-near-bent function and ∈ F 2 .

If D 1 f = then Îf = tr + Remark: According to the definition of Îf this lemma means that if D 1 f = then f (x) = 0 if and only if tr(x) + = 1.

Proof. D 1 f = means that f (x + 1) = f (x) + . The transform τ :

x → x + 1 is a permutation of F 2 2t-1 and then preserves the weight of every (2t -1)-Boolean function. Thus:

{x | f (x) + tr(ux) = 1} = {x | f (x + 1) + tr(u(x + 1)) = 1}. (E) {x | f (x) + tr(ux) = 1} = {x | f (x) + + tr(ux) + tr(u) = 1}.
If tr(u) + = 1 the right hand member of (E) is:

{x | f (x) + tr(ux) = 0} = 2 2t-1 -{x | f (x) + tr(ux) = 1} Hence (E) becomes: {x | f (x) + tr(ux) = 1} = 2 2t-1 -{x | f (x) + tr(ux) = 1} In other words w(f + t u ) = 2 2t-1 -w(f + t u ) and thus: If tr(u) + = 1 then w(f + t u ) = 2 2t-2 which is equivalent to f (u) = 0.
For every the number of u such that tr(u) + = 1 is 2 2t-2 and this is also the number of u such that f (u) = 0 (see (P 6 ) ). Then, immediately: f (u) = 0 if and only if tr(u) + = 1. This means Îf 0 = tr + 3.2.2. Results on bent functions. The next theorem is a key point of this work.

Theorem 10. Let F = [f 0 , f 1 ] be a (2t)-bent function and let F = [ f0 , f1 ] be its dual function. Then f0 + f1 = Îf 0 Proof. From the definitions of bent and near-bent functions, every a in F 2 2t-1 belongs to one of the following sets:

A 1 = {a ∈ F 2 2t-1 | f0 (a) = -2 t and f1 (a) = 0} A 2 = {a ∈ F 2 2t-1 | f0 (a) = 0 and f1 (a) = -2 t } A 3 = {a ∈ F 2 2t-1 | f0 (a) = 2 t and f1 (a) = 0} A 4 = {a ∈ F 2 2t-1 | f0 (a) = 0 and f1 (a) = 2 t }
The definition of the dual of F induces that (a, η) is in the support of F if and only if F (a, η) = -2 t . From Lemma 5: F (u, 0) = f0 (u) + f1 (u) and F (u, 1) = f0 (u) -f1 (u).

Therefore: F (a, 0) = -2 t if and only if a ∈ A 1 or a ∈ A 2 , F (a, 1) = -2 t if and only if a ∈ A 1 or a ∈ A 4 , We deduce that, (a, 0) is in the support of F if and only if a ∈ A 1 ∪ A 2 . and (a, 1) is in the support of F if and only if a ∈ A 1 ∪ A 4 . In other words the support of f0 is A 1 ∪ A 2 and the support of f1 is

A 1 ∪ A 4 . It follows that the support of f0 + f1 is A 2 ∪ A 4 which is nothing but the set {x ∈ F 2 2t-1 | f0 (x) = 0}. A bent function F = [f 0 , f 1 ] is self-dual if F = F where F is the dual of F . Corollary 11. If a (2t)-bent function F = [f 0 , f 1 ] is self-dual then f 1 = f 0 + Îf 0 .
Proof. By Theorem 10, f0 + f1 = Îf 0 . Since F is self-dual then f0 + f1 = f 0 + f 1 which gives the result.

The converse of Theorem 7 is not true. In other words, it is not true that if f 0 + f 1 = tr then D 1 f 0 is a constant function, as it will be seen with several examples. However, we have a pseudo-reciprocoal theorem. The next Theorem is an improvement of a Theorem of [START_REF]Wolfmann Special Bent and Near-Bent Functions in Advances in Mathematics of Communication[END_REF] Theorem 12. Let F = [f 0 , f 1 ] be a bent function and let F = [ f0 , f1 ] be its dual function. Let be in F 2 .

f 0 + f 1 = tr + if and only if D 1 f0 = .

Proof.

Step 1: assume

f 0 + f 1 = tr + with ∈ F 2 .
From (P 1 ): f0 (u) = 2 2t-1 -2w(f 0 + t u ). Thus:

(R) f0 (u + 1) = 2 2t-1 -2w(f 0 + t u+1 ). • If = 0: Since f 0 + f 1 = tr then: f1 (u) = 2 2t-1 -2w(f 0 + t 1 + t u ) = 2 2t-1 -2w(f 0 + t u+1 ).
By using (R) we have f1 (u) = f0 (u + 1). We deduce from lemma 5: F (u, 0) = f0 (u) + f0 (u + 1) and F (u + 1, 0) = f0 (u + 1) + f0 (u) whence F (u, 0) = F (u + 1, 0). According to (P 4 ): f0 (u) = f0 (u + 1)) and this means D 1 f0 = 0.

• If = 1. Since f 0 + f 1 = tr + 1 then: f1 (u) = 2 2t-1 -2w(f 0 + t 1 + t u + 1) = 2 2t-1 -2w(f 0 + t u+1 + 1). Observe that w(f 0 + t u+1 + 1) = 2 2t-1 -w(f 0 + t u+1 ). This implies: f1 (u) = 2 2t-1 -2[2 2t-1 -w(f 0 + t u+1 )] = -2 2t-1 + 2w(f 0 + t u+1 ).
With (R) it follows f1 (u) = -f0 (u + 1). From lemma 5: F (u, 0) = f0 (u) -f0 (u + 1) whence F (u + 1, 0) = f0 (u + 1) -f0 (u) and finally F (u, 0) = -F (u + 1, 0). This means that F (u, 0) and F (u + 1, 0) are not equal to 2 t in the same time. From (P 4 ): f0 (u) = f0 (u + 1)) which yields D 1 f0 = 1.

Conclusion: if

f 0 + f 1 = tr + then D 1 f0 = .
Step 2: Conversely, assume D 1 f0 = . Lemma 9 shows that Î f0 = tr + . Hence by Theorem 10:

f 0 + f 1 = tr + . Corollary 13. If a bent function F = [f 0 , f 1 ] is self-dual then f 1 = f 0 + tr + if and only if D 1 f 0 = . Proof. This a direct consequence of Theorem 12 since if F is self-dual then f0 + f1 = f 0 + f 1 .
Another result on Îf is given in [START_REF] Leander | Construction of Bent Functions from Near-Bent Functions Journal of Combinatorial Theory[END_REF].

Theorem 14. (McGuire and Leander)

Let f be a near-bent function.

D 1 ( Îf ) = 1 if and only if [f, f + tr] is a bent-function.
Proof. see [START_REF] Leander | Construction of Bent Functions from Near-Bent Functions Journal of Combinatorial Theory[END_REF].

Corollary 15.

Let F = [f 0 , f 1 ] be a bent function and let F = [ f0 , f1 ] be its dual.

D 1 (f 0 + f 1 ) = 1 if and only if [ f0 , f0 + tr] is a bent function.
Proof. Theorem 10 says that Îf 0 = f0 + f1 and Theorem 14 proves that D 1 ( f0 + f1 ) = 1 if and only if [f 0 , f 0 + tr] is bent. The expected result is obtained by interchanging the roles of F and F .

Example: (obtained by computing ) t = 4, f 0 = tr(x 29 + x 27 + x 23 + x 21 + x 5 + x), f 1 = tr(x 29 + x 27 + x 23 + x 21 + x 9 )

f 0 + f 1 = tr(x 9 + x 5 + x) then D 1 (f 0 + f 1 ) = 1.
The dual of [f 0 , f 1 ] is [tr(x 13 + x 7 ), tr(x 19 + x 11 + x 7 )] then [tr(x 13 + x 7 ), tr(x 13 + x 7 + x)] is bent. The next theorem is the generalisation of a result of [START_REF]Wolfmann Special Bent and Near-Bent Functions in Advances in Mathematics of Communication[END_REF].

Theorem 17. Define two conditions on a (2t)-Boolean function G = [g 0 , g 1 ]: (T ) : g 0 + g 1 = tr + ξ with ξ ∈ F 2 and (C) : D 1 g 0 = µ with µ ∈ F 2
Let F be a bent function. If F meets condition (T ) with ξ = then:

A) The pseudo-duals F0 and F1 are bent functions.

F0 meets (T ) with ξ = 0 and (C) with µ = .

F1 meets (T ) with ξ = 0 and (C) with µ = + 1. Remark: We deduce from the previous result that if D 1 f 0 is not a constant function then F0 and F are distinct.

The following examples are obtained with computer assist.

Example 1: t = 4, D 1 f 0 = 0. F : f 0 (x) = tr(x + x 3 + x 7 + x 11 + x 19 + x 21 ), . f 1 (x) = f 0 (x) + tr(x) F : f0 (x) = tr(x 7 + x 11 + x 19 + x 21 ) f1 (x) = f0 (x) + tr(x) F0 = F and F0 = F F1 = F + [tr, tr].
Example 2: t = 4, D 1 f 0 = 0, 1.

F : f 0 (x) = tr(x 7 + x 13 ) . f 1 (x) = f 0 (x) + tr(x) F : f0 (x) = tr(x 5 + x 7 + x 9 + x 13 + x 19 + x 21 ) f1 (x) = f0 (x) + tr(x + x 5 + x 9 ) F0 : f (0) 0 (x) = f0 (x), f (0) 1 (x) = f0 (x) + tr(x) F1 : f (1) 0 (x) = f1 (x), f (1) 1 (x
) = f1 (x) + tr(x) F0 : f (0) 0 (x) = tr(x + x 7 + x 9 + x 13 + x 19 + x 21 ), If f 0 meets one of the two following cases and if F = [f 0 , f 0 + tr] then F, F , F0 , F1 , F0 , F1 are bent functions and they satisfy condition (T ).

f (0) 1 (x) = f (0) 0 (x) + tr(x) F1 : f (1) 0 (x) = tr(x + x 3 + x 7 + x 13 + x 19 + x 21 ) f (1) 1 (x) = f (
(1) Kasami-Welch case f 0 is a (2t -1)-Boolean function such that:

f 0 (x) = tr(x 4 s -2 s +1
) with 2t -1 ≡ 0 mod 3, 3s ≡ ±1 mod (2t -1), s < t.

(2) Quadratic case f 0 is a (2t -1)-near-bent function such that: f 0 (x) = t-1 j=1 c j tr(x 2 j +1 ) with c j ∈ F 2 .

Proof. In order to apply Theorem 17, we have to check in every case that F is bent and meets (T ) with ξ = 0.

• In case (1), the derivative D 1 f 0 is not a constant function but it was proved (McGuire, Leander, [START_REF] Leander | Construction of Bent Functions from Near-Bent Functions Journal of Combinatorial Theory[END_REF]) that [f 0 , f 0 + tr] is bent.

• In case [START_REF] Canteaut | Decomposing Bent Functions[END_REF] it is easy to prove that D 1 f 0 is a constant function and then [f 0 , f 0 + tr] is bent.

Remark: Several examples of near-bent functions in case (2) are given in [START_REF] Khoo | Stinson A new characterisation of semi-bent and bent functions on finite fields Des[END_REF].

Adding new functions.

Let F = [f 0 , f 1 ] be bent with f 0 + f 1 = tr.

Question: find g such that [f 0 + g, f 1 + g] is bent.

Remark : In this case, f 1 + g = f 0 + g + tr.

Recall that if u ∈ F 2 2t-1 then t u is the (2t -1)-function defined by t u (x) = tr(ux).

Theorem 20.

Let v be in

F 2 2t-1 .
If one of the following conditions is satisfied:

(i) f 0 is a (2t -1)-near-bent function such that the derivative D 1 f 0 is a constant function (ii) F = [f 0 , f 1 ] is a bent function such that f 0 + f 1 = tr. then the Boolean function F † = [f † 0 , f † 1 ] with f † 0 = f 0 + t 1 t v and f † 1 = f † 0 + tr is bent.
In order to prove Theorem 20 we need the two following lemmas. The second one is due to Canteaut and Charpin ( [START_REF] Canteaut | Decomposing Bent Functions[END_REF], Theorem 8).

Recall that the inner product <, > which is used for the calculation of the Fourier coefficients is now defined by: T (a,η) (x, ν) =< (a, η), (x, ν) >= tr(ax) + ην.

Lemma 21.

(L 1 ): The indicator of (u, 0), (v, 0) ⊥ is:

[t u t v , t u t v ] + [t u + t v + 1, t u + t v + 1]. (L 2 ): If G = [g 0 , g 1 ] then: D (b,0) D (a,0) G = [D b D a g 0 , D b D a g 1 ].
(L 3 ): If the derivative D 1 g is a constant function then for every u:

D 1 D u g = D u+1 D u g = 0
Proof. Proof of (L 1 ): Let I be the indicator (u, 0), (v, 0) ⊥ . (x, ν) is othogonal to (u, 0) if and only if < (u, 0), (x, ν) >= tr(ux) = 0. In other words the indicator of (u, 0) ⊥ is t u +1. Similarly, the indicator of (v, 0) ⊥ is t v + 1. This means that I(x, ν) = 1 if and only if: (tr(ux) + 1)(tr(vx) + 1) = 1. This result is independent of ν and

f 0 (u) = I(u, 0) = f 1 (u) = I(u, 1) = (t u + 1)(t v + 1) that is I = [(t u +1)(t v +1), (t u +1)(t v +1)]. Since (t u +1)(t v +1) = t u t v +t u +t v +1 then I = [t u t v , t u t v ] + [t u + t v + 1, t u + t v + 1]. Proof of (L 2 ): G(x, y) = (y + 1)g 0 (x) + yg 1 (x). D (a,0) G(x, y) = (y + 1)g 0 (x + a) + yg 1 (x + a) + (y + 1)g 0 (x) + yg 1 (x). = (y + 1)[g 0 (x + a) + g 0 (x)] + y[g 1 (x + a) + g 1 (x)] = (y + 1)D a g 0 (x) + yD a g 1 (x) Then D (a,0) G = [D a g 0 , D a g 1 ]. By replacing G by D (a,0) G and a by b: D (b,0) D (a,0) G = [D b D a g 0 , D b D a g 1 ]. Proof of (L 3 ): From D u g(x) = g(x + u) + g(x): D 1 D u g(x) = g(x + 1 + u) + g(x + 1) + g(x + u) + g(x), D u+1 D u g(x) = g(x + 1 + u) + g(x + 1) + g(x + u) + g(x) then D 1 D u g = D u+1 D u g. Note from the above calculation that D 1 D u g(x) = D 1 g(x+u)+D 1 g(x). Since D 1 g is a constant function then D 1 D u g = 0.
E = [t 1 t v , t 1 t v ] + [t 1 + t v + 1, t 1 + t v + 1].
Hence, in view of Lemma 22:

F = F + [t 1 t v , t 1 t v ] + [t 1 + t v + 1, t 1 + t v + 1] is bent. Since L = [t 1 + t v + 1, t 1 + t v + 1
] is an affine linear form then: F = F + L = [f 0 , f 1 ] with f 0 = f 0 + t 1 t v and f 1 = f 1 + t 1 t v also is bent.

Remark: Instead of t 1 t v is it possible to add some other products of two forms. This give rise to bent functions equivalent to those of the kind [f 0 , f 1 ] by means of the addition of a linear form. Example:

F † † = [f † † 0 , f † † 1 ]
with f † † 0 = f 0 + t u t u+1 and f † † 1 = f † † 0 + tr. It is easy to check that t u t u+1 = t 1 t u + t u . Therefore: F † † = [f 0 +t 1 t u +t u , f 0 +t 1 t u +t u +tr] = [f 0 +t 1 t u , f 0 +t 1 t u +tr]+[t u , t u ] Corollary 23. If (i) or (ii) of Theorem 20 is satified then f † 0 = f 0 + t 1 t v is near-bent. Proof. Since F † = [f † 0 , f † 1 ] is bent then according to Proposition 3, (a): f † 0 and f † 1 are near-bent.

Conclusion

We have introduced a way to construct new bent functions starting from a near-bent functions f such that D 1 f is a constant function or from a bent function such that the sum of the two components is a Boolean function of degree 1. An open question now is to describe explicitely the near-bent functions of the first kind, for example by means of the trace function. Another question is to generalise the study to the case where the degree of the sum of the two components is greater than 1.

Proposition 3 .

 3 F is a bent function if and only if: (a) f 0 and f 1 are near-bent. (b) ∀a ∈ F 2 2t-1 | f0 (a) | + | f1 (a) |= 2 t Remark: (b) means that one of | f0 (a) | and | f1 (a) | is equal to 2 t and the other one is equal to 0.

3. 3 .

 3 Pseudo duality. The results of the previous theorems lead to introduce a new definition. Definition 16. Let G = [g 0 , g 1 ] be a (2t)-bent function and let G = [g 0 , g1 ] be its dual function. The pseudo-duals of G are the two (2t)-Boolean functions: Ḡ0 = [g 0 , g0 + tr] and Ḡ1 = [g 1 , g1 + tr].

2 )

 2 If D 1 f 0 = 0 then 1) shows that F0 = F whence f0 + f1 = tr. Thus F1 = [ f1 , f1 + tr] = [ f0 + tr, f0 + tr + tr] = [ f0 + tr, f0 ] = [ f0 , f0 + tr] + [tr, tr]. Finally F1 = F0 + [tr, tr] = F + [tr, tr].

3 . 4 .

 34 1) 0 (x) + tr(x + 1) New infinite families of bent functions. Proposition 19.

  Lemma 22. ([2],Theorem 8) Let F be a (2t)-bent function. Let A and B two distinct non-zero elements of F 2 2t and E = A, B ⊥ . Let Φ E be the indicator of E. Then the function F + Φ E is bent if and only if the dual function F satisfies D B D A F = 0. Proof. of Theorem 20. If v = 0 or v = 1 the result is trivial. If v = 0 and v = 1: From Theorem 7, in case (i) the Boolean function [f 0 , f 0 + tr] is bent. Consequently, in both cases (i) and (ii) we have to consider the bent function F = [f 0 , f 0 + tr]. Let F = [ f0 , f1 ] be the dual function of F . Theorem 12 proves that D 1 f0 = 0 and D 1 f1 = 1. Therefore, from (L 3 ): D 1 D v f0 = D 1 D v f1 = 0. It follows from (L 2 ) that D (1,0) D (v,0) F = [0, 0]. Now, using the notations of Lemma 22 with A = (1, 0), B = (v, 0), v = 0, v = 1 then from (L 1 ): Φ

B) The dual F0 of F0 meets (T ) with ξ = and (C) with µ = 0. C) The dual F1 of F1 meets (T ) with ξ = + 1 and (C) with µ = 1.

Proof.

Notations:

1 ] Proof of A): Since f 0 + f 1 = tr + then, from Theorem 12 and Lemma 6: ( ) D 1 f0 = and D 1 f1 = + 1. By the definition F0 = [ f0 , f0 + tr]. That is F0 meets (T ) with ξ = 0. Furthermore, we deduce from Theorem 7 that F0 is a bent function and meets (C) with µ = . Similarly, again from Theorem 12 and from Theorem 7, F1 is a bent function and meets (T ) with ξ = 0 and (C) with µ = + 1.

Proof of B) and C):

From the definition of the duals: Starting from a bent function F = [f 0 , f 1 ] such that f 0 + f 1 = tr + , Theorem 17 provides new bent functions. The question now is: are F, F , F0 , F1 , F0 , F1 distinct functions or not ? A first answer is given by the next proposition.


Lukáš Ručka

Igor Peterlík

FAST RECONSTRUCTION OF IMAGE DEFORMATION FIELD USING RADIAL BASIS FUNCTION

Keywords: image warping, image interpolation, medical imaging, non-rigid registration

Fast and accurate registration of image data is a key component of computer-aided medical image analysis. Instead of performing the registration directly on the input images, many algorithms compute the transformation using a sparse representation extracted from the original data. However, in order to apply the resulting transformation onto the original images, a dense deformation field has to be reconstructed using a suitable inter-/extra-polation technique.

In this paper, we employ the radial basis function (RBF) to reconstruct the dense deformation field from a sparse transformation computed by a model-based registration. Various kernels are tested using different scenario. The dense deformation field is used to warp the source image and compare it quantitatively to the target image using two different metrics. Moreover, the influence of the number and distribution of the control points required by the RBF is studied via two different scenarios. Beside the accuracy, the performance of the method accelerated using a GPU is reported.

INTRODUCTION

Image registration is the process of overlaying multiple images of the same scene, acquired by different modalities, from various viewpoints and/or at different times [START_REF] Zitová | Image Registration Methods: A Survey[END_REF].

The goal of registration is thus to find the optimal transformation W , mapping coordinates from source image S onto target image T [START_REF] Glasbey | A review of imagewarping methods[END_REF]. From mathematical perspective, image registration an is optimization of objective function given as

M (W (S), T ) + R M (W ) (1)
which consists of similarity measure M , typically complemented with a regularization term R M , which introduces a priori knowledge into the optimization process. It is often used to limit the set of possible transformations to those which are physically feasible [START_REF] Sotiras | Deformable Medical Image Registration: A Survey[END_REF]. Instead of using the dense representation typically represented by the images T and S, the registration process can be designed as optimization of objective function evaluated on a sparse representations t and s, typically geometric primitives such as points or convexes, extracted from the original data. In this case, a sparse transformation w is computed and the metric M is usually replaced by a different criterion m such as Euclidean or Hausdorff distance. The sparse registration is then summarized as follows: SR1 Sparse representation s and t such as set of features or a mesh is extracted from the source and target images.

SR2 Sparse transformation w is computed by optimization of the criteria m applied to t and w(s). The process can be again regularized by a priori information introduced via a suitable term.

SR3 Dense transformation W is constructed from w via suitable interpolation or approximation technique [START_REF] Wang | Non-Rigid Multi-Modal Image Registration Using Cross-Cumulative Residual Entropy[END_REF].

SR4 Similarity measure M is evaluated over target image T and transformed source image W (S).

The approach described above has been employed in a large number of algorithms. A finite element model (FE model) is core component of ITK filters for physics-based non-rigid registration presented in [START_REF] Clatz | Robust Nonrigid Registration to Capture Brain Shift from Intraoperative MRI[END_REF][START_REF] Liu | An ITK Implementation of a Physics-Based Non-Rigid Registration Method for Brain Deformation in Image-Guided Neurosurgery[END_REF]. They use the data from image block matching as input for a finite element modeling of the underlying tissue, which yields sparse transformation. The deformation field W is then obtained as a result of hybrid iterative approximation-interpolation algorithm with outlier block rejection step. Non-linear finite element modeling has been also deployed for brain registration [START_REF] Mostayed | Biomechanical Model as a Registration Tool for Image-Guided Neurosurgery: Evaluation Against BSpline Registration[END_REF], interpolating the sparse transformation with algorithm proposed in [START_REF] Joldes | Cortical Surface Motion Estimation for Brain Shift Prediction[END_REF]. The interpolation is based on B-spline deformation of Gaussian curvature.

A compound dense transformation is used in [START_REF] Mattes | PET-CT Image Registration in the Chest Using Free-Form Deformations[END_REF], consisting of free-form deformation with homogeneous coarse grid and B-spline convolution kernel for interpolation inbetween the latices. This approach has been adapted in [START_REF] Wang | Non-Rigid Multi-Modal Image Registration Using Cross-Cumulative Residual Entropy[END_REF], improving the results with a cross-cumulative residual entropy metric. The dense deformation field, obtained by the CLI-TPS algorithm [START_REF] Johnson | Consistent Landmark and Intensity-Based Image Registration[END_REF] is subjected to constraint-based smoothing in lung CT registration [START_REF] Li | Pulmonary CT Image Registration and Warping for Tracking Tissue Deformation During the Respiratory Cycle Through 3D Consistent Image Registration[END_REF]. The constraints presented ensure that the dense field preserves the topological properties of the original structure.

The sparse transformation obtained by the registration process often provides values in positions which are not organized in a regular structure. Therefore, an interpolation from scatter data such as radial basis functions (RBF) interpolation is used to compute the dense transformation [START_REF] Amidror | Scattered Data Interpolation Methods for Electronic Imaging Systems: A Survey[END_REF]. A series of criteria for sparse deformation field interpolation via RBF was given in [START_REF] Arad | Image Warping Using Few Anchor Points and Radial Functions[END_REF]. Furthermore, computational efficiency of RBF interpolation has been addressed in [START_REF] Livne | Fast Multilevel Evaluation of Smooth Radial Basis Function Expansions[END_REF] and later evaluated [START_REF] Lapeer | An Optimised Radial Basis Function Algorithm for Fast Non-Rigid Registration of Medical Images[END_REF]. Remarkably, in this work, the uncertainty at input was greater than the error introduced by approximation of RBF solution.

In this paper, we focus mainly on steps SR3 and SR4 as given in the list above. Both SR1 and SR2 are performed by an external registration process. We have chosen a meshbased registration presented in [START_REF] Peterlik | Model-Based Identification of Anatomical Boundary Conditions in Living Tissues[END_REF]. The method presented therein is regularized with co-rotational formulation of linear elasticity, therefore the resulting sparse elastic transformation is physically admissible and capable of handling very large displacements. We apply the RBF to convert the sparse transformation and study the influence of the RBF method for both the domain covered by the FE mesh (interpolating mode), as well as in the areas outside of the mesh (RBF in extrapolating mode). Method is assessed using two different image metrics and visual inspection for different RBF splines, resolutions of the image data and different inputs provided by the sparse registration.

METHOD

In this section, we first shortly describe the RBF interpolation method. Then, we couple this technique with a sparse FEbased registration in order to obtain dense transformation that can be used to warp 3D images.

Radial Basis Function Interpolation

In general, the RBF interpolation in point x is given as a weighted sum of radial base function ϕ(r), centered in control points X i :

f (x) = A • x + i w i • ϕ( X i -x ) ( 2
)
where A is an affine transformation term. The interpolation is computed in two steps. First, given known values of the function f in m control points X j , j = 1 . . . m, the vector w of weights is found by substituting f (X j ) for f (x) on the left-hand side, and X j for x on the right-hand side of the Eq. 2 and solving the resulting system of linear equations. Second, for any new position x, the value f (x) is computed directly by Eq. 2. By its nature, RBF interpolation is unidimensional; interpolation of multidimensional values is considered as combination of independent unidimensional transformations. As for the choice of the base function, in this paper, we consider following options:

ϕ(r) = r 2 • log(r)
Thin plate splines (TPS) -solution of bending energy functional, minimizing second order bending energy potential in 2D [START_REF] Lapeer | An Optimised Radial Basis Function Algorithm for Fast Non-Rigid Registration of Medical Images[END_REF]. ϕ(r) = r 2 (biharmonic RBF spline) -solution of bending energy functional, minimizing second order bending energy potential in 3D [START_REF] Lapeer | An Optimised Radial Basis Function Algorithm for Fast Non-Rigid Registration of Medical Images[END_REF]. ϕ(r) = r 3 (triharmonic RBF spline) -solution of bending energy functional, minimizing third order bending energy potential in 3D [START_REF] Lapeer | An Optimised Radial Basis Function Algorithm for Fast Non-Rigid Registration of Medical Images[END_REF].

ϕ(r) = r 3 • log(r), ϕ(r) = r 2 • log(log(r)), ϕ(r) = r • log (log(r))
-test base functions, illustrating different aspects of dense transformations they yield 1 . In the following, these base functions are referred to as excited triharmonic (ETRI), damped thin-plate spline (DTPS) and excited linear with damping (EDL) respectively.

From Sparse Registration to Dense Transformation

In the scenario considered in this paper, we employ the surfacematching registration presented in [START_REF] Peterlik | Model-Based Identification of Anatomical Boundary Conditions in Living Tissues[END_REF] which performs a meshto-mesh registration based on constrained biomechanical simulation. Indeed, the method does not work directly with the images, but it is necessary to extract the sparse representation given by the FE meshes of the object of interest in source and target configurations, M s and M t , respectively. During the registration process regularized by the co-rotational formulation of the linear elasticity, the M s is deformed due to the application of sliding constraints which are used to match the shape of M s to the shape of M t . The procedure results in the registered mesh Ms : although being deformed, this mesh has the same topological structure as M s . The vectors of nodal positions of the two meshes constitute the sparse transformation computed by the mesh-to-mesh registration. Now, to obtain the dense transformation, RBF is employed as follows. Using the nomenclature introduced in the previous section, the nodes of the deformed mesh Ms are used as control points and f (ñ s ) = n s -ñs is the 3D function which is used to compute three weight vectors w x , w y , w z and affine transformation A. Then, the value of transformation is interpolated over the original image, i.e., the Eq. 2 is evaluated in each voxel of the image.

RESULTS

Data and Sparse Registration

The evaluation was performed used a set of CT scans of a female pig liver (target organ) and kidneys (auxiliary organs). Two abdominal CT scans were acquired in flank and supine positions with SOMATOM R Definition AS 128 device. A bounding box containing the liver was extracted from the original The FE meshes extracted from the original images were composed of 843 nodes (liver) and 165 and 207 nodes for kidneys. The sparse registration was performed using the method presented in [START_REF] Peterlik | Model-Based Identification of Anatomical Boundary Conditions in Living Tissues[END_REF]: the supine configurations of the liver and kidney meshes were registered to the corresponding flank configurations; the registration of each organ was performed separately resulting in three complementary sparse transformations. While the kidneys displayed important translation and rotation, the liver underwent large deformations resulting in large displacements of lateral lobes (up to 6 cm).

Testing Scenarios

The dense transformation was constructed for each base function from section 2.1. Further, two different scenarios were considered. Since we are interested in registration of liver and its surroundings, the dense transformation was first constructed using only the control points provided by the liver mesh. In the second scenario, the control points obtained by the registration of kidneys as auxiliary organs were also added to the liver mesh nodes to augment the set of control nodes. Since the three organs perform relatively independent movements, the goal was to test the capabilities of the RBF method.

All the scenarios were evaluated for both the low-and highresolution images using two metrics: the interpolation mode was evaluated by the Dice metric [START_REF] Zijdenbos | Morphometric Analysis of White Matter Lesions in MR Images: Method and Validation[END_REF] computed using the binary maps of the liver. Further, Mattes mutual information (MMI) [START_REF] Mattes | Nonrigid multimodality image registration[END_REF] was computed using the target image and the source image warped by the dense transformation. Hence, the Dice metric was used to assess the interpolating mode locally, while MMI was employed as a global metric quantifying the quality of warping that was obtained by combining the physicsbased approach (the sparse registration) and geometric RBF method.

Evaluation Using Warped Images

Dense transformations were constructed using the RBF method using a PC equipped with Intel R Core TM i7-3770 CPU running 3.40 GHz, 32 GB RAM and NVIDIA GeForce GTX 560. The second phase of the method, where the RBF is evaluated for each voxel, was accelerated using the GPU. The computing times are reported in Tab. The quantitative results showing the MMI and Dice metrics are reported in Tab. 1. Beside the values of metrics, the relative differences are shown on the right-hand side. The best results are displayed in colors. Moreover, an illustration of results is shown in Fig. 1. We have deliberately chosen an axial slice located below the liver, which represents the area covered with the control points, in order to illustrate the behavior of the extrapolated dense deformation field for different choices of RBF spline.

Both TPS and triharmonic splines perform well in terms of Dice metric. Nonetheless, although the liver keeps its shape, the near vicinity suffers from vast drifts. While triharmonic spline renders the segmented liver smoother, TPS creates artificial concavities. Furthermore, the results of EDL, biharmonic and DTPS splines display visual similarity, without any noticeable issues. However, none of the splines allows for approaching the deformation of the flank image outside of the liver region. In the scenario using also the registered kidneys as auxiliary sources of control points, the Dice metric performs best with both Thin plate splines and DTPS. However, TPS achieves inferior score in terms of Mattes MI metric. On the other hand, the performance of other spline is acceptable for Mattes mutual information (MMI) metric as well. The best results for MMI are achieved using the biharmonic spline. Upon visual inspection of the segmented liver, all three results look naturally. Nevertheless, the TPS seems to emphasize bulges on liver surface. Furthermore, visual inspection shows slight image deterioration along the borders, clearly visible for TPS warp. In addition, the DTPS exhibits stable results throughout distinct scenarios.

Out of all RBF splines tested, the ETRI performs the worst in terms of both Dice and MMI metrics. Upon visual inspection, the surface of the liver seems to suffer from material drift. Inspection of the deformation field reveals large displacements in image parts uncovered by RBF control points, rendering the intensity image completely useless. Correspondingly, adding displacement contribution from kidneys lead to complete deformation field breakdown, thus degenerating segmented liver to point cloud.

DISCUSSION AND CONCLUSIONS

While including auxiliary organs slightly corrupts the Dice metric, it improves the MMI metric in all cases where single organ warping does not already exhibit deterioration. The improvement is, however, only decimal order higher than metric sensitivity, computed as the difference of MMI values for noninterpolated warps. Thus, as the changes in Dice metric are negligible (< 0.0003), including auxiliary organs improves the overall MMI result. However, neither the improvement alone, nor MMI score guarantee that the result is acceptable. For example, the triharmonic spline achieves reasonably-high score in MMI, yet the image is severely degenerated. The only significant deterioration of Dice result in combination with improvement in MMI suggest degenerated image.

It is equally important to understand the nature of all the splines considered in this paper: the further the inter-/extrapolated voxel from a control point, the stronger the effect of given control point. The data used also displays strong rotational changes for kidneys, as pig lacks ligaments which are present in human body. This observations combined probably explain the poor performance of the triharmonic and higher splines. As noted in [START_REF] Arad | Image Warping Using Few Anchor Points and Radial Functions[END_REF], the distance property can introduce significant artifacts.

Beside testing other intensity-based metrics, we plan to compare this warping technique to other global interpolation methods such as natural neighbor interpolation. However, exploring different interpolation frameworks or base functions with better locality properties (e.g. Gaussian) is beyond the scope of this paper and will be subject of the future work.

Fig. 1 :

 1 Fig. 1: Comparison of the warped images obtained via different RBF splines. (a): source image, (b): target image On the left, scenarios using control points only from liver mesh. On the right, scenarios using control points from liver and kidney meshes. Types of the RBF spline: (c),(d): biharmonic, (e),(f): DTPS, (g): EDL, (h): triharmonic.

Table 1 :

 1 Evaluation of target organ driven warping with respect to Dice and Mattes MI metrics. Blue and green † colors denote the best and the second best result achieved for given metric (maximum for Dice and minimum for Mattes MI).scans in two versions: first, the original resolution was preserved resulting in high-resolution images (HR) composed of 480 × 452 × 388 voxels with spacing 0.5 × 0.5 × 0.5 mm 3 . The same volume was resampled to lower resolution (LR): 240 × 226 × 194 voxels, spacing 1 × 1 × 1 mm 3 .

			Liver only			Liver + kidneys			Cross-scenario improvement
	Interpolation		Dice		MMI		Dice	MMI	Dice	MMI
		LR	HR	LR	HR	LR	HR	LR	HR	LR	HR	LR	HR
	No interpolation	0.5123	0.5148	-0.1032	-0.098	0.5123	0.5148	-0.1032	-0.098	-	-	-	-
	EDL	0.9333	0.9495	-0.2395	-0.2237	0.9333	0.9496	-0.2586	-0.2457	0	0.0001	-0.0191	-0.0219
	Biharmonic	0.9341	0.9519	-0.2306	-0.2301	0.9339	0.9519	-0.2635	-0.2607	-0.0002	0	-0.0329 †	-0.0306 †
	DTPS	0.9358	0.9526	-0.239 †	-0.2247 †	0.9357 †	0.9526	-0.2608 †	-0.2543 †	-0.0001	0	-0.0218	-0.0296
	Thinplate	0.9363 †	0.9528 †	-0.2236	-0.2187	0.936	0.9526 †	-0.2533	-0.2435	-0.0003	-0.0002	-0.0297	-0.0247
	Triharmonic	0.9371	0.9531	-0.179	-0.171	0.8531	0.8621	-0.2495	-0.2335	-0.084	-0.091	-0.0705	-0.0625
	ETRI	0.9286	0.9446	-0.1703	-0.1622	0.1453	0.1528	-0.0547	-0.0377	-0.7833	-0.7918	0.1156	0.1245

Table 2 :

 2 2. Processing

	Scenario	GPU L [s] LR HR	GPU L+K [s] LR HR	CPU L [s] LR HR
	EDL	2.64 21.13 4.76 37.88 338.9 2 695
	Biharmonic	1.70 13.45 3.23 25.64 132.2 1 052
	DTPS	2.59 20.57 4.69 38.17 333.8 2 645
	Thinplate	2.19 17.39 4.10 32.61 231.1 1 833
	Triharmonic 1.72 13.62 3.24 25.79 133.6 1 067
	ETRI	2.17 17.20 4.06 34.22 253.2 2 012
		times of individual scenarios -liver (L),
	liver + kidneys (L+K) -retrieved via code instrumentation.
	CPU times for liver registration (CPU L) are provided for
	illustration.

As noted in[START_REF] Amidror | Scattered Data Interpolation Methods for Electronic Imaging Systems: A Survey[END_REF], certain base functions can yield good results even though they are not mathematically justified.