Anne-Marie Barthe-Delanoë
email: anne-marie.barthe@mines-albi.fr

Matthieu Lauras
email: matthieu.lauras@mines-albi.fr

Sébastien Truptil
email: sebastien.truptil@mines-albi.fr

Frédérick Bénaben
email: frederick.benaben@mines-albi.fr

Hervé Pingaud
email: herve.pingaud@univ-jfc.fr

A Platform for Event-Driven Agility of Processes: A Delivery Context Use-Case

Keywords: event-driven architecture, agility, complex event processing, detection, adaptation

The French SocEDA project aims at providing a Service Oriented Architecture (SOA) platform for dynamic and complex event-driven interactions between large highly distributed and heterogeneous service systems. This platform should be able to combine a publish/subscribe mechanism (to collect events coming from heterogeneous and distributed services) with complex event processing of these collected events in order to detect interesting situations a service should react on. The SocEDA platform should also offer suggestions about relevant changes (adaptation) to do in running processes and services. The context of city logistics is chosen to show some of the basic abilities of the SocEDA platform and the relevance of event management in everyday situations.

Introduction

In the literature of the last decade, many authors have proposed the idea of "city logistics" to solve the issues of urban freight transport like congestion and negative impacts on the environment (such as air pollution and noise) [START_REF] Ruske | City-Logistics -Solutions for urban commercial transport by cooperative operations management[END_REF] [START_REF] Taniguchi | City Logistics: Network Modelling and Intelligent Transport Systems[END_REF].

But, as it is widely recognized that companies are nowadays operating in a complex economical environment where markets are more open and globalized, the freight transport companies have also to deal now with uncertainties and short-term changes in demand or supply. They need to have agile processes to respond quickly to these changes, whatever their nature (traffic environment, congestion, safety and energy savings, uncertainties on supply or demand).

This article is organized as follows: Section 2 presents the architecture of the SocEDA platform. Section 3 develops the agility service, including the detection of the relevant moment where change is needed, and the proposed adaptation of running processes and services. Finally, Section 4 shows a use case inspired by the delivery of drugs to French pharmacies, with the associated events and business rules.

SocEDA Platform

As the SocEDA platform is a global (Internet level scale) structure to combine events from many sources with the goal of connecting and orchestrating services, things and people, it can contributes to support the coordination of shippers, carriers and movements in an integrated logistics system. Through a federated middleware, the platform provides a transparent and distributed infrastructure to allow users and services to be connected in a publish/subscribe mechanism through the whole SocEDA platform. We distinguish the main components that form the federated middleware and the components that are plugged in and using this middleware.

Fig. 1. Overview of SocEDA platform's components

As shown in Figure 1, the main components of the SocEDA platform are:

• Event-Cloud stores and forwards events between the bus and DiCEPE (Distributed Complex Event Processing), which are filtered by the social filter. Every event in the platform must be forwarded to the Event-Cloud to be stored and to be handled asynchronously (pub/sub mechanism) or synchronously by simple query retrieval. • The Distributed Service Bus (DSB) is in charge of connecting Web Services consumers and providers to the platform (SOA layer). It includes an Event

Broker to integrate events. The DSB implements the Web Service Notification (WS-N) standard [START_REF]OASIS: Web Services Base Notification 1.3 OASIS Standard[END_REF] in order to provide the core EDA feature. • The Agility Service has to ensure the agility of the dynamic of the collaboration, in accordance with the evolutions of the context. Its objective can be divided into two sub-objectives: (i) detect the evolution of the context and (ii) adapt the dynamic of the collaboration to the new context. • The Governance component to govern the services and events.

• SeaCloud allows to subscribe to a specific event producer, to add Complex Event Processing rules or deploy process. The components of the platform that connects to the federated middleware are:

• The DiCEPE component: reasons on incoming events to generate complex events to the Event-Cloud to be sent to subscribers. It contains a Complex Event Processing (CEP) engine. It detects relevant events or combinations of events by applying event patterns (as described in [START_REF] Etzion | Event Processing in Action[END_REF]) on these events. • Event Sources, i.e. business services in charge of sending/receiving events to/from the federated middleware (partners' devices, sensors, reports, etc.).

Agility Mechanism

The delivery processes defined in the city logistics area are subject to issues related to the dynamic adjustments of vehicle, as underlined by [START_REF] Taniguchi | Intelligent transportation system based dynamic vehicle routing and scheduling with variable travel times[END_REF]. For a decade, several commercial products and research projects have been attempting to provide agility to collaborative workflows. On the commercial side, we can cite Bonita (a suite of tools to design, execute and monitor processes) and the Architecture of Integrated Information Systems (ARIS) tools. ARIS has the ability to combine determined process fragments according to received events. In a way, the ARIS approach manages workflow adaptation (but in a determinist manner) [START_REF] Scheer | ARIS Architecture and Reference Models for Business Process Management[END_REF]. On the research works side, several projects like the WORKPAD project [START_REF] Catarci | WORKPAD : process management and geo-collaboration help disaster response[END_REF] and the CRISIS project [START_REF]Crisis Project website[END_REF] focus more on recovering the disconnecting nodes through specific tasks (WORKPAD project) or supporting collaboration into crisis situation and on exploring decision-making under conditions of uncertainty (CRISIS project). Another research project, the European PLAY project proposes an adaptation recommender service [START_REF] Verginadis | Service Adaptation Recommender in the Event Marketplace: Conceptual View[END_REF] which allows to adapt the ongoing processes on pre-determined milestones, through the addition of relevant pieces of processes (extracted from a knowledge database).

The Agility Service in the SocEDA project is a tool designed to adapt the dynamic of the situation, i.e. processes, to the evolution of the context. In this section, we will detail the mechanisms of the Agility Service. The adopted approach of agility in our research work is given by the formula below (which is close to the definition given by [START_REF] Charles | A model to define and assess the agility of supply chains: building on humanitarian experience[END_REF]):

Agility = (Detection + Adaptation)*Reactivity . (1)
• Detection: this step consists in detecting an evolution of the situation that could not be solved by the ongoing processes or making the running processes not relevant to the current situation, • Adaptation: when an evolution is detected, the adaptation step is executed to modify the ongoing processes in order to make them relevant to the current situation, • Reactivity: detection and adaptation steps have to be done as fast as possible (and real-time if possible).

Agility allows on one hand to detect if the ongoing processes meet the requirements of the current situation, on the other hand to adapt the ongoing processes if necessary.

Detection

The approach to detect a divergence is based on three major steps:

1. Define the situation model to obtain expected situation model (picture of the expected situation) and field situation model (picture of the ongoing situation), 2. Comparing expected situation model and field situation model, 3. Calculate the difference between these two models to check the adequacy of the running processes with the field situation model. A situation model is an instant capture of the running collaborative processes, the collaboration itself (all the actors and their services, their objectives, etc.) and the environment in which the processes are running (risks, opportunities, goods, etc.). We can use these events to track the changes inside the collaborative situation model by this method (illustrated by Figure 2):

• First, we duplicate the initial model of the collaborative situation (i.e. model at time 0), • Then, both models are automatically updated with the received events (the Agility Service subscribes to all the events emitted by the DiCEPE). We obtain two models through this update:

Fig. 2. Principles of detection of the divergence

The measure of ∂ is given by the following formula:

(

Where:

• ∂ i is a difference identified between the expected model and the field model. E.g.: a partner is not available, a new risk has appeared, a resource is free, • w i is the weight applied to ∂ i according the type of element concerned by the identified difference (e.g. partner, risk, resource, activity, etc.) and the kind of difference, called operation here (added, deleted, updated). This weight is used to qualify each detected difference, as each difference has not the same impact on the relevancy of the processes. For example, the addition of an order (the deletion of a partner) has more negative impact on the processes than the deletion of an order (the addition of a new partner) (see Table 1). The matrix {element type × operation} is a parameter defined by the partners of the collaboration among the execution of the Agility Service. Among the difference detection, a threshold (named ∂ threshold ) is also defined by the partners. If the calculated ∂ is equal or over this threshold, the Agility Service will execute the adaptation step.

Adaptation

We propose combining two kinds of system adaptations (to the most appropriate conformation): (i) the ability to evolve in a predetermined closed geometry and/or (ii) the ability to redesign a new structure fitting the situation.

The first point refers more to a "Design-time" agility, based, for instance, on risk studies and leading to the building of models including a number of conditional branches to optimize coverage of the possibilities. The second point refers more to a "Run-time" agility where the (re)building of the best possible conformation has to be improvised, at the convenient moment.

To decide what kind of adaptation should be proposed to the user, the adaptation part is based on the study of the difference details gathered by the ∂ calculus in the detection step. First of all, we defined a set of business rules describing how the adaptation step can advice users on the loop to choose (partial or total redefinition, etc.), based on the details of each detected difference and on the weight of each concept participating into the calculated ∂.

Then, when all theses rules are run, the Agility Service is able to indicate the best solution(s) for adaptation to the users, considering all the detected differences. It is very important to note that the final choice is let to the users. The Agility Service is a decision support system in order to choose a relevant adaptation of the processes, not the system that takes the decision of the adaptation.

City Logistics And Transportation Use-Case

Our use-case is based on the process of French drug deliveries to pharmacies. It is inspired by the one described in [START_REF] Lauras | An event-driven platform to manage agility: Behavior adaptation in delivery context[END_REF]. In France, drugs and Over-the-counter (OTC) medicines are produced by various pharmaceutical manufacturers and mainly distributed through a two-levels network composed of wholesalers and pharmacies (about 80% of the national sales of drugs and OTC medicines [START_REF] Hénaff | Médicaments, une distribution fébrile[END_REF]) as shown on Figure 3. As the French law requires wholesalers to have a minimum stock of each product in their warehouses, it is considered in this use-case that there is no shortage at wholesaler level. The focus is put on the products distributed from wholesalers to pharmacies. A wholesaler supplies several pharmacies, contained in a given geographical area, about three to four times per day (in average) [START_REF] Hénaff | Médicaments, une distribution fébrile[END_REF]. The average delay between the reception of the order by the wholesaler, and the delivery to the pharmacy is 2h15. Usually, the deliveries are organized as a set of scheduled rounds.

Considering the fact that the deliveries are dependent on traffic conditions, lorries conditions and availability, weather conditions, it is crucial to (i) gather relevant information in real time and (ii) to manage it very quickly. Thus any changes, any evolution, any information that could challenge the rounds and the quality of service of the delivery process have to be managed in order to be able to adapt the routes of the wholesaler's lorries in real time. According to [START_REF] Etzion | Event Processing in Action[END_REF], [START_REF] Chandy | Event Processing: Designing It Systems for Agile Companies[END_REF] and [START_REF] Luckham | Event Processing Glossary -Version 1.1[END_REF], these elements that happened and that embedded data can be considered and managed as events.

We have listed four macro types of events covering the whole kinds of events that can happen in various contexts as delivery management, as shown in Table 2. These events help us to gather knowledge about the evolution of the situation (i.e. process execution and environment of execution). In the context of the delivery use-case, these events can embed at least (i) the traffic information provided by the city hall to assess the traffic situation and potential new routes (situational event), (ii) the pharmacies' stock to assess the priorities in terms of supply (resource event), (iii) the patients' urgent orders to assess the emergencies in term of demand (consequence event), (iv) the ongoing preparations at wholesaler's warehouse to adapt the time delivery schedule (activities event), etc.

Based on these events, we defined some business rules in order to be able to adapt the behavior of the delivery system within the SocEDA platform. Below is a sample of the rules used in this use-case:

• If there is an important traffic jam 1 km around pharmacy N then go to pharmacy N+1 , • If the planned route is modified then inform all the remaining pharmacies for the delay, • If there is an urgent request for a specific drug then proceed to the delivery of this drug first.

The considered scenario concerns one vehicle, which has to deliver five pharmacies (P) located in various places of a given area. The scheduled round is: So the Agility Service of the SocEDA platform modifies its round schedule and proposes an adapted route to the driver, in order to avoid the traffic jam and to deliver P#34 instead of P#22. The SocEDA platform also informs P#22 that the delivery should be effective around 10:45 am. Finally, he arrives at P#34 at 9:18 am and leaves at 9:28 am.

The new schedule indicates to the driver the route to go to P#11. While he is driving to P#11, he receives an alert about an urgent demand received by P#14 and for which the medicines are in the delivery order he is transporting. The Agility Service of the SocEDA platform adds an Order element into the field situation model (and nothing into the expected situation model). According the matrix {element type × operation}, ∂ = 2 and ∂ = ∂ threshold . Consequently, the platform suggests him to change his route to deliver first P#14 instead P#11. And the pharmacist of P#11 receives a notification from the platform at 9:38 am about the move of its delivery (from 11:00 am to 10:25 am). Finally the driver reaches P#14 at 9:58 am and leaves at 10:13 am. The medicine is delivered without any delay to the patient. The driver continues his route normally and arrives à P#11 at 10:24 am and leaves at 10:35 am. The driver goes to P#22 following the route proposed by the platform: he reaches P#22 at 10:45 am and leaves at 10:55 am.

Finally the driver goes back to the wholesaler's warehouse at 11:05 am, with only a delay of 15 min on the planned schedule. This example shows how the SocEDA platform allows detecting, analyzing and reacting on events and can help people avoiding losing time and means by adapting processes.

Conclusion

The SocEDA platform provides an environment to manage processes in an eventdriven perspective, allowing the users to be dynamically connected to events coming from the environment and the processes themselves. Events are gathered through a publish-subscribe mechanism based on topic and content subscription. This kind of subscriptions is very interesting because it allows the user to get the relevant information at the right time, without any necessity of knowing all the event sources providing this information. Events gathered by the SocEDA platform are filtered and aggregated to detect interesting situations and, through the SocEDA's Agility Service, to propose relevant changes to existing and running processes and services. Regarding the overall proposition, two reproaches might be done: (i) what about the robustness of such a system? And (ii) how to deal with access to events and unreliable data?

Concerning the first point, it is clear that as we are based on organizations' Information Systems, we are thus technically dependent on the network. If the network goes down, a part or, in the worst case, the entire system fails too. Hardware security measures should be taken to protect the physical network, in addition to the security measures taken to protect the data network against acts of piracy.

Concerning the second point, we agree on that events can only be used if they are trustable: a possible solution to tackle this point may be the definition of a social filter to evaluate the event providers and the events themselves [15] [16]. This social filter operates on a social network of services to compute the strength of the relationships between them through a trust interference algorithm. The social filter is one of the components of the SocEDA platform and will be integrated on the CEP. Another complementary solution could be the use of a governance tool to manage subscriptions to event types.

Fig. 3 .

 3 Fig. 3.French drugs and OTC medicines distribution network (according[START_REF] Lauras | An event-driven platform to manage agility: Behavior adaptation in delivery context[END_REF]).

•

  The expected situation model: the planned and expected situation model at time t (i.e. what we expect to obtain when we apply the collaborative workflows to the collaborative situation). It is obtained by updating the initial model with monitoring events, • The field situation model: the real situation model of the collaboration at time t, whatever the applied collaborative workflows are (i.e. the "what actually happened" situation at time t). It is obtained by updating the initial model with events coming from the field. • At any time t, a measure of the divergence ∂ is made between the expected model and the field model.

Table 1 .

 1 Sample of an element types × operation matrix.

	Element	Added	Deleted Updated
	Partner	1	3	2
	Traffic Jam Constraint 4	3	2
	Order	2	1	1

Table 2 .

 2 Macro event types.

	Event type	Description
	Situational	Used for information about the situation (e.g. information coming from
		sensors or human made description)
	Consequence	Used to transmit the result of one activity/a sequence of activities
	Activity	Used to give information about the state of a service: waiting, work in
		progress, done, failed
	Resource	Used for information about the resources (i.e. information such as
		reports, or objects or human means)

The Agility Service of the SocEDA platform adds a Traffic Jam Constraint element into the field situation model (and nothing into the expected situation model). According the matrix {element type × operation}, ∂ = 4 and ∂ > ∂ threshold .

  Go back to warehouse at 10:50 am. The wholesaler has chosen a ∂ threshold equals to 2. At 8:25 am the vehicle is ready to start. Its round begins at 8:30 am. At 8:43 am, the driver arrives to P#53 and leaves at 8:51 am.While he is driving to P#22, he received an alert about a traffic jam due to a street protest in the street just near P#22.

	1. Start from wholesaler's warehouse at 8:30 am,
	2. Arrive to P#53 at 8:45 am and leave at 8:55 am,
	3. Arrive to P#22 at 9:15 am and leave at 9:25 am,
	4. Arrive to P#34 at 9:30 am and leave at 9:40 am,
	5. Arrive to P#11 at 10:00 am and leave at 10:10 am,
	6. Arrive to P#14 at 10:20 am and leave at 10:30 am,
	7.

Acknowledgments

This work has been partially funded by the French Research Agency (ANR) regarding the research project SocEDA (SOCial Event Driven Architecture) (Grant ANR-10-SEGI-013). This project aims to provide dynamic and adaptive workflows to collaborative situations through EDA and CEP. The authors would like to thank the project partners for their advices and comments regarding this work.


L Lecoq Abcd

D Flick

A Plana-Fattori Bcd

O Laguerre

CHARACTERIZATION OF EVAPORATION RATE, TEMPERATURE, VELOCITY AND HUMIDITY FIELDS IN A COLD CHAMBER

Keywords: food, cooling chamber, humidity, temperature, air flow velocity, food safety

In a food factory some pathogenic bacteria can grow even at low temperature, especially Listeria Monocytogenes which resists to temperature around 0°C. Currently, two ways are used to prevent bacterial growth; product treatment (pasteurization, pH decrease) and cleaning. However, this is not sufficient to eliminate all microorganisms in practice. This work takes part of a four-year research project called EcoSec (2013EcoSec ( -2016) ) in France. The practical objective is to propose an ecological and economical way to destroy Listeria Monocytogenes by humidity reduction. The principle idea is based on the dehumidification of the supply air of a plant. To achieve this objective, a method to predict the evaporation rate on wet surfaces in a refrigerated room was developed using experimental and numerical approaches (CFD, Computational Fluid Dynamics). The study was performed inside a cold room (volume 29m 3 ) in which the temperature and velocity of air were controlled. The results show that the air humidity and its local velocity near the wet surface play a major role.

NOMENCLATURE

Cp air Specific heat of air at 5°C (J.kg -1 .K -1 ) C pl Specific heat of the plate (stainless steel) (J.kg -

INTRODUCTION 1.Background

The study of water evaporation inside a production plant is a topic of major interest to prevent Listeria Monocytogenes growth in the product. These bacteria survive especially thanks to water and high relative humidity inside the room. The main objective of this work is to propose optimal operating conditions to reduce the humidity in a production plant and technical solutions to control food contamination. More precisely, the objective of this paper is to develop a numerical approach to predict the water evaporation rate in a cold chamber. This rate varies due to the position inside a cold room because of the heterogeneity of temperature, velocity and humidity inside the chamber. An evaporation model was developed and the predicted results were compared with the experimental ones. After validation, the model will be used to predict the evaporation rate in other conditions unexplored experimentally. This can help to propose operating conditions to completely dry the production plant as fast as possible after cleaning.

Literature review

Many studies have been performed concerning the prediction of the velocity, temperature and relative humidity fields in refrigerated equipment. Laguerre and al. (2009a [1], 2009b [START_REF] Laguerre | Airflow, heat and moisture transfers by natural convection in a refrigerating cavity[END_REF], 2010a [START_REF] Laguerre | Numerical simulation of simultaneous heat and moisture transfer in domestic refrigerator[END_REF], 2010b [START_REF] Laguerre | Study of water evaporation and condensation in a domestic refrigerator loaded by wet product[END_REF], 2012 [START_REF] Laguerre | Heat transfer modelling in a refrigerated display cabinet: the influence of operating conditions[END_REF]) and [START_REF] Hoang | Heat transfer modelling in a ventilated cavity loaded with food product: Application to a refrigerated vehicle[END_REF] [START_REF] Hoang | Heat transfer modelling in a ventilated cavity loaded with food product: Application to a refrigerated vehicle[END_REF] worked into the prediction of this parameters into refrigerated display cabinets, domestic refrigerators and refrigerated trucks. For the numerical approach, the conservation equations of mass, momentum and energy are solved allowing the knowledge of the temperature, humidity, velocity fields. Sometimes evaporation/condensation at the product surface are also taken into account. [START_REF] Delele | Optimization of the humidification of cold stores by pressurized water atomizers based on a multiscale CFD model[END_REF] [START_REF] Delele | Optimization of the humidification of cold stores by pressurized water atomizers based on a multiscale CFD model[END_REF] studied, for example, the influence of a humidification system to avoid fruits to be dried. A multiphase model was developed and the droplets evaporation was taken into account as well as the product and system characteristics (humidification rate, on and off cycle of the cold production system). The ambient conditions used in this study were: velocity (around 0.5m.s -1 ), temperature (20°C, 10°C, 1°C) and humidity (50%, 83.6% and 96%). This study shows a possibility of a CFD approach to analyze the influence of the humidification parameters and the possibility to define the optimal conditions. The coupling of heat and mass transfer equations has been studied in different conditions. [START_REF] Kaya | Numerical modeling of heat and mass transfer during forced convection drying of rectangular moist objects[END_REF] [START_REF] Kaya | Numerical modeling of heat and mass transfer during forced convection drying of rectangular moist objects[END_REF] studied, for example, the heat and mass transfer during a wet object drying. Several geometries were used and the object was placed inside a working area (0.6x0.2 m²) in which a low velocity (0.33m/s) was maintained. 2D simulations (commercial software Fluent) have been achieved to predict the heat and mass transfer coefficients. This study also shows the boundary layer and temperature field around this object. [START_REF] Bower | A study of the Sherwood-Rayleigh relation for water undergoing natural convection-driven evaporation[END_REF] [START_REF] Bower | A study of the Sherwood-Rayleigh relation for water undergoing natural convection-driven evaporation[END_REF] used an electronic balance to measure the evaporation rate of a bulk of water (different widths; 15.2cm, 30.5cm, 45.7cm and 60.9cm) under natural convection with different ambient conditions (22°C to 32°C and 51% to 58% relative humidity). These authors propose a correlation between the mass transfer and airflow ( with a, b and c constants determined experimentally). An important point in the evaporation is also the radius of the droplets or the water film thickness. Kelly-Zion et al. (2013) [START_REF] Kelly-Zion | Correlation for the convective and diffusive evaporation of a sessile drop[END_REF] have shown that the smaller the radius of the droplet is, the smaller the evaporation rate will be. D'Agaro et al. ( 2006) [START_REF] D'agaro | Numerical simulation of glass doors fogging and defogging in refrigerated display cabinets[END_REF] studied the evolution of the shape of a droplet during evaporation. Two periods were observed: one with constant wet surface and decrease of contact angle and one with decrease of wet surface and constant contact angle.

MATERIALS AND METHODS

Cold chamber

The experiments were carried out inside a cold chamber at Irstea were the temperature and velocity can be controlled. The cold chamber is 3.4m long, 3.4m wide and 2.5m high which represent a total volume of 29m 3 (Figure 1). The cooling unit is constituted of two axial fans of 30cm diameter which allow the flow rate to vary from 1000 m 3 .h -1 to 2450 m 3 .h -1 , and the average inlet air velocity from 2 to 4.8 m.s -1 (same order of magnitude as in productions plants). The evolution of temperature and relative humidity inside the chamber was followed using thermocouples and hygrometers. The water was sprayed on a stainless steel plate (Width x Length x Thickness 0.9m x 0.9m x 0.001m) with a 10 cm thickness thermal insulation (extruded polystyrene) underneath and placed on an electronic balance located on the floor. The evolution of water evaporation rate could be followed by the recorded weight. The experimental evaporation rate at different inlet air velocities was compared with the one calculated from the developed model (presented in 2.4.).

Temperature and hygrometry measurement point (a) (b) In order to prevent the influence of the air blowing inside the room from having an impact of the value measured by the balance, the ventilators were stopped every 20 minutes during 2 minutes and the weight was taken during this time.

Measurement devices

Five thermocouples (T type, 1 mm diameter, 0.2°C precision) and hygrometers (3% RH precision) were used to measure air temperature and relative humidity at different positions (Figure 1a). The thermocouples were calibrated at 5 different temperatures (-10°C, 0°C, 10°C, 20°C, 30°C) and the hygrometers at 3 relatives humidity (11.3%, 75.5%, 97.7%). The velocity measurement was performed with a hot wire anemometer (0.03 m.s -1 precision) in front of the fans at different positions.

CFD simulation

In order to understand the phenomena, particularly the airflow pattern inside the cold chamber, simulation of the air flow was performed on Ansys-Fluent 14.5 with the standard k-ε turbulence model and the enhanced wall treatment. The mesh was based mainly on tetrahedral elements with refined prismatic cells nears the wall. It consisted in approximately 9.1*10 5 elements. The inlet air temperature was fixed at 5°C and the upper part of the plate temperature at 10°C, with the thermal insulation underneath of 10 cm thickness. The convective heat transfer coefficient (h) was identified by the simulation and this value was used to develop the evaporation rate model (part 2.4).

Water evaporation rate model

The stainless steel plate was insulated at the bottom using extruded polystyrene. Thus, it can be considered that the exchanges by convection and evaporation with the air take place only on the top part of the wet plate (Figure 2).

Figure 2 -Physical phenomenon on the plate

The equation 1 describes the variation of the water weight due to the evaporation. The equation 2 describes how heat accumulation in the plate is related to the difference between the heat exchange by convection and water evaporation.

With:

As a first approach to predict the evaporation rate, the thermal inertia of the plate is considered negligible. Equation 2 becomes then: By replacing the equation 3 in the equation 1, the evaporation rate (dm/dt) can be calculated knowing the plate temperature.

At the beginning of evaporation, it can be assumed that the droplets or the puddles occupy almost the whole surface (S m =S). Therefore, the plate temperature is supposed to be constant and equal to the wet bulb temperature. The evolution of the water weight can then be predicted by the equations 4 or 5:

It is to be emphasized that for 0°C<T pl <5°C, C sat(Tpl) can be calculated approximately by the following expression (obtained by linear regression): It is to be reminded that the mass transfer coefficient (k) can also be estimated from the convective heat transfer coefficient (h) using the Lewis analogy (equation 7):

Using the equations 3, 6 and 7, the wet-bulb temperature can be approximated by following expression:

So, using together the equations 4 and 8, it is easy to predict the evolution of the mass of water at the beginning of evaporation. It can be observed that, the higher the relative humidity, the lower the evaporation rate. Indeed, for a high humidity the wet-bulb temperature (or plate temperature) tends toward the air temperature. So that, according to equation 4, the mass of water decreases slower for higher air humidity.

RESULTS AND DISCUSSION

Velocity field in the chamber

The velocity field inside the chamber obtained from CFD simulation (3D) is presented in figures 3 and 4 for an average inlet velocity of 4.8 m.s -1 . These figures show the directions of the airflow and the velocity magnitude. It is then possible to examine the position of high or low velocity zones. This is an important point because the velocity plays a major role for the determination of the convective heat transfer coefficient: the higher the local velocity, the higher the transfer coefficient. For the investigated position, the simulation showed that the convective heat transfer coefficient between the plate and the air is 7.6 W.m -2 .K -1 .

Evaporation rate considering water droplets

The value of heat transfer coefficient obtained from the CFD simulation was applied to the model of water evaporation (eq. 4, constant wet surface) for 3 air inlet conditions: 3.0 m.s -1 , 5.4°C , 90% RH -4.0 m.s -1 , 5.4°C, 90% RH and 4.8 m.s -1 , 5.1°C, 85% RH. The comparisons between the experimental and calculated data are shown in figure 5. The first part of the evaporation seems to be predicted correctly while the second part is under estimated by the model because the wet surface area is then decreasing. The influence of air velocity on the evaporation rate is well shown in this figure: the higher the inlet velocity, the faster the evaporation. This can be explained by the fact that a higher air velocity leads to a higher mass transfer coefficient.

CONCLUSION AND PERSPECTIVES

A CFD simulation was carried out to better understanding some phenomena inside a cold chamber, especially for the airflow direction and velocity field. It allowed also determining the convective heat transfer coefficient which was then used for the prediction of water evaporation rate in the chamber. The model equations show the importance of the air relative humidity inside the room on the evaporation rate: the higher the air humidity, the lower the evaporation rate. The numerical and experimental results also showed the influence of airflow rate: the higher the air velocity of the fans, the higher the evaporation rate.

The first part of the evaporation was well predicted while the second part was under estimated by the model. To improve the model, it is necessary to understand the mechanism of droplet surface evolution which requires the development of an experimental setup. This will be carried out in the future. Currently, the inlet air temperature and relative humidity were constants during experiments. These parameters will be varied allowing studying their impact on the evaporation rate at different positions in the chamber.

The result of this study leads to conclude that to avoid humidity stagnation in a production plant, one can increase the airflow rate or decrease the air relative humidity. The later solution is more comfortable for workers compared to the one of high air velocity.

Figure 1 -

 1 Figure 1 -Cold chamber (a) -Side view (b) -Top view

  Conference on Sustainability and the Cold Chain, ICCC 2014 23/06/2014 -25/06/2014, London UK pp. 267-274

Figure 3 -Figure 4 -

 34 Figure 3 -Velocity vector Side view (symmetry plane)

Figure 5 -

 5 Figure 5 -Evolution of the water weight -Comparison with the developed model

3rd IIR International Conference on Sustainability and the Cold Chain, ICCC 2014 23/06/2014 -25/06/2014, London UK pp. 267-274

ACKNOWLEDGEMENT

The research leading to this result has received funding from the French National Research Agency (EcoSec Project, ANR-12-ALID-0005-04)


Rémi Wolf

Josselin Duchateau

Philippe Cinquin

Sandrine Voros

3D tracking of laparoscopic instruments using statistical and geometric modeling

Keywords: laparoscopic surgery, image-based localization of surgical instruments, Condensation algorithm

During a laparoscopic surgery, the endoscope can be manipulated by an assistant or a robot. Several teams have worked on the tracking of surgical instruments, based on methods ranging from the development of specific devices to image processing methods. We propose to exploit the instruments' insertion points, which are fixed on the patients abdominal cavity, as a geometric constraint for the localization of the instruments. A simple geometric model of a laparoscopic instrument is described, as well as a parametrization that exploits a spherical geometric grid, which offers attracting homogeneity and isotropy properties. The general architecture of our proposed approach is based on the probabilistic Condensation algorithm.

Introduction

Minimally invasive surgery (MIS) has become more and more popular with surgeons and the public in the recent years, although its superiority on open surgery has not been systematically proven in terms of per and post-operative complications ( [START_REF] Touijer | Laparoscopic radical prostatectomy: A critical analysis of surgical quality[END_REF]). Furthermore, the mini-invasive approach is more expensive and more challenging technically for the surgeons than open surgery (loss of depth information, loss of tactile information, limited field of view, inverted movements of the instruments due to their insertion in the patient's abdominal cavity through "fixed" points on the abdominal wall).

These limitations highlight the need to develop further new methodologies dedicated to MIS. These methodologies can be divided in two categories: skills (learning, practicing and evaluating) and tools to perform these skills (endeffectors, endoscopic cameras). New approaches to learn MIS have already been developed, especially using virtual reality simulators [START_REF] Sherman | Assessing the learning curve for the acquisition of laparoscopic skills on a virtual reality simulator[END_REF] and video recordings of surgeries for the evaluation of the surgeon's performance [START_REF] Aggarwal | Toward feasible, valid, and reliable video-based assessments of technical surgical skills in the operating room[END_REF]. Extensive research has also been performed on the optimization of laparoscopic tools, the most outstanding being the development of robots, either holding the endoscope inside the abdominal cavity (AESOP R , ViKY R ) or offering more degrees of freedom than the standard instruments (DaVinci R ). Here we will concentrate on the enhancement of the endoscope manipulator. It is bothersome and distracting for the surgeons to control the positioning of the robot manually or by voice. Therefore there has been some research to visually track the tip of the surgeon's instruments and adapt the position of the endoscope to follow them [START_REF] Doignon | Real-time segmentation of surgical instruments inside the abdominal cavity using a joint hue saturation color feature[END_REF][START_REF] Uecker | Automated instrument tracking in robotically assisted laparoscopic surgery[END_REF].

In this paper, we present a new method to localize a laparoscopic instrument from video images in real time, using only a priori knowledge on the geometry of the instrument, the position of its insertion point and a statistical model for the displacement of the instrument based on the Condensation algorithm. Compared to other methods described until today, this algorithm allows us to automatically retrieve the 3D orientation and tip position of the instrument in the camera referential, offering new opportunities such as analysing the motion of the instruments to classify the operating task.

Materials and methods

General framework

The goal of our method is to track the motion of a laparoscopic instrument inside the abdominal cavity without modifying it (ie. without a physical tracking system such as infrared or magnetic localizers), using only information from the endoscopic image and a geometric model of the tool.

After an initial calibration procedure, we determine the insertion point of each instrument in 3D using a priori knowledge about its geometry. We then build a geode centered on the insertion point describing the different orientations that can be taken by the tool. A propagation algorithm using particles evolving over the geode's elements according to a determinist and a stochastic law is then computed (the Conditional Density Propagation or Condensation algorithm [START_REF] Isard | Condensation-conditional density propagation for visual tracking[END_REF]), with measures performed at each time step to guess the most likely orientation of the tool, after reconstruction and projection of the instrument model in the 2D image. Finally, we determine the tip depth along this orientation with an Otsu-based [START_REF] Otsu | A threshold selection method from gray-level histograms[END_REF] segmentation method over a sliding window.

System calibration

Our method requires system calibration: in particular it must be possible to calculate the 3D line in the fixed field of reference (FOR) that corresponds to a 2D image point on the screen. This requires an intrinsic calibration and distortion correction performed using Zhang's procedure [START_REF] Zhang | A flexible new technique for camera calibration[END_REF] with a planar chessboard, which is validated in the medical imaging community.

If the endoscope is displaced by a robot, it is also necessary to know the camera's displacement in the robot's referential, which requires a "hand-eye" calibration of the robotic system. In this paper, we choose to work with a fixed camera, since we concentrate on the presentation of our localization method.

Instrument model

Laparoscopic tools have differents appearances, according to their color or to their tip shape. However, they all have in common a cylindrical shape and a standardized radius, in order to fit in trocars. It is thus possible to model a tool as a finite cylinder of known diameter [START_REF] Cano | Laparoscopic tool tracking method for augmented reality surgical applications[END_REF]. Our search, given a predefined insertion point (see 2.4), is thus limited to two angular parameters (pan and tilt) that define a direction axis, depicted by a 3D unit vector e 1 , and a translational parameter, γ, representing the instrument's depth along its axis (see 2.5).

The knowledge of the axis of the tool and of its insertion point allows us to reconstruct its borders as seen by the camera (see fig. 1, left): the vector CT is the translation vector representing the instrument insertion point T in the camera's FOR (where C is the camera's center of projection). We can construct a referential linked to the tool by defining two new vectors as follow:

e2 = T⊗e 1 T⊗e 1 , e3 = e1 ⊗ e2
The instrument axis can be defined by the line (T ,e1) but also by the line (P ,e1), where CP = T ⊗ e1 e3, as shown in fig. 1, left. The plane (P,e2,e3) is the only right section of the cylinder representing the tool that contains the camera's center of projection C. Working in this right section (fig. 1, right), we define two points S1 and S2 that belong to the tangent planes to the instrument running through the camera center. S1 and S2 satisfy PS1,2 = ρ, where ρ is the instrument's radius and CS1,2 •PS1,2 = 0, i.e. the following system:

αe2

+ βe3 = ρ (αe2 + ( T ⊗ e1 + β)e3) • (αe2 + βe3) = 0 (α, β) ∈ R 2
This system admits two solutions (one for each tangent plane) :

(α, β) = ±ρ 1 - 1 T ⊗ e1 , - ρ 2 T ⊗ e1
L = T ± |α|e2 + βe3 + γe1 γ ∈ R

Determination of the insertion point

Several methods have already been developed in order to retrieve the lines representing the borders of a laparoscopic instrument in 2D images [START_REF] Uecker | Automated instrument tracking in robotically assisted laparoscopic surgery[END_REF][START_REF] Doignon | Real-time segmentation of surgical instruments inside the abdominal cavity using a joint hue saturation color feature[END_REF]. Here we use manual selection of these lines in a sequence of images with one camera position and a moving instrument, which allows us to find (thanks to the calibration and the geometric model of the tool) a beam of 3D lines each representing an instrument's axis. The intersection of these lines, which corresponds to the insertion point, is computed by solving a rank-3 system with a SVD decomposition. We will see in the discussion that this method could easily be automatized.

Instrument localization

Our tracking system is based on the Condensation algorithm [START_REF] Isard | Condensation-conditional density propagation for visual tracking[END_REF]: it uses factored sampling in which randomly generated sets represent the probability distribution of a certain feature. The random set is propagated using a dynamic model and visual observations, which is of particuliar interest in our situation in which background clutter can often cause direct observation to fail. We therefore need to build a specific dynamic model for our process, along with the corresponding observation model. Given our model construction, it is preferable to detect the instrument 3D position in two steps, determining the angular parameters before the translational parameter, since this parameter only influences the length of the apparent segments, and not their orientation.

Geodesic geometry Two parameters are needed to describe the angular orientation of the tool in a 3D space. The typical representation of these parameters would be standard latitude and longitude coordinates (pan and tilt), but these are not well suited to our problem as they are not isotropic (non-uniform distribution) and computational stability near the poles is non-trivial. A solution developed in the field of climate modeling is the use of quasi-uniform spherical geodesic grids [START_REF] Randall | Climate modeling with spherical geodesic grids[END_REF] among which icosahedral grids give almost homogeneous and quasi-isotropic coverage of the sphere. These grids are obtained by the successive subdivisions of an isocahedron (fig. 2a). Parametrization is easily obtained by dividing the geodesic sphere into 5 panels (fig. 2b) composed almost exclusively of hexagonal cells. Each panel, once rotated and twisted, can be represented as a rectangular array (fig. 2c).

Each 3D unit vector that runs from the center of the sphere to a cell on its surface corresponds to a unique cell on one of the 5 rectangular arrays, as illustrated by the hexagonal cell marked in blue on the sphere in fig. 2a, which corresponds to an unique cell of the rectangular array in fig. 2c. This vector can easily be computed given the 2D position of the cell in the rectangular array, and in our algorithm correspondences between 2D cells and 3D vectors are computed and stored in correspondence tables in order to increase speed for the rest of the execution.

Particle model In our model, each particle is described at a given time by five variables: the current and previous 2D coordinates of the particle st = [xt, yt, xt-1, yt-1] expressed in the panel variable pt.

st = Ast-1 + Bwt, where A =     0 0 1 0 0 0 0 1 -0.8 0 1.8 0 0 -0.8 0 1.8     (1)
Matrix A represents the deterministic evolution of the particle, wt is a vector of standard normal variates and BB T is the process noise covariance. The coefficients of A have been set empirically for now, describing an exponentially decreasing speed. The determination of B is achieved by noting that probabilities of random motion are isotropic: the correponding probability density function is described in Eq. 2, where i, j, k represent the three directions on the hexagonal grid (see fig. 2d). Fitting the hexagonal grid on the rectangular array parametrized by x = i and y = k coordinates, the probability density function can be expressed as Eq. 3.

p(i, j, k) ∝ exp(-(i 2 + j 2 + k 2 )) (2) p(x, y) ∝ exp(-(x 2 + y 2 + (x + y) 2 )) (3)
which has a covariance matrix proportional to :

BB T = 1 1 2 1 2 1
Particle dynamics are therefore almost defined for a single panel. The scale factor for the covariance matrix is however unknown and is determined empirically.

Measurements At each time step, the different cell values are extracted using the instrument model. For each particle, the 3D vector v that corresponds to the 2D cell it occupies is extracted from the precomputed correspondence table. The segments that correspond to the projection in the 2D image of the borders of a cylinder of direction vector v passing through the insertion point T are then extracted according to the instrument model previously described.

The image measurements we use are based on contour detection: image derivatives along the x and y axis are approximated by convoluting the image with corresponding Sobel kernels. The two components being extracted, the gradient vector G i is available in each image point. We define a score, calculated for each segment : σj = m j i=1 (G i •n j ), where j = (1, 2), n j is the unit vector normal to the segment j and mj is the number of points in segment j. Finally the score π = |σ1 -σ2| is given to the cell and the particle (n1 and n2 have opposite directions). The instrument axis is then computed as the average of all the particle vectors weighted by their score.

Determination of the instrument depth Once the instrument axis has been computed, we determine the tip position along this axis by looking for a transition point between the instrument and the background. We developed an algorithm using a sliding window moving along the projection of the instrument's axis in the 2D image. For each position of the sliding window, we find the optimal threshold value dividing the pixels into two distinct classes according to Otsu's method [START_REF] Otsu | A threshold selection method from gray-level histograms[END_REF]. We consider that the transition between the intrument's body and its tip is the point where the interclass variance will be maximal. Therefore we store the position of the sliding window for which the inter-class variance is maximal,the position of the middle of this window returning the position of the tip of the instrument.

Results

Our first evaluation parameter is the angular error between the theorical instrument's axis and the computed one, in 3D and in the projected image. First, we tested the validity of our method with simulated data. We built a fictive sequence using an abdominal background image on which we pasted a black mask, corresponding to known 3D positions and orientations of the instrument. Its dynamic law was governed by the same equation as the particles (see Eq. 1).

We then tested our method on a testbench modeling a laparoscopic surgery, and compared our results with 1) manual segmentation of the tool in the image followed by 3D reconstruction in the same fashion as in 2.4 and 2) 3D results obtained using an Optical localizer. Our testbench consisted of a fixed endoscope imaging a surgical instrument inserted through a 5mm trocart, with a background image of an abdominal cavity containing another tool (the blue tool in fig. 4). The angular error for each setup according to the number of particles and geode resolution can be found fig. 3. Fig. 4 shows typical image results that we obtained on our testbench. Working with 768x576 images acquired by a clinical endoscope (no downsampling was performed), a geode resolution of 163842 cells and 500 or 2000 particles, our algorithm computes at the speed of respectively 16 and 8.6Hz with an Intel Xeon 2.67GHz, 3.48Go RAM PC.

Our second evaluation parameter was the 2D tip position in the image. It was compared to the tip position manually clicked, with the error expressed in pixels. Our mean error is evaluated at 27.8 pixels for the same dataset as fig. 3. Our tracking system allows a surgeon to track the motion of his tools during an intervention, in order to either automatically control the movements of the robot, or to evaluate the quality of the surgeon's gesture. As shown by the results, our precision in the 2D orientation of the instrument inside the image is precise, robust and quick enough to allow the first anticipated use. The precision in the 3D orientation and the depth evaluation are still imperfect and should call for several improvements in our algorithm.

Increasing the number of cells in the geode (i.e. its resolution) does not slow down the algorithm as the correspondence table between the cells and 3D orientations is computed only once. However a limit is reached when the order is superior to 128 (163842 cells), as the surface of a cell at a distance of 35cm (the instrument's length) from the insertion point will be 3,13mm 2 , which is approximatively the surface of the instrument tip.

Simulated data highlights the intrinsic error of the method, linked to the statistical and dynamical model, as the insertion point remains fixed, when testbench data shows the error due to the motion of the insertion point as well.

We now plan several improvements of our method to improve the 3D orientation detection and ensure its compliance with clinical conditions: first, we plan to allow the insertion point to move freely around its computed position , as it is the case in real surgery (in an amplitude range described in [START_REF] Voros | Automatic detection of instruments in laparoscopic images: A first step towards high-level command of robotic endoscopic holders[END_REF]), following a white noise or an evolution law which needs to be determined. Secondly, we plan to detect automatically the insertion point at the beginning of the procedure by extracting lines using the Hough transformation (not adapted to real-time tracking) and the same pivot procedure as presented in section 2.4.

Our next step will be to test our method on anatomic specimens mimicking true surgical interventions, in order to improve the deterministic and stochastic matrices describing the evolution of the particles -currently empirically evaluated -and to study

Fig. 1 :

 1 Fig. 1: Description of the tool model

Fig. 2 :

 2 Fig. 2: Representation of the geodesic grid: (a) The spherical geodesic grid (b) Cutting the geodesic sphere in 5 panels. (c) Twisted and rotated panels can be represented in rectangular arrays. (d) A particle evolves freely in three independent directions in the hexagonal tiling. Note that j = i + k. (images source http://kiwi.atmos.colostate.edu/BUGS/geodesic)

Fig. 3 :Fig. 4 :

 34 Fig. 3: Results for 1050 consecutive frames, using 2000 particles and a geode resolution of 163842 cells covering the 3D sphere

their influence on the results. Improvements will also be brought to the tip detection in order to include an evolution model similar to that used for the axis orientation, thus smoothing aberrant measurements which pollute our mean result.

This work has been supported by French National Research Agency (ANR) through TecSan program (project DEPORRA ANR-09-TECS-006)


Sandrine Voros
email: sandrine.voros@imag.fr

Jean-Alexandre Long

Philippe Cinquin
email: philippe.cinquin]@imag.fr

Automatic localization of laparoscopic instruments for the visual servoing of an endoscopic camera holder

The use of a robotized camera holder in laparoscopic surgery allows a surgeon to control the endoscope without the intervention of an assistant. Today, the orders that the surgeon can give to robotized camera holders remain limited. In order to provide higher level interactions between the surgeon and a robotized camera holder, we have developed a new method for the automatic tracking of laparoscopic instruments which works in near real-time. The method is based on the measurement of the 3D positions of the insertion points of the instruments in the abdominal cavity and a simple shape model of the laparoscopic instruments. We present the results of our first experimentation on a cadaver.

Introduction

Laparoscopic surgery offers a lot of benefits for the patients compared to open surgery: less blood loss, less transfusions, a smaller consumption of analgesia and a shorter hospitalisation time [START_REF] Makhoul | Laparoscopic radical nephrectomy for T1 renal cancer: the gold standard? A comparison of laparoscopic vs open nephrectomy[END_REF]. However, because of the limited field of view, the specific instruments and the difficult coordination between the assistant manipulating the camera and the surgeon, laparoscopy is a more complex procedure for the surgeon than open surgery. Robotized systems were developed to assist the surgeon in performing laparoscopy, including robotized camera holders: they stabilize the image [START_REF] Kavoussi | Comparison of robotic versus human laparoscopic camera control[END_REF], reduce the staining of the endoscope and enable solo surgery. Several robotized endoscopic holders have been commercialised, including Aesop (Computer Motion) and EndoAssist (Armstrong Healthcare) [START_REF] Ballester | Comparison of task performance of robotic camera holders EndoAssist and Aesop[END_REF].

Today, the applications allowing the surgeon to control robotized endoscopic holders only integrate simple commands. Providing higher level commands of the robot would be a major clinical added value: it would allow the surgeon to concentrate on his surgical procedure, rather than on the control of the robot. In this objective, the automatic tracking of the instruments is an intuitive approach, since the moving instruments must stay in the field of view of the camera for security reasons.

Several approaches were proposed to detect the instruments in laparoscopic images: [START_REF] Krupa | Autonomous 3-D positioning of surgical instruments in robotized laparoscopic surgery using visual servoing[END_REF] designed a robotized instrument holder equipped with a laser beam. The laser beam is detected in the images in order to track the tip of the instrument. This approach raises the questions of cost issues and sterilisation. [START_REF] Wei | Real-Time Visual Servoing for Laparoscopic Surgery. Controlling Robot Motion with Color Image Segmentation[END_REF], [START_REF] Zhang | Application of Visual Tracking for Robotic-Assisted Laparoscopic Surgery[END_REF]and [START_REF] Tonet | Image analysis-based approach for localization of endoscopic tools[END_REF] taped a color mark on the instrument. The main limit of this approach is the sterilisation of the mark. These four approaches put constraints and add complexity to the procedure. [START_REF] Wang | A new framework for vision enabled and robotically assisted minimally invasive surgery[END_REF] used color classification to identify colors that are likely to correspond to instruments, but the computation time is not indicated. [START_REF] Doignon | Detection of grey regions in color images : application to the segmentation of a surgical instrument in robotized laparoscopy[END_REF] proposed to detect the grey regions in the laparoscopic images, which are likely to correspond to the metallic tip of the instrument. A good computation time is achieved. However, some laparoscopic instruments do not have a metallic tip and the tip might be hidden behind an organ. Finally, [START_REF] Climent | Automatic instrument localization in laparoscopic surgery[END_REF] proposed an approach based on the detection of the lines in the image, but some parameters have to be tuned by the user. None of these approaches give solutions to identify the instrument that has to be tracked among the other instruments.

Our approach is based on the measurement of the 3D positions of the insertion points of the instruments. With a calibrated robot and a calibrated camera, we can project the insertion point of the instrument we wish to track on the image plane of the camera (even if the insertion point is not visible on the image). This information provides an important constraint on the position of the instrument on the image. It is used to automatically detect the axis, edges and tip of the instrument. The knowledge of the 2D position of the tip allows to perform a visual servoing of the orientation of the camera. For this purpose, we do not require subpixel precision but rather need to roughly determine the position of the tip. As suggested in [START_REF] Tonet | Image analysis-based approach for localization of endoscopic tools[END_REF], the orientations of the edges can be used to find the depth of the instrument and include the zoom in the visual servoing loop. [START_REF] Voros | Automatic detection of instruments in laparoscopic images: a first step towards high level command of robotized endoscopic holders[END_REF] showed encouraging results on laparoscopic images and presented our visual servoing on a test bench. In this paper, we focus on the results of our first experimentation on a cadaver.

Method

Calibration of the robotized laparoscopic holder and the camera

Our robotized camera holder [START_REF] Berkelman | Development of a Compact Cable-Driven Laparoscopic Endoscope Manipulator[END_REF] was calibrated, as well as the endoscopic camera [START_REF] Zhang | A Flexible New Technique for Camera Calibration[END_REF]: the 3D position of the camera is known in a fixed frame linked to the robot (fig. 1). The 3D coordinates of a point can be measured by stereovision: two views of this point for different positions of the robot are sufficient.

Measurement of the 3D positions of the insertion points of the instruments

The 3D positions of the insertion points of the instrument are measured at the beginning of the intervention. Since the surgeon gives orders to the robot with a vocal command, we have developed a 'vocal mouse' which allows to move a cursor on the endoscopic image. For two different positions of the camera in which the insertion point is visible, the surgeon selects the position of the insertion point in the image with the vocal mouse. This initialisation step could be easily integrated to the surgical protocol since the surgeon creates the insertion points at the beginning of the intervention under visual control. We show in section 3 our measurements of the 3D position of an insertion point: the hypothesis that they are relatively fixed during an intervention is validated.

Model of the laparoscopic instruments in the images

We define P = proj(T ) as the projection of an insertion point T in the image plane of the camera. In a laparoscopic image, an instrument can be represented as (fig. 2) :

three lines: a symmetry axis ∆, and two edges which are symmetrical compared to ∆. The positions of these three lines are constrained by the position of P : C is a circle centered on P of radius R such that the symmetry axis and the two edges intersect C ; a point S, which represents the tip of the instrument. S must belong to ∆. To detect the instrument inserted through T, we search for a symmetry axis and two edges which intersect C of radius R, and we search for the tip S of the instrument on the symmetry axis.

Main steps of the detection algorithm

We suppose that the surgeon has selected an instrument to track (by selecting an insertion point T ). The projection P of the insertion point on the image plane is computed.

Step 1: segmentation:

In this step, we compute the edge points in the image by a simple gradient method. The orientation of the edges corresponding to the tool are constrained by the position of P . After this step, we obtain a cloud of edge points, most of which correspond to the instrument.

Step 2: search for the symmetry axis of the instrument: We use the edge points computed in the previous step to find the symmetry axis of the instrument. To do this, we use a Hough method [START_REF] Duda | Use of the Hough Transformation To Detect Lines and Curves in Pictures[END_REF] centered on P . This allows to consider only the lines that have a satisfying orientation, and the computation time is very reasonable compared to a classical Hough method.

Step 3: search for the edges of the instrument:

All the edge points which voted for the symmetry axis were stored in the previous step. The symmetric of the cloud of edge points compared to the symmetry axis is computed to obtain more points, and linear regressions are performed for all the points above the axis and below the axis.

Step 4: search for tip of the instrument:

The tip of the instrument is searched on the symmetry axis. The points of the axis are ordered, according to the position of the insertion point. An Otsu threshold [START_REF] Otsu | A threshold selection method from gray level histograms[END_REF] is used to separate the points of the axis in two classes. The pixels of the axis belonging to the same class are grouped into zones. The difficulty of this step is to determine which of the two classes corresponds to the instrument, and to group the zones correctly in order to find the tip. We developped two approach, one for the first image of the tracking, the other for the next images.

For the first image, we search for the longest zone of each class and select the class corresponding to the instrument according to its postion with regard to the tip. The average color of the zone corresponding to the tool and the average color of the zone corresponding to the instrument are stored for the next images. In the next images, we use this color information to determine which zone corresponds to the instrument and to the background.

Once the zone corresponding to the instrument has been found, we look if it can be concatenated with other zones corresponding to the tool according to color and length considerations.

Results

The method was tested on anatomical pieces. Our first step was to validate our hypothesis that the insertion points are fixed during the intervention. We measured along time the position of the insertion point of an instrument (section 2.2). To check that the first computation of the insertion point was satisfying, we projected the insertion point on the image plane (fig. 4). The distance between the measured points and the reference was inferior to 5mm except for measure n • 3, but this was due to leak in of gas in the abdominal cavity. As soon as the leak was stopped the results were satisfying. We were able to perform the tracking of an instrument when the specular reflexions were limited. The computation time is around 100 ms for a 200x100 image. It allows a smooth tracking when the instrument is moved at a normal speed. It could still be reduced by optimizing the code allowing to take into account abrupt and large movements. A resolution of 200x100 pixels is enough for the tracking of a laparoscopic instrument, since we are not looking for subpixel precision to control the camera. We considered that the tool was centered as long as the distance between the detected tip and the previous tip was less than 11 pixels. Fig. 5 shows the distance in pixels between the tip found by the method and the tip manually selected. In 70% of the images, the error is inferior to 5 pixels, and in 87 % the error is less than 11 pixels. We analysed the images for which the detection was wrong (fig. 6): fig. 6 a) corresponds to the first error peak (frame 9). The detection of the axis is correct, but the tip is wrongly positioned: the red segments correspond to the pixels labelled as pixels belonging to the instrument and the yellow segments as pixels labelled as the background. Since the background close to the tip is very dark, it is wrongly labelled by the Otsu threshold as a region belonging to the instrument. This is also the cause of the error for frames 35, 36, 61, 63 and 64.

Fig 6 b) corresponds to the second error peak (frame 21). Again, the axis is correctly detected, but not the tip. This time, the error is due to specular reflexions, which are labelled as background. This is also causing errors for frames 67 and 69. The important error for frame 37 is caused by a false detection of the axis of the tool (fig. 6 c)): the edges of the instrument are barely visible in the image because of the lack of contrast, and the organ on the far hand side of the image provides a lot of edges, which deteriorates the axis detection.

When the specular reflecions are important, they do not deteriorate the detection of the axis, since they are along the axis of the instrument. However, they cause the failure of tip detection, because the wrong classes are attributed to the instrument and background.

Discussion

This experiment allowed us to validate our approach by showing that the insertion point of the instruments are relatively fixed during the intervention. We were able to robustly detect the axis of an instrument, except in cases for which the tool is in a border of the image with a low contrast, or very close to the insertion point. To enhance further this detection, we could give more credit to radial lines than lines with other orientations, since the instruments are usually directed towards the center of the image, as suggested by [START_REF] Climent | Automatic instrument localization in laparoscopic surgery[END_REF]. The detection of the edges is also quite robust, and we plan to integrate the zoom in the visual servoing of the robot.

The detection of the tip works well on images with little specular reflexions. It must be improved to deal with them. We are currently investigating a solution in which the specular reflexions along the axis of the instrument would be detected and removed. Since we detect correctly the axis, we could perform this detection only on a region of interest around the axis, with a reasonable computation time.

To further enhance the speed of the detection, we also plan to take into account the position of the tip of the instrument in the previous images. However, considering the very strong magnification of an endoscope, the instrument might very well move a lot in the image between two frames, so the region of interest (ROI) used must be carefully chosen. Moreover, if the ROI is too small, only a tiny portion of the instrument will be visible, and our shape model will not be accurate, but we consider making a more complete model of the tool.

Finally, we must avoid as much as possible false detections. Our method deals partly with them: if the amount of edge points which voted for the axis is too low, or if the length of the instrument is too small an error is raised. If the tool was correctly detected in the first image, its color and the background color are stored for the next images. Thus, if the detection algorithm finds a tool which color is very close to the color of the background, an error is raised. We must now find a method to reject a false detection when the image is the first image of the tracking. One solution could be to gather a few images of the abdominal cavity without an instrument at the beginning of the intervention to obtain information about the color of the background.

Conclusion

We developed a novel method to detect the instruments in laparoscopic images. The method is automatic and works almost in real-time. Our first cadaver test showed encouraging results and allowed us to find the difficult cases still to solve. Several high level control modes based on this method could be implemented to supervise surgical tasks such as the insertion/removal of an instrument or suturing tasks.

Fig. 1 .

 1 Fig. 1. Calibrated system: if the coordinates of P are known in R0, its projection on the image plane can be computed thanks to the calibration matrixes T (geometric model of the robot) and C (pinhole model of the camera).

Fig. 2 .

 2 Fig. 2. Representation of the instruments in a laparoscopic image. The position of each instrument in the image is constrained by the position of the projection of its insertion point in the abdominal cavity.

Fig. 4 .

 4 Fig. 4. Variation of the position of the insertion point along time. On the right image, the insertion point was projected on the image plane.

Fig. 5 .Fig. 6 .

 56 Fig.5. Error in pixels between the tip found by the method and the tip selected manually and an example of correct detection (yellow dots: candidate points, red points: points corresponding to the axis, blue dot on the axis: tip of the instrument).

Acknowledgement: the authors would like to thank Bruno Thibaut from the anatomy laboratory for his help during these experiments.


Thomas Durieux

Youssef Hamadi

Martin Monperrus

Production-Driven Patch Generation

We present an original concept for patch generation: we propose to do it directly in production. Our idea is to generate patches on-the-fly based on automated analysis of the failure context. By doing this in production, the repair process has complete access to the system state at the point of failure. We propose to perform live regression testing of the generated patches directly on the production traffic, by feeding a sandboxed version of the application with a copy of the production traffic, the "shadow traffic". Our concept widens the applicability of program repair, because it removes the requirements of having a failing test case.

I. INTRODUCTION

Program repair requires the presence of a failing test case to reproduce a failure that has happened in production. Writing such a failing test case is a really hard task, because the developer in charge of reproducing a failure has little access to the system state at the point of failure (she basically only has logs). The difficulty of reproducing production failures has a direct impact on applicability of program repair: with no failing test, there is no patch generation. We aim at weakening the requirements of program repair by removing the mandatory presence of a failing test case.

Our intuition is to perform program repair directly in production, so that the repair process has a direct access to the system state at the point of failure. This paper presents an architecture, called Itzal, it generates patches without requiring a failing test case. The process of Itzal is as follows. First, Itzal uses production assertions or runtime exceptions to detect failures. Second, right after the failure is detected in production, a patch is searched in a sandboxed environment that mimics the production one. If a patch fixes the failure, it is a "candidate patch". Third, the patches are tested for regression, directly in production, based on traffic that is an exact copy of the production traffic -we call it shadow traffic. Itzal has been realized in a prototype implementation for Java which focuses on generating source code patches for null dereferences.

This is a new line of research in automatic repair. Compared to classical test-suite based patch generation (e.g. [START_REF] Le Goues | Genprog: A generic method for automatic software repair[END_REF]), Itzal does patch generation online, i.e. as soon as the failure happens, with no need for reproducing the failure. Yet, Itzal is not a classical runtime repair technique either (e.g. [START_REF] Rinard | Enhancing server availability and security through failureoblivious computing[END_REF]): while the patches are generated online in production, the system state is never altered. The Itzal patches are applied later, once the developer has validated them.

To sum up, our contributions are:

Figure 1: The Architecture of Itzal. The key idea is to duplicate production traffic via a "shadower", the duplicated traffic is used to search for patches and to validate candidate patches.

• Itzal, an architecture for patch generation in production that does not require a failing test case. • The use of shadow production systems and shadow traffic in the context of automatic repair to generate patches in production.

• The design and implementation of a Java implementation of this vision for null pointer exceptions.

This paper is based on content from Arxiv's document #1609.06848 [START_REF] Durieux | Production-driven patch generation and validation[END_REF] and is structured as follows. Section II presents Itzal. Section III presents the related works and Section IV concludes.

II. ITZAL

We now present Itzal, a novel software repair technique for generating patches without requiring failing test cases, directly in the production environment. Itzal is a "testless" patch generation approach.

A. Intuition

The intuition behind Itzal is twofold. First, one can use production runtime contracts to drive the generation of source code patches. This includes classical pre-and post-conditions and implicit contracts such that an accessed variable must not be null. The latter is important because the violations of those implicit contracts come for free in any modern runtime, usually under the form of runtime exceptions.

The second intuition is that one can use the diversity of the production inputs to perform in-the-field regression testing on the synthesized patches. This has the advantage that the regression exactly corresponds to what actually matters.

B. Architecture

The Itzal architecture is composed of seven components, as shown in Figure 1.

1) The Unmodified Application (see Section II-B1) is the application onto which automatic patch generation is plugged. Algorithm 1 shows the interactions between each component of the Itzal. Itzal receives the request from the client (line 1). Then it redirects the request to the Unmodified Application (line 2). Once the request has been handled by the Unmodified Application, the response is sent back to the client, yet it is additionally sent to the Request Oracle Service (arrow a in Figure 1) which verifies the viability of the output (line 4). If the Request Oracle Service determines that there is a failure, the request is sent to the Patch Synthesis Service (arrow b in Figure 1 and line 5). The patches generated by Patch Synthesis Service which pass the Request Oracle Service (i.e. fix the failure at hand) are sent to Regression Assessment Service (line 6). If the request has succeeded (no failure on the original application), the request is also sent to the Regression Assessment Service (line 8) where all the previously generated patches are being validated on-the-fly against the new request. When the Regression Assessment Service has identified valid patches with no regression, it sends them to the Patch Reporting Service.

To sum up, Itzal does patch generation online, i.e. as soon as the failure happens, directly in production. However, while the patches are generated online, they are applied later, once the developer has validated them. The side effects of patch search or regression testing on the production state are completely sandboxed, with no interference with the production environment.

1) The Unmodified Application: Itzal augments a production application with automatic patch generation capabilities. The requirement to deploy Itzal is that the application must use requests, ie. must have a message-driven architecture. A web The type of request may vary between applications, for example a request in a web application will be the request sent by a user's browser to a web-server, in a micro-service application, the request will typically be a REST message, in a mobile application, a request would be a touch event triggered when a user touches a mobile device's screen.

2) Request Oracle Service: The responsibility of the Request Oracle Service is to verify whether the application has succeeded to answer the request. For instance, in a web-server, the Request Oracle Service can check the HTTP request return code ("assert response code != 500 (internal server error)"), of check the presence or not of an exception. Itzal works with generic oracles such as checking the absence of exceptions (e.g. in a web request container or in a thread monitor), and it can also work with domain-specific oracles written by software engineers on top of domain concepts and data (e.g. the returned XML must comply with a specific schema). When possible, the Request Oracle Service provide some information about the failure to help the Patch Synthesis Service to search for a patch. The information provided by default is the stack trace when the failure is based on an exception.

The Itzal does not require a perfect Request Oracle Service, i.e. the Request Oracle Service may miss some failures (false negatives). In the case of false negatives, when the Request Oracle Service misses the failure detection, Itzal simply does not generate patches: this is unfortunate but it impacts neither the original application nor the patches generated for the other failures. In the case of false positives, when the Request Oracle Service detects a failure when there is no failure in reality, the Patch Synthesis Service would generate patches, yet they would be lilkely benign if they pass regression testing done by the Regression Assessment Service.

3) Patch Synthesis Service: The Patch Synthesis Service is the service that synthesizes patches that fix a failing request. Itzal can work with any patch synthesis approach compatible with the Request Oracle Service. The patches are applied to the failure point, for instance at the line where an exception has occurred. For a given failure point, the Patch Synthesis Service performs an exhaustive application of all possible patches.

For each tentative patch, Patch Synthesis Service calls the Request Oracle Service (arrow c in Figure 1) to verify that the request has been correctly handled by the patch template under consideration (the failure has been fixed). Because the Patch Synthesis Service generates the patches only based on one request (the failing one), the patches may break the behavior of the application for other requests, in other word, the may introduce a regression. Thus, if the patch is successful on the failing request, the corresponding patch is transferred to the Regression Assessment Service (arrow d in Figure 1) that will further validates its correctness based on other requests.

The application and execution of candidate patches can change the state of the application in runtime. Consequently, each execution is done in a sandboxed environment, this nullifies the potential side effects of the request or of the patch templates. The sandboxed environment contains a shadow state of the application, which is regularly synchronized with the production one. Since the space of patches is sometimes large, Itzal uses a time budget. It explores the patch alternatives sequentially until they are all explored or until the time budget is consumed.

Itzal can work with any patch model, whether domainspecific (such as out-of-bounds exception) or generic (à la Genprog [START_REF] Le Goues | Genprog: A generic method for automatic software repair[END_REF]). Similarly, Itzal can be applied to binary code if the patch synthesis technique supports it. Our current prototype generates patches for null dereferences. If the patch model generates too much patches, i.e. the search space is too large, this would be a problem because it would represent a huge computation effort on the Patch Synthesis Service and much more importantly on the Regression Assessment Service.

4) Regression Assessment Service:

The patches generated by the Patch Synthesis Service can introduce regressions because their generation only involves one request (the failing one). The Regression Assessment Service has the responsibly to check the behavior of the application when the generated patches are injected against other requests. It detects these regressions by comparing the output of the Unmodified Application against the output of the patch-augmented application. If the output is different, it means that the patch has introduced a regression, and the patch is consequently marked as invalid. This comparison is done on-the-fly, directly on production traffic. Doing regression testing "live" has the advantage that there is no need to record the potentially enormous amount of production data.

5) Regression Oracle:

The Regression Oracle compares the output of the Unmodified Application (arrow g in Figure 1) and the output of a patched version in the Regression Assessment Service for the same request. If the outputs are different, the Regression Oracle marks the current patch as invalid. For example, a regression oracle for a web server compares the HTML text of both versions. The comparison is not necessarily a byte-to-byte one, it can include heuristics to discard transient information such as time, cookie identifiers, etc.

6) Shadower: The role of the Shadower is to create shadow traffic from actual end-user traffic coming into the application. The "shadow traffic" is made of production requests that are duplicated one or several times and sent to sandboxed shadow applications. In our case, the shadow applications are the Patch Synthesis Service and Regression Assessment Service.

In Itzal, the Shadower receives the requests from the clients duplicates them and sends one duplicate to each service of the architecture (arrows a, b, and f in Figure 1). The response is also shadowed for the regression oracle service (arrow g in Figure 1).

In the context of web applications, the concept of running multiple instances of an application is well known and heavily used: this is done for load balancing and rolling deployment. The difference between a load balancer and a Shadower is twofold: first, a load balancer does not duplicate the traffic; second, a load balancer does not send requests to sandboxed "sinks" as Itzal does.

Since Itzal is a production technique, it must have a reasonable impact on the performance of the application. In order to minimize the impact on the Unmodified Application, Itzal computes the Regression Assessment Service and the Patch Synthesis Service asynchronously. Indeed, the goal of Itzal is to perform patch generation, not automatic error recovery system. Hence, the Shadower directly sends the output as soon as the Unmodified Application has handled a request (even if there is a failure). Itzal does not have to wait for the end of the patch search or the end regression testing for sending the response back to the client. The Shadower is thus the only component that impacts the performance of the Unmodified Application. The Shadower requires to 1) copy and reroute requests on the fly and 2) maintain an appropriate shadow state of the system under consideration. In a typical HTTPbased setup, the cost of the former is similar to that of classical web proxies and load balancers, which are extensively used in production systems. The latter point is more an open question: very few works study production state shadowing, neither in academia nor in industry. It may be a costly operation if databases are naively copied for instance. However, we envision piggy-backing on the latest advances in efficient online backups and copies of file systems.

7) Patch Reporting Service: The Patch Reporting Service is the service that communicates the results of Itzal to the developers (arrow i in Figure 1).

It happens that multiple patches (corresponding to multiple patch templates) successfully pass the regression test over production traffic. Consequently the Patch Reporting Service has to sort the patches in order to first propose the most useful ones to the developers. To sort, the Patch Reporting Service uses the number of execution of the patched line in the Regression Assessment Service (the number of requests that execute the patch). The idea is that the more a patch has been executed by the Regression Assessment Service, the less likely it is to introduce a regression.

We now discuss the reporting medium to the human developer. There are several types of communication that can be used in the Patch Reporting Service. In the current prototype, we have a dashboard where the developers follow in real time the failures, the generated patches and the progression of the patch validation.

C. Prototype Implementation for Java

We have implemented a prototype of Itzal for Java in a tool named Itzal4j, dedicated to reactive applications based on HTTP. Itzal4j generates patches for null dereference failures. In Itzal4j, the Request Oracle Service is based on exceptions. Any uncaught exception happening during the processing of a request is considered as a failure. The Patch Synthesis Service is dedicated to null pointers and uses the NPEFix technique [START_REF] Durieux | Dynamic patch generation for null pointer exceptions using metaprogramming[END_REF] for searching the space of possible patches for null dereferences. Sandboxing of patch search is achieved using Docker, a major software containerization platform which provides powerful sandboxing (both disk and IO based). In our implementation, the Patch Synthesis Service sends all candidate patches to the Regression Assessment Service using an HTTP-based protocol. For the Regression Oracle, we compare the body of the HTTP response of the Unmodified Application against the output produced by the patched application (e.g. the HTML body text). If the outputs match, the patch is considered validated for the current request otherwise the patch is permanently marked as invalid. The Shadower is implemented on top of a HTTP proxy implementation in Java called "Jetty Proxy". The Patch Reporting Service of Itzal4j is a web dashboard, where the developers can access in real time the current patches of Itzal: the ones that have fixed at least one failure and the patch that are under regression testing. For each patch, they can visualize the number of failures of the system detected by the Request Oracle Service, see the actual patch code, and the patch validation metrics such as the number of executions done by the Regression Assessment Service.

III. RELATED WORK

Our work is much inspired by the classical work on runtime repair. Rinard et al. [START_REF] Rinard | Enhancing server availability and security through failureoblivious computing[END_REF] present a technique called "failure oblivious computing" to avoid illegal memory accesses by adding additional code around each memory operation during the compilation process. Assure [START_REF] Sidiroglou | Assure: automatic software self-healing using rescue points[END_REF] is a self-healing system based on error-virtualization. Long et al. [START_REF] Long | Automatic runtime error repair and containment via recovery shepherding[END_REF] proposes the concept of "recovery shepherding" in a system called RCV. Those techniques do not produce patches and do not perform regression testing in production.

Gu et al. [START_REF] Gu | Automatic runtime recovery via error handler synthesis[END_REF] presents Ares a runtime error recovery for Java exceptions using JavaPathFinder (JPF). The two majors differences with Itzal4j are: first Itzal4j is safer, it does not modify the production state of the application as Ares does, and secondly, while Ares performs runtime repair, Itzal4j produces source code patches that are then communicated to the developers.

Perkins et al. [START_REF] Perkins | Automatically patching errors in deployed software[END_REF] propose ClearView, a system for automatically repairing errors in production. Itzal and ClearView both perform repair in production, yet they are very different: 1) ClearView does not produce source code patches while Itzal does; 2) ClearView modifies the production state, while Itzal only modifies the sandboxed shadow requests and state (this means that ClearView can mess up the application while Itzal never does so); 3) ClearView works with learned invariantbased oracles, while Itzal uses human designed request oracles.

The concept of shadow traffic is related to the execution of multiple versions of the same software in parallel, called in the literature "multi-version execution" [START_REF] Hosek | Safe software updates via multi-version execution[END_REF], or "parallel execution" [START_REF] Trachsel | Variant-based competitive parallel execution of sequential programs[END_REF]. However, none of the related work uses shadow traffic to generate patches.

IV. CONCLUSION

In this paper, we have presented Itzal, an approach for generating patches in production. The failure detection that triggers the patch search is achieved with runtime assertions, and the regression assessment is done on live production traffic. In Itzal, the patch search is done in a fully sandboxed environment, with no interference with the production data. Our future work now consists of devising an approach to fully automatically store a shadow of the production state, and to efficiently synchronize the shadow state with the actual production state.

  Algorithm 1 The main Itzal algorithm Input: A: the Unmodified Application Input: G: the Patch Synthesis Service Input: V: the Regression Assessment Service Input: O: Request Oracle Service 1: while new request r client from Client do

	2:	out put = A(r client )
	3:	send out put to Client
	4:	if O(out put) is failure then
	5:	patches = send r client to G
	6:	push patches to V
	7:	else
	8:	send r client to V for regression
	9:	end if
	10:	if ∃ validated patches ∈ V then
	11:	p ← order the patches
	12:	report p to developers
	13:	end if
	14: end while
	application, or a web REST service are examples of message-
	driven applications.