E Coupeté
email: eva.coupete@mines-paristech.fr

S Manitsaris

F Moutarde

Real-time recognition of human gestures for collaborative robots on assembly-line

Keywords: Technical gestures recognition, collaborative robotics (cobotics), factory assembly-line

We present a framework and preliminary experimental results for real-time recognition of human operator actions. The goal is, for a collaborative industrial robot operating on same assembly-line as workers, to allow adaptation of its behavior and speed for smooth human-robot cooperation. To this end, it is necessary for the robot to monitor and understand behavior of humans around it. The real-time motion capture is performed using a "MoCap suit" of 12 inertial sensors estimating joint angles of upper-half of human body (neck, wrists, elbows, shoulders, etc…). In our experiment, we consider one particular assembly operation on car doors, which we have further subdivided into 4 successive steps: removing the adhesive protection from the waterproofing sheet, positioning the waterproofing sheet on the door, pre-sticking the sheet on the door, and finally installing the window "sealing strip". The gesture recognition is achieved continuously in real-time, using a technique combining an automatic time-rescaling similar to Dynamic Time Warp (DTW), and Hidden Markov Model (HMM) for estimating respective probabilities of the 4 learnt actions. Preliminary evaluation, conducted in realworld on an experimental assembly cell of car manufacturer PSA, shows a very promising action correct recognition rate of 96% on several repetitions of the same assembly operation by a single operator. Ongoing work aims at evaluating our framework for same actions recognition but on more executions by a larger pool of different human operators, and also to estimate false recognition rates on unrelated gestures. Another interesting potential perspective is the use of workers' motion capture in order to estimate effort and stress, for helping prevention of physical causes of some musculoskeletal disorders.

Introduction

Recent advances in motion capture technologies make it now easier to perform real-time monitoring of human activities, for various purposes. Meanwhile, a current trend in manufacturing robotics is the development of collaborative robotics (cobotics), in which robots perform activities jointly, or at least side-by-side, with human operators. As already highlighted long ago by [START_REF] Inagaki | Behavior-based intention inference for intelligent robots cooperating with human[END_REF], for humans and robots to have a common goal and work cooperatively, human intention inference by robots is required. This in turn necessitates robots to be able to recognize human actions. Prototyping a framework (from sensors to dataprocessing) for real-time human technical gestures recognition in the context of automotive assembly line is precisely one of our research goals within the "Chaire PSA Peugeot-Citroën on Robotics and Virtual Reality". Technical gestures are invented by a specialist minority for use strictly within the limits of their particular activity. These gestures are meaningless to anyone outside the specialization, and make sense only in their narrow operating field. In automotive assembly lines and more precisely in the experimental cell of PSA Peugeot Citroen, technical gestures can be "to remove the sealing sheet", "to fit the windows sealing sheet", "to screw up without tightening", "to hammer" etc. Recognition of human gestures and activities has become an important research area with numerous potential applications including sign language interpretation, automated surveillance and human robot collaboration. In order to capture human motions, several systems can be used, as illustrated on Figure 1. The oldest and most commonly used is the video from RGB cameras. More recently the apparition of real-time depth cameras, like the Microsoft Kinect, brought new possibilities for human motion capture and monitoring. Thirdly, inertial motion sensors (IMS) like MotionPod of Movea (see http://www.movea.com/), or the motion capture suit developed by AnimaZoo (see http://www.animazoo.com/), enable the direct acquisition of motion information in real time. When using information provided by RGB camera, in order to be able to perform gesture recognition, some features have to be extracted. In [START_REF] Laptev | Space-time interest points[END_REF] and [START_REF] Schuldt | Recognizing Human Actions: A Local SVM Approach[END_REF], authors introduce space time interest points, features computed in a spatio-temporal domain. More specifically to human, [START_REF] Gorelick | Actions as space-time shapes[END_REF] use body motion as a space time feature. In [START_REF] Laptev | Retrieving actions in movies[END_REF] authors use optical flow histograms and spatial gradient histograms. In [START_REF] Lv | Single View Human Action Recognition using Key Pose Matching and Viterbi Path Searching[END_REF] authors compare silhouettes extracted in the frames with key poses, using an extension of the PMK (Pyramid Match Kernel) described in [START_REF] Grauman | The Pyramid Match Kernel: Discriminative Classification with Sets of Image Features[END_REF]. Similarly in [START_REF] Sullivan | Recognizing and Tracking Human Action[END_REF] authors compare key frames to the frames computing shape equivalences. With inertial motion sensors, the output is the angles of joints, which can directly be used as features for gesture recognition. However, there are several possible choices for rotational representation: for instance Euler angles as in [START_REF] Piovan | On Coordinate-Free Rotation Decomposition: Euler Angles About Arbitrary Axes[END_REF], quaternions as in [START_REF] Bachmann | Orientation Tracking for Humans and Robots Using Inertial Sensors A Quaternion Attitude Filter[END_REF] or Direction Cosine Matrix (DCM) as in [START_REF] Bar-Itzhack | Direction cosine matrix estimation from vector observations using a matrix Kalman filter[END_REF]. As for depth cameras, a skeleton posture has to be estimated first, before extracting the joint angles. For example, in case of the Microsoft Kinect, a default skeletonization is performed by the provided SDK. It enables real-time pose estimation, restricted however to human standing up and facing the sensor. This limitation is because this skeletonization relies on a randomized decision forest [START_REF] Lepetit | Randomized Trees for Real-Time Keypoint Recognition[END_REF] trained on a database of depth images of human poses, very large but limited to standing-up and facing postures [START_REF] Shotton | Real-time human pose recognition in parts from single depth images[END_REF]. Various approaches have been proposed to handle dynamic gestures recognition by exploiting the features previously calculated. Most of the problems have been solved with statistical techniques: Hidden Markov Models (HMM) as in [START_REF] Lovell | Model Structure Selection & Training Algorithms for an HMM Gesture Recognition System[END_REF], [START_REF] Kellokumpu | Human Activity Recognition Using Sequences of Postures[END_REF] or [START_REF] Oka | Real-time tracking of multiple fingertips and gesture recognition for augmented desk interface systems[END_REF] and/or Principal Component Analysis (PCA) as in [START_REF] Kim | Three Dimensional Gesture Recognition Using PCA of Stereo Images and Modified matching Algorithm[END_REF]. SVM (Support Vector Machine) are also becoming a popular way for visual spatio-temporal pattern recognition as in [START_REF] Schuldt | Recognizing Human Actions: A Local SVM Approach[END_REF]. Templatematching is also sometimes used for matching with temporal templates, as in [START_REF] Bobick | The Recognition of Human Movement Using Temporal Templates[END_REF]. This article presents our research work on development of a technical gesture recognition system for human-robot collaboration in an experimental cell of PSA Peugeot Citroen. Section 2 presents our approach, sensors and algorithms used, and our experimental protocol. Section 3 then provides and analyzes results of our experiments. Finally, section 4 contains our conclusions and perspectives for this work.

Methodology

Technical gestures are strongly related to the movement of parts of the body, what is called "motion". But, what we mean by the word "motion" and how can we analyze it? Motion analysis intends to describe the body movements, but in most applications not all the details of the human body is required. The representation and structure of human motion should be simplified according to the application. In the case study of the Human-Robot Collaboration (HRC), a good simplification of the complexity of the human body is to estimate only position of joints, or angle of segments, since they can provide a sufficient representation of the human posture. Other information, such as clothing, tendons, muscles etc, are absolutely unnecessary for the technical gesture recognition for HRC in general, and in particular for our specific case study of PSA. In order to avoid self-occlusions and scene occlusions of the gestures, we have chosen to capture rotational information of body segments using inertial sensors.

Inertial motion sensors for technical gesture recognition

Two different types of IMS have been tested with virtually performed gestures: a) 3 individual and non-hierarchical inertial sensors from Movea (MotionPod); b) 12 hierarchical inertial sensors for the upper part of the body from AnimaZoo (IGS120+). Both of these technologies are unaffected by occlusions (contrary to vision-based techniques), and they provide a rotational representation of the gestures.

For using Movea inertial sensors, two of them have been mounted on the wrists of the worker and the third one on his torso. The latter has been used as a reference point for the other two sensors. The acquired data are angular accelerations for the X, Y, Z axes for the wrists, and they have been normalised a posteriori to Euler angles in order to train the HMMs. The IGS120+ is a product designed for an industrial use, with robustly packaged inertial sensors and specific sensor data correction to reduce potential disruptions due to possible magnetic fields that may exist in the industrial environment. On the other hand, the price of a single sensor of the IGS120+ is twice the price of a single Movea sensor. Additionally, the IGS120+ is more appropriate for an industrial use for security reasons. More precisely, the IGS120+ gives a complete representation of the upper part of the worker's body, which is an important information required in order to be able to detect body postures in case of accidents, or other abnormal situations. The IGS120+ provides exactly the same type of data as those from Movea, but with 12 sensors placed on different body segments (see Figure 2).

Gestures recognition technique

We perform recognition of technical gestures using the processing pipeline illustrated on Figure 3. Gestures execution speed, as well as details of movements can vary significantly from one realization to another. We therefore use, for learning and recognition, a hybrid approach combining Hidden Markov Model (HMM) with a time-rescaling similar to Dynamic Time Warping (DTW); this technique was designed and developed by [START_REF] Bevilacqua | Wireless sensor interface and gesturefollower for music pedagogy[END_REF][START_REF] Bevilacqua | Gesture in Embodied Communication and Human-Computer Interaction[END_REF], and we use the implementation they provide in the "Gesture Follower" tool 1 . This is a template-based method which allows us to use a single gesture to define a 1 http://imtr.ircam.fr/imtr/Gesture_Follower gesture class, by training a HMM. One of the advantages of this method is the DTW-like online time-alignment, illustrated on Figure 4, which permits good robustness of recognition to variation in execution of same technical gesture. A HMM is a statistical modeling method which can be used to establish a model for a time-series with spatial and temporal variability. As described by (Rabiner,1989), a HMM is composed of N hidden states {S 1 , S 2 ,…, S N }, each associated with M possible observations. We denote the actual state at time t as q t and the observation at time t as O t . A HMM is defined by its structure and its parameters. The whole is denoted by λ. The structure of a HMM is the number of states and the non-zero transition probabilities. The parameters are:

-The starting state probabilities: , with . -The state transition probabilities: with -The observation probabilities: with .

To simplify the HMM learning procedure, each learning example is defined by a left-to-right Markov chain. The left-to-right model does not have backward paths, the state transition probabilities allow to stay in the same state or to go in the upper state. This model is appropriate to represent a temporal system. The learning example is down-sampled, with a constant sampling rate, and, as in DTW, each sample is associated to a state. The observation probability for each state is a Gaussian model centered on the sample value. The recognition procedure can be done in real time using the well-known forward procedure which computes the likelihoods for the different models [START_REF] Bevilacqua | Gesture in Embodied Communication and Human-Computer Interaction[END_REF].

Figure 4 : The red curve represents the learning example and the blue curve an observation. The red arrows illustrate the time alignment between the reference gesture and the performed gesture. Each sample of the reference gesture is used to generate a state of a left-toright HMM.

Experimental protocol

The gestures that have been analysed and recognized with our system are the following:  G1: "to remove protective paper"  G2: "to fit on door the waterproofing-sheet"  G3: "pre-stick the waterproofing-sheet"  G4: "to fit the window sealing strip".

The durations of all 4 gestures are about 8-10 seconds. It has been asked to the worker to perform 5 times these 4 gestures. Recognition performances of the system have been evaluated offline on pre-recorded isolated gestures. In our use case, a dataset contains observations of all the 4 gestures (co-articulated, i.e. one gesture after the other, without interruption). In total, 5 observations for each gesture have been recorded and different learning and recognition databases have been used in 5 iterations. We have first evaluated our gesture recognition approach using the « jackknife » method. Here, jackknifing means estimation of the precision of the recognition accuracies for isolated gestures by using subsets of the available gestural data. The basic idea behind the jackknife variance estimator lies in systematically recomputing the statistical estimate leaving out one or more observations at a time from the sample set. For each iteration, one dataset (4 observations) is left out to be used as the learning database and train the HMM of each gesture, until all the datasets are used once. The rest of the datasets are used as a test database.

Results

Gesture recognition performances

The usual precision and recall metrics have been used to evaluate gesture recognition performances of our approach.  Precision is defined by: where is the number of test sequences for which the maximum likelihood from the HMM corresponds to the correct actually performed gesture, while is the number of test sequences for which HMM output is wrong. Precision thus measures the proportion of correct results among all examples for which output by the system corresponds to the considered gesture.

 Recall is defined by: where is the number of test sequences for which the HMM does not output maximum likelihood for the actual gesture. Recall therefore estimates the proportion of correct results among all examples of considered class in the recognition database.

Table 1 shows the results for the 5 iterations of the jackknifing for motionPods data, as well as the Precision and Recall per gesture. The average recognition performances obtained are:

 precision  86%  recall  86%

Gesture G3 ("pre-sticking the waterproofingsheet") is perfectly recognized, with 100% recall. Conversely, recognition rate of gesture G4 ("to fit the window sealing strip") is the worst one, with only 75% recall. One reason that may explain the very good results for gesture G3 can be the fact that the worker makes circular movements with his/her right hand without walking in his/her workspace, in contrast with the other 3 gestures. In other words, from a stochastic point of view, G3 is easier to

Output (maximum likelihood) G1 G2 G3 G4 Recall Observation (Gesture) G1 18 1 1 - 90% G2 1 16 1 2 80% G3 - - 20 - 100% G4 1 3 1 15 75%
Precision 90% 80% 87% 88% 86% discriminate, and present less variability. Regarding precision, gesture G2 ("to fit on door the waterproofing-sheet") is the one for which performance is lowest (80%); the significant rate of confusion between gestures G4 and G2 can be explained by actual similarity of postures in those two gestures.

Table 2 shows the recognition results using AnimaZoo data. The average recognition performances are therefore the following:  precision  96%  recall  96%.

These rates are significantly higher than the 86% obtained with 3 MotionPods, which is not surprising, as the 12 inertial sensors of AnimaZoo half-suit provide much more information on upperbody movements. Also note that per-gesture results are coherent with those obtained with 3 motionPods: gesture G4 ("to fit the window sealing strip") has the lowest recall (85%), and gesture G2 has lowest precision (87%), with significant confusion of gesture G4 with G2.

Test of real-time early recognition

While the above performance evaluation has been done offline on isolated pre-recorded gestures, it is essential for our application to be able to recognize gestures online in real-time, as early as possible.

For testing this, we analyzed the probabilities returned on every timestep for each learnt gesture, on a continuous sequence of actions. We used the AnimaZoo™ for the learning and the testing.

As we can see on top of Figure 5, there is almost always a gesture with an estimated probability around 1 issued by one of the trained HMM. However, if we look at the gesture with the highest probability in real-time (see blue curve on bottom of Figure 5), there are some fluctuations at each gesture transition due to overlapping movements that occur during the gesture sequence (gesture coarticulation). For short time periods, the gesture model with the highest probability may not correspond to the actually performed gesture. Indeed, the system needs a minimum amount of data to establish which gesture is currently being executed. These fluctuations are longer for the fourth gesture G4 because the beginning of the gesture is similar with the G2 gesture, and only the end of G4 is really characteristic of this gesture.

To prevent false recognition and add robustness to the output of our system, we process data in realtime to establish which gesture is being performed.

For each timestep, we analyze instantaneous likelihoods in a time-window with a length equal to 1/3 of duration of the shortest gesture. In our experiment the duration of this time-window is 1.64 second.

Figure 5: On top, gesture probabilities output in real-time by our recognition system. On bottom, the gesture with the highest probability in real-time (blue), the final result (estimate of performed gesture) obtained after time-window filtering (red), and the ground truth (black). The four performed gestures are successively G1, G2, G3 and G4. Their time intervals are represented with color background: blue for G1, red for G2, green for G3 and yellow for G4.

Output (maximum likelihood) G1 G2 G3 G4 Recall Observation (Gesture) G1 20 - - - 100% G2 - 20 - - 100% G3 - - 20 - 100% G4 - 3 - 17 85%
Precision 100% 87% 100% 100% 96%

We calculate which gesture has the largest number of highest probabilities in this window, taking into account only the probabilities above 0.7. If this gesture has the highest probabilities for at least 75% of the timesteps in the window, the system returns this gesture, otherwise the system returns 0. The '0' output means that no gesture is recognized with enough certainty. As we can see on bottom of Figure 5 (red curve), our time-window filtering prevents most of the false recognitions. At the beginning of gesture G2 and gesture G3 and in the middle of gesture G3, the system returns 0. There is still an error at the beginning of the fourth gesture.

The confusion between gesture G2 and gesture G4 is too long to be ignored. A solution could be to extend the window but we will take the risk to not recognize gesture G4. Note that this is totally coherent with the lower recall and precision observed for gesture G4 on tables 1 and 2. Finally, it can be noted that the average time between the actual gesture beginning and its recognition by our system is 3.4 seconds, which is an acceptable earlyrecognition for our application.

Conclusions & perspectives

In this paper, we have presented a framework and preliminary experimental results for real-time recognition of human operator actions on factory assembly-line. The goal is to enable a collaborative industrial robot, to understand behavior of humans around it, so that it can adapt its own actions and execution speed accordingly.

Our experiments shows that it is possible to obtain very good gesture recognition rates (96% precision and recall on several repetitions of the same assembly operations by a single operator), using real-time motion capture with a "MoCap suit" of 12 inertial sensors (AnimaZoo) estimating joint angles of upper-half of human body (neck, wrists, elbows, shoulders, etc…). The great advantage of the chosen motion capture technology is total avoidance of any occlusion problems, contrary to visual-based motion capture. Also, we find that somewhat lower but nevertheless rather good recognition rates (86% precision and recall) can be obtained with only 3 inertial sensors placed respectively on the torso and the two wrists of the operator. Finally, we have tested real-time recognition capabilities of our framework using MoCap suit, during continuous sequences of actions performed by a human operator. These tests show that our system is able to perform online early recognition, after an average delay of less than 4 seconds after gesture beginning, with very little erroneous output.

Ongoing work aims at evaluating our framework for same actions recognition but on more executions by a larger pool of different human operators, and also to estimate false recognition rates on unrelated gestures. Another interesting potential perspective is the use of workers' motion capture in order to estimate effort and stress, which could be helpful for prevention of physical causes contributing to initiate or worsen some musculoskeletal disorders.

Figure 1 :

 1 Figure 1 : Framework for human motion capture and gesture recognition.

Figure 2 :

 2 Figure 2 : On the left, a worker fitting the windows sealing-sheet on a door in the experimental cell of PSA, while he wears the AnimaZoo IGS120+; on the right 3 inertial sensors from IGS120+.

Figure 3 :

 3 Figure 3 : Gesture recognition pipeline based on the AnimaZoo suit of wireless motion sensors for the upper-part of the body, and one-shot learning of HMM + time-alignment for robust recognition.

Table 1 :

 1 Precision and Recall per gesture obtained with motionPods capture, estimated by jackknifing.

Table 2 :

 2 Precision and Recall per gesture obtained with AnimaZoo inertial "MoCap suit" IGS-120+, estimated with jackknifing.

Acknowledgements

This study was carried out thanks to the financial support of "Chaire PSA Peugeot-Citroën on Robotics and Virtual Reality".

Authors wish to thank Pascal Ligot, from PSA, for his important contribution to the experimentations conducted in PSA experimental cell.


Denis Hatebur
email: denis.hatebur@uni-duisburg-essen.de

Maritta Heisel
email: maritta.heisel@uni-duisburg-essen.de

Jeanine Souquières
email: jeanine.souquieres@loria.fr

A Method for Component-Based Software and System Development

Keywords: Component based approach, interoperability, formal method

We propose a method for component-based software and system development, where the interoperability between the different components is given special consideration. The method uses existing notations and languages with their associated tools: context diagrams for analyzing and structuring the problem, composite structure diagrams for describing the overall system in terms of components and interfaces, sequence diagrams to describe the behavior of each component, and the formal method B for specifying the interfaces of the different components and for proving their interoperability. The method proposes to integrate these different notations; at the end of the process, the interoperability is guaranted by the use of the B method with its underlying concept of refinement, and its powerful tool support, the B prover.

I. INTRODUCTION

The idea underlying the paradigm of component orientation [START_REF] Heineman | Component-Based Software Engineering[END_REF], [START_REF] Szyperski | Component Software[END_REF] is to develop software systems not from scratch but by assembling pre-fabricated parts, as is common in other engineering disciplines. As in object orientation, components are encapsulated, and their services are accessible only via interfaces and their operations. To really exploit the idea of component orientation, it must be possible to acquire components developed by third parties and assemble them in such a way that the desired behavior of the system to be implemented is achieved. These considerations lead to the following requirements on how to describe components in such a way that they can be assembled into systems in a systematic way:

• The specification of a component must contain sufficient information to decide whether or not to acquire it for integration in a new system. This requirement concerns the access to the component's source code that may not be granted in order to protect the component producer's interests. Moreover, component consumers should not be obliged to read the source code of a component to decide if it is useful for their purposes or not. Hence, the source code should not be considered to belong to the component specification. • It does not suffice to describe the provided interfaces of a given component. Often, components need other components to provide their full functionality. Hence, also the required interfaces must be part of a component specification. • For different components to interoperate, they must agree on the format of the data to be exchanged between them. Hence, each interface of a component must be equipped with an interface data model describing the format of the data accepted and produced by the component. • It does not suffice to give only the signature of interface operations (e.g., operation foo takes two integers and yields an integer as its result) as is common in current interface description languages. It is also necessary to describe what effect an interface operation has (e.g., operation foo takes two integers and yields their sum as a result). • The interface operations provided or required by some component usually cannot be called in an arbitrary order. Instead, certain communication protocols must be adhered to. These protocols also must be provided in a comprehensive component description. The method for component-based software and system development we present in this paper was designed to fulfill these requirements. It is based on existing languages and notations with their associated tools to appropriately describe existing components or components to be constructed, and to check if two components can be connected via given interfaces or not. In particular, we use the following languages:

• Jackson's context diagrams [START_REF] Jackson | Problem Frames. Analyzing and structuring software development problems[END_REF] serve to analyze and structure the problem in terms of domains and shared phenomena (which will be organized into interfaces later). • UML 2.0 [START_REF]OMG Unified Modeling Language: Superstructure[END_REF] composite structure diagrams serve to express the overall architecture of the system. • UML 2.0 sequence diagrams serve to express the visible behavior of components. They show communication sequences with other connected components. • The formal method B [START_REF] Souquières | The B Book[END_REF] serves to specify interfaces.

An interface specification consists of a data model and the specification of the different operations provided or required by the interface. Since the B concepts of a machine and of refinement fit well with components and their interoperability, we use a B tool for proving component interoperability. (Using for example the object constraint language OCL [START_REF] Warmer | The Object Constraint Language: Precise Modeling with UML[END_REF], and generating verification conditions from scratch would be much more tedious.) In the following, we present an overview of the B method in Section II. We then describe our method in Section III and illustrate it with the case study of an access control system in Section IV. In Section V, we discuss related work. The paper finishes with some concluding remarks in Section VI.

II. THE FORMAL METHOD B

The B method [START_REF] Souquières | The B Book[END_REF] is a formal software development method based on set theory. Because of its rigor and powerful tool support, it is often used to develop software for critical systems. The B method supports an incremental development process, using refinement. A development begins with the definition of an abstract specification, which can be refined step by step until an implementation is reached. The refinement of models is a key feature for incrementally developing models from textual descriptions, preserving correctness in each step. Thus, B follows the proof-based development paradigm [START_REF]Obligations de preuve: Manuel de référence[END_REF]. The method has been successfully applied in the development of several complex real-life applications, such as the METEOR project [START_REF] Behm | METEOR: A Successful Application of B in a Large Project[END_REF]. It is one of the few formal methods which has robust and commercially available support tools for the entire development life-cycle from specification down to code generation [START_REF] Bert | Adaptable Translator of B Specifications to Embedded C Programs[END_REF]. B specifications consist of abstract machines, which are very close to notions well-known in programming under the names of modules, classes or abstract data types. Each abstract machine consists of a set of variables, invariant properties of those variables, and operations. The state of the machine, i.e. the set of variable values, is modifiable by operations, which must preserve its invariant. An example of a B machine is given Fig. 1. The invariant clause characterizes the sensible states that are permitted for the machine. The machine should never arrive at a state in which some part of the invariant clause is false. The B method provides structuring primitives that allow one to compose machines in various ways. Large systems can be specified in a modular way and in an object-based manner [START_REF] Meyer | A systematic approach to transform OMT diagrams to a B specification[END_REF], [START_REF] Ledang | Modeling class operations in B: application to UML behavioral diagrams[END_REF]. Proofs of invariance and refinement are part of each development. The proof obligations are generated automatically by support tools such as AtelierB [START_REF]Obligations de preuve: Manuel de référence[END_REF] or B4free [START_REF] Clearsy | B4free[END_REF], an academic version of AtelierB. Checking proof obligations with B support tools (either through automatic or interactive proofs) [START_REF] Abrial | Click'n'Prove : Interactive Proofs Within Set Theory[END_REF], is an efficient and practical way to detect errors introduced during development.

III. A METHOD FOR COMPONENT-BASED

DEVELOPMENT Our goal is to provide a method for component-based software and system development that pays special attention to the question how the interoperability between different components can be guaranteed. Components are specified as black boxes, so that component consumers can deploy them without knowing their internal details. Hence, component interface specifications play an important role, as interfaces are the only access points to a component. Our method consists of the following steps.

A. Set up a context diagram

Context diagrams as introduced by Jackson [START_REF] Jackson | Problem Frames. Analyzing and structuring software development problems[END_REF] serve to analyze and structure a given development problem. They consist of rectangles that are connected by lines (for an example, see Fig. 2). The rectangles denote domains, and a connecting line between two (or more) domains indicates that there are shared phenomena between the connected domains. Shared phenomena may be events, operation calls, messages, and the like. They are observable by at least two domains, but controlled by only one domain. For example, if a user pushes a turnstile, this is a phenomenon shared by the user and the turnstile, which is controlled by the user. For each connecting line, the corresponding shared phenomena are given in a context diagram. In the original version of context diagrams, each such diagram contains exactly one machine domain, denoted by a rectangle with a double vertical stripe. This is the domain that has to be constructed to solve the problem. All other domains are given domains that already exist1 . Their properties (called domain knowledge) have to be taken into account in the problem solving process. In our method, domains correspond to components or actors, and connecting lines and their shared phenomena correspond to interfaces between components or components and actors.

In contrast to Jackson, we allow several or even zero machine domains in a context diagram. If a component has to be constructed, it is expressed as a machine domain. If a component exists, it is expressed as a given domain. Thus, context diagrams with no machine domain express the fact that a system shall be assembled from existing components only, and that we have to check the interoperability of these components. If a context diagram contains at least one machine domain, we have more flexibility in the development process. For example, we might change the interfaces of a component to be constructed if this is necessary in order to interoperate with the given components. Hence, the domains contained in our context diagrams may play one of the following roles:

• A machine domain corresponds to a software component to be constructed, (e.g., Controller in Fig. 2). • A given domain can be -an existing software component to be used (e.g., Data Base in Fig. 2) -an existing hardware component to be used (e.g., Turnstile Entry in Fig. 2) -an actor communicating with the system (e.g., User in Fig. 2) Note that for existing components we need not distinguish between hardware and software components, as both are described in the same way. Thus, our method not only covers component-based software development, but also component-based system development.

A second generalization of Jackson's context diagrams concerns the control of shared phenomena. In Jackson's method, this information is only added in a later step. In contrast, we add it to the context diagram. The notation U!{push} means that the User component controls the shared phenomenon push. Setting up a (generalized) context diagram, we structure the given development problem by identifying the different actors and the components that are available or must be constructed. Specifying the shared phenomena between the components is a preparation for the next step of the method, where all provided and required interfaces must be identified.

B. Construct the system-or software architecture

Based on the context diagram developed in the first step of the method, we now decide on the provided and required interfaces of all the identified components. For this purpose, we have to inspect the shared phenomena expressed in the context diagram, and consider the following rules:

• A shared phenomenon observable but not controlled by component C corresponds to an operation or signal in a provided interface of C, because C must in some way be notified of the occurrence of the shared phenomenon. • Conversely, a phenomenon controlled by component C will often correspond to an operation or signal in a required interface of C, because controlling a phenomenon often corresponds to sending a signal or message to another component (thus invoking an operation of a provided interface of the receiver component). However, there are exceptions to this rule. They concern passive components that are mere data structures. Such components do not send messages to their environment without an external stimulus. Hence, the phenomena controlled by passive domains correspond to return values of operations of a provided interface. An example is the data base in Fig. 2. Although it controls the phenomenon authorized, it only has a provided interface, where the authorized phenomenon is the return value of a data base query issued by the controller component. The software or system architecture is expressed as a UML 2.0 composite structure diagram. Such diagrams contain named rectangles, called parts. These parts are the components of the system. Parts may have ports, denoted by small rectangles, and ports may have interfaces associated to them. Interfaces may be required or provided. Provided interfaces are denoted using the "lollipop" notation, and required interfaces using the "socket" notation. For an example, see Fig. 3. The procedure for setting up the architectural diagram for the system is as follows. Each connecting line between two domains of the context diagram must be transformed into one or two connections between two components. Such a connection consists of a provided interface for one of the components and a corresponding required interface of the other component. For a bi-directional communication, we need two connections (e.g., between Turnstile Entry and Controller in Fig. 3); for a one-directional communication, we only need one connection (e.g., between Data Base and Controller in Fig. 3). The decision on provided and required interfaces must be taken according to the rules given above. Finally, we delete all domains corresponding to actors from the diagram.

Lights

C. Specify components

For each component of the architecture, a specification must be set up containing:

• Sequence diagrams describing the visible behavior of the specified component; the sequence diagrams must contain all operations or signals of all interfaces and all other components the specified component is connected with. Messages received by the specified component from its environment must correspond to an operation or signal of some provided interface, while messages sent out by the specified component must correspond to an operation or signal of some required interface. • A B machine for each provided and each required interface. For all interfaces that connect the specified component with the same outside component, the interface data models (IDM, see Section I) must be the same, and the IDM must be encoded in the B machine by specifying:

-the types used in the interface -a data state as far as necessary to express the effects of operations -invariants on that data state When a component manipulates data, it is possible to use a UML class diagram to express the interface data model for reasons of readability. This class diagram can then automatically be transformed into a B specification [START_REF] Meyer | A systematic approach to transform OMT diagrams to a B specification[END_REF]. Each machine specifies the operations belonging to its corresponding interface. An operation specification consists of its signature (i.e., the types of its input and output parameters), its precondition expressing under which circumstances the operation may be invoked, and its postcondition expressing the effect of the operation. Both pre-and postcondition will refer to the interface data model. The behavioral specification and the operation specifications must be coherent. To ensure coherence, the sequence diagrams should be annotated with states, which are also used to express the pre-and postconditions of the B operations. In this way, we make sure that the communication sequences required by the sequence diagrams do not violate any preconditions of interface operations. Such a specification contains all the information that is needed to decide if a given component can be used in a given context or not.

D. Prove interoperability

In component-based development, the components must be connected in an appropriate way. To guarantee interoperability of components, we must consider each connection of a provided and a required interface contained in the architecture and try to show that -after some syntactic transformations -the provided interface is a B refinement of the required interface. This means that the provided interface constitutes an implementation of the required interface, and we can conclude that the two components can be connected as shown in the architectural diagram. The process of proving interoperability between components is described in [START_REF] Chouali | Proving Component Interoperability with B Refinement[END_REF]. If we cannot demonstrate the interoperability of an intended connection of the architectural diagram, we either can try to change the specification of a component to be developed (corresponding to a machine domain in the context diagram), or we can try to develop an adapter. In both cases, the architectural diagram and the component specifications will have to be adjusted accordingly. If the latter approach is not possible either, we have to conclude that the components cannot be connected as intended.

In the case where all components are given, the result of our method is a proof that the system can be assembled as specified in the architectural diagram. In the case where some components do not yet exist but must be constructed, our method yields detailed specifications of all the components to be developed.

IV. CASE STUDY: ACCESS CONTROL SYSTEM

We illustrate our method with the case study of a simple access control system. The access to a building is to be controlled. Persons (called users) who are authorized to enter the building have access cards with some identification stored on it. There are two turnstiles, one at the entrance to the building, and one at the exit. At the entrance, there is also a card reader as well as a red and a green light. A user who wants to enter the building inserts his or her card into the card reader. The information on the card is read, and a data base is queried to decide if the user is granted access or not. If access is granted, the green light is turned on for some time, the card is ejected, and the entry turnstile -which is normally blocked -is unblocked. The entry turnstile is re-blocked either after the user has entered or after some timeout. If the user does not take the card within some time limit, the card is retracted and kept. If access is denied, a red light is turned on for some time, and the entry turnstile remains blocked. The card is ejected. Again, if the user does not take the card within some time limit, the card is retracted and kept. The number of persons present in the building must be counted. Therefore, there is also an exit turnstile which is never blocked, but just serves to observe when a person leaves the building.

A. Context diagram

A context diagram for this access control system is shown in Fig. 2. We have one component to be developed, namely Controller, and one domain corresponding to an actor, namely User. The other domains correspond to existing components.

As an example, we consider the interfaces c and g: the controller can give commands to the entry turnstile to block or unblock it. Conversely, the entry turnstile sends a signal entered to the controller. This signal is sent as a reaction of a push event caused by the user (interface g). Since the user can observe if the turnstile is locked or not, the interface contains two phenomena locked and unlocked controlled by the entry turnstile.

B. Architecture of the system

The architecture of the system is given in Fig. 3. We have used the following naming conventions for interfaces. Each interface name has the form X Y Z, where

• X is an abbreviation of the name of the component the interface belongs to, e.g., its first letter; • Y is either "P" for provided or "R" for required;

• Z is the abbreviation of the name of the component the given component is connected to. As an example, C P TE is the provided interface of the Controller component which is connected to the Turnstile Entry component. Note that the User is no longer contained in the diagram, and that the data base only has a provided interface, because it is a passive component (even though it controls the phenomenon authorized). In all other interfaces, controlled phenomena correspond to operations in required interfaces, and observed phenomena correspond to operations in provided interfaces.

C. Component Specifications

We want to use two copies of the same turnstile component available on the market to be installed as an entry and an exit turnstile, respectively. For reasons of space, we only present the specification of that turnstile component, and parts of the controller specification.

1) Specification of the turnstile component: The available turnstile component has the following interfaces:

• T P C, its provided interface to a controller component, with two operations lock c and unlock c; • T R C, its requested interface to the controller component, with an operation pushed; • T P U, its provided interface to an external user, with an operation push. There is no required interface to the user because the user finds out if the turnstile is locked by trying to push it and not be being sent a message from the turnstile. The turnstile component has two states, namely locked and unlocked, and its initial state is locked. It has no invariant property. The externally visible behavior, i.e., its usage protocol, is given in Fig. 4. The first sequence diagram expresses that initially the turnstile is locked. After it receives an unlock c command from the controller, it changes its state to unlocked. In that state, an arbitrary number of push events can be received from the user. When the turnstile component receives a lock c command from the controller, it re-enters the state locked. The second sequence diagram expresses the fact that in state locked, no push events can be received from the user. The neg construct has been used to express this forbidden scenario. Note that the given turnstile component is specified using different names for the states and the operations than the ones used in Fig. 2, because it cannot be assumed that component consumers and component producers choose the same names independently of each other.

We now give B specifications for two of the tree interfaces. The interface data models of interfaces that connect the same two components must be the same. Hence, T P C and T R C have the same IDM, including the state of the turnstile with the two possible values, locked and unlocked, and the initial state locked. The interface specifications are shown in Figs. 1 and5. The operations lock c and unlock c change the state of the turnstile, whereas the operation pushed corresponds to sending a signal to the controller and does not change the state of the turnstile. As specified in the sequence diagram, it can only be invoked if the turnstile is in state unlocked. 2) Specification of the controller component: In the same way, we have to specify the controller component. Since the controller is connected with five other components, its behavioral specification is much more complex than the behavioral specification of the turnstile. Fig. 6 shows the behavior of the controller when a user enters the building. It captures exactly the behavior sketched in the informal description of the access control system. Figure 7 shows the behavior of the controller component when a user leaves the building. As an example of an interface specification, we give the B specification of the required interface of the controller with the entry turnstile, C R TE in Fig. 8. It looks quite similar to the specification of the interface T P C, and in fact, we will show that the two interfaces can be connected. The specification of interface of the controller with the exit turnstile is given in Fig. 9. Note that here we require that the turnstile is always unblocked, as specified in the behavioral specification, see Fig. 7.

MACHINE

D. Proving Interoperability

We now must prove that the controller can be connected with two of the given turnstile components as specified Fig. 3. For this purpose, we must carry out three refinement proofs. As example, we show that a transformed version of the interface T P C is a refinement of the interface C R TE, see Fig. 10. The refinement machine NewT P C is a transformation of the machine T P C (Fig. 1), where unlock c is renamed to unlock and lock c is renamed to We might now try to perform a similar refinement proof for the interfaces T R C and C P TX, i.e., to show that the controller and the exit turnstile can be connected as intended, see Fig. 11. However, this proof fails, because the initializations of the two machines are incompatible. The problem is that the given turnstile component comes in state locked, but to be usable as an exit turnstile in the access control system, it must be in state unlocked. Hence, we must change the specification of the controller. The controller must get a new interface C R TX which it can used to initialize the turnstile component. As an initialization operation, the operation unlock c provided by the turnstile component can be used. Figures 12 and13 show the new interface of the controller and its interoperability with the interface T P C of the turnstile component. These changes have to be propagated to the other specifications as follows:

REFINEMENT

• The initialization of the machine C P TX must be changed, so that it also allows an initial state locked. • A connection between the controller and the exit turnstile must be added to the architectural diagram. It connects the interfaces C R TX and TPC. • The initialization of the exit turnstile must be achieved by the controller. The corresponding sequence diagram is shown in Fig. 14. Thus, we have shown how our method supports developers in assembling systems from components, always guaranteeing that the different components are able to interoperate in the intended way.

V. RELATED WORK

In an earlier paper, we have investigated the role of component models in component specification [START_REF] Heisel | Toward a formal model of software components[END_REF]. The specification of a component model makes it possible to obtain more concise specifications of individual components. In this paper, we investigate the necessary ingredients a component specification must have in order to be useful for assembly of a software system out of components. These ingredients are independent of concrete component models. Several proposals for component specification have already been made. They have in common that they have no counterpart of our interface data model and that they do not consider interoperability issues, but only the specification of single components. A working group of the German "Gesellschaft für Informatik" (GI) has defined a specification structure for business components [START_REF] Turowski | Standardized Specification of Business Components[END_REF]. That structure comprises seven levels, namely marketing, task, terminology, quality, coordination, behavioral, and interface. Our specification structure covers the layers terminology, coordination, behavioral, and interface by proposing concrete ways of specifying each of those levels. The other layers of the GI proposal have to do with non-functional aspects of components. Beugnard et al. [START_REF] Beugnard | Making components contract aware[END_REF] propose to define contracts for components. They distinguish four levels of contracts: syntactic, behavioral, synchronization, and quality of service. However, they do not introduce data models for their interfaces. Hence, it cannot easily be checked if two components can be combined. The component specification approach of Lau and Ornaghi [START_REF] Lau | A formal approach to software component specification[END_REF] is closer to ours, because there, each component has a context that corresponds to our interface data model. A context is an algebraic specification, consisting of a signature, axioms, and constraints. In contrast, we deem it more appropriate to allow for an object-oriented specification of the data model of a component interface. This makes it possible to take side effects of operations into account and to use inheritance, concepts that are frequently used in practice. Cheesman and Daniels [START_REF] Cheesman | UML Components -A Simple Process for Specifying Component-Based Software[END_REF] propose a process to specify component-based software. This process starts with an informal requirements description and produces an architecture showing the components to be developed or reused, their in-terfaces and their dependencies. For each interface operation, a specification is developed, consisting of a precondition, a postcondition and possibly an invariant. Our specification of component interfaces is inspired by Cheesman and Daniels' work because that work clearly shows that for each interface, a data model is necessary. However, Cheesman and Daniels do not consider the case that already existing components with possibly different data models have to be combined, and hence they do not define a notion of interoperability. Canal et al. [START_REF] Canal | Extending CORBA interfaces with protocols[END_REF] use a subset of the polyadic π-calculus to deal with component interoperability only at the protocol level. The π-calculus is well suited for describing component interactions. The limitation of this approach is the lowlevel description of the used language and its minimalistic semantics. Bastide et al. [START_REF] Bastide | Formal specification and prototyping of CORBA systems[END_REF] use Petri nets to specify the behavior of CORBA objects, including operation semantics and protocols. The difference with our approach is that we take into account the invariants of the interface specifications. Zaremski and Wing [START_REF] Zaremski | Specification matching of software components[END_REF] propose an interesting approach to compare two software components. It is determined whether one component can be substituted for another. They use formal specifications to model the behavior of components and the Larch prover to prove the specification matching of components. Others [START_REF] Han | A comprehensive interface definition framework for software components[END_REF], [START_REF] Vallacillo | Object interoperability[END_REF] have also proposed to enrich component interface specifications by providing information at signature, semantic and protocol levels. Despite these enhancements, we believe that in addition, a data model is necessary to perform a formal verification of interface compatibility. The idea to define component interfaces using B has been introduced in an earlier paper [START_REF] Chouali | Verifying the compatibility of component interfaces using the B formal method[END_REF]. The use of the B refinement to prove that two components are compatible at the signature and semantics levels has been explored in [START_REF] Chouali | Proving Component Interoperability with B Refinement[END_REF].

VI. CONCLUSION

We have presented a method for component-based software and system development. In this method, components are considered as black boxes. They are only described by their visible behavior and by their interfaces. This approach makes it possible to describe hardware and software components as well as existing components and components to be constructed in the same way. Interoperability between component is defined rigorously and can be checked with tool support. Furthermore, methodological guidance is given to developers, as our method consists of four well-defined steps. Thus, the ideas underlying the concept of componentbased development can be fully exploited, and the way of constructing systems becomes more similar to other engineering disciplines. To construct a working system out of components, adapters have to be defined that implement the transformation of required interface data into provided interface data and vice versa. In some cases, we can relax compatibility conditions and use more liberal versions of specification matching, e.g. plug-in-post matching [START_REF] Zaremski | Specification matching of software components[END_REF]. In this case, an adapter must check if the precondition of the provided operation holds. If not, it has to take appropriate actions other than calling the provided operation. Thus, the construction of adapters becomes a program synthesis problem. This problem becomes more complex for weaker versions of specification matching. We are currently working on alternative versions of compatibility and their mappings to refinement in B, and to give patterns for the corresponding adapters.

Fig. 1 .

 1 Fig. 1. An example of a B Specification

Fig. 2 .

 2 Fig. 2. Context diagram for the access control system

Fig. 3 .

 3 Fig. 3. Architecture of the global access control system

Fig. 4 .

 4 Fig. 4. Behavioral specification of the turnstile component

Fig. 5 .Fig. 6 .

 56 Fig.5. Specification of T R C

Fig. 7 .

 7 Fig. 7. Sequence diagram for the controller component (exit)

Fig. 8 .

 8 Fig. 8. Specification of C R TE

Fig. 10 .

 10 Fig. 10. Specification of New T P C

Fig. 12 .

 12 Fig. 12. Specification of C R TX

Fig. 13 .Fig. 14 .

 1314 Fig. 13. Refinement InitT P C of C R TX

  This renaming is necessary because the definition of refinement in B requires that the refining machine defines operations with the same names as the refined machine. The linking invariant, state = locked ⇔ e state = blocked, relates the states of the different machines.

	MACHINE
	C P TX
	SETS
	Exit States Turnstile = {blocked, unblocked}
	VARIABLES
	x state
	INVARIANT
	x state ∈ Exit States Turnstile
	INITIALISATION
	x state := unblocked
	OPERATIONS
	left	=
		PRE x state = unblocked
		THEN SKIP END
	END
		Fig. 9. Specification of C P TX
	lock. REFINEMENT
	New T P C
	REFINES
	C R TE
	SETS
	Turnstile States = {locked, unlocked}
	VARIABLES
	state
	INVARIANT
	state ∈ Turnstile States ∧
	(state = locked ⇔ e state = blocked)
	INITIALISATION
	state := locked
	OPERATIONS
	unlock =/* unlock c */
		PRE state = locked
		THEN state := unlocked END;
	lock = /* lock c */
		PRE state = unlocked
		THEN state := locked END
	END

There are also designed domains, which are data structures to be developed. We do not use designed domains in this work; hence, they will not be mentioned further.


Alice Galdi

P Orgiani

L Maritato

Laurence Méchin

come

Correlation between structural properties and resistivity critical behavior in SrRuO3 thin films


Note that the total variation is not differentiable, so it cannot be exactly minimized by simple gradient descent, although it is convex.

• In this section, we discuss discretization and the most straight-forward exact solver I know of.

Solving the TV-L 2 model x 0

• The TV-L 2 model argmin

u ∇u 1 + 1 2λ u -f 2 2 .
is a building block for several more general algorithms. It is therefore important to fully understand it.

Grayscale images

x 0

In the continuous setting, an image is a function or scalar field u : Ω → R with Ω ⊂ R 2 . For simplicity, we assume in this course that the domain is always two-dimensional.

Discretization of an image

• the domain Ω is a regular grid of size M × N

• an image u is a matrix u ∈ R M×N • the gray value or intensity at a pixel (i, j) is written as u i,j

Bastian Goldluecke Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"

Color images

x 0

Color images have several components, usually red, green and blue for the primary colors. Thus, they are functions u : Ω → R 3 . In the discretization, each component u 1 is represented by its own matrix of intensity values.

Vector fields

x 0

A vector field is a map p : Ω → R 2 . One can imagine it as an assignment of a little arrow to each point:

Discretization of a vector field

• a vector field p consists of two components p 1 , p 2 , each of which is a matrix of size M × N

• the vector at a pixel (i, j) is written as

p 1 i,j p 2 i,j
.

Bastian Goldluecke Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"

The gradient of a scalar field

x 0

In the continuous setting, one of the basic operators acting on a function is the derivative. The vector field of the two partial derivatives of a function u : Ω → R is called the gradient of u,

∇u = ∂ x u ∂ y u .

Discretization of the gradient with forward differences

For i < N and j < M:

(∇u) i,j = u i+1,ju i,j u i,j+1u i,j .

Neumann boundary conditions ∂u ∂n = 0 are used for functions, i.e.

(∇u) 1 M,j = 0 and (∇u) 2 i,N = 0.

Bastian Goldluecke Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"

The divergence of a vector field

x 0

In the continuous setting, the divergence of a vector field p is the scalar field or function div(p) : Ω → R defined as

div(p) = ∂ x p 1 + ∂ y p 2 .

Discretization of the divergence with backward differences

Derivatives on vector fields use Dirichlet boundary conditions p| ∂Ω = 0. Thus,

(div(p)) i,j = p 1 i,j -p 1 i-1,j if 1 < i ≤ M p 1 i,j if i = 1 + p 2 i,j -p 2 i,j-1 if 1 < j ≤ N p 2 i,j if j = 1 Bastian Goldluecke
Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"

Discrete total variation

x 0

Let u ∈ R M×N be a matrix, then its total variation is defined as

TV(u) := M i=1 N j=1 (∇u) i,j 2 ,
where the gradient is computed with forward differences and Neumann boundary conditions.

Note:

Actually computing the discrete TV requires to compute the sum over all elements in a matrix, which is a not fully parallelizable type of operation called a reduction.

Dual formulation

x 0

The Euclidean norm in R n can be written for any vector v as

|v | 2 = sup p∈R n ,|p| 2 ≤1 (v , p),
where (v , p) denotes the usual inner product.

We can use this fact to rewrite the total variation of a matrix as follows:

TV(u) = max p∈K i,j (∇u) i,j , p i,j ,
where the (convex) set K consists of all vector fields such that their vectors in each point have length at most one,

K := p ∈ R M×N × R M×N : p i,j 2 ≤ 1 for all i, j .
Instead of finding a minimizer for ROF, we can look for a saddle point (u, p) of min

u∈R M×N max p∈K i,j (∇u) i,j , p i,j + 1 2λ (u i,j -f i,j ) 2 .
• Note that the energy is now differentiable, so as a simple method to find a saddle point we can iterate gradient descent in u and ascent in p starting from an initial point (u 0 , p 0 ).

• We will formulate a slightly faster method, which is based on the observation that for given p, one can solve the minimization in u explicitly.

Gradient ascent in the dual variable

x 0

The energy of the saddle point problem is

E(u, p) = i,j (∇u) i,j , p i,j + 1 2λ (u i,j -f i,j ) 2 .
Its derivative with respect to p is just

d dp E(u, p) = ∇u,
so a (discrete) gradient ascent step for p is given by setting

p n+1 = p n + τ ∇u n ,
where (p n , u n ) is the nth iterate for the saddle point and τ > 0 is a small enough step size (to be determined later).

Projection onto K

x 0

The allowed range for p is the set K . After the update step, p n might lie outside K , so we have to back-project it onto K . This can be done by projecting each single vector of the field p back onto the unit disk,

Π K (p n ) i,j = p n i,j max 1, p n i,j 2
.

Bastian Goldluecke Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"

Updating the primal variable

x 0

The reason the finite difference operators introduced in the first lecture were defined in this way is that they need to satisfy i,j

(∇u) i,j , p i,j =i,j u i,j (div(p)) i,j

like their continuous counterparts.

This means we can rewrite the saddle point energy as

E(u, p) = i,j -u i,j (div(p)) i,j + 1 2λ (u i,j -f i,j ) 2 .
The gradient with respect to u is

d du E(u, p) = -div(p) + 1 λ (u -f ),
At a minimum for u, this must be zero. This leads to the update equation

u n+1 = f + λdiv(p n+1 )
which computes the exact solution for u n+1 explicitly, given the previously computed p n+1 .

Complete algorithm

x 0

Algorithm to solve ROF (Bermudez-Moreno 1981, Chambolle 2005)

• Given an image f to be denoised and a smoothing parameter λ > 0.

• Set the step size to τ := 1 2λ ∇ = 1 8λ .

• Start with initial u 0 , p 0 , where u 0 is a scalar field and p 0 is a vector field. Initial values can be arbitrary, choose e.g. zero for everything.

• Iterate until convergence:

p n+1 = Π K (p n + τ ∇u n ) u n+1 = f + λdiv(p n+1 )

Introduction to variational methods

Convex optimization primer: proximal gradient methods

Variational inverse problems

Segmentation and labeling problems 3D reconstruction

Summary

Setting and notation

x 0

• A Hilbert space V with

• inner product (•, •) and

• norm • = (•, •),
in this tutorial, this will be V = R N or V = L 2 (Ω).

• The dual space V * -in our setting, you can substitute V * = V if you are not familiar with the concept

• We consider "proper" functionals f : V → R ∪ {∞}, i.e. they can take the value infinity, but we exclude the (boring) functional f = ∞ for technical reasons.

The epigraph of a functional

x 0

The epigraph epi(f ) of a functional f : V → R ∪ {∞} is the set "above the graph", i.e.

epi(f

) := {(x, µ) : x ∈ V and µ ≥ f (x)}. epi(f ) f V dom(f ) R Bastian Goldluecke
Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"

Convex functionals

x 0

Definition • A functional f : V → R ∪ {∞} is called convex if epi(f ) is a convex set.
• The set of all proper and convex functionals on V is denoted conv(V).

We choose the geometric definition of a convex function here because it is more intuitive, the usual algebraic property is a simple consequence.

Bastian Goldluecke Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics" non-convex energy convex energy

Convex energies can be globally minimizedfor non-convex energies this is usually impossible.

Bastian Goldluecke

Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"

The subdifferential

x 0 Definition • A vector ϕ ∈ V * is a subgradient of f at x ∈ V if f (y ) ≥ f (x) + y -x, ϕ for all y ∈ V.
• The set of all subgradients of f at x is the subdifferential ∂f (x).

Geometrically: the function on the right hand side is an affine function

h(y ) = f (x) + y -x, ϕ
with slope ϕ, which touches the epigraph of f in x. The condition says that h needs to lie below f .

Example: the subdifferential of f : x → |x| in 0 is

∂f (0) = [-1, 1].
Bastian Goldluecke Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics" 29 Subdifferential and Fr échet derivative

x 0
The subdifferential is a generalization of the Fr échet derivative (or the gradient in finite dimension), in the following sense.

Theorem

Let f ∈ conv(V) be Fr échet differentiable at x ∈ V. Then

∂f (x) = {df (x)}.
The proof of the theorem is surprisingly involved -it requires to relate the subdifferential to one-sided directional derivatives. We will not explore this here, see e.g. [Rockafeller 1970].

Euler-Lagrange equation for convex functionals

x 0

Theorem

Let f ∈ conv(V). Then x is a global minimum of f if and only if

0 ∈ ∂f (x).
In particular, a convex functional does not have minima which are local but not global.

Geometrically intuitive, but the proof is also surprisingly simple. Rewriting the subgradient definition, one sees that 0 ∈ ∂f (x) just means that

f (y) -y , 0 ≥ f (x) -x, 0 for all y ∈ V.
This already implies the equivalence.

Proximity operator

x 0

• The proximity operator for a closed convex functional H and step size parameter τ > 0 is defined as

prox τ H (x) = argmin y 1 2τ y -x 2 + H(y) . • At the minimizer x = prox τ H (x), 0 ∈ 1 τ (x -x) + ∂H(x) ⇔ x ∈ x -τ ∂H(x).
The proximity operator computes an implicit subgradient descent for the convex functional H with step size τ .

• Let H = TV , then

prox λTV (x) = argmin y TV (y) + 1 2τ x -y 2
is an instance of ROF-denoising.

• Solve e.g. with Bermudez-Moreno.

Bastian Goldluecke Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"

• Let C ⊂ V be a convex set, and

δ C (x) := 0 if x ∈ C, ∞ otherwise be its indicator function. • Then prox δ C (x) = argmin y ∈V 1 2τ y -x 2 + δ C (y ) = argmin y ∈C y -x 2
is the projection onto C.

• Adding an indicator function to a functional is equivalent to an optimization constraint.

• Let H = • 1 (sparsity inducing norm), e.g. interesting in compressive sensing or dictionary learning, next talk.

• Then (component-wise)

prox τ • 1 (x) i = (|x i | -τ )sgn(x i ),
which is called the shrinkage operator.

Bastian Goldluecke Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"

Applicability: Minimization problems of the form

argmin u∈V {G(u) + F (u)}
such that

• G is convex

• F is convex and differentiable

• dF is Lipschitz continuous, i.e. there exist L > 0 such that

dF (u) -dF (v ) ≤ L u -v for all u, v ∈ V.
• The proximation for G is (comparatively) easy to evaluate.

Algorithm: Initialize with u 0 = 0, then compute alternating "forward-backward steps":

• gradient descent in F :

v n+1 = u n - 1 L dF (u n ).
• implicit subgradient descent in G:

u n+1 = prox 1 L G (v n+1 ).
only O(1/n) convergence (slow). Fast ISTA (FISTA)

x 0

Algorithm: Initialize with u 0 = 0, ū0 = 0, R ∋ t 0 = 1, then iterate

• gradient descent in F :

v n+1 = ūn - 1 L dF (ū n ).
• implicit subgradient descent in G:

u n+1 = prox 1 L G (v n+1 ). • extragradient step: t n+1 = 1 2 (1 + 1 + 4t 2 n ), ūn+1 = u n + t n -1 t n+1 (u n+1 -u n ).

Beck and Teboulle SIIMS 2009

Bastian Goldluecke Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"

• Leads to O(1/n 2 ) algorithm (optimal for a first-order method)

• Good general-purpose solver, easy to implement (in particular if proximation is already available)

Example: constrained differentiable problems When applying FISTA, δ C is the non-differentiable part taken care of via the prox operator, which equals the projection onto C.

• Alternating projection (projection onto convex sets)

• Douglas-Rachford splitting, alternating direction method of multipliers ADMM (popular e.g. in learning)

• Split-Bregman (basis pursuit problems, i.e. constrained L 1 -minimization)

• Chambolle-Pock (swiss army knife if both F and G are non-differentiable, but convex with simple prox-operators).

All of these use the same basic building blocks as the ones shown in the previous section.

Introduction to variational methods

Convex optimization primer: proximal gradient methods

Variational inverse problems

Segmentation and labeling problems

3D reconstruction

The convolution b * u with a kernel b of total mass 1 can be interpreted as a blurring operation.

TV deblurring model

Recover a deblurred image u as the minimizer of

argmin u ∇u 1 + 1 2λ b * u -f Example x 0

Original Blurred + Noisy Solution

FISTA for the general inverse problem

x 0

• The differentiable part is

F (u) = 1 2λ Au -f 2 2 ,
which has functional derivative

dF (u) = 1 λ A T (Au -f ),
which is Lipschitz-continuous with L = 1 λ A T A 2 .

• The non-differentiable convex part is the regularizer, e.g. TV, with the proximation given by solving an instance of ROF.

Inpainting problem

x 0

• Restoration problem: recover missing regions of an image • Input:

• damaged region Γ ⊂ Ω • partial image f : Ω \ Γ → R • TV inpainting: Solve min u ∇u 1 such that u = f on Ω \ Γ.

Damaged image f Recovered image u

TV inpainting results ("textured" images)

x 0

Original Damaged Inpainted

• TV inpainting is unconvincing for highly textured images if the missing regions are larger. The reason is that no structure is inferred from surrounding regions, and only boundary values of Γ are taken into account.

• A better variational model for inpainting can be found e.g. in papers on non-local TV by Osher et al., or check methods based on dictionary learning.

Object removal

x 0

Once we have an inpainting algorithm, we can employ it to remove unwanted regions in an image by marking them as damaged.

Uncaging a bird

Optic flow x 0

Idea: introduce auxiliary variable v

• Decouples regularizer from the data term.

• Let θ > 0 be a coupling parameter, define family of energies in the variables u, v ,

E θ (u, v ) = J (u) + 1 2θ u -v 2 2 + Ω ρ(x, v (x)) dx.
• For θ → 0, the quadratic term forces u to be close to v , so the solution of E θ approximates the solution of E.

x 0

Algorithm

Set θ > 0, start with initial u 0 , v 0 . Then iterate the following steps until convergence.

• Solve an instance of TV-L 2 :

u k+1 = argmin u E θ (u, v k ).
• Solve the point-wise problem

v k+1 = argmin v E θ (u k+1 , v ).
Bastian Goldluecke Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics" 53

Real-time depth Example: histograms of user input (naive, illustration only)

x 0

Foreground histograms (RGB)

Background histograms (RGB)

P fg (x) =
number of marked foreground pixels with color f (x) total number of marked foreground pixels P bg (x) = number of marked background pixels with color f (x) total number of marked background pixels Straightforward:

• Since F (u) = (c, u), the derivative is constant,

dF (u) = a.
• It is Lipschitz-continuous with arbitrary Lipschitz constant L > 0, one can e.g. just choose L = 1.

• Thus, the gradient descent step in F just subtracts c from the solution.

• Proximation in the TV-regularizer is an instance of ROF. TV (best result)

Length does not work well with elongated structures (curvature regularity is better in that case, but complex)

Regularization x 0 Σ γ 2 γ 1
The regularization penalty is proportional to the label distance times the length of the interface.

In this example d(γ 1 , γ 2 ) • L(Σ)
Euclidean representation of the label distance:

• Each label γ is represented by a point a γ ∈ R k . • Label distance d(γ, µ) = |a γ -a µ | 2 .

Important special cases

Different Multilabel Regularizers

x 0

A comparison using the "triple-junction" problem Zach, Gallup, Frahm, Niethammer '08

J 1 (u) = 1 2 γ Ω |Du γ |
Lellmann, Kappes, Yuan, Becker, Schn örr '08

J 2 (u) = Ω γ Du γ 2 A , v A = √ v T A T Av
Chambolle, Cremers, Pock '08

J 3 (u) = Ω Ψ (
Labeling regions should have certain spatial relationships, i.e. heaven is always above ground.

Direction-Aware New Regularizer

x 0

Strekalovskiy and Cremers, ICCV 2011

The labeling penalty may depend also on the normal n of the interface between two regions,

d : Γ × Γ × S → R + .
New direction-aware regularizer: Find surface Σ ⊂ R 3 with texture T : Σ → R 3 which optimally matches the input images.

J (u) = sup p∈C γ Ω p γ , ∇u γ dx, with C = {(p γ : Ω → R n ) : p µ -p γ , n ≤ d(γ, µ, n) ∀γ,
Given n images I i : Ω i → R 3 with projections π i : R 3 → R 2

Find surface Σ ⊂ R 3 with texture T : Σ → R 3 which optimally matches the input images.

  V be a convex constraint set, i.e. the solution shall be restricted to lie in C. Thus, we want to solve argmin x∈V {F (x)} subject to x ∈ C, which could be e.g. a constrained least squares problem as in the first talk. This is equivalent so solving argmin x∈V {F (x) + δ C (x)} .

  We observe the image f , which was created from an unknown original viaf = b * u + Gaussian noise.More general, f = Au + noise (linear inverse problem).

  time t I t+1 at time t + 1 flow field u : Ω → R 2 TV-L 1 optical flow (Zach et al. DAGM 2007) Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"

•

  An observed image f on a domain Ω, from which one computes (e.g. via a trained classifier) • a foreground probability distribution P fg (x) = Probability that the point x is in the foreground • A background probability distribution P bg (x) = Probability that the point x is in the background Wanted: • The MAP estimate for a binary segmentation u : Ω → {0, 1}, where u(x) = 1 means a point is in the foreground, and u(x) = 0 means it is in the background: argmax u:Ω→{0,1} P(u | f ) = argmax u:Ω→{0,1} P(f | u)P(u). Bastian Goldluecke Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics" 57

•

  Can be solved globally with functional lifting [Pock, Sch önemann, Graber, Bischof, Cremers '08] • Continuous version of [Ishikawa '03] a γ = e γ ∈ R N Potts model • Example: segmentation • No globally optimal solution possible if N > 2 • Continuous version of [Potts '52]

  Du), Ψ(q) = sup p { q, p : |p γp µ | ≤ d γµ } Bastian Goldluecke Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics" • Stereo assignment cost c γ (x) = |I left (x) -I right (x + γ)|• Linear discontinuity penalty ⇒ globally optimal solution.Images from UCSD lightfield repository Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"

  I i : Ω i → R 3 with projections π i : R 3 → R 2

•

  Find a surface Σ ⊂ R 3 which minimizes the photo-consistency error, argmin Σ Σ ρ(s) ds.• Point on surface• projections look similar in different views • small value of ρ Bastian Goldluecke Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"Classical variational formulation (Faugeras and Keriven, 1998) Find a surface Σ ⊂ R 3 which minimizes the photo-consistency error, argmin Σ Σ ρ(s) ds. • Point not on surface • projections look different in different views • large value of ρ Bastian Goldluecke Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics" Convex functional minimizes photo-consistency error: argmin u:Ω→{0,1} Ω ρ |Du| Silhouette constraints to avoid constant solutions A ray through the silhouette must intersect the surface Kolev, Klodt, Brox, Cremers IJCV'09 Bastian Goldluecke Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics" Convex functional minimizes photo-consistency error: argmin u:Ω→{0,1} Ω ρ |Du| Silhouette constraints to avoid constant solutions A ray through the background must miss the surface Kolev, Klodt, Brox, Cremers IJCV'09 i E i ) • φ i ) = 0 on T φ i = π i • τ 30 cameras, input image resolution 768 × 576 • Initial geometry: Kolev and Cremers, ECCV 2008 • Texture resolution 2048 × 2048

Segmentation example: elongated structures

  Input with User marks log(p bg /p fg ) (normalized)

			x
			0
	Input image Input image	User marks User marks	log(p bg /p fg ) log(p bg /p fg )	Result Result
	... but is not very robust against noise Simple model works pretty well with easy, clean input ... (better classifiers requires, i.e. random forests)

Bastian Goldluecke

Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"

Bastian Goldluecke Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"

Bastian GoldlueckeVariational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"

x 0

2 ,which is again the MAP estimate for Gaussian noise with TV prior.Bastian GoldlueckeVariational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"

x 0Wedel et al. ECCV 2008Bastian GoldlueckeVariational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"

Here, J(u) is the prior. A typical prior is the length of the interface.

Bastian Goldluecke

Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics" 59

The binary case • In T : energy is convex Idea: assume the projection parameters π are unknowns in the superresolution energy,

Optimize alternatingly for texture and projection. In a way, this can be thought of as a continuous version of bundle adjustment.

Goldluecke, Aubry, Kolev, Cremers, IJCV 2013. • Variational methods: a powerful tool for realistic modeling

• Non-convex energies: gradient descent, local minimum, usually requires good initialization

• Convex energies: no local minima, powerful optimization tools (e.g. proximal gradient methods)

• For argmin u {F (u) + G(u)} with convex G and differentiable and convex F , FISTA is a good general-purpose method with optimal convergence rate. • Numerous appliations: variational inverse problems, segmentation and labeling, 3D reconstruction ...

Bastian Goldluecke

Variational Methods, Eurographics 2014 Tutorial "Optimization Techniques in Computer Graphics"