Michel Taïx
email: taix@laas.fr

David Flavigné
email: flavigne@isir.upmc.fr

Etienne Ferré

Human interaction with Motion Planning Algorithm

published or not. The documents may come   L'archive ouverte pluridisciplinaire

Introduction

Manufacturing industries (automotive, aerospace, shipbuilding) take benefits from numerical simulations for several decades. The computer is used almost in all the stages of the product life cycle: design, manufacturing, maintenance. One of the most recent applications of the numerical simulation is the digital mockup for assembly path definition. Thanks to the numerical model of the product, design engineers are able to build the manufacturing processes at the early stages of the product definition. This simulation avoids the building of a physical mockup, which is costly, and time consuming. Moreover, the time needed to build a physical mockup makes the prototype always outdated as soon as it is available. Thanks to the numerical model, the engineers always have access to an up-to-date representation of the current product. Thus, by checking the assembly motion of the product while the design phase, the project team can minimize the risk induced by a non-feasible product. Indeed, an error in the feasibility of assembling a product leads to several months of delay. The engineers use various types of tools to define assembly paths. The goal of the study is the creation of a collision free motion of one part from a position outside the assembly to its mounted position. The assembly problem inside the digital mockup can be solved generally by two approaches : manually or by algorithms. The user can manually create such a motion thanks to highly efficient collision detector warning of a clash while he/she moves the target component. Another solution is the automatic path planner. These algorithms use the numerical model to automatically calculate the collision free trajectory.

When the assembly is done manually, the user can be immersed in a virtual reality environment. Virtual reality allows him to have a better world perception and can use a haptic arm in order to interact with CAD software [START_REF] Bayart | Evaluation of an haptic step-guidance algorithm through a teaching 2d/3d path scenario[END_REF]. Haptics are a recent enhancement that provides an additional perceptual modality in virtual environments [START_REF] Burdea | Haptics issues in virtual environments[END_REF] and can be used for learning . A force feedback device allows the user to test collision-free paths in assembly tasks during the process of industrial design or benchmark maintenance [START_REF] Galeano | Artificial and natural force constraints in haptic-aided path planning[END_REF]. The contribution of haptics can improve the use of virtual reality [START_REF] Ye | A comparative study of assembly planning in traditional and virtual environments[END_REF][START_REF] Lecuyer | Can haptic feedback improve the perception of self-motion in virtual reality?[END_REF][START_REF] Chabal | Emm-3d : a virtual environment for evaluating maintenability from cad models[END_REF]. The use of a space mouse for interaction, is less expensive than a haptic device and can give very good results in many design processes with digital mockup [START_REF] Chen | A 3d mouse for interacting with virtual objects[END_REF][START_REF] Su | Virtual assembly platform based on pc[END_REF]. Manufactured goods using such technology are already available in the market and are used in car and aeronautical industries. The main applications concern the operations of disassembly for maintenance checks as well as during final phase assembly [START_REF] Lecuyer | A haptic prototype for the simulations of aeronautics mounting/unmounting operations[END_REF].

But in complex scenes the user may get lost and need help in finding the solution easily. It is very easy to move a ball in a 3D environment with a haptic device because in that case the orientation of the object has no importance. But this is not the case with a real non-convex object in a complex scene.

Elsewhere, in robotics, a lot of work in motion planning [START_REF] Latombe | Robot Motion Planning[END_REF] has been done to compute free paths in digital models for mechanical systems. With the recent results in random planning algorithm [START_REF] Lavalle | Planning Algorithms[END_REF] it is possible to solve automatically problems for systems with many degrees of freedom. The algorithm computes a collision-free roadmap in a configuration space (CS) where the object is reduced to a point. This point represents the robot's model in the environment and the dimension of CS is equal to the number of degrees of freedom of the mechanical system. The algorithm searches a free path for a point in the roadmap. One example of applying motion planning techniques in real studies is in the maintenance and operation of nuclear power stations [START_REF] Siméon | Computer aided motion: Move3d within molog[END_REF][START_REF] Siméon | Computer aided motion for logistics in nuclear plants[END_REF]. Now, companies use automatic path planning to solve PLM applications [START_REF] Ferré | Progresses in assembly path planning[END_REF][START_REF] Laumond | Kineo cam: a success story of motion planning algorithms[END_REF].

There are mainly two families of methods for building a roadmap, the Probabilistic Roadmap (PRM) [START_REF] Kavraki | Probabilistic roadmaps for path planning in high-dimensional configuration spaces[END_REF] and the Rapidly-exploring Random Tree (RRT) [START_REF] Lavalle | Rapidly-exploring random trees: Progress and prospects[END_REF]. The RRT approach is more interesting to our study because it is faster in the single query case, when the roadmap computation must be limited in time. For really constrained environments [START_REF] Ferré | An iterative diffusion algorithm for part disassembly[END_REF], there are very effective methods to solve problem. But this is not always the case, especially when the object must pass through a narrow passage. It is difficult for the algorithm to find the passage entrance (see for example [START_REF] Boor | The Gaussian sampling strategy for probabilistic roadmap planners[END_REF][START_REF] Wilmarth | Maprm: A probabilistic roadmap planner with sampling on the medial axis of the free space[END_REF][START_REF] Amato | OBPRM: An obstaclebased PRM for 3D workspaces[END_REF][START_REF] Sun | Narrow passage sampling for probabilistic roadmap planning[END_REF][START_REF] Strandberg | Augmenting RRT-planners with local trees[END_REF][START_REF] Zhang | An efficient retraction-based rrt planner[END_REF][START_REF] Flavigné | Improving motion planning in weakly connected configuration spaces[END_REF][START_REF] Rantanen | A connectivity-based method for enhancing sampling in probabilistic roadmap planners[END_REF]). To do so it will have to investigate randomly all the CS dimensions. On the other hand when the algorithm finds the narrow passage entrance in the CS, it will quickly progress to find a solution. Random algorithms progress easily in narrow passages but have no global vision of the solution.

On the contrary, in the case of a free or low-constrained environment, a user will easily find a solution (for example, we quickly steer the key towards the keyhole). On the other hand, a user will have much more difficulties to progress in narrow passages. Thus, the user has a good global vision of the problem, but finds it difficult to make fine movements.

It is noticeable that the algorithm and a user with an interaction device, can both work in a complementary fashion. The complementarity between user and algorithm is the essential idea of this work. The main contribution of the method is to allow simultaneous cooperation between the loop of roadmap search and the loop of interaction for the user. As the roadmap tree iteratively enlarges, it uses the directions of motion of the user device to favor the directions of search extensions.

In this paper, we choose to integrate the human into the planning loop via an interaction device, which can be a space mouse or a haptic arm. The user can thus act on the search, showing areas to explore or to go through. On the other side, the algorithm can gather information while exploring the scene (independently from the user's input) and display them visually and/or haptically. The purpose is not to provide an optimal solution ( [START_REF] Jaillet | Sampling-based path planning on configuration-space costmaps[END_REF][START_REF] Karaman | Sampling-based algorithms for optimal motion planning[END_REF]) but to simultaneously accelerate the search of a solution, and to improve the operator guidance during a motion task in a virtual environment with the help of an automatic path planning algorithm. These results are an extension of [START_REF] Flavigné | Interactive motion planning for assembly tasks[END_REF].

This algorithm can also propose other possible solutions to the user. For example, during the insertion of window glass in a car's door, it may be interesting to test the various solutions to insert the window. For this purpose, human and computer must cooperate in an efficient way along the process.

The main asset of our method is the deep integration of Human into the planning algorithm. This allows a narrow collaboration between human and algorithm (where other human interaction methods have a two-step decomposition that separates the human and the algorithm, our method proposes to make them work in the same loop). This collaboration facilitates better consideration of human constraints that cannot be easily converted to mathematical or logical expressions (preferred solutions because of easier manipulation for the human, or because of flexible obstacles that are not modeled for example). Also, a better perception of the search and solution by the operator : it helps him to decide if the computed solution is applicable (regarding the constraints) in practice and to learn the solution for in real life application.

The paper is organized as follow. We begin by describing previous work in interactive motion planning. The following section describes the basic RRT algorithm and extensions that have been made to answer our problem. The interactive motion planning algorithm, I-RRT, and the main steps are presented in section 4. Simulation results with an industrial application are presented in section 5. Finally, we analyse our results and concludes the paper.

Related works

We can find works where the main problem was to move a robot arm in a virtual environment for off-line programming industrial robot by means of interactive path planning method. In [START_REF] Hea | Haptic-aided robot path planning based on virtual tele-operation[END_REF] the user selects critical robot configurations with a haptic interface to facilitate the automatic path planning research with PRM approach. The advantages of human intuition are exploited to facilitate the robot path planning process. A virtual tele-operation system provides a convenient tool for manipulating a virtual robot arm. They defined a semi-automatic path planning without a real feedback loop between user's motion and motion planning method. [START_REF] Chen | Haptic-based interactive path planning for a virtual robot arm[END_REF] used potential field method to guide user's control on the virtual robot through a haptic device. This approach provides little information to the user because of the many local minima induced by the force fields in complex environment. In these two works the motion planning issue is not the important point of the approach.

On the other hand, idea of taking into account user input has been studied by various authors specialized in motion planning. In [START_REF] Bayazit | Enhancing randomized motion planners: Exploring with haptic hints[END_REF] the authors use user input to work with an automatic motion planner. They show that randomized techniques are quite useful to transform an approximate user input path with collisions into a collision-free path. The idea is to push an approximate user path to free space in CS. The use of a path planning technique based on harmonic functions to generate guiding forces that aid a user in a virtual environment is introduced in [START_REF] Vazquez | Haptic guidance based on harmonic functions for the execution of teleoperated assembly tasks[END_REF], [START_REF] Rosell | Motion planning for haptic guidance. Journal of Intelligents[END_REF]. The idea is to compute a solution channel by cell decomposition that connects the initial and the final configuration and then two harmonic functions are computed over the CS to find a guiding path. RRT algorithms including heuristics based on a study of the workspace are presented in [START_REF] Ladezeve | Interactive motion planning in virtual reality environments[END_REF], [START_REF] Ladeveze | Haptic assembly and disassembly task assistance using interactive path planning[END_REF]. The idea is to discretize the workspace using an unbalanced octree decomposition and generate a free continuous volume between initial and final configuration using A * . In a second time a free path for the object is computed by using RRT approach.

In all these studies, it is necessary to note that the interaction between the automatic search by motion planning algorithm and the user is simplified by a decomposition in two stages. The stage that includes the user can also include new kinds of peripheral devices that allow users to interact in a more natural way with machines. These various devices such as space mice or haptic arms (figure 1) can be used to improve motion planning algorithm with a natural user input, thus speeding-up path computation. A space mouse is a kind of mouse that can move in three dimensions, and even with orientation, the drawback is that the user has no feedback from the machine. This is done by haptic arms, which can send back forces and torques to the user. The Virtuose 6D Desktop includes three articulated branches whose end is fixed on the sphere of seizure. This haptic interface has a kinematics with 6 degrees of freedom (dof), and the return of effort is exerted on the 6 degrees.

RRT Path Planning

We present the principles of the RRT that we modify thereafter to take into account the interaction with the user. The Rapidly-exploring Random Tree (RRT) was introduced in [START_REF] Lavalle | Rapidly-exploring random trees: A new tool for path planning[END_REF] as an efficient sampling method to quickly search in high dimensional spaces for single query motion planning.

The RRT algorithm is shown in Algorithm 1 and the development is illustrated in Fig. 2. Starting at a given initial configuration, RRTs incrementally search the configuration space for a path connecting the initial and the goal configuration. At each iteration a new configuration, q rand , is sampled (RANDOM CONFIG) and the extension from the nearest node, q near , in the tree (NEAREST NODE) toward this sample is attempted. If the extension succeeds (CONNECT ) a new node, q new , and an edge are created in the roadmap.

Algorithm 1

The RRT-CONNECT algorithm T (q init ) for i = 0 to N do q rand ← RANDOM CONFIG() q near ← NEAREST NODE(q rand , T ) if CONNECT(T, q rand ,q near ,q new ) then Add Vertex(T , q new ) Add Edge(T , q near ,q new ) end if end for By definition the RRT tends to grow rapidly in the unexplored regions of the CS and the probability to find a path (if it exist) approaches one as the number of RRT nodes tends to infinity. Many different planners exploit the exploration properties of the basic RRT and try to increase the rapidly exploring property. RRT planner performs well in many applications but gives unsatisfactory results when narrow passages are far from both initial and goal configurations. To improve performances it is necessary to guide exploration by leveraging available information. In our case it is necessary to take into account the interaction with the user.

The main modifications of the basic RRT are to be made on the following stages:

• By definition, the state of the node is binary. If the node is free and the connection path is also collision free, the node/edge is created in the roadmap, otherwise, any information is memorized. It is interesting to memorise information for the user when the CONNECT function return false.

• The sampling area in RANDOM CONFIG is a crucial point to extend efficiently the roadmap. User's interaction must be taken into account in the definition of the sampling area.

• The NEAREST NODE function has not only to take into account the Euclidean distance but also the most interesting nodes with regard to the task of the user.

• At each iteration, the new peripheral device position must be taken into account to compute these different steps

To take into account these various points we propose a new algorithm called Interactive Rapidly-exploring Random Tree, I-RRT, which is presented in the next section.

Interactive Motion Planning System

In this section, the first part gives an overview of the method and the second part details each step of the I-RRT algorithm.

We note that interactive device can be a space mouse or a haptic device. If the interactive device is a haptic device the method uses force through the haptic feedback, in the case of a space mouse, the motion of the mouse is transformed to be used.

Approach Overview

The idea behind the I-RRT algorithm is to take advantage of both a human operator and the RRT algorithm capacities to solve a motion planning problem. Motion planning methods that include a human in the loop are focusing on two steps approaches that separate user input from algorithm execution. Our main contribution is to gather these two steps in a single general loop, shown on Fig. 3.

The goal is to improve the user guidance using motion planning algorithms while these algorithms will benefit from the user's experience. One of the key problem in motion planning is the narrow passage problem. It illustrates the purpose of our method. Traditionnal sampling methods have difficulties to find and go through narrow passages in the configuration space, but a human can easily find the entrance of narrow passages and guide the search. For example, steering a key in a keyhole might be a difficult problem for an automatic algorithm but is easy for a human. This was one of the motivation of integrating a human into the planning loop. The operator will deal with narrow passages issues while the algorithm can help in cases that are difficult for him.

The principle is to let the user move an object (the avatar) in the virtual scene to guide the RRT exploration by introducing a bias toward the avatar during the planning phase. The operator uses an interactive device to move the virtual object. Such algorithm has to take two constraints into account :

• The operator movement has to lead the roadmap extension.

• The planning algorithm has to provide useful information to the user concerning environment, obstacles and exploration status. The proposed planning method is based on a RRT algorithm and can be divided into two distinct parts, represented on Fig. 3 by two loops running in parallel :

1. In the first loop ("interactive algorithm"), the planner explores the configuration space in an automatic way, looking for a collision-free path between start and goal configurations. The sampling and the extension steps are different from the original RRT algorithm, in order to take the user input into account. The user interaction is represented by the pseudo-force F u on Fig. 3. We use pseudo-force instead of force because data are not a force in the physical sense, it is a vector proportional to a force. In this article, we use the same name to represent both a real force and this kind of vector. Using this pseudo-force, the planner extends its roadmap and gathers data on the virtual scene. This data is sent back to the user through another pseudo-force called F a . This pseudo-force is also used in the next iteration of the general loop, along with F u , to extend the roadmap.

2. In the second loop ("interactive device"), the operator moves the avatar in the virtual scene. A pseudo-force F s is applied to the interactive device if it is possible (haptic device) or directly on the avatar otherwise. This pseudo-force may be interpreted as a disturbance by the operator but is intended to act as a guide towards a better area to explore in the configuration space. The device interaction with the virtual scene is interpreted as a pseudo-force F d , which can include F s in the case of a non-haptic device.

These two loops are interacting through the virtual scene. First, the user uses the interactive device to move the avatar (F d ). The avatar movement is used to generate F u (which includes F s either directly or applied to the interactive device if possible). Then, the algorithm builds a RRT roadmap using data gathered during previous iterations (F a ) and F u . This roadmap is used to gather data for next iterations and give feedback to the user (F s ).

The different nature of the two parts of the algorithm allows us to distinguish three different modes of working for the I-RRT algorithm :

• The user is moving the interactive device faster than the algorithm builds its roadmap because he knows where to go to get closer to the goal configuration. In this case, the algorithms tracks the avatar position in the virtual scene. The user's intention is considered predominant. If the solution is obvious to the user, the algorithm follows his movements.

• The user does not move the device. The algorithm keeps exploring the configuration space to gather data on the environment and to give visual and/or (depending on the device) haptic hints to the user. This happens if the user is stuck and does not know how to move the avatar. This mode is close to the original RRT algorithm.

• The user slowly moves the device. In this case, the algorithm extends the roadmap around the user's position while still exploring other part of the configuration space, without the user's influence, to gather more data on the environment.

These three different modes of working are not three separate states of the algorithm, but rather two extremities and the middle of a behaviour line (Fig. 4). Each loop is running at its own speed. The interactive device loop is faster than the algorithm loop (which can vary due to the increasing number of nodes, thus increasing the number of required collision tests). This is not a problem since they are nonblocking. Only pseudo-forces are shared between each loop. In the next section we detail each step of the algorithm that is different from the original RRT.

Description

This section is divided into two main parts, each of them describing one of the two loops involved in the interactive motion planning system.

The interactive device loop

This loop is simple as its only purpose is to deal with the interactive device data. It computes necessary pseudo-forces for feedback from and to the user.

Algorithm 2 Interactive device algorithm loop F s ← ALGORITHM FEEDBACK() if USING HAPTIC() then APPLY FEEDBACK(F s ) Q device ← RETRIEVE DEVICE STATE() F d ← COMPUTE DEVICE PSEUDO FORCE(Q device ) else Q device ← RETRIEVE DEVICE CONFIG() F d ← COMPUTE DEVICE PSEUDO FORCE(Q device ,F s ) end if end loop
The algorithm 2 shows the principle of the interactive device loop. First, it starts by retrieving information from the algorithm to guide the user in the virtual scene (ALGORITHM FEEDBACK). This information is represented by the pseudo-force F s , which points to areas to explore according to the algorithm. Then, it handles two different cases:

1. The interactive device is a haptic device. F s is sent to the device to provide the user with a feedback from the algorithm (APPLY FEEDBACK). The state of the device resulting from this feedback and the user's action is retrieved (RETRIEVE DEVICE STATE) and used to compute a pseudo-force F d (COM-PUTE DEVICE PSEUDO FORCE) that is used for moving the avatar in the virtual scene, movement which is transformed into F u .

2. The interactive device does not have haptic features. In this case the state of the device Q device is mixed with F s to compute F d (COMPUTE DEVICE PSEUDO FORCE).

The interactive algorithm loop

The algorithm loop is based on the RRT motion planning algorithm. It is shown on algorithm 3. The main steps of interactive RRT (I-RRT) are :

• Configuration sampling using computed pseudo-forces (F u , F a ). We are looking for an efficient way of combining user movements with the algorithm automatic search.

• Choosing the nearest neighbour in the roadmap.

• Extending the roadmap.

• Updating node data given the result of the roadmap extension.

• Computing a pseudo-force F a aimed at attracting the user to relevant areas to explore. It is computed using data contained in the roadmap nodes.

• Computing a pseudo-force F u from the motion of the avatar in the virtual scene. This motion is generated by the user interaction with the device (F d ) and the virtual scene (obstacles may stop the movement).

• Sending back some useful information to the user by using haptic feedback or visual hints.

Algorithm 3 I-RRT algorithm T ← q init ,q goal while !PATH(q init ,q goal ) do q rand ← SAMPLING(F a , F u ) q near ← NEAREST NODE(q rand , T ) if CONNECT(T, q rand , q near , q new ) then T ← q new ,EDGE(q near , q new ) end if GATHER DATA(T , q near , q new ) F a ,F u ← COMPUTE PSEUDO FORCES(T , q user ) end while Algorithm 4 COMPUTE PSEUDO FORCES(T , q user ) L near ← NEAREST NEIGHBORS(q user , T ) B ← Compute Barycenter(L near ) F a ← Compute FAlgo(B, q user );

F u ← Compute FUser(q user ); F r ← Compute FResult(F a ,F u );
In a first time, the pseudo-forces computation is described because they are the basis of interaction in the method and are used in other parts of the algorithm. Then, the parts of the original RRT algorithm that were modified for our algorithm is explained. Finally, the data gathering and storage and the haptic feedback are described.

Attractive pseudo-force computation

Interaction being the basis of the I-RRT, the contribution of each interacting part (the user and the algorithm) has to be defined. These contributions have to be flexible and have to reflect the planning constraints (obstacles, goal, ...). This means that they vary continuously during the planning, according to the virtual scene. This is achieved through computation of pseudo-forces which are aimed at representing each interacting part influence over the path searching process.

The contribution of the algorithm part is defined by the computation of the attractive pseudo-force F a (Compute FAlgo in algorithm 4). This computation is illustrated on Fig. 5. First, the p nearest neighbours of the user position q user are retrieved from the RRT roadmap (NEAREST NEIGHBORS) and Fig. 5-a). The barycenter B of these nodes is then computed (Compute Barycenter). The data contained in the nodes are used as weights during this computation (Fig. 5-b and5-c). The computation of weights is detailed below in the paragraph "Node Labelling". F a is the vector starting from q user and directed towards the barycenter (Fig. 5-d).

User pseudo-force computation

The other interacting part of the system is the user. Its contribution is defined by the pseudo-force F u (Compute FUser). F u results from the avatar movements in the virtual scene and its interaction with obstacles. As the avatar movements are controlled by the user and the pseudo-force F d , F u can be deduced from F d most of the time. However, this is different if the avatar is not allowed to collide with obstacles during the path search, when there is a contact with an obstacle. The user will apply a force on the device (F d = 0), but the avatar will not move if it implies a collision with the obstacle (F u = 0).

Resulting pseudo-force computation

As the path search has to reflect both part of the interaction, the two previously computed pseudo-forces are mixed into a resulting pseudo-force called F r (Equ. 1 and Fig. 6-a). F r takes into account the user's intention and the algorithm suggestions for the path search. We can retrieve the three modes of working illustrated on Fig. 4 by studying the variations of F r :

(a) (b) (c) (d)
Figure 5: Computing F a pseudo-force.

• The user moves fast : F u > F a . In this case, F a is negligible compared to F u and then F r is almost equivalent to F u . The user is predominant and leads the path search.

• The user doesn't move : F u = 0 and F a > 0. In this case, F r = F a . The algorithm is predominant.

• The user moves slowly : The user and the algorithm are working together.

F r = α • F u + (1 -α) • F a (1)
In equ. 1, α can be used to modify manually the influence of each part. If α = 0, the algorithm is equivalent to an RRT algorithm. If α = 1, the algorithm is only tracking the user position.

Sampling method

The resulting pseudo-force F r is computed in order to guide the sampling of the RRT algorithm. The sampling is the heart of probabilistic methods. If the sampling can be controlled, the path planning is controlled as well. The F r pseudo-force computed in the previous step is used by the sampling method to guide the path search (Fig. 6). The sampling method is a gaussian sampling along F r direction. If the user does not move or move slowly, the sampling is equivalent to a gaussian sampling around the user position but if he is moving fast the sampling anticipates the user's movements. This sampling method can be seen as a position tracking method. The standard deviation σ of the gaussian sampling is used to tune the influence of the algorithm on the search.

But this sampling method has some drawbacks. If the user moves too fast, the RRT algorithm may not be able to link newly sampled nodes to the existing roadmap due to obstacles. To deal with this problem, three other sampling methods were investigated.

• The first method allows creation of several local trees rooted at nodes that cannot be linked directly to the main tree. It is useful in the case of I-RRT method to keep track of user's movements even if the connection with the main tree was lost (this is the drawback of the original sampling method). The missing connection is made by the algorithm part of the interactive motion planner.

• The second method stores the path traveled by the user and samples around this path (gaussian sampling). The advantage is that the connection to the main tree is never lost, but it slows down the resolution. In this method, the influence of the algorithm can also be tuned using σ .

• The third method adds a uniform sampling in the configuration space. This method allows a global exploration of the virtual scene, independently of the user's position, and then the data gathering is also more global compared to the original sampling method. To preserve the interactivity, the uniform sampling should be used with one of the other methods. A parameter β indicates the ratio between the usage of the two used methods.

All sampling-based algorithm are based on random sampling which is dense with probability one to ensure probabilistic completeness. In the worst case, without human intervention, the algorithm will proceed as an RRT and inherits the property of their completeness. With only a sampling around F r it is not possible to guarantee this property. It is necessary to add a uniform sampling to the main search direction to ensure probabilistic completeness in case the user does not move the interface. This is why the uniform sampling method is always used in combination with one of the other methods described above.

Nearest neighbour

The choice of the nearest neighbour q near for extension is made using the distance d (Equ. 2) between the newly sampled node and each node of the tree (NEAR-EST NODE). This distance is also used for the choice of the nearest neighbours of the user's position (NEAREST NEIGHBORS) during barycenter computation.

d = n ∑ i=1 ω i • (do f i sample -do f i node ) 2
(2)

Extend method

The extension step of the I-RRT is the same as a RRT-Connect algorithm : the nearest neighbour in the tree is extended towards the new sample until it reaches this sample or collide with an obstacle. Given the result of this extension, the extend method can return three different states :

• Collision : the extension failed because q near was near an obstacle. No new node is added to the tree.

• Obstacle : the extension was partially done because an obstacle was between q rand and q near . A new node q new is added to the tree.

• Success : q near and q rand were successfully connected. q rand is added to the tree (q new = q rand ).

For the original RRT algorithm, this is the last step of an iteration. For the I-RRT algorithm, we need to define the feedback for the user. This feedback is computed as hints that can be either haptic hints or visual hints (or both). With this feedback, the interaction loop is closed. The steps needed for this are detailed in the following paragraphs.

Node Labelling

The first type of feedback is visual. The algorithm uses data stored in nodes to display visual hints for the user. These visual hints are modeled as a coloring of nodes and edges according to labels computed from return states of the extend step. These labels are the data stored in the nodes. In the original RRT algorithm, configuration are sampled randomly and the only piece of information returned is their collision state with the scene obstacles. In the I-RRT, these collision tests are used to label roadmap nodes. Two types of labels were defined for visual feedback :

• The distance label is a function to guide towards the goal. It is the distance from the current node to the goal node.

• The collision label is a function of the covered distance during the extension, the number of extension of the current node and the number of collisions resulting from extension attempts. It is computed using the return states of the extend method. It reflects the constraint of workspace to extend this node.

However, as the roadmap grows the screen can be overloaded and the user will not benefit from this feedback. In practice, the graph display can be switched off to ease the search and overcome this drawback.

Haptic feedback computation

The second type of feedback is the haptic feedback. A pseudo-force is inserted in the interaction loop to give a more intuitive guidance for the user and to avoid cluttering the screen with visual artifacts. To compute this pseudo-force F s , the previously computed attractive force is used :

F s = k × F a ( 3
)
In equ. 3, k is a scale factor set by the user before algorithm execution.

In the case of the operator using a haptic device, F s is sent directly to the device. The user constantly feels an effort from the device which guides him towards areas to explore, according to the algorithm. In other cases, the haptic feedback is applied to the avatar movements in the virtual scene. The user does not feel the feedback directly but will be able to see the perturbations due to F s on the movements of the avatar.

F s also allows the user to stay within range of the algorithm roadmap and thus limit the problems mentioned in paragraph Sampling method.

Results

System implementation

The algorithm was developed in C++ using the HPP (Humanoid Path Planner) platform developed at LAAS and based on the KineoWorks software 1 . Experiments were made on a dual-core PC, 2,1GHz and 2GB RAM.

In practice, using a local trees sampling method mixed with a uniform sampling method minimizes the user's influence on the path search. The main idea behind the I-RRT method is to make user and algorithm interact equally, for this reason these two methods will not be used together in most of the cases. In the following experiments, the uniform sampling method will be used with the sampling method centered on the user's position (along F r ) and with the sampling method along the user's path. In each case, β = 3. The pseudo-forces F u and F a have the same influence on the search (α = 1 2 ). The nearest neighbours search is made by comparing each node of the tree. The number of nodes p used to compute F a is five. For this computation, all the neighbours should remain in the same connected component.

The required time to solve the path planning is dependent on the user's ability to manipulate the interactive device. If he does not move or is stuck in a dead-end, the planning time can increase and the algorithm will not benefit from the user's guidance. In the following case, we decided to let the user go through obstacles (contact forces are not considered) so that user's movements will not be limited, and he will still be able to guide the algorithm. The F s pseudo-force replaces the contact forces to inform the user of obstacles limiting the exploration capacity of the algorithm in the current area. With this, even if the user path is in collision with the environment, the algorithm can still find a collision-free solution.

The RRT generated path is usually a broken line. In the context of a physical application, this path is not feasible. To make it physically feasible, this solution has to be deformed to a smooth path. One of the most common solution is to use the final path nodes to build a bezier curve. The drawback is that the deformed path might not be collision free, and making it so is not trivial. Here it is not a problem because the movement of the object is achieved by the user who performs an implicit path smoothing.

Simple 2D case

The main goal of this first case is to illustrate the principle of the interactive algorithm, and the different modes that can be observed during a planning query. In this case, the user moves the avatar with a space mouse in the virtual scene shown in Fig. 7, 8 and9. The coloring of nodes was not activated for this example.

The user has to reach q goal from q init . The avatar (the green block) has three degrees-of-freedom : x, y, θ . The first part of the environment is simple and user and algorithm can both easily find the path. The user being faster, he leads the search during this part. Then, the environment becomes tighter and there is a crossroad (Fig 7-a). The way is obvious for the user, but in a non-interactive case the algorithm would have difficulties to find the narrow passage. The user is predominant and leads the search into the narrow passage.

After the narrow passage, the environment has two rooms with no exit. These rooms would have required a large amount of time to be explored in a non-interactive case. Instead, the user avoids them and leads the search in the right direction (Fig. 7-b). The algorithm is not predominant but can initiate a slow exploration towards regions left out by the user (small branches are developed towards the rooms). On Fig. 8-a, the exploration made by the algorithm part is more visible. Meanwhile, the user saw a possible passage and tries it. The passage is too narrow for the user to find a solution here. While the user is trying to develop nodes in this passage, the algorithm continues exploration in other parts of the environment (uniform sampling), as shown on Fig. 8-b. The algorithm finds a solution left out by the user because it seemed to be too narrow for the robot (red obstacle on Fig. 8-b). The user mistake is due to his point of view of the virtual scene. The red obstacle is placed below the plan in which the robot can move. This trap may seem artificial but it illustrates perfectly problems that can be met in complex 3D environments.

(a) The I-RRT roadmap after the solution was found (b) A Bi-RRT roadmap Figure 9: Comparison between the roadmap obtained after using the I-RRT algorithm and after using a Bi-directional RRT algorithm.

On Fig. 9, we can compare the roadmap built by the I-RRT algorithm with the roadmap of a bi-directional RRT algorithm. The interaction between the user and the algorithm allows a lower number of nodes and a faster exploration (few minutes) than a RRT algorithm. The RRT algorithm made a dense exploration and covered the maximum of the space. This kind of exploration requires a lot of time and resources (> 1 hour).

Simple 6D case

In this example, shown on Fig. 10-a, the user has a haptic arm (Virtuose 6D Desktop -Haption) to manipulate a S-shaped object in order to cross a wall through a hole in its center. The wall dimensions are 700 × 700 × 20 and the hole's are 80 × 80. The object is composed of three sticks assembled along three edges of a virtual cube. The nodes of the roadmap represent the position of the center of this cube (and thus are not a point on the object itself). Each stick is 20 × 20 × 100. The environment has only one difficult passage which requires the help of the user. The sampling method used are user tracking (using F r ) and local trees are allowed. The node coloring is activated. The configuration space is bounded to force the user to go through the hole, and the object has six degrees-of-freedom.

Solving this planning query requires that the user has a good ability to manipulate the interactive device, the passage being really narrow. The user has to slowly guide the search through the hole in order to improve the search. On Fig. 10-b and Fig. 11, the developed roadmap is displayed with nodes colored according to the labelling. Nodes that are closer to the goal are displayed in a darker color than others. The algorithm part covered both part of the wall (due to the local trees) so that the user has only to guide the part in the narrow passage. The roadmap density can testify to the dexterity of the user. Here, the roadmap is not dense nor sparse, which shows that the user has an average dexterity. Figure 11: The roadmap is developed on both sides of the wall.

Industrial case

The following example is based on a real automotive case. The goal is to find a disassembly path for a car silencer. In this example, the configuration space is bounded. The user is using the same haptic arm as in previous example (Fig. 13) to move the silencer and extract it from its location. The extension of a tree from the goal configu-ration is allowed, because it does not modify the interaction user-algorithm : the goal configuration is in a large free space, easily accessible after the silencer is removed. The silencer has to be remove from bottom to be placed on the top.

The query is solved in few minutes by the user (Fig. 12). The main difficulty is to extract the silencer from the start area, under the green part. This disassembled query can be solved more quickly by specific automatic algorithms [START_REF] Ferré | An iterative diffusion algorithm for part disassembly[END_REF]. This problem does not favor interaction between user and algorithm but, in the context of guiding the user, the I-RRT method is still valid because the interactivity allows to control the solution and to take into account constraints imposed by the user. On Fig. 12-b, we can see that the roadmap around the goal configuration is developed mainly by the algorithm part and the roadmap around the start configuration by the user part.

Analysis

Due to its dual nature, the interactive RRT performances are strongly dependent on several parameters. The first parameters, on the user side, are the user capacities to visualize the virtual scene and the potential solutions and his skills at manipulating the interactive device. Then, on the algorithm side, the choice of the sampling method to use is determinant. Last, the parameters needed for exchanges between each part of the I-RRT are also to take into account (k for the haptic feedback and α for computing F r ). The standard deviation used for gaussian sampling methods is also important for exchanges.

The settings of the algorithm parameters will influence its performances and the aspect of the resulting path. They can be used to adapt the algorithm to the user skills.

In this section, we choose the three parameters that seemed to be the most representative of each type of parameters. We analyze the influence of these parameters on the algorithm:

• The user dexterity for the user part.

• The sampling method for the algorithm part.

• The standard deviation for the exchange part.

For testing purposes, we choose a simple environment with square block robot in a three-dimensional configuration space (X, Y , θ ). For the analysis α = 1/2 and k = 100.This environment is shown on Fig. 14. The user dexterity is a parameter difficult to model. In the tests, it is represented by the user speed along a predetermined path path p (Fig. 14). For each iteration of the algorithm part (each time a new configuration is sampled), the user's position is moved by a step ε v along this path. ε v is used as equivalent of the user speed. The tests were run for four different values of ε v .

Four different sampling methods were chosen for the tests:

• POS: Gaussian sampling around the user's position.

• POS +CS: Gaussian sampling around the user's position and uniform sampling in the configuration space.

• PAT H: Gaussian sampling around the user's current path (path part along tra j p between start position and user's position).

• PAT H +CS: Gaussian sampling around the user's current path and uniform sampling in the configuration space.

The sampling method defines how to use the user input to gain information, explore the environment and find the solution. All the methods are based on gaussian sampling as it allows randomness to be guided (by user input in our case). The standard deviation σ is the parameter that defines the influence of the guide on the gaussian sampling. Four different values of standard deviation were chosen for the tests.

Visual and haptic feedback are not considered during the tests. The results shown are average over 50 runs. The time limit for each test was fixed to 900 seconds.

Sampling

Fig. 15 shows the computation time given the user speed for each values of the standard deviation. Curves (a) to (d) shows the results for each particular sampling method.

On each curve, we can observe that the computation time is higher (almost reaching the time limit) for lowest user speed values. This is due to the fact that in this case, the user is not fast enough to find a solution before the time limit. When a uniform sampling method is used, however, the problem can be solved within time limit.

The curve (a) shows the results for the POS sampling method. In this case, high values of σ (red and magenta on Fig. 15-a) shows the worst performances because it makes the user position tracking less efficient, and the narrow passages of the environment reduce the probability of linking a new sample to the existing roadmap for this method. On the contrary, if the sampling is precise enough when following the user (low values of σ , blue and green), the increase of the user speed will benefit to the algorithm.

The curve (b) shows the results for the POS +CS sampling method. In this case, a uniform sampling in the configuration space is added to the gaussian sampling around the user position. (one sampling of three is chosen with the uniform method). We observe on the curve that it introduces a latency in computation time for low σ values, but speeds up the resolution for high values. The addition of the uniform method reduces significantly the dependence of the algorithm on the standard deviation.

The curve (c) shows the results for the PAT H sampling method. In this case, σ has no influence over the path planning process. The resolution time always benefits from the increasing user speed. The memorization of the user path allows to deal with narrow passages and problems of link between the samples and the existing tree.

The curve (d) shows the results for the PAT H + CS sampling method. The computation time increases a little bit compared to the PAT H sampling method, and the differences due to σ are visible, but the addition of the uniform sampling has no significant influence on the planning. This figure allows us to compare each sampling method directly given a single value of standard deviation. On curves (a) and (b) (low standard deviation values), we can observe that the addition of a uniform component to the sampling method increases the computation time (green and magenta). But while the PAT H + CS method follow the same kind of path as PAT H and POS (due to the memorization of the user path), the method POS +CS has worse results as the user speed increases. This is due to the slow down caused by the uniform sampling and the linking problem when approaching narrow passages. The POS method gives the best performances: this method benefits from the user speed more than others.

Standard deviation

On curves (c) and (d), we can see that the two POS based methods are giving worse results because of the increase of the standard deviation. We can observe on curve (d) that, without the supply of a uniform component, the POS method cannot resolve the problem within the time limit (blue). On curve (a), we can observe that computation time are high for every sampling method. As said earlier, this is due to the user speed being to low to reach the goal before the time limit. Sampling methods without a uniform component cannot resolve the path planning problem in this case, but we can see that other methods can do it before the time limit.

User dexterity

On curve (b), we can observe that no matter the σ value and the sampling method, the performances are similar when ε v = 1.0.

On curves (c) and (d), we see that PAT H based sampling methods are not affected either by the user speed or the standard deviation. The POS + CS method is affected by the user speed but not by the standard deviation (green). This is due principally to the uniform sampling because if we look at the POS method, the computation time increases with the standard deviation (narrow passages and linking problem). However, for low values of σ and high user speed, this method is the most efficient.

In this section, three parameters of the I-RRT algorithm were analyzed. The results showed that each of these parameters are important and modifies the performances of the algorithm. They must be adapted to the user capacities and the problem environment.

Sampling methods without a uniform component (POS,PAT H) are more adapted for an experienced user that can move the avatar in the virtual scene fast and with precision. They can take full benefit from the user guidance. Sampling methods with a uniform component (POS +CS,PAT H +CS) are more efficient when the user is inexperienced (slow user speed) or cannot find the solution. They ease the problem resolution, but require more time. Path tracking methods (PAT H,PAT H + CS) are very little affected by user speed or standard deviation. They constitute a good compromise for unknown user capacities or environment.

Conclusion

In this paper, we presented a new method allowing a human operator and a path planning algorithm to cooperate for finding a path in an efficient way. The main idea is to gather the algorithm and the user input into the same planning loop, and make them exchange information about the search using visual hints and by exchanging pseudoforces. The method is made of two separate loops one including the user interaction and the other the planning algorithm. These two loops are communicating through the virtual scene using the pseudo-forces. These pseudo-forces are generated by the user with the interactive device and by the algorithm which uses collected data to compute them. The principal advantage is that user and algorithm are not decoupled in a two step algorithm but work together at the same time, with a varying influence over the path search. Our approach based on RRT algorithm generates non-optimal paths. [START_REF] Karaman | Incremental sampling-based algorithms for optimal motion planning[END_REF] shows that RRT does not converges to an optimal solution almost surely as the number of samples increases. In our problem, the path optimality is less a priority than finding a path while taking into account the operator constraints (constraints that are not mathematically described easily but are known by the user : relevant work experience, business processes). Different cases are presented to show the validity of the method : a simple two-dimensional case which illustrates the principle, a difficult six-dimensional case, and an industrial case from the automotive industry. The studied cases show that the search benefits from each part : the user speeds up the search by going through difficult passages and the algorithm gives global information to help the user and improve approximate solutions. The usage of both a haptic arm and space mouse in these cases shows that the method can be easily adapted to different devices. Finally, the analysis shows that different parameters can be tuned to adapt to different situations.

In this paper, the contact forces between the object and the obstacles are not considered because they can cause difficulties to the user. They can be useful in the path search if slip is considered and if the user has the possibility to easily deactivate these contact forces when he needs it. In several manipulation tasks slip can guide the user to the solution (ex: steering a key in its keyhole). But if the user get stuck in a dead-end or in a very narrow passage, he must be able to deactivate these contacts to get out of this situation easily and resume the path search in another part of the virtual scene.

It could be interesting to analyze the influence of the different parameters on the algorithm performances for finding paths (number of nodes and number of iterations) and on the shape of the solution path (is it a broken line or a smooth curve). The dynamic tuning of the parameters during the path search phase could be done using this analysis. This allows the algorithm to automatically adapt to the topology of the different parts of the virtual environment or to changing environments. The comparison between solutions obtained by a search leaded by the user (the user is more influential) with solutions obtained by a search leaded by the algorithm can show the respect of human-defined different criteria (is the path natural? intuitive? human feasible?). As the interactive device is an important part of the method, it could be interesting to study its influence on the search and the difference of performances between haptic and non-haptic devices. The visual feedback can also be improved by changing colors (according to labels) in a continuous gradation in a background image instead of displaying each node and edge. This would lighten the display. The cases studied in this paper focus on free-flying objects and assembly problems, but the method can be applied to different situations such as moving animated characters or multi-user path search for very complex problems (using several interactive devices). The algorithm can theoretically be ran on embedded systems but, as is, it will not be suitable for real-time applications as it needs computing power.

acknowledgements

Thanks to E Ferré to make available the Haption haptic Device in KineoCAM Company.

Figure 1 :

 1 Figure 1: The Virtuose6d and the space mouse.

Figure 2 :

 2 Figure 2: RRT Algorithm.

Figure 3 :

 3 Figure 3: General interaction loop between the algorithm and the user (represented by the interactive device) through the virtual scene.

Figure 4 :

 4 Figure 4: Behaviour line of the algorithm. The algorithm state can be anywhere between fully automatic behaviour (left) and fully manual behaviour (right).

Figure 6 :

 6 Figure 6: Sampling a new configuration using user and algorithm interaction.

  (a) Few nodes are created (b) The user has a dominant influence on the search

Figure 7 :

 7 Figure 7: The user is leading the search and the algorithm is tracking his position.

Figure 8 :

 8 Figure 8: The user is stuck in a dead-end, the algorithm is searching other possibilities to reach the goal.

Figure 10 :

 10 Figure 10: A case with a 6-dimensional configuration space.

  Figure12: The silencer has to be remove from bottom to be placed on the top.

Figure 13 :

 13 Figure 13: Extracting the silencer using an interactive device and the I-RRT algorithm.

Figure 14 :

 14 Figure 14: Testing environment. The user's path is displayed in black.

Figure 15 :

 15 Figure 15: Comparing the standard deviation effect on resolution time given the user speed for different sampling methods.

Fig. 16 Figure 16 :

 1616 Fig.16shows the computation time given the user speed for each sampling method. Curves (a) to (d) shows the results for each value of σ .

Fig. 17 Figure 17 :

 1717 Fig.17shows the computation time given the standard deviation for each sampling method. Curves (a) to (d) shows the results for each value of ε v .

KineoWorks is the path planning dedicated Software Development Kit developed by KineoCAM


R Meftah

D Duhamel

Modèle de contact dynamique pneumatique/chaussée

Dans ce travail, nous proposons une démarche de modélisation de contact dynamique entre pneumatique et chaussée. La démarche est basée sur les techniques de convolutions. Les fonctions de Green du pneumatique sont supposées connues. La démarche consite d'abord à décomposer ces fonctions dans une base modale, ensuite les paramètres identifiés sont utilisés pour construire une convolution rapide. Le déplacement est calculé par une convolution rapide entre forces de contact et fonctions de Green. Le problème de contact est traité en ajoutant une condition cinématique sur les vitesses. La démarche est appliquée sur un modèle simplifié de pneumatique sur une chaussée sinusoïdale.

Mots clefs : contact dynamique ; Convolution ; pneumatique-chaussée 1 Introduction

Le contact pneumatique/chaussée représente la source principale du bruit de trafic routier, en particulier pour des vitesses supérieures à 50 km/h. L'une des difficultés majeures de la modélisation de ce phénomène est la prise en compte de l'effort de contact dans le comportement dynamique du pneumatique. La complexité de la géométrie du couple chaussée/pneumatique rend le problème difficile. Ce genre de problème est souvent traité à l'aide de l'intégration temporelle en utilisant des approches basées soit sur la méthode de pénalisation [START_REF] Kanto | A dynamic contact buckling analysis by the penalty finite element method[END_REF] soit sur la méthode de multiplicateurs de Lagrange [START_REF] Bathe | A solution method for planar and axisymmetric contact problems[END_REF]. La plupart de ces méthodes sont souvent difficiles à mettre en oeuvre ou bien mènent à des calculs lourds. Les techniques de convolution peuvent être utilisées pour traiter ce genre de problème. La réponse dynamique est calculée par un produit de convolution entre les fonctions de Green du système dynamique et les forces de contact. Cette approche est largement utilisée dans la littérature [START_REF] Nordborg | Wheel/rail noise generation due to nonlinear effects and parametric excitation[END_REF][START_REF] Mcintyre | On the oscillations of musical instruments[END_REF][START_REF] Wang | New analysis method for a thin beam impacting against a stop based on the full continuous mode[END_REF][START_REF] Kropp | A model to describe the rolling noise of a smooth radial tire on a rough road surface[END_REF][START_REF] Larsson | Modelling of dynamic contact-exemplified on tire/road interaction[END_REF][START_REF] Andersson | Time domain contact model for tyre/road interaction including nonlinear contact stiffness due to small-scale roughness[END_REF]. La démarche que nous proposons dans ce travail est basée sur la décomposition de la fonction de Green dans une base modale. Les paramètres modaux identifiés sont ensuite utilisés pour construire une convolution plus rapide. Le contact est géré en utilisant les conditions de contact classique, puis en ajouttant une condition cinématique sur les vitesses de déplacements.

Convolution rapide

La fonction de Green est supposée connue. L'idée est d'approximer cette fonction par une combinaison linéaire de la contribution de chaque mode. Les paramètres modaux sont identifiés à l'aide de l'algorithme LSCE.

G(ω) ≈ k=N m ∑ k=1 A k -ω 2 + 2 jξ k ωω k + ω 2 k (1)
Après avoir déterminé les paramètres modaux A k , ω k et ξ k et tronqué la décomposition à un ordre N, le déplacement dans le domaine temporel peut être obtenu par une transformation de Fourier inverse de l'équation (1)

g(t) = k=2N m ∑ k=1 A k ω k d e -ξ k ω k t sin(ω d k t) (2) Avec ω k d = ω k 1 -ξ 2 k
Le déplacement peut être calculé par un produit de convolution

u(t) = Z t 0 g(τ)q(t -τ)dτ = Z t 0 g(t -τ)q(τ)dτ (3)
On remplace l'expression de g(t) dans l'équation (3), puis on sépare les variables t et τ , nous pouvons écrire le déplacement sous la forme

u(t) = k=2N m ∑ k=1 A k ω d k e -ξ k ω k t sin ω d k t α k (t) -cos ω d k t β k (t) (4)
où α k (t) et β k (t) sont calculés par les formules

α k (t) = Z t 0 e ξ k ω k τ cos(ω d k τ)q(τ)dτ (5)
β k (t) = Z t 0 e ξ k ω k τ sin(ω d k τ)q(τ)dτ (6)
3 Problème de contact

Démarche générale

Les conditions classiques de contact sont données par les équation suivantes

u(t) = u r (t) ; f c (t) > 0 (7) u(t) > u r (t) ; f c (t) = 0 (8)
Où u(t) est le déplacement d'un point du système, u r (t) est l'amplitude du profile de la chaussée au même point et f c (t) la force de contact au même point. Le principe de calcul est détaillé ci-dessous :

1. Premièrement, l'historique du déplacement u h (n∆t) est calculé en supposant que f c (n∆t) = 0.

u h (n∆t) = k=N ∑ k=1 A k ω d k e -ξ k ω k n∆t sin ω d k n∆t α k ((n -1)∆t) -cos ω d k n∆t β k ((n -1)∆t) (9)
2. L'historique du déplacement est comparé avec l'amplitude du profile de la chaussée.

∆x(n∆t) = u r (n∆t) -u h (n∆t) (10)
3. Si l'hypothèse de non contact (∆x < 0) est vérifiée, le déplacement est calculé par :

u(n∆t) = u h (n∆t) (11)
Ensuite, la même procédure est répétée pour le prochain pas de temps ( f c ((n + 1)∆t) = 0).

4. Si ∆x ≥ 0, la force de contact est calculée par le principe de déconvolution détaillé dans la prochaine section.

5. La même procédure est répétée pour le prochain pas de temps ( f c ((n + 1)∆t) = 0).

Déconvolution

la force et la fonction de Green dans le domaine temporel, le déplacement peut être calculé par un produit de convolution :

u(t) = Z t-∆t 0 g(t -τ)q(τ)dτ + Z t t-∆t g(t -τ)q(τ)dτ = Z t-∆t 0 g(t -τ)q(τ)dτ u h (t) + Z ∆t 0 g(τ)q(t -τ)dτ = u h (t) + Z ∆t 0 g(τ)q(t -τ)dτ (12)
Si on dérive cette équation (12) on trouve

v(t) = u ′ (t) = Z t-∆t 0 g ′ (t -τ)q(τ)dτ + Z t t-∆t g ′ (t -τ)q(τ)dτ = Z t-∆t 0 g ′ (t -τ)q(τ)dτ v h (t) + Z ∆t 0 g ′ (τ)q(t -τ)dτ = v h (t) + Z ∆t 0 g ′ (τ)q(t -τ)dτ (13) Notant Y = [u(t) v(t)] T et Y h = [u h (t) v h (t)] T , on peut écrire Y = Y h + Ψ(q) (14)
Où Ψ est un opérateur intégral qui traduit l'influence de la force de contact à l'instant t.

Dans le cas où, il y a contact, le déplacement est imposé par le profile de la chaussée qui est supposé rigide. On suppose aussi que le point de contact suit la chaussée, ce qui peut être traduit par les conditions suivantes :

u(t) = u r (t) ; v(t) = du r (t) dt ( 15
)
Sous la forme vectorielle

Y = Y r = [u r (t) du r (t) dt ] T (16)
En utilisant la décomposition modale de la fonction de Green, on peut écrire

u(t) = k=N ∑ k=1 A k ω d k e -ξ k ω k t sin ω d k t α k (t) -cos ω d k t β k (t) (17) et v(t) = - k=N ∑ k=1 A k ξ k ω k ω d k e -ξ k ω k t sin ω d k α k (t) -cos ω d k t β k (t) + k=N ∑ k=1 A k e -ξ k ω k t cos ω d k t α k (t) + sin ω d k t β k (t) (18)
Les équations ( 17) et (18), peuvent s'écrire sous la forme matricielle

Y = k=N ∑ k=1 B k x k (19) où x k = [α k (t) β k (t)] T (20) et B k = A k e -ξ k ω k t     sin(ω d k t) ω d k - cos(ω d k t) ω d k cos ω d k t -ξ k ω k ω d k sin ω d k t sin ω d k t + ξ k ω k ω d k cos ω d k t     (21)
Les conditions de contact

Y r = Y h + Ψ(q) (22)
L'objectif est de calculer la force de contact qui vérifie à la fois ces deux conditions. La force de contact est reliée avec le déplacement par l'opérateur Ψ.

∆Y = Y r -Y h =    Z ∆t 0 g(t -τ)q(τ)dτ Z ∆t 0 g ′ (t -τ)q(τ)dτ    (23)
Ces intégrales peuvent être calculées en utilisant deux points de Gauss. Les valeurs de la force de contact aux deux points de Gauss sont calculées en inversant l'opérateur Ψ.

q = q 1 q 2 = Ψ -1 (t 1 ,t 2 )∆Y (24)
Cet opérateur est donné par

Ψ =   g(t -t 1 ) g(t -t 2 ) g ′ (t -t 1 ) g ′ (t -t 2 )   ∆t (25)
où,

g(t) = k=N ∑ k=1 A k ω d k e -ξ k ω k t sin ω d k t (26)
Dans l'intervalle de temps [t t + ∆t], la force de contact est la moyenne des valeurs de la force aux deux points de Gauss t 1 et t 2 qui sont définis par

t 1 = t + 1 - 1 (3) ∆t 2 ; t 2 = t + 1 + 1 (3) ∆t 2 (27)
Connaissant les valeurs de la force de contact aux deux point de Gauss t 1 and t 2 , on peut calculer les paramètres α k (t + ∆t) et β k (t + ∆t) par les équations de récurrence suivantes

α k (t + ∆t) = α k h (t) + e ξ k ω k t 1 cos(ω d k t 1 )q 1 + e ξ k ω k t 2 cos(ω d k t 2 )q 2 2 ∆t ( 28
)
β k (t + ∆t) = β k h (t) + e ξ k ω k t 1 sin(ω d k t 1 )q 1 + e ξ k ω k t 2 sin(ω d k t 2 )q 2 2 ∆t (29)

Conclusion

Une nouvelle formulation de contact dynamique entre pneumatique et chaussée a été présentée. Le formalisme a deux avantages : premièrement sa rapidité en utilisant une convolution basée sur une décomposition modale. Deuxièmement, l'introduction de la condition cinématique de contact facilite le calcul des efforts de contact. La méthode a été appliquée dans le cas d'un modèle simplifié de pneumatique sur une chaussée sinusoïdale, elle peut être encore plus avantageuse en terme de temps de calcul dans le cas d'un pneumatique réel sur une chaussée réelle.

FIG. 2 -

 2 FIG. 2 -Forces de contact aux différents points de contact

4 Application

 4 La modélisation d'un pneumatique réel avec tous les détails est très complexe. Dans la littérature, le modèle d'anneau circulaire sous fondation élastique est souvent adopté pour modéliser un pneumatique dans les basses fréquence [0 500Hz] . Dans ce modèle la bande de roulement est modélisée par une poutre de Bernoulli, les propriétés élastiques du flanc et de la jante sont modélisées par des ressorts[START_REF] Meftah | Efficient computation of tire/road contact forces using an ARMA model of the Green function[END_REF].On suppose que la ligne de contact contient cinq points. La matrice des fonctions de Green est alors de taille 5x5, les figures 1 et 2 montrent respectivement les déplacements et les forces de contact de ces cinq points de contact.

		8	x 10 -3
		7
		6
	Displacement [m]	1 2 3 4 5
		0
	-1
	-2	0	0.01	0.02	0.03	0.04 Time [s] 0.05 0.06	0.07	0.08	0.09	0.1
	FIG. 1 -Déplacement des différents points de contact


A Gütlein

C Ciemniak

F Von Feilitzsch

N Haag

M Hofmann

C Isaila

T Lachenmaier

J.-C Lanfranchi

L Oberauer

S Pfister

W Potzel

S Roth

M Von Sivers

R Strauss

A Zöller

Solar and Atmospheric Neutrinos: Background Sources for the Direct Dark Matter Searches

Keywords: Coherent neutrino nucleus scattering, direct dark matter search, solar neutrinos, atmospheric neutrinos

Direct dark matter search

There are many hints for the existence of non-baryonic cold dark matter [1], [2], [3], [4]. Until now, however, it has not been possible to directly detect dark matter. At present there are several experiments searching for dark matter particles, in particular for WIMPs motivated by supersymmetric theories [1], [4]. Due to the very low expected count rates these experiments need a highly effective background suppression. Most background events originate from the interaction of incident particles (e.g. gammas, betas) with the electrons of the target material (electron recoils), while WIMP interactions are expected to induce nuclear recoils in the target material. Most of the present detectors for direct dark matter searches can distinguish between electron recoils and nuclear recoils, see, e.g., ArDM [5], CDMS [6], CRESST [7], EDELWEISS [8], WARP [9], XENON [10], XENON100 [START_REF] Aprile | XENON100 collaboration[END_REF]. In this way, such background events can be discriminated. Neutrinos, however, could mimic a WIMP event due to coherent neutrino nucleus scattering.

The present best limit for the WIMP-nucleon cross section is ∼ 3•10 -8 pb =3 • 10 -44 cm 2 for a WIMP mass in the region of 50 -70 GeV [6], [START_REF] Aprile | XENON100 collaboration[END_REF]. In [6] this sensitivity was achieved with a raw exposure of 612 kg-days (germanium) showing 2 events which cannot be excluded as signal but which could also be background with 23% probability [6]. Thus, for the next-generation ton-scale experiments, in addition to the already challenging background sources like neutrons also solar neutrinos could become a dangerous background with a rate of several events per ton-year.

Coherent neutrino nucleus scattering as a background for dark matter searches was first suggested in [START_REF] Cabrera | [END_REF] and [13], and has been calculated in [14], [15], [16]. The calculations in the present paper extend beyond previously published work by considering calcium tungstate (CaWO 4 ) and sodium iodide (NaI) as target materials and include all parts of the solar neutrino spectrum (see sections 3 and 4). Some of our results are in disagreement with those reported in [15] and [16].

Coherent Neutrino Nucleus Scattering (CNNS)

CNNS is a neutral-current interaction and thus independent of the neutrino flavour. A neutrino scatters off a nucleus by exchanging a virtual Z 0 boson. For low transferred momenta, the wavelength of the Z 0 boson is in the same order of magnitude as the radius of the nucleus. Thus, the neutrino scatters off all nucleons coherently leading to an enhanced cross section [17]:

dσ d cos θ = G 2 F 8π Z 4 sin 2 θ W -1 + N 2 E 2 ν (1 + cos θ) |f (q)| 2 (1)
dσ dE rec = G 2 F 8π Z 4 sin 2 θ W -1 + N 2 M 2 - E rec M E 2 ν |f (q)| 2 (2)
σ tot = G 2 F 4π Z 4 sin 2 θ W -1 + N 2 E 2 ν |f (q)| 2 (3)
where G F is the Fermi constant, Z the proton number of the target nucleus, N the neutron number of the target nucleus, θ W the Weinberg angle, E ν the neutrino energy, θ the scattering angle in the laboratory frame, M the mass of the target nucleus and E rec its recoil energy. For larger transferred momenta q = √ 2M E rec the Z 0 wavelength is smaller than the target nucleus. Thus, the cross section is modified by a form factor |f (q)| 2 . For this work we used the Helm form factor [18]:

|f (q)| 2 = 3 j 1 (qR 0 ) qR 0 2 e -q 2 s 2 (4)
where j 1 is the spherical Bessel function of the first kind and order 1, R 0 = 1.14A 1 3 fm [19] the effective radius of the target nucleus, and s = 0.9 fm [19] the skin thickness of the nucleus. For the calculation of this form factor the distribution of the charge is assumed to be a convolution of a constant inner part and a gaussian distribution modeling the skin of the nucleus. We used the Helm form factor to compare our results with those of [14].

Since sin 2 θ W = 0.23 [20], the total cross section σ tot is approximately

σ tot ≈ G 2 F 4π N 2 E 2 ν |f (q)| 2 = 4.2 • 10 -45 N 2 E ν 1 MeV 2 |f (q)| 2 cm 2 .
(5)

Count rate calculation

The spectra for solar [21], [22] and atmospheric neutrinos [23] are depicted in figure 1. The count rate for neutrinos is given by where N t is the total number of target nuclei, Φ(E ν ) the neutrino flux, and

R = N t ∞ 0 dE ν Φ(E ν ) 2E 2 ν M 0 dE rec dσ(E ν , E rec ) dE rec (6)
2E 2 ν
M the maximal recoil energy (i.e., back scattering of the neutrino at θ = π).

For the simple case of monoenergetic solar pep neutrinos the flux is given by [21]:

Φ(E ν ) pep = δ(E ν -1.442 MeV) • 1.41 • 10 8 1 cm 2 sMeV . (7)
In the case of monoenergetic pep neutrinos the second integral in equation ( 6) can be calculated leading to the expression for the total cross section in equation (3). With the approximation (5) for the total cross section, the count rate for pep neutrinos is

R pep = N t • 1.41 • 10 8 1 cm 2 s • 4.2 • 10 -45 N 2 • 1.442 MeV 1MeV 2 cm 2 , ( 8
)
where |f (q)| 2 ≈ 1 in good approximation (the corrections are a few percent only). For a germanium target, for example, pep neutrinos produce a count rate of R pep ∼ 539 events per ton-year (for zero energy threshold, see below).

For the calculation of count rates and recoil spectra of other parts of the solar neutrino spectrum, the two integrals in equation ( 6) can be swapped.

Then, the integration limits have to be changed too (the lower limit is determined by the minimum neutrino energy E ν leading to the recoil energy E rec ):

R = N t ∞ 0 dE rec ∞ √ ErecM 2 dE ν Φ(E ν ) dσ(E ν , E rec ) dE rec . (9)
The recoil spectrum dR(Erec) dErec is given by the second integral in equation ( 9)

dR (E rec ) dE rec = N t ∞ √ ErecM 2 dE ν Φ(E ν ) dσ(E ν , E rec ) dE rec . ( 10
)
In order to evaluate the number of events due to CNNS in a detector, the number of events above a given energy threshold E th is calculated. It is important to state what energy threshold means in this context. The energy threshold which is applied in the following calculations is the minimal nuclear recoil energy used for the dark matter search. The count rate above an energy threshold E th is given by:

R th = ∞ E th dE rec dR (E rec ) dE rec = N t ∞ E th dE rec ∞ √ ErecM 2 dE ν Φ(E ν ) dσ(E ν , E rec ) dE rec . (11)

Count rates for solar neutrinos

The solar neutrino fluxes are taken from [21] and the spectra (see figure 1) of the different solar neutrinos from [22]. Figure 2 depicts the recoil spectra of solar neutrinos for different target materials. In figure 3, the integrated count rates above the recoil-energy threshold (horizontal axis) are plotted. In table 1 the count rates for various target materials and different energy thresholds are summarized. These count rates are calculated per ton of target material and year of exposure for an energy window between the minimum energy thresholds as given in table 1 and a maximum energy threshold of 100 keV.

For example, the current energy threshold for xenon is 4.5 keV [10] leading to 1.16 • 10 -3 events per ton-year. For germanium and calcium tungstate (CaWO 4 ) the energy threshold used is ∼ 10 keV [6], [7] producing no events (∼ 10 -8 ) in germanium and ∼ 1 event per ton-year in calcium tungstate mainly due to oxygen recoils. All these events are only due to solar 8 B and hep neutrinos. The other parts of the solar neutrino spectrum become important only for energy thresholds below 1 keV. The results of table 1 are in good agreement with the results of [14] for germanium, xenon and argon and they are also comparable to the simple approximation applied to the pep neutrino count rate in equation (8). However, our results are lower than those obtained in [15] by a factor of ∼ 3 and are higher by a factor of up to ∼ 10 3 than those of [16]. Since the changes due to the form factor |f (q)| 2 (see equation ( 4)) are only a few percent, the disagreement with the results of [15] and [16] must have a different origin. In the present paper we include also results for calcium tungstate and sodium iodide.

Discussion

Solar neutrinos

In table 2 the maximum nuclear recoil energies due to solar neutrinos for different target nuclei are summarized. For energy thresholds above these maximum recoil energies, solar neutrinos are not a background source for direct dark matter searches.

Calcium tungstate (CaWO 4 ) (CRESST)

Figure 4 shows the recoil spectra of solar neutrinos for the different nuclei of calcium tungstate. The maximum recoil energy of tungsten 4.09 keV (see table 2). Thus, for recoil energies above ∼ 4 keV neutrinos only scatter off calcium or oxygen.

the neutrino) is 4.1 keV for tungsten as target nucleus. Thus, the events above an energy threshold of ∼ 4 keV are only due to scattering off calcium or oxygen. If a calcium tungstate detector is able to distinguish between tungsten recoils (i.e., the incident particle scatters off tungsten as expected for WIMPs) and calcium or oxygen recoils [24], only tungsten recoils can be used for the WIMP search and the solar neutrino events above 4 keV can be rejected. However, without this discrimination solar neutrinos lead to ∼ 1 event per ton-year for the present energy threshold of ∼ 10 keV. Thus, the discrimination between tungsten recoils and recoils of oxygen or calcium is crucial for the use of calcium tungstate in future ton-scale experiments.

In order to perform a background-free measurement the count rates of solar neutrinos should be 0.1 counts per ton-year, so that the probability for no events is 90.5% due to the Poisson distribution. For a given WIMP mass of ∼ 60 GeV a background-free experiment with an exposure of one ton-year can reject all WIMP-nucleon cross sections that predict a count rate higher than 2.3 counts per ton-year. To estimate the number of WIMP events, here and in the following examples we have assumed spin-independent elastic WIMP-nucleon scattering [19] and a WIMP mass of ∼ 60 GeV.

If the discrimination between tungsten recoils and oxygen or calcium recoils can be achieved [24] the energy threshold for calcium tungstate should be above 2.05 keV so that solar neutrinos are leading to ∼ 0.1 tungsten recoil events per ton-year. With this energy threshold a background-free experiment using calcium tungstate as target material can reach a sensitivity of ∼ 1.1 • 10 -11 pb for the WIMP-nucleon cross section with a confidence level of 90%.

Without the discriminetion between tungsten recoils and oxygen or calcium recoils the energy threshold should be 16.31 keV to achieve 0.1 counts per ton-year for solar neutrinos. With this energy threshold the reachable sensitivity for the WIMP-nucleon cross section is ∼ 6.8 • 10 -11 pb.

The same conclusion we reach for zinc tungstate (ZnWO 4 ) which is presently also investigated as target material for direct dark matter searches [25].

Germanium (CDMS, EDELWEISS)

The present energy threshold for dark matter experiments with germanium is ∼ 10 keV [6] which is in agreement with the maximum recoil energy of solar neutrinos (see table 2). For recoil energies below 4.91 keV solar neutrinos are leading to ∼ 0.1 events per ton-year. For this energy threshold a background-free experiment with germanium as target material can achieve a sensitivity of ∼ 1, 7 • 10 -11 pb for the WIMP-nucleon cross section with a confidence level of 90%.

The calculated solar neutrino count rates for germanium are in agreement with [14]. However, these count rates are smaller by a factor of ∼ 3 than those in [15].

Xenon (XENON)

For the present energy threshold (4.5 keV) of dark matter experiments using xenon as target material solar neutrinos would produce about 1.2•10 -3 events per ton-year. Already for the slightly lower energy threshold of 2.89 keV, solar neutrinos are leading to ∼ 0.1 events per ton-year. For this energy threshold experiments using xenon as target material can reach a sensitivity of ∼ 8.6•10 -12 pb for the WIMP-nucleon cross section with a confidence level of 90%.

Our count rates for xenon are in agreement with [14], but our count rates are smaller by a factor of ∼ 3 than those of [15] and higher by a factor of ∼ 1000 than those of [16].

Argon (ArDM, WARP)

For argon as target material a recoil energy threshold of 8.62 keV would lead to ∼ 0.1 solar neutrino events per ton-year. For this energy threshold experiments using argon can reach a sensitivity of ∼ 5.0 • 10 -11 pb for the WIMP-nucleon cross section with a confidence level of 90%. The present recoil energy threshold of these experiments is between 30 keV and 50 keV [5], [9], well above 8.62 keV.

The calculated solar neutrino count rates for argon are in agreement with [14] but our count rates are smaller by a factor of ∼ 3 than those of [15].

Atmospheric neutrinos

As shown in the previous section, the solar-neutrino background will not become a problem, if the energy thresholds applied for the dark-matter search (see section 3) are selected with care. However, atmospheric neutrinos have higher energies leading to higher recoil energies. Figure 5 shows, for comparison, the recoil spectra of atmospheric and solar neutrinos for different target materials. The spectra of atmospheric neutrinos were calculated in an analogous way as for solar neutrinos described in section 3. The flux and energy spectra ( 100 MeV) of atmospheric neutrinos (ν µ , ν µ , ν e , ν e ) were taken from [23]. For neutrino energies above 100 MeV the form factor becomes very small. Thus, for this calculation we took only atmospheric neutrinos with energies below 100 MeV into account.

At energies 100 MeV, the flux of atmospheric neutrinos has an uncertainty of ∼ 20% [14]. The spectra in figure 5 were calculated with the Helm form factor of equation (4). For the present energy thresholds of the dark matter experiments, atmospheric neutrinos lead to ∼ 10 -2 events per tonyear [14], [15]. Our results are in agreement with those of [14]. In [15] only calculations for carbon have been included. They derive a value which is a factor of ∼ 100 smaller. This is due to the much lower number of neutrons in carbon as compared to xenon.

As a consequence, atmospheric neutrinos will be a serious background for WIMP sensitivities better than ∼ 10 -12 pb. This sensitivity can be reached with an exposure of ∼ 10 ton-years assuming the energy thresholds mentoined in sections 5.1.1 -5.1.4. To acquire even better sensitivities the recoil spectra of atmospheric and solar neutrinos would have to be fitted to the measured data. Assuming that for a reasonable fit at least ∼ 100 events due to neutrinos are necessary, an exposure of a few 10 3 ton-years would be required to reach enough statistical data for a spectral fit. To achieve such high statistics is very challenging. Thus, if the next generation of darkmatter experiments fails to find WIMPs, the direct detection of WIMPs might become impossible with the present detection techniques and ideas. The DAMA/LIBRA collaboration has been taking a different -in this respect more favourable -route by looking for a characteristic WIMP signature, an annular modulation [START_REF] Bernabei | DAMA[END_REF] using a much larger target mass than other direct dark matter experiments.

Conclusions

Naturally occurring neutrinos are a dangerous background for direct dark matter searches. To reach sensitivities better than ∼ 10 -10 pb for the WIMPnucleon cross section, the solar neutrino events have to be rejected by cuts with nuclear recoil-energy thresholds of 2 keV for calcium tungstate (if tungsten recoils can be singled out, 16 keV otherwise), 5 keV for germanium, 3 keV for xenon, and 9 keV for argon as target material. Thus, next-generation experiments like EURECA [START_REF] Kraus | [END_REF], SuperCDMS [28], ArDM [5] and XENON1T [29] should not only focus on a reduction of the present energy thresholds but mainly on the production of large absorber masses.

Atmospheric neutrinos become a dangerous background for sensitivities better than ∼ 10 -12 pb for the WIMP-nucleon cross section. Atmospheric neutrinos can not be rejected because it is not possible to select a proper recoil-energy window where their contribution can be neglected without a drastic reduction of the sensitivity for the WIMP search. Thus, for sensitivities better than 10 -12 pb direct dark matter searches would suffer from neutrinos mimicking WIMPs. To reach better sensitivities it would be necessary to fit the recoil spectra of solar and atmospheric neutrinos to the measured data. Such a spectral fit would need high statistics (we assume at least about 100 events, see section 5.2) and thus an exposure of a few 10 3 ton-years. To reach such high exposures with the present techniques in direct dark matter experiments seems too demanding. However, detectors with, e.g., directional sensitivity (see, e.g., [30]) could reject solar neutrinos also for low-energy thresholds. Applying such a technique the sensitivity for the WIMP-nucleon interaction could be increased for the same exposure. A promising approach is the confirmation of an annular modulation of the WIMP intensity due to the motion of the earth relative to the WIMP stream [START_REF] Bernabei | DAMA[END_REF]. The DAMA/LIBRA collaboration observes a modulation with the same phase as the expected modulation of the WIMP intensity. According to our calculations, with an exposure of 1.17 ton-years [START_REF] Bernabei | DAMA[END_REF], there should be ∼ 15 solar neutrino events (if we assume 100% efficiency for the detection of solar neutrinos and a recoil energy between 2 and 6 keV without quenching) within the complete measuring period which corresponds to ∼ 10 -5 events per keV, kg, and day. The amplitude of the measured modulation is ∼ 10 -2 events per keV, kg and day. In addition, and even more important, the phase of the expected modulation of the solar neutrino flux, due to the elliptic orbit of the earth around the sun, is shifted by about half a year compared to the modulation observed by DAMA/LIBRA. Thus, neutrinos can not be the reason for the measured modulation.

Acknowledgment

Figure 1 :

 1 Figure 1: Energy spectra of solar (solid lines) and atmospheric (dashed lines) neutrinos (see section 5.2).

Figure 2 :

 2 Figure 2: Nuclear recoil spectra from solar neutrino scattering for different target materials. For the energy thresholds of the present dark matter experiments only 8 B and hep neutrinos contribute to the count rate. The dashed lines indicate the contribution of the CNO cycle to the neutrino flux. The kinks in the spectra for CaWO 4 and NaI arise from the different contributions of the elements in these materials.

Figure 3 :

 3 Figure 3: Integrated count rates from solar neutrino scattering above recoil-energy threshold for different target materials. For other details, see caption of figure 2.

Figure 4 :

 4 Figure 4: Recoil spectra of solar neutrinos for the different nuclei of calcium tungstate. The maximum recoil energy of tungsten 4.09 keV (see table2). Thus, for recoil energies above ∼ 4 keV neutrinos only scatter off calcium or oxygen.

Figure 5 :

 5 Figure 5: Recoil spectra of atmospheric (ν µ , ν µ , ν e , ν e ) and solar neutrinos for different target materials. The feature in the recoil spectrum of atmospheric neutrinos for CaWO 4 at a recoil energy of ∼ 50 keV is due to the form factor (see equation (4)).

  This work was supported by funds of the Deutsche Forschungsgemeinschaft DFG (Transregio 27: Neutrinos and Beyond), the Munich Cluster of Excellence (Origin and Structure of the Universe), and the Maier-Leibnitz-Laboratorium (Garching).

Table 1 :

 1 Nuclear recoil count rates from solar neutrino scattering per ton-year for various target materials and different energy thresholds. The count rates for the different nuclei of calcium tungstate are given separately.

	Material 0 keV 1 keV	2 keV	3 keV	5 keV 10 keV
	Ar	8.0•10 3 1.3•10 2	6.5•10 1	3.2•10 1	6.2 1.5•10 -2
	Ge	1.6•10 4 1.4•10 2	3.6•10 1	7.2 7.9•10 -2 2.5•10 -8
	NaI	2.4•10 4 9.6•10 1	1.4•10 1	6.3	2.7 2.4•10 -1
	Xe	3.0•10 4 9.5•10 1	4.7 6.2•10 -2 9.0•10 -5	-
	ZnWO 4	2.7•10 4 6.6•10 1	1.8•10 1	9.7	4.4	1.0
	CaWO 4	2.8•10 4 5.8•10 1	1.8•10 1	1.2•10 1	5.4	1.1
	O	5.8•10 2 1.4•10 1	1.1•10 1	8.1	4.6	1.1
	Ca	9.0•10 2 1.4•10 1	7.2	3.5 6.9•10 -1 1.7•10 -3
	W	2.6•10 4 2.9•10 1 1.4•10 -1 2.3•10 -3	-	-

  The maximum recoil energy (i.e. backscattering of

	Nucleus Max. recoil energy
	O	47.1 keV
	Na	37.6 keV
	Ar	18.8 keV
	Ca	18.8 keV
	Zn	11.6 keV
	Ge	10.3 keV
	I	5.9 keV
	Xe	5.7 keV
	W	4.1 keV

Table 2 :

 2

Maximum nuclear recoil energies for different target materials for an energy of solar neutrinos of 18.79 MeV (i.e. the maximum hep neutrino energy).