Aude Bolopion

Barthélemy Cagneau

Stephane Redon

Stéphane Régnier

Haptic molecular simulation based on force control

In this paper, force control is proposed to connect a molecular simulator to a haptic device. Most of the works dealing with this kind of simulators use position control to manipulate the molecules, with major stability concerns. These two control modes are compared in terms of adequacy with the molecular simulator. Stability with respect to the scaling coefficients introduced to connect the macro and the nanoworlds is also considered. The theoretical results and the experiments carried out confirm that position control is sensitive to the gain tuning. Force control enables to get stable force feedback for varying gains, and is thus a promising coupling to perform manipulations on complex molecular systems.

I. INTRODUCTION

Molecular simulation is one of the fields in which haptic devices greatly improves the user knowledge as well as its ability to perform complex operations, such as prototyping bio-nanorobots [START_REF] Hamdi | Prototyping bio-nanorobots using molecular dynamics simulation and virtual reality[END_REF]. The works can be classified in three main topics: evaluation of the benefits of using haptic, how to interact with large virtual environments using haptic devices, and the choice of the coupling.

The evaluation of haptic for molecular simulators has demonstrated that it helps operators to understand nanoscale phenomena. Its use is recommended for educational purposes [START_REF] Persson | Designing and evaluating a haptic system for biomolecular education[END_REF], [START_REF] Murayama | An educational environment for chemical contents with haptic interaction[END_REF], but is not limited to academic courses. Chemists and biologists also benefit from haptic to find specific locations in complex molecular systems, such as docking sites [START_REF] Brooks | Project GROPE -Haptic displays for scientific visualization[END_REF].

The benefits of haptic feedback depends on the coupling used. In particular, how to reach the entire virtual environment using a haptic device with a limited workspace is a key issue to get an interactive system. Several techniques are proposed, from the concept of clutching (freezing the displacement of the virtual object while enabling the user to modify the position of the haptic handle), to the Bubble technique (combination of position and rate control) [START_REF] Dominjon | The "Bubble" technique: interacting with large virtual environments using haptic devices with limited workspace[END_REF]. As an alternative the concept of Active Haptic Workspace is considered [START_REF] Subasi | A new haptic interaction and visualization approach for rigid molecular docking in virtual environments[END_REF]. In all cases, when the user is close to the point of interest, a classical position coupling scheme is used to control the object and send the interaction forces to the user. Thus, even if complex techniques are used, the choice of the haptic coupling used when dealing with a precise area of the workspace remains unchanged.

In most of the works dealing with molecular systems A. Bolopion stephane.redon@inria.fr position control is used (the user sets the position of the molecule, and feels the interaction force through the haptic interface). However, stability is difficult to ensure due to long computation times, scaling factors used to link the macro and nano worlds, and the high variation of the forces.

Stability is usually guaranteed at the expense of the fidelity of the force feedback. Either the accuracy of the molecular interactions computed is decreased by using simpler models [START_REF] Lee | Smoothing haptic interaction using molecular force calculations[END_REF], or the damping added to the coupling deteriorates the transparency [START_REF] Daunay | Stable six degrees of freedom haptic feedback for flexible ligand-protein docking[END_REF]. This manipulation mode is also not adapted to molecular dynamics simulators, since setting the position of the molecule leads to potentially physically unacceptable positions before shifting to the next step of the simulation, as we will see in this paper. It seems that none of the proposed systems (simulator and haptic coupling) gives satisfying results in terms of the trade off between the accuracy of the computed interaction forces, the stability and the transparency of the force feedback.

In this paper, we present a molecular simulator SAMSON, which enables to simulate complex systems in a few hundreds of milliseconds, and takes into account the flexibility of all the considered molecules. To avoid stability concerns due to the commonly used position control, we propose to use force control. Force control is well known for macro scale teleoperations, but is not used for molecular simulations. The molecular dynamics equation solved by the simulator, on which an additional force can be included, enables to implement this control mode. Some of the potential applications of our system (conception of new molecules, analysis of molecular properties, ...) are described in [START_REF] Bolopion | Haptic feedback for molecular simulation[END_REF]. We now present a comparison of position and force control. A detailed description of the molecular simulator and the control schemes, as well as how they are connected, is made to compare the adequacy of the two control modes to molecular dynamics simulators. The stability sensitivity of the coupling schemes, measured as the variation of the magnitude of the roots of the control schemes' transfer functions with respect to given parameters, is considered for the specific application of molecular interactions.

The rest of the paper is organized as follows. In Section II, an overview of the simulator and the coupling is given. Based on the available inputs and outputs, the two control modes (position and force control) are detailed and compared in Sections III and IV in terms of adequacy with the simulator. The stability sensitivity is analyzed in Section V. Experimental results are given in Section VI.

II. COUPLING BETWEEN THE HAPTIC INTERFACE AND THE SIMULATION SOFTWARE

A. Molecular simulator SAMSON SAMSON (System for Adaptative Modeling and Simulation Of Nano-Objects) is based on a quasi-statics method to simulate the motion of the molecular system (Figure 1) [START_REF] Rossi | Adaptive torsion-angle quasi-statics: a general simulation method with applications to protein structure analysis and design[END_REF], [START_REF] Grudinin | Practical modeling of molecular systems with symmetries[END_REF]. The force field used is derived from a well known molecular mechanics force field, CHARMM [START_REF] Brooks | CHARMM: A program for macromolecular energy, minimization, and dynamics calculations[END_REF], which models interactions through van der Waals, electrostatic and dihedral contributions. Our simulator has several specific properties:

• as the simulator solves the equation of motion, the value of the force applied to each atom is directly known (whereas many simulators compute the energy of the system [START_REF] Daunay | Stable six degrees of freedom haptic feedback for flexible ligand-protein docking[END_REF])

• both external and internal efforts are computed • the flexibility of the molecules is simulated, and systems with thousands of degrees of freedom can be handled

• it enables fast computation based on a tree representation of the molecules. The simulation period is only a few hundreds of milliseconds for complex systems. This simulator enables to use two different manipulation modes based on the equation of motion solved for each rigid body. Using spatial notations, this equation is [START_REF] Featherstone | Robot dynamics algorithms[END_REF]:

Ia = F ext + F int -v × Iv (1)
where a (resp. v) is the spatial acceleration (resp. velocity) of the molecule and I is its inertia. F ext are the forces applied by other molecules and F int is the sum of the internal forces. Considering (1), the user can either manipulate the molecule by setting the position of its center of mass, or by applying a force F i on it. The first mode is thus the position control in which the instantaneous position taken into account while solving (1) is the position that the user sets as the desired position for the molecule. For force control, the user applies an additional force to the molecule to be manipulated, so that the equation of motion is modified to:

Ia = F i + F ext + F int -v × Iv (2)
Equations ( 1) and ( 2) are integrated with the molecular time T m . This parameter is introduced since the molecular dynamics do not allow to use the same time basis as the one we use. A single simulation step (moving the molecule, reconfiguring the system, sending the interaction forces to the user) takes around T s = 100ms (simulation loop period), and this corresponds to a period of around T m = 2fs for the evolution of the molecular system [START_REF] Stone | A system for interactive molecular dynamics simulation[END_REF].

The force fed back to the user is F m :

F m = F ext + F int (3)
Thus the user can choose between feeling the internal forces by manipulating a single atom of a molecule, or the total force applied to a given molecule by its environment (while dealing with a whole molecule the internal forces cancel out).

Since no modification on our system are necessary the user can easily switch from one mode to the other.

B. Haptic coupling schemes

The quality of the force rendering depends on the control scheme used. The choice of the coupling is thus of utmost importance. When connecting a haptic device to a molecular simulator the delays due to long computations of complex interactions must be taken into account while analyzing the stability. The transfer functions are thus represented using the discrete time variable z. Since the variables coming from the simulator are delayed by one sampling period (period of the simulator), an explicit one step delay is introduced between the simulation and the coupling. This is modeled by the delay operator z -1 . The discrete time transfer function H(z) (which input, resp. output, is the force applied to the haptic handle, resp. its position) is computed from the continuous time domain using the Z-transform function [START_REF] Ogata | Discrete-time control systems -second edition[END_REF], and taking into account the effects of the sampler and zero order hold:

H(z) = 1 B h T s (z -δ) -(1 -δ)(z -1) M h B h (z -1)(z -δ) (4)
where δ = e

- B h Ts M h
and T s is the simulation's sampling period. M h and B h are respectively the master's inertia and viscosity, as given in [START_REF] Diolaiti | Stability of haptic rendering: Discretization, quantization, time delay, and coulomb effects[END_REF].

Two scaling factors are used to connect the haptic interface and the molecule. The force (resp. position) scaling factor is denoted by A f (resp. A d ).

Based on the two manipulation modes described above two control schemes with different inputs and outputs can be used, and are described in next sections. Note that, although we have used SAMSON simulator in our implementation, the control schemes we will present may be connected with any simulator which provides access to the required variables (positions, forces, etc.). However, the accuracy of the obtained results (reconfiguration of the molecule, force computed, ...) depend on the accuracy of the simulator.

III. POSITION CONTROL

A. Haptic coupling

The classical impedance display or Direct Force Feedback (DFF) coupling, is based on position control. It is the simplest structure to connect the haptic device to the simulation

F m F h x i F op H(z) A -1 d A f z -1 + - x h Simulation Coupling
Master device Delay Fig. 2. DFF control scheme using the scaled position of the haptic handle x i as the input of the molecular dynamics simulator (Figure 2). Only two scaling factors are necessary to control the molecule.

B. Details of the DFF algorithm

To compare the adequacy of position and force control with the molecular simulator it is necessary to give some details about how the haptic coupling schemes and the simulation are connected. This is done in this section for the DFF coupling, and in Section IV-B for force control. Figure 3 gives more details about the notations used.

F l e = F l int + F l ext R 1 m ...R j m ...R n m R 1 e ...R l e ...R q e
Center of mass position xm

Rigid bodies of the manipulated molecule

Rigid bodies of the environment

The desired position of the center of mass xi is set using the haptic device Fig. 3. Notations and principle of position control

As mentioned in [START_REF] Bolopion | Haptic feedback for molecular simulation[END_REF], the calculations of positions and forces performed by the simulation software are made on each rigid body. However, for ease of manipulation, the user manipulates the molecule as a whole, controlling the position of the center of mass x i . As explained in Section II, the instantaneous position of the molecule considered when solving the equation of motion (1) is the desired position of the rigid body. Internal and external forces are computed and used to update the positions of the molecules in the environment and their internal degrees of freedom (to simulate their flexibility). This reconfiguration step involves a modification of the rigid bodies positions, and thus the actual position of the center of mass of the manipulated molecule x m might be different from the one set by the user using the haptic interface x i (see Algorithm 1 1 ).

Algorithm 1 Position control algorithm 1) Compute the desired position of the center of mass x i (k + 1) set by the operator using the control scheme (position of the haptic handle scaled by A d ). 2) Set the position of the rigid bodies to x j i (k + 1), computed as the relative positions with respect to x i .

3) Update all interatomic forces F 1 m (k +1), ..., F n m (k +1) (resp. F 1 e (k + 1), ..., F q e (k + 1)) applied to the rigid bodies which are inside (resp. outside) the manipulated molecule using the molecular simulator. x l e ) of the rigid bodies which are inside (resp. outside) the manipulated molecule using the quasi-statics simulator. Note that in the following equations x j i (k + 1) (resp. F j m (k + 1) and F l e (k + 1)) are known since they have previously been computed: Step 2 (resp. 3) of the algorithm.

x j m (k + 1) = x j i (k + 1) + T 2 m m j F j m (k + 1)
x l e (k + 1) = x l e (k) +

T 2 m m l F l e (k + 1)
The DFF coupling has one main drawback: the position the user sets using the haptic interface might be physically unacceptable (i.e., with high potential energy due to atomic clashes). Even if the next step of the simulation corrects the position, this can lead to instabilities, and the forces sent to the user may misrepresent the correct molecular interactions.

IV. FORCE CONTROL

A. Haptic coupling

In the Force-Force (FF) coupling scheme (Figure 4), the input of the simulation is the force applied by the user to the molecule. In addition to the scaling factors A f and

A d F m x m F h x h F i F op H(z) z -1 z -1 A -1

B. Details of the FF algorithm

Contrary to the DFF algorithm, the simulator's input is the force F i that the user wants to apply to the entire molecule. This force is set by controlling the distance between the haptic handle and the molecule positions. The corresponding forces applied on each rigid body F j i of the manipulated molecule as well as the environment forces are used to update the position of all the rigid bodies (inside and outside the manipulated molecule) using (2) (see Algorithm 2). Algorithm 2 Force control algorithm 1) Compute the force applied by the operator to the molecule using the control scheme:

F i (k + 1) = k i x h (k+1) A d -x m (k)
2) Compute the force F j i (k + 1) = F i (k + 1)/n, applied to all the j th rigid bodies of the manipulated molecule.

3) Update the positions x j m (resp. x l e ) of the rigid bodies which are inside (resp. outside) the manipulated molecule using the quasi-statics simulator. Note that in the following equations, F j i (k + 1) is known since it has previously been computed in Step 2:

x j m (k + 1) = x j m (k) + T 2 m m j (F j m (k) + F j i (k + 1))
x l e (k + 1) = x l e (k) +

T 2 m m l F l e (k)
4) Update all interatomic forces F 1 m (k +1), ..., F n m (k +1) (resp. F 1 e (k + 1), ..., F q e (k + 1)) applied to the rigid bodies which are inside (resp. outside) the manipulated molecule using the molecular simulator. Using force control enables to integrate the user's input while solving the mechanical equation [START_REF] Persson | Designing and evaluating a haptic system for biomolecular education[END_REF]. Contrary to the DFF manipulation mode, the position of the molecule is always physically acceptable, which avoids instabilities due to atomic clashes. Thus, this manipulation mode is more adapted to the molecular simulator than position control.

V. COMPARISON OF THE CONTROL SCHEMES STABILITY

In this section, the sensitivity of the coupling schemes' stability is compared. Since passivity is conservative (it ensures stability for any operator and any environment as long as they are themselves passive) [START_REF] Khalil | Nonlinear systems[END_REF], it is not used in this work. To compare the haptic coupling schemes, the action of the operator is not taken into account, and a simplified model of the environment is considered (an equivalent spring constant k e ). This is obviously not enough to model the complex interactions between molecules, but comparisons on the control schemes' stability can be made. The sensitivity analyzed is the variation of the magnitude of the roots of the transfer functions' characteristic equations with respect to given parameters. We choose to analyze here the influence of the scaling factors.

A. Characteristic equations

The discrete time transfer function of the DFF control scheme is derived from the coupling presented in Figure 2:

Z x h (s) F op (s) = X h (z) F op (z) (5) = T s z-δ z-1 -(1 -δ) M h B h B h (z -δ) + T s z-δ z-1 -(1 -δ) M h B h A f A d k e z -1 (6)
where Z is the Z-transform function, and δ is the parameter defined in (4). The corresponding characteristic equation is:

a 3 z 3 + a 2 z 2 + a 1 z + a 0 = 0 (7)
where:

a 3 = B h , a 2 = -B h (1 + δ) a 1 = δB h + T s -(1 -δ) M h B h A f A d k e a 0 = -δT s + (1 -δ) M h B h A f A d k e
The transfer function of the FF coupling is computed according to Figure 4:

X h (z) F op (z) = H(z) 1 + H(z)G(z) A f A d (8)
where:

G(z) = Fm(z) X h (z) = kikeT 2 m mz(z-1)+kiT 2 m +keT 2 m
and m is the total mass of the molecule. The characteristic equation is:

b 4 z 4 + b 3 z 3 + b 2 z 2 + b 1 z + b 0 = 0 (9)
where:

b 4 = mB h , b 3 = -(2 + δ)mB h b 2 = (k i + k e )T 2 m B h + mB h (1 + 2δ) b 1 = -(1 + δ)(k i + k e )T 2 m B h -δmB h + T s -(1 -δ) M h B h k i k e T 2 m A f A d b 0 = B h δ(k i + k e )T 2 m + -T s δ + (1 -δ) M h B h k i k e T 2 m A f A d
Analytical stability criteria can be obtained for both of these control schemes using the Jury criterion [START_REF] Ogata | Discrete-time control systems -second edition[END_REF], as it is done in [START_REF] Gil | Stability boundary for haptic rendering: Influence of damping and delay[END_REF] for a control scheme similar to the DFF one. However, the relations derived are complex and the influence of the control scheme parameters on stability is not highlighted. Thus, in the next section a numerical analysis is performed.

B. Numerical comparison of stability

For both the DFF and FF control schemes two parameters can be tuned to ensure stability: the force and the displacement scaling factors. The force scaling factor A f is chosen so that the user can distinctly feel the interaction forces between and/or within the molecules. Thus, only the displacement coefficient can be modified. A numerical analysis of the influence of A d on each control scheme is performed: the variations of the highest magnitude of the roots (absolute number of the largest pole) of the characteristic equations ( 7) and ( 9 The variations of the magnitude of the roots are of particular interest. For both of the control schemes, small values of A d and/or stiff environments may lead to an unstable system. The magnitude of the roots of the DFF control scheme quickly increases when A d is decreased, or k e is increased. On the contrary, the roots obtained with the FF control scheme are almost constant. Thus, the DFF control scheme stability is more sensitive than the FF one to the variations of A d for a given environment stiffness k e . This result is confirmed experimentally in the next section.

VI. EXPERIMENTAL COMPARISON OF THE CONTROL

SCHEMES

An experiment involving a classical blocker of potassium channels (KcsA), the tetraethylamonium (TEA) is performed. The docking of the TEA into KcsA is studied in particular in the context of mutagenesis studies to determine the influence of aromatic residues on the external blockade of the TEA. The experiment consists in approaching the TEA from the potassium channel and then moving it away. The complexity of this example (5 molecules involved, for a total of 1428 degrees of freedom) is representative of the molecular systems biologists and chemists consider. This experiment performed for both the DFF and FF control schemes for varying values of A d is presented to support the conclusions of Section V-B. More examples can be found in [START_REF] Bolopion | Haptic feedback for molecular simulation[END_REF].

A. Characteristic parameters

The parameters of this experiment are: m = 1.8 • 10 -25 kg (mass of the TEA : the molar mass is 110g/mol in the CHARMM19 representation in which the 20 hydrogen atoms are not repesented), T s = 138 • 10 -3 s (simulation period, as it varies depending on the number of atoms involved in molecular interactions, the largest one is given), T m = 2 • 10 -15 s (physical integration period), d = 8 Å (cutoff distance used to limit the computational time due to the molecular interactions that can be felt at infinite distance), A f = 0.3.10 9 (force amplification), k i = 24N.m -1 (for the FF coupling). The 124 dof the most important for the accuracy of the molecular interaction computation are activated by our adaptive simulator [START_REF] Rossi | Adaptive torsion-angle quasi-statics: a general simulation method with applications to protein structure analysis and design[END_REF]. The position scaling factor depends on the user needs (precise or fast manipulation). To test the sensitivity of the control schemes' stability with respect to this parameter, three experiments are performed, with A d = 0.50 • 10 9 (Experiment A), A d = 0.25 • 10 9 (Experiment B) and A d = 0.125 • 10 9 (Experiment C). The haptic interface used is a Phantom Omni device from Sensable (Figure 1), with an inertia M h = 0.072kg and a viscosity B h = 0.005N.s.m -1 [START_REF] Diolaiti | Stability of haptic rendering: Discretization, quantization, time delay, and coulomb effects[END_REF].

B. Experimental results

For each experiment, the tetraethylamonium is approached and retracted from the potassium channel. The movement is directed along one axis. The forces returned to the operator are depicted in Figures 6, 7 and 8 2 .

For both of the coupling schemes, the user distinctly feels the interaction forces since they are rendered with a force in the order of the newton. As predicted in Section V-B, as A d is decreased, the DFF control scheme becomes unstable. In experiment C, atom clashes produce high magnitude forces and instabilities. On the contrary, the FF coupling scheme remains stable. The maximum force varies depending on the exact displacement performed by the user. The variations on the force profile of the DFF control scheme (see circles in Figures 6 and7) are due to the manipulation mode. As the user freezes the movement of the molecule to reach a more comfortable position with the haptic interface the rest of the system applies forces to the molecules, which moves according to (1). This is not the case for the FF control scheme since the user does apply a force on the molecule even while clutching its position.

Both the conclusions of Sections III-B and IV-B, the numerical comparison of stability (Section V-B), and the experiments highlight that force control is well suited for molecular manipulation and enables a stable coupling. More experiments were made on different molecules (HIV virus, water molecules, ...), and similar conclusions were observed. The results are not displayed here since they do not provide additional information. More examples can be found in [START_REF] Bolopion | Haptic feedback for molecular simulation[END_REF].

Adding damping to the DFF coupling would improve the stability, but transparency would be deteriorated. Since molecular interactions involve complex variations of the forces and efforts of small magnitude (in particular for attractive interactions), transparency is a major concern to provide to the user a high fidelity force feedback. It can be noted that for both the DFF and FF control schemes the user feels the interaction forces scaled by the force coefficient. Thus, high fidelity force feedback is ensured.

VII. CONCLUSION

Position control, which is commonly used to connect a haptic device to a molecular simulator has two major drawbacks: due to the equation of motion solved, the position set by the user can be physically unacceptable before being updated, which can produce instabilities. The second drawback is that the stability of the control scheme is very sensitive to the scaling factors used, which can limit its potential applications. On the contrary, force control offers both a good adequacy with the simulation software and the stability is not sensitive to the chosen gain parameters. Since the molecular simulator used is able to simulate the flexibility of complex molecules composed of thousands of degrees of freedom, our system improves the user's understanding of molecular interactions.

Besides the two criterion used in this paper to compare the control schemes, the user feeling should also be taken into account, as well as the intuitiveness of the manipulation mode. User based tests should thus be carried out in future works to perform a detailed comparison. Solutions to automatically determine the most appropriated scaling factors and gains for the control schemes should also be considered. Adaptive gains might be a solution to improve the DFF control scheme performances.

Fig. 1 .

 1 Fig. 1. Manipulation of a tetraethylamonium (TEA) around a potassium channel (KcsA), 1428 degrees of freedom.

  Position of the l th rigid body

4 )

 4 Compute the sum F m (k + 1) = n j=1 F j m (k + 1) of the forces applied to the molecule, and send it to the user. 5) Update the positions x j m (resp.

Fig. 4 .

 4 Fig. 4. FF control scheme already considered on the DFF coupling, a proportional gain k i is introduced. This gain adjusts the magnitude of the force applied to molecule, which is based on the difference between the molecule and the haptic handle positions.

5 )

 5 Compute the sum F m (k + 1) = n j=1 F j m (k + 1) of the forces applied to the molecule, and send it to the operator.

Fig. 5 .

 5 Fig.5. Variations of the magnitude of the roots against A d , for different environment stiffness ke. For clarity, only the highest magnitude is plotted for each control scheme given the stiffness.

Fig. 6 .Fig. 7 .

 67 Fig. 6. Haptic force during Experiment A (124 activated dof )

Fig. 8 .

 8 Fig. 8. Haptic force during Experiment C (124 activated dof )

  and S. Régnier are with Institut des Systèmes Intelligents et de Robotique, UPMC -CNRS UMR 7222, 4 Place Jussieu, 75005 Paris, France. {bolopion, regnier}@isir.upmc.fr B. Cagneau is with Laboratoire d'Ingénierie des Systèmes de Versailles, EA 4048, Université de Versailles Saint Quentin, 45 avenue des Etats Unis, 78035 Versailles, France. barthelemy.cagneau@uvsq.fr

	S.	Redon	is	with	NANO-D	-	INRIA	Rhône-Alpes,
	655	avenue	de	l'Europe,	38334	Saint-Ismier,	France.

For the sake of clarity only the translations are considered in this algorithm. Rotations are implemented using the same methodology. This remark is also valid for Algorithm

To avoid damages on the haptic interface forces are limited to

3N

VIII. ACKNOWLEDGMENT This work was supported by the French National Agency of Research, through the PACMAN project.


Aude Bolopion

Barthélemy Cagneau

Stéphane Redon

Stéphane Régnier

Haptic feedback for molecular simulation

In this paper, a new tool dedicated to the analysis and the conception of molecules is presented. It is composed of an adaptive simulation software and a haptic device used to interact with molecules while feeling either the forces applied by the environment or the internal forces. The adaptive articulated body algorithm allows fast simulations of complex flexible molecules. To handle the coupling with the force feedback device, two different control schemes designed for nanoscale applications and providing high transparency rendering are proposed and compared. The system we propose is highly flexible since either a single rigid body or the entire molecule can be manipulated via the haptic device. The user can choose between setting a desired position/orientation of the molecule, or apply forces/torques to manipulate it. It allows the operator to control each stage of the design process of new molecular structures. The validity of this tool is demonstrated through examples of haptic interaction between the HIV protease and its inhibitors, and unfolding one of these drugs.

I. INTRODUCTION

Molecular mechanics simulators are efficient tools to study molecular interactions and docking, with applications in many fields, including pharmacology, cosmetics, nanosciences and surface physics. In particular, interactive molecular simulations may be a significant help to understand molecular interactions. For these simulations, a user can control the motion of one or more atoms to guide the evolution of the molecular system.

Most of the work deal with molecular docking. In [START_REF] Ouh-Young | Using a manipulator for force display in molecular docking[END_REF], the authors explain how the user's intervention may facilitate the docking, by controlling the position of the ligand and feeling the interaction forces. Other studies have then been carried out. The simulator and the bilateral couplings proposed differ in many ways. For example, in the simulation, the force field can be computed thanks to the derivation of the potential energy [START_REF] Daunay | Energy-field reconstruction for haptic-based molecular docking using energy minimization processes[END_REF]. Even if some simulators take into account the flexibility of the ligand and the receptor [START_REF] Lai-Yuen | Energy-field optimization and hapticbased molecular docking and assembly search system for computeraided molecular design (camd)[END_REF], most of them do not [START_REF] Teodoro | Molecular docking: a problem with thousands of degrees of freedom[END_REF]. Concerning the coupling, wave variables are also used to guarantee the passivity of the system even when time delays increase [START_REF] Daunay | 6 dof haptic feedback for molecular docking using wave variables[END_REF]. However, this method may lead to a lack of transparency.

Besides molecular docking, other emerging fields of research concern molecular interactions. The molecules can be used as a channel to transmit an information between two nanomachines ( [START_REF] Moore | Interfacing with nanomachines through molecular communication[END_REF] and [START_REF] Atakan | An information theoretical approach for molecular communication[END_REF]). Another research field is the A. Bolopion stephane.redon@inria.fr study of molecules in order to develop biomicroelectromechanical system devices (bio-MEMS) ( [START_REF] Rossi | Molecular and supramolecular systems for sensing and actuation[END_REF] and [START_REF] Grove | Packaging actomyosinbased biomolecular motor-driven devices for nanoactuator applications[END_REF]). Tools to deal with specific applications have been developed. However, emerging applications as the conception of new materials and alloys, or molecules with specific properties require more flexibility. Both the internal structure of the molecule and the interactions with its environment must be understood. The manipulation of the molecule must be done in real time, and be intuitive for users with no training. Thus, we propose a system, composed of a molecular simulation software and a haptic interface, able to deal with complex flexible molecules in real time (figure 1). The simulation software is based on adaptive simulation of articulated body, so that the computation time is decreased dramatically. The algorithm used has already been exposed in [START_REF] Morin | A force-feedback algorithm for adaptive articulated-body dynamics simulation[END_REF] and [START_REF] Rossi | Adaptive torsion-angle quasi-statics: a general simulation method with applications to protein structure analysis and design[END_REF]. Both internal forces and interactions between the molecule and its environment are computed. To interact with the simulation, a Phantom Omni haptic device is used to manipulate the molecule, and render forces to the user. To handle the coupling between the haptic device and the simulation, two different control schemes are implemented. The user can manipulate a single rigid body, or an entire molecule. He or she can either set the desired position/orientation, or apply force/torque to control the object. We performed several experiments using this system. Interactions between the HIV protease and its inhibitors are used to compare the two control schemes. The feeling of internal forces is demonstrated by unfolding one of its inhibitors.

This paper is organized as follows. In section II, the molecular simulator is briefly introduced. The two control schemes designed for nanoscale applications are presented in section III. In section IV, examples of manipulations that can be handled by our tool are presented. The advantages and drawbacks of each control scheme are also discussed. Finally, concluding remarks are exposed in section V.

II. MOLECULAR MECHANICS

In this paper, we use an adaptive molecular quasi-statics method to simulate the motion of the molecular system [START_REF] Rossi | Adaptive torsion-angle quasi-statics: a general simulation method with applications to protein structure analysis and design[END_REF]. The adaptive method relies on a recursive, divide-andconquer representation of the molecular system, where any non-rigid group of atoms is considered to be the union of two groups of atoms. Thus, any molecular system is associated to a binary assembly tree, in which the leaf nodes represent the user-defined rigid bodies, the internal nodes represent the sub-assemblies, and the root node represents the complete molecular system (figure 2). The distinctive feature of the adaptive algorithm is that it allows a user to finely trade between precision and computational cost. Precisely, the user may arbitrarily choose the number of degrees of freedom that are simulated at each time step, and the adaptive method automatically determines which degrees of freedom are activated at each time step (a sub-tree of the binary assembly tree), so as to approximate as faithfully as possible the motion which would have been obtained by a non-adaptive simulator. In practice, the adaptive algorithm focuses the computational resources on the joints with the largest accelerations, that is, the most flexible regions. Using this approach, complex molecular systems may be simulated on cheap personal computers while retaining molecular flexibility. In docking applications, the adaptive simulator automatically focuses the motions computations on the docking site.

Besides adaptivity, the simulator follows traditional approaches: at each time step, the interatomic forces are first computed, then the acceleration of the system is determined and used to update its state.

A user may influence the state of the molecular system by either directly controlling the position and orientation of a selected molecule, or by applying a force and torque to a group of atoms. In the latter case, the selected group of atoms may be either a complete molecule (e.g. an inhibitor of an HIV protease model), or a rigid body from a molecule.

The force fed back to the user is the total force applied by the molecular system to the group of atoms selected by the user. These forces directly depend on the force field being used in the simulation to model interatomic interactions. In our current implementation, the force field is derived from a well known molecular mechanics force field, CHARMM [START_REF] Rossi | Adaptive torsion-angle quasi-statics: a general simulation method with applications to protein structure analysis and design[END_REF], which models interactions through van der Waals, electrostatic and dihedral contributions.

• Molecule. Assume the user manipulates an entire molecule M, composed of rigid bodies R 1 , ..., R n . In this case, we determine the molecular forces F 1 , ... F n applied on each rigid body by the other rigid bodies in the environment (i.e. van der Waals and electrostatic contributions from surrounding atoms, including those resulting from other rigid bodies in the selected molecule), and return to the user the sum F = n i=1 F i of these forces. Note that, in this sum, the contributions from rigid bodies within the selected molecule cancel out.

• Rigid Body. When only one rigid body R i from a molecule is selected, we compute the molecular force F i applied to the rigid body by the environment as before, but add to this force the kinematic constraint forces, i.e. the forces which maintain the kinematic constraints between the selected rigid body and its topological neighbors [START_REF] Rossi | Adaptive torsion-angle quasi-statics: a general simulation method with applications to protein structure analysis and design[END_REF]. This flexibility allows us to propose two control schemes to handle the coupling between the simulation and the haptic device, that we detail in the next section.

III. BILATERAL COUPLING

In order for our tool to be used by untrained users, the control schemes must fulfill several requirements. In particular, the manipulation of the molecule via the master arm must be intuitive. The forces felt by the user must reflect the interaction forces computed by the simulation. The coupling has therefore to be transparent, so that even small modifications on the force profile are rendered to the operator.

Regarding the specificities of the software described in the previous section, different bilateral couplings can be investigated, based indistinctly on impedance or admittance control.

A. Direct Force Feedback control

The Direct Force Feedback control scheme (DFF) appears as the most intuitive when the simulation's input is the position of the rigid body being manipulated and forces sent to the user are the ones computed by the simulation. It is represented in figure 3, using the continuous time Laplace variable s. The different parameters are:

F m q s q h F op H(s) A -1 d A f Simulation Coupling Master device
• A d = A d,T 0 0 A d,R
where A d,T and A d,R are respectively used to scale down the positions and the rotations of the haptic handle. A f is used to enhance the forces sent to the user. These coefficients are necessary to adjust the scale difference between macro-and nanoworlds.

• q s = [x s , θ s ], q h = [x h , θ h ] are respectively the position vectors of the rigid body's center of mass and of the haptic handle.

• H(s) is the haptic device transfer function. It is composed of a mass and viscous friction [START_REF] Diolaiti | Stability of haptic rendering: Discretization, quantization, time delay, and coulomb effects[END_REF].

The operator applies a force F op on the master device. The position of the haptic device handle is scaled by A d , and used as the desired position of the rigid body's center of mass. The orientation of the rigid body is computed with respect to the rotations of the haptic device handle. For a given time step, the user moves the rigid body to a desired position q s . The simulator modifies this position using a criterion based on energy minimization and taking into account the influence of other atoms. This new position will be used as the reference position for the next time step. Thanks to this method, all the displacements make sense and represent the molecular dynamics. If the user wants to manipulate multiple rigid bodies as an unique entity (for example to manipulate a flexible molecule at its center of mass), the same displacement and orientation are applied to all the rigid bodies. In that case, the simulation also updates the configuration of the rigid bodies based on internal and environmental interactions, and determines the new sub-tree of rigid bodies as described in section III.

The simulation remains stable if the displacements are small enough for a given time step, so that no atoms interpenetrate each other. The matrix A d is used to adjust the displacements of the user so that the simulation can deal with them. The molecular forces F m are sent to the user after being increased by a factor A f .

B. Force-Force control

The second approach, Force-Force control (FF), consists in using only forces for both manipulation of the rigid body and haptic feedback. This bilateral coupling is represented on figure 4. This control scheme is similar to the DFF as the feedback path is also composed of molecular forces F m amplified by A f . The difference is that the user applies a force F s and a moment M s on the rigid body to manipulate it. These values are computed based on the position of the rigid body's center of mass q m and the position of the haptic handle q h :

F m q m q h F s , M s
F s = k T (x h /A d,T -x m ) M s = k R (θ h /A d,R -θ m ) (1)
During a simulation step, the software takes into account the forces applied by the user on the rigid body, as well as the forces due to the interactions between atoms. These forces are used to derive the new position/orientation of the rigid body. Compared to the DFF, both the forces from the user and from the other atoms are considered simultaneously. The position does not need to be corrected afterwards, contrary to section III-A, and therefore is less subject to oscillations. When the user manipulates a molecule composed of several rigid bodies, the force applied on each of them is computed so that their accelerations are the same. This is done using adapted values of k, and considering the mass of each rigid body. As previously, for flexible molecules, the configuration of these rigid bodies are updated at each time step using a criterion based on energy minimization.

Compared to the DFF formulation, the ergonomy of the system is modified when a force is used as the simulator's input. The user does not control the position but the dynamics of the molecule. Depending on the application, users can choose between these two modes.

IV. SIMULATION RESULTS

Experiments are performed on a 3GHz processor Intel PC, using Windows operating system. The available haptic device is a Phantom Omni haptic interface from Sensable1 . This master arm is a 6 degrees of freedom (dof ) movement and a 3 dof force feedback device. It is controlled using the OpenHaptics development toolkit.

A. Manipulation of molecules

The tool presented offers the possibility to the user to interact in real time during the conception process of new molecular structures. This functionality is of primary importance to design new materials and alloys with specific properties. A first experiment is performed using the HIV protease and its inhibitors (pdb code: 2AZ8) to demonstrate that our tool is adapted to such kind of manipulations. The system considered is composed of 848 degrees of freedom. To simulate it in real time so that haptic feedback can be used, the adaptive feature of the algorithm is used, and only 432 dof are activated. This allows to simulate complex systems in a few dozens of milliseconds.

The manipulation performed consists in extracting an inhibitor of the molecule, and then reinserting it (see figure 5). It is repeated twice, using both of the bilateral couplings presented in section III. The molecule is manipulated as a whole, so that it would be equivalent for the user to manipulate its center of mass in the case of a rigid molecule. However, the molecule is indeed flexible, and can be deformed by the simulation (interactions with other molecules). Using this manipulation mode, even an untrained user can move precisely the HIV's inhibitor.

The movement performed by the user is mainly directed along the y axis (see figure 5), so that it can easily be reproduced in order to compare the two control schemes. However, rotations are allowed.

The forces along the x, y and z axes are depicted in figures 6, 7 and 8. They are represented with respect to the position along the y axis (direction of the main displacement).

The parameters of the controllers are chosen so that the forces delivered by the simulator (in the order of dozens of nanonewtons) are amplified to be easily felt by the user. To select the position scaling factors, a trade off between stability of the system and ease of manipulation due to the limited workspace of the haptic device must be considered. These issues have been partially addressed in our previous work [START_REF] Bolopion | Tuning the gains of haptic couplings to improve force feedback stability in nanorobotics[END_REF]. In this experiment, the force scaling factor is the same for both the DFF and FF control schemes:

A f = 0.3 • 10 9 (2)
Regarding the velocity scaling factors, they are set to:

A DFF d = 0.25 • 10 9 0 0 5 ; A FF d = 0.025 • 10 9 0 0 5 (3)
Forces obtained with the DFF and FF couplings are similar. Three phases can be clearly distinguished (corresponding to the three states of the molecule depicted in figures 6, 7 and 8). The first phase is the extraction of one of the inhibitor, when it is situated at less than 15 Å of its initial position. The second phase consists in removing the inhibitor from the rest of the structure. When the inhibitor is further than 15 Å, no force is applied on it. The reinsertion into the system occurs when the position is again less than 15 Å. For both of the couplings, the magnitude of the forces is in the order of the newton so that the interaction forces can clearly be felt by the user. When the DFF control scheme is used, some oscillations can be felt. As seen in section III, this is due to the structure of the simulation since the rigid body first takes the position set by the user. This position is then modified according to the interactions with other atoms. Adding damping and changing the structure of the coupling has been considered in previous work, but decreases the transparency [START_REF] Bolopion | Tuning the gains of haptic couplings to improve force feedback stability in nanorobotics[END_REF].

Even if the forces felt using the DFF and FF are similar (especially on figures 7 and 8), some differences can be noted (e.g. dotted circle in figure 6). As stated earlier, rotations are free. For a given position along the y axis, the orientation can be different. The interactions between the manipulated inhibitor and the HIV protease are thus not the same which explains the differences on the force fed back. As demonstrated in this example, using both of the couplings, the system remains stable and the forces felt by the user allow him or her to clearly distinguish the interactions between the molecules. Compared to [START_REF] Daunay | 6 dof haptic feedback for molecular docking using wave variables[END_REF], no damping is added in our system which allows a good transparency. Even untrained users did not experience any difficulty in using either position or force control mode. Therefore, our system can indeed be used during a conception process of molecular structures, or to test docking properties of a ligand.

B. Internal forces

The tool that we present in this paper is not limited to molecular interactions. It can also be useful to manipulate a single rigid body of a flexible molecule, and feeling the internal forces and the inertial effects. This is useful in particular to synthesize a new molecule by increasing progressively the number of atoms. Therefore, it is possible to change the topology of the molecule and to control its specificities. The experiment presented in this paragraph consists in manipulating a single rigid body of a flexible molecule, and feeling internal forces and inertia. No other molecules are present. To perform this experiment, the user selects a rigid body of an inhibitor of the HIV virus. Using the haptic device and the FF coupling, he or she applies a force on this inhibitor to unfold it. The other extremity is fixed. When the molecule is stretched, the user moves back the haptic device to the initial position. Doing so, the molecule is refolded. Since the molecule is composed of only 17 degrees of freedom, the adaptive feature of the algorithm is not useful to obtain a fast simulation of the whole molecule (the simulation's sampling period is 14ms with all the dof The same position scaling factor A FF d than in section IV-A is used. The force scaling factor is A f = 0.1 • 10 9 .

The figure 9 represents the norm of the force felt by the user against the norm of the position of the center of mass. First, the force increases as the molecule is stretched (part 1). When the molecule is deployed, the user applies an additional force and the internal forces increase (part 2). Then, the user moves back the rigid body and he or she releases the constraint applied on it. The forces quickly decrease because of the reduction of internal forces (part 3). Finally, the molecule is moved toward its initial position. The interaction forces due to the other rigid bodies are felt (part 4).

These results show that our system is able to render small variations of the forces to the user (in the order of dozens of nanonewtons). It demonstrates that the system is well suited to achieve high transparency. The new possibility to manipulate a single rigid body of a molecule has been asked by Sanofi-Aventis. With this possibility, one can expect to evaluate in real-time the influence of the molecular structure on its dynamics.

V. CONCLUSION

We have presented a new tool to perform molecular simulations with haptic feedback. Using the force feedback device, the user can manipulate the entire HIV inhibitor. The simulator allows for real-time computation of the forces thanks to an adaptive method. To provide flexibility, two different manipulation modes are proposed (position and force) thanks to two different bilateral couplings. Both of them allow the operator to feel interaction forces. All the stages of the design process of specific molecular structures can be handled by our new system. The internal forces of a molecule can be felt by the user. This possibility may be relevant for studying the possible conformational changes of a molecule and its properties. The experiment we performed, and which has been tested by Sanofi-Aventis engineers, shows that the simulator and the coupling achieve high performances. Indeed, the low values of the internal forces are felt by the user in a transparent and stable way.

This paper demonstrates the first results obtained with this system. User-based tests must be carried out to evaluate the effectiveness of the proposed system compared to conventional devices such as a mouse or a joystick. Tests should also be carried out to compare complex manipulation tasks performed with or without haptic feedback.

Fig. 1 .

 1 Fig. 1. Manipulation of a tetraethylamonium (TEA) around a potassium channel (KcsA). A haptic device is used to control the molecule and feel the interactions between the ligand and the protein derived thanks to the adaptive simulation

Fig. 2 .

 2 Fig. 2. The binary assembly tree of a tetra-alanine.

Fig. 3 .

 3 Fig. 3. DFF control scheme

Fig. 4 .

 4 Fig. 4. FF control scheme The parameter k = k T 0 0 k R is a stiffness matrix used to tune the forces and the moments applied on the rigid body.

Fig. 5 .

 5 Fig. 5. Manipulation of one of the HIV protease's inhibitors using the haptic device

3Fig. 6 . 4 5Fig. 7 .

 647 Fig. 6. Force felt by the user along the x axis of the haptic interface while manipulating an inhibitor of the HIV molecule (interaction forces between molecules)

Fig. 8 .

 8 Fig. 8. Force felt by the user along the z axis of the haptic interface while manipulating an inhibitor of the HIV molecule (interaction forces between molecules)

Fig. 9 .

 9 Fig. 9. Force felt by the user while unfolding an inhibitor of the HIV virus (internal forces). Non adaptive simulation

  , B. Cagneau, S. Régnier are with Institut des Systèmes Intelligents et de Robotique, Université Pierre et Marie Curie, CNRS UMR 7222, 4 Place Jussieu, 75005 Paris, France. {bolopion, cagneau, regnier}@isir.fr

	S. Redon is with NANO-D -INRIA Rhône-Alpes, 655
	avenue	de	l'Europe,	38334	Saint-Ismier	Cedex,	France

http://www.sensable.com

VI. ACKNOWLEDGMENT

This work was supported by the French National Agency of Research, through the PACMAN project.


Julien Totems
email: julien.totems@onera.fr

Véronique Jolivet
email: veronique.jolivet@onera.fr

Jean-Philippe Ovarlez
email: jean-philippe.ovarlez@onera.fr

Nadine Martin
email: nadine.martin@gipsa-lab.inpg.fr

Jean-Philippe Ovarlez

Nadine Martin Signal

J.-P Ovarlez

Signal Processing Methods and Poly-Pulse Waveforms for Laser Vibrometry in Pulsed Mode

come   L'archive ouverte pluridisciplinaire

Introduction

Coherent laser radars can measure the small amplitude vibrations of remote surfaces targeted with the laser, by determining the Doppler shift of the back-scattered wave. The interference between this scattered wave and an artificially shifted local oscillator, on a detector, is the source of a heterodyne signal around this artificial frequency shift, with an additional shift which is proportional to the scattering surface's radial velocity. This technique has been successfully applied to the detection of structural damage in buildings [START_REF] Jolivet | On the full scale dynamic behaviour of RC-buildings using coherent laser radar vibrometer[END_REF] , and target recognition for military purposes [START_REF] Kranz | Target Classification by Laser Vibration Sensing[END_REF] . In both cases, the aim is to identify the vibration's modal frequencies, after frequency demodulating the signal.

Today, one of the major goals is to lengthen the range of these instruments. To this end, a pulsed mode has been developed [START_REF] Hill | Fiber-based 1.5 µm lidar vibrometer in pulsed and continuous modes[END_REF] , which allows the signal's peak power to be increased, for the same mean laser power. But until now, the performance of this operation mode, which usually works with pulse pairs, was limited by its sensitivity to phase noise and Doppler ambiguities that appear for poorly chosen pulse separations [START_REF] Gatt | Poly-Pulse Waveforms for Coherent Lidar Measurements[END_REF] . Gatt et al. [START_REF] Kachelmyer | Spectrogram Processing of Laser Vibration Data[END_REF] have shown the potential of poly-pulse waveforms to counter these problems, but did not describe the best suited signal processing. We compare various nonparametrical processing strategies for these waveforms, and evaluate their interest in regard to the classical CW and pulsepair operating modes, in the context of long range identification of moving targets.

Signal model

This study is based on the following expression of the signal, in a complex form:

) ( )) ( exp( ). ( . ). ( ) ( ) ( ). ( ) ( 0 t i t j t m i t t i t i t t i b vib b het S        (1)
where µ(t) is an amplitude modulation applied to the laser emission (µ(t) = 1 in CW mode), i het (t) is the heterodyne current, with a mean amplitude i 0 , m(t) is a complex multiplicative noise (speckle noise), gaussian valued, circular and centered, with a variance set to 1, φ vib (t) = 4π.x vib (t)/λ is the phase modulation caused by the targeted surface vibration x vib (t) for the laser wavelength λ, and i b (t) is an additive complex noise (detector / photon noise), white, gaussian valued, circular and centered, with variance σ b 2 . This simplified model is valid only if we assume that any global Doppler shift has previously been removed, and if we neglect the phase effects of atmospheric turbulence and laser frequency noise (which can be mitigated by decreasing optical path difference between local oscillator arm and measuring arm) in regard to speckle noise, which is mostly due to the coherent scattered figure, from the target's rough surface, moving in front of the reception aperture when the target rotates. Thus the multiplicative noise is directly linked to the characteristics and dynamics of so-called speckles.

The two main noise parameters are B m , the inverse of the multiplicative noise's coherence time, so that its autocorrelation function is: Γ m (τ) = exp(-B m 2 τ 2 ), and the Carrier to Noise Ratio (CNR), defined by <|i het | 2 > / <|i b | 2 > = i 0 2 / 2σ b 2 . B m will be chosen of a few kiloHertz, close to the modulation bandwidth created by the vibration (which corresponds to mm/s order vibration velocities). Such values are consistent with the actual parameters expected for long range moving vehicle identification (previous studies 3 considered speckle noise with lower bandwidth, few hundred Hertz). We remind that the goal is to determine v vib (t) = dx vib (t)/dt = λ/2 f inst (t), which is equivalent to a frequency demodulation in the presence of strong additive and complex multiplicative noises.

Pulsed waveforms are created by modulating the amplitude with µ(t), which takes values 0 or µ. Their parameters are defined in figure 1. Our comparison between various waveforms assumes in particular an equal mean laser power; the peak value µ is thus √1/(N p t p PRF). We also suppose that we have previous knowledge of the target's distance, and know precisely which samples contain signal. In practice, the radar that was used to detect the target can provide us with telemetry accurate enough for this purpose.

Instantaneous frequency estimators

We present various IF estimators that we developed, adapted, or found in laser radar or signal processing literature:

Table 1 : Summary of studied IF estimators CW mode Autocorrelation First Lag                K k s s inst t k t i t k t i t t f 1 * ) ). 1 ( ( ) . ( arg 2 1 ) ( ˆ SpectroGram Centroid                 B inst df B f j f t SG B t f   2 exp ). , ( arg 2 ) ( ˆ . Lee's Spectral Matching   ) , ( ) ( max arg ) ( ˆf t SG f S t f ref f inst   . PP mode Pulse-Pair   * 1 , 2 , arg 2
1 ) / ( ˆk s k s S inst i i T PRF k t f        . Poly-Pulse mode Poly Pulse-Pair                  p N l l k s l k s S inst i i T PRF k t f 2 * 1 , , arg 2
1 ) / ( ˆ Correlation FT             d f j f PRF k t f k i f inst s ) 2 exp( ). ( ). ( max arg ) / ( ˆ, Poly-Pulse Spectrum           k Multiplet ) 2 exp( ). ( ). ( max arg ) / ( ˆdt ft j t t i PRF k t f s f inst   IF Likelihood   ) .( . ) ( max arg ) / ( ˆ2 , 1 2 , ft j k s H ft j k s f inst e i Q e i PRF k t f         (Q = Q m + Q ib ,

sum of covariance matrices of multiplicative and additive noises)

For CW mode, we have implemented three popular estimators in literature: the short term « coherent average » of phase difference between consecutive samples (Autocorrelation First Lag 3 , AFL), the centroid of a spectrogram's columns (with a circular mean, in order to avoid bias from the non-zero-mean noisy background) (SpectroGram Centroid 5 , SGC), and a short-term spectrum matching, also on the spectrogram (Lee's Spectral Matching 6 , LSM).

Pulse Pair mode is processed by calculating the phase difference between each pair of pulses (Pulse-Pair, PP). This method is also implemented in Poly-Pulse mode: the phase difference between consecutive pulses is coherently averaged here again (Poly Pulse-Pair, PPP). But various works like Youmans' 6 insist on the benefit of differentiating phase over non-consecutive pulses, in which case the problem boils down to linear regression of the autocorrelation's phase, i.e. the search of its Fourier Transform's (FT) maximum. That is why we consider a second class of estimators for Poly-Pulse mode, which we call "spectral". The first one is based on the FT of the autocorrelation function, after a proper windowing to select relevant samples and lags (Correlation Fourier Transform, CFT). Another one, theoretically equivalent to the latter, uses the spectrum of the received waveform multiplied by the emitted waveform (Poly-Pulse Spectrum, PPS). This corresponds to radar matched filter, which is optimal when dealing with additive white noise. In case of multiplicative noise, though, PPS is not optimal, and we propose a new estimator based on the likelihood function of the poly-pulse's instantaneous frequency (Instantaneous Frequency Likelihood, IFL, using a variable change from Ghogho et al. [START_REF] Ghogho | Cramér-Rao Bounds and Maximum Likelihood Estimation for Random Amplitude Phase-Modulated signals[END_REF] ). It does require prior knowledge of noise characteristics (CNR and B m , which can be evaluated on the signal), and remains sensitive to deviations from the given signal model, but it is theoretically the best estimator for this model. For these spectral estimators, we also propose to build a pseudo Time-Frequency Representation (TFR), by juxtaposing in columns the function we need to maximize for each poly-pulse. In case of strong additive noise, peak tracking between columnsbetween poly-pulses-can help discriminate the peak at the true signal IF from the transitory peaks due to strong noise. We implemented a simple temporal smoothing on the TFR in order to show the potential of this option (we then call the processing methods: PPS-s, CFT-s and IFL-s). The same smoothing is applied to CW mode TFRs for a fair comparison (SGC-s, LSM-s).

Poly-pulse waveforms theoretical performance

Before applying these estimators to vibrometry, we study the error on the velocity estimate over a single waveform, in the case of pulse-pair and poly-pulses. We evaluated numerically the theoretical Cramér-Rao Lower Bound (CRLB) on this error (Q only on non-zero samples):

2 / 1 2 1 )] ) ( [tr(      Q f Q v  with I e e i u t Q b ft j u t B m 2 2 ) .( 0 2 . ) , ( 2 2        (2)
We had to add the effect of CNR variation in time (so-called fading) to better predict velocity error. On figure 2, we see that theory and simulation agree on better precision with waveforms containing six pulses rather than only two, as predicted by Gatt et al. [START_REF] Kachelmyer | Spectrogram Processing of Laser Vibration Data[END_REF] . Also, we notice that our IFL estimator almost reaches optimum precision when CNR is high, whereas PP and PPP estimators are limited by multiplicative noise. On the other hand, the CRLB cannot account for the degradation seen at low CNR; only simulation can help qualify waveforms and estimators in case of strong additive noise.

Simulation results

Using our heterodyne signal simulator, we have qualified all presented estimators, on the criterion of Signal-to-Noise Ratio (SNR), which is evaluated on the vibration spectrum, as the ratio of the power spectral density at the modal frequency over the power spectral density of the noise floor. The simulated vibration had 5 modal frequencies between 20 and 120 Hz, with maximum velocities around 1 mm/s. The retained SNR value is the average of the 5 individual SNRs. The multiplicative noise bandwidth B m is set to 5 kHz. For each given result, we average the SNR of 200 runs. In pulsed mode, with a strong additive noise (CNR = -25 dB in 1 MHz, cf. figure 3), our IFL estimator with a temporal smoothing of the RTF performs a little better than other processing methods. CW mode obtains equivalent SNR, unless smoothing is applied to the spectrogram, in which case SNR increases by about 5 dB, thanks to background noise mitigation.

When multiplicative noise is predominant (CNR = 20 dB in 1 MHz, cf. figure 4), CW mode also yields much higher SNR, thanks to an efficient phase noise averaging allowed by numerous measurements. In pulsed mode, the IFL estimator allows a 5 dB gain compared with estimators that are only adapted to additive noise, and a 10 dB gain compared with pulsepair mode. Our objective to improve the robustness of pulsed mode against speckle noise is met, with a computing time comparable to the one in CW mode.

For the previous results, it was assumed that the analysis bandwidth, defined by the filter used in CW mode or 1/Ts in pulsed mode, was closely adjusted to the frequency modulation bandwidth due to target vibration : B = 1/T S = 4 v vib,max / λ. Yet in practice, without prior knowledge of vibration maximum velocity, a larger analysis bandwidth has to be chosen. We show that, in case of strong additive noise, temporal smoothing allows TFR-based processing to be much less sensitive to bandwidth enlargement (cf. figure 5). With the chosen PRF, CW mode still yields higher SNR, but if the PRF is increased in pulsed mode, so as to better average phase noise, poly-pulse waveform with IFL processing allows over 3 dB gain, despite the loss in peak CNR due to the greater number of pulses. An optimum is found at PRF = 1500 Hz, which depends on relative strength of multiplicative and additive noises.

Conclusion

We have simulated and compared performance of various IF estimators for signal processing for coherent laser radar vibrometry, in a long range, moving target identification context. Our results indicate that the estimator we proposed, which is based on a pseudo time-frequency representation of the IF likelihood, is better than the other estimators available today for poly-pulse waveforms in pulsed operation mode. We confirm the advantage of this waveform over the classical pulsepair, particularly for its robustness against the complex multiplicative noise encountered. CW mode remains generally more efficient, when dealing with large bandwidth vibration-induced frequency modulation and comparably large bandwidth multiplicative noise. Yet, pulsed mode as we implemented it can be better in practice, without prior knowledge of the target's vibration speed. It is also intrinsically relevant in applications that require monostatic instruments and work with slow targets, since neither optical path isolation nor difference in carrier frequency can shield reception from noise due to signal parasitic reflections. Our future works will involve development of peak tracking techniques for optimal extraction of the vibration from the time-frequency representation, study of the effects of other perturbations like high frequency motor vibrations, and laboratory experimentation of poly-pulse vibrometry.

(Figure 1 :

 1 Figure 1: Poly-pulse waveform parameters

  vs  = B / (2 V max / ) @ CNR cw ,1MHz = -20 dB, B m = 5 kHz CW / SGC CW / SGC-s Pulse-Pair / PP Poly-Pulse (6) / IFL Poly-Pulse (6) / IFL-s Poly-Pulse (6) / IFL-s, PRF = 1 kHz

Figure 3 :

 3 Figure 3: SNR results of estimators under strong additive noise

Figure 4 :

 4 Figure 4: SNR results of estimators under predominant multiplicative noise

Figure 5 :Figure 2 : 6 )

 526 Figure 5: SNR results in case of analysis bandwidth mismatch

Acknowledgments

We thank Chris Hill from Qinetiq and Matthieu Valla from Onera for useful discussion.


Abdelhay Aboulaich

Malgorzata Geszke

Lavinia Balan

Jaafar Ghanbaja

Ghouti Medjahdi

Raphaël Schneider
email: raphael.schneider@ensic.inpl-nancy.fr

Water-based route to colloidal Mn-doped ZnSe and core/shell ZnSe/ZnS quantum dots

Relatively monodisperse and highly luminescent Mn 2+ -doped zinc blende ZnSe nanocrystals were synthesized in aqueous solution at 100°C using the nucleation-doping strategy. The effects of the experimental conditions and of the ligand on the synthesis of nanocrystals were investigated systematically. It was found that there were significant effects of molar ratio of precursors and heating time on the optical properties of ZnSe:Mn nanocrystals. Using 3-mercaptopropionic acid as capping ligand afforded 3.1 nm wide ZnSe:Mn QDs with very low surface defect density and which exhibited the Mn 2+ -related orange luminescence. The post-preparative introduction of a ZnS shell at the surface of the Mn 2+ -doped ZnSe QDs improved their photoluminescence properties, resulting in stronger emission.

A 2.5-fold increase in photoluminescence quantum yield (from 3.5 to 9%) and of Mn 2+ ion emission lifetime (from 0.62 to 1.39 ms) have been observed after surface passivation. The size and the structure of these QDs were also corroborated by using transmission electron microscopy, energy dispersive spectroscopy and X-ray powder diffraction.

Introduction

Semiconductor nanocrystals or quantum dots (QDs) have become one of the most attractive areas of current research because of their tremendous potential applications in photonic devices, biological detection and labelling and imaging. [1][2][3][4][5][6] The wide interest in QDs arises from their unique optical and electrical properties, including broad absorption and narrow emission spectra, large extinction coefficients, resistance to photobleaching, long fluorescence lifetime, and size-tunable emission. [7][8][9] Moreover, QDs with different emission colors can be simultaneously excited with a single light source, with minimal spectral overlap, thus providing significant advantages for multiplexed detection of molecular targets. [10][11][12] During the last ten years, cadmium-based QDs, and especially core/shell CdSe/ZnS, have been the workhorse of the QD emitters. Because of the ultimate elimination of highly toxic class A elements like Cd, Hg and Pb, numerous researches are currently conducted to develop QDs without any heavy metal.

Transition element-doped nanocrystals (d-dots) can be an alternative to Cd-based ones for numerous applications. [13][14] Over the various non-cadmium II-VI semiconductors that constitute hosts for dopants, zinc selenide ZnSe with a room temperature bulk band gap of 2.80 eV (442 nm) is especially interesting because it is widely used for various applications like light-emitting diodes (LEDs), photovoltaic solar cells, photocatalysts, and biological sensors. [15][16][17][18] In d-dots, the dopants lead to phenomena not found in the bulk state because their electronic states are confined in a very small volume. Mn 2+ is one of the most common dopant for many II-VI semiconductors because (i) its location inside the crystal lattice of the host material can be confirmed by electron paramagnetic resonance (EPR), and (ii) it exhibits a typical photoluminescence from the Mn 2+ 4 T 1 → 6 A 1 transition centered at ca. 595 nm. 19 Over the last ten years, a few groups have reported the synthesis of ZnSe [20][21][22][23][24][25] and of the related Mndoped ZnSe nanocrystals [26][27][28][29][30][31][32][33] via organometallic synthesis. This approach generally utilizes pyrolysis of organometallic complexes in coordinating solvents such as tri-n-octylphosphine (TOP) or tri-noctylphosphine oxide (TOPO). The exchange of the hydrophobic TOP or TOPO ligands with hydrophilic ligands, such as thioacids, and subsequent transfer of QDs from oil to aqueous solution can significantly reduce the photoluminescence quantum yield (PL QY) of the nanocrystals. 2 Direct synthesis of thiol-capped QDs in water is a promising alternative route to organometallic reactions and offers the following advantages: (1) it uses less toxic reagents, (2) it is cheaper and simpler, (3) it uses lower reaction temperatures (100°C) to obtain QDs with comparable PL QY and size-tunable fluorescence, (4) the surface functionalization with water-soluble ligands occurs during the synthesis, and (5) the QDs obtained have a small diameter, which is suitable for biological applications. Recent studies have demonstrated that good quality ZnSe QDs can be prepared in water either under hydrothermal conditions [34][35][36][37][38][39] or under microwave irradiation 40 by reaction of a Zn(+2) salt with NaHSe, generally in the presence of a thioalkyl acid stabilizer. Using a mixture of the dopant and of the host precursor during the nucleation step (also called nucleation-doping strategy), Wang et al. succeeded recently in the preparation of Mn 2+ -doped ZnSe QDs with a PL QY of 2.4%, which represents the first direct synthesis of fluorescent Mn 2+ -doped ZnSe QDs in aqueous solution. 41 In the present study, we report an alternative aqueous route to prepare ZnSe:Mn QDs stabilized by 3mercaptopropionic acid (MPA). The trap state emission at 480 nm of the d-dots obtained was found to be very weak and could completely be suppressed by passivating the Mn-doped ZnSe core material with the wider band gap shell material ZnS. Using this strategy, MPA-capped core/shell ZnSe:Mn/ZnS nanocrystals with PL QYs up to 9% in water at neutral pH could be readily prepared. , 99%), mercaptopropionic acid (MPA, 99%), selenium powder (99.5%, 100 mesh), sodium borohydride (NaBH 4 , 98%), and ethanol (HPLC grade) were used as received without additional purification. All solutions were prepared using Milli-Q water (18.2 MΩ.cm -1 , Millipore) as the solvent.

Experimental Section

Synthesis of ZnSe:Mn

QDs. The preparation of NaHSe was performed according to Klayman et al. 42 with some modifications. Under an argon atmosphere, 76 mg (2.009 mmol) of NaBH 4 was added to a small flask cooled with ice containing 1 mL of ultrapure water. 79 mg (1.004 mmol) of selenium powder was then added and a small outlet was connected to the flask to discharge the hydrogen pressure generated by the reduction of Se into NaHSe. After 3 h at 4°C, the black selenium powder disappeared and a white NaHSe solution was obtained. The solution was diluted with 20 mL of argon-saturated water and the concentration of the final NaHSe solution was 0.05 M.

Typical procedure for the preparation of MPA-capped ZnSe:Mn QDs: Solutions of 0.1 M ZnSO 4 . 7H 2 O (5 mL), 0.01 M Mn(OAc) 2 . 4H 2 O (2 mL) for a Mn 2+ doping of 4% relative to Zn 2+ and 0.5 M MPA (20 mL) were mixed and the pH of the mixture was adjusted to 10.3 by dropwise addition of a 2 M NaOH solution. The solution was placed in a three-necked flask fitted with a septum and valves, and further stirred under argon for 1 h. 9 mL of a freshly prepared NaHSe solution (0.05 M) were then injected through a syringe into the mixture at room temperature. The molar ratio Zn 2+ /Se 2-/MPA in the solution was 1/0.9/20. The growth of the MPA-capped ZnSe:Mn nanocrystals proceeded on refluxing at 100°C for 24 h under argon flow with a condenser attached. After cooling to room temperature, the MPA-capped ZnSe:Mn nanocrystals were precipitated by ethanol, the precipitate was centrifuged, washed 3 times with ethanol, and then dried in vacuum (30 mmHg, 2 h, rt).

Thioglycolic acid (TGA), 6-mercaptohexanoic acid (MHA), mercaptosuccinic acid (MSA), and cysteine(Cys)-capped ZnSe:Mn QDs were prepared at pH = 10.3 using a molar ratio of Zn 2+ /Se 2-/ligand of 1/0.9/20, similarly to MPA-capped QDs.

Synthesis of MPA-capped ZnSe:

Mn/ZnS core/shell nanocrystals. For the ZnS shell growth, 10 mL of 0.2 M Zn(OAc) 2 . 2H 2 O solution and 0.7 mL of MPA were mixed together and the solution was diluted to 88 mL with water. The pH was adjusted to 10.3 with 4 M NaOH and the solution was saturated with argon by bubbling for 1h. The ZnSe:Mn core solution was prepared by dispersing 20 mg of ZnSe:Mn QDs in 130 mL of water, transferred into a three-necked flask fitted with a septum and valves and purged by argon bubbling for 1h. Then, 20 mL of the Zn 2+ -MPA complex were added dropwise to the ZnSe:Mn solution and the mixture was heated at 100°C for 10 h. After cooling to room temperature, the reaction mixture was concentrated down to approximately 15 mL, nanocrystals were precipitated with ethanol and dried in vacuum at room temperature for 12 h before measurements.

Instruments.

Transmission electron microscopy (TEM) images were taken by placing a drop of the particles in water onto a carbon film-supported copper grid. Samples were studied using a Philips CM20 instrument operating at 200 kV equipped with Energy dispersive X-ray Spectrometer (EDX). The X-ray powder diffraction (XRD) diagrams of all samples were measured using Panalytical X'Pert Pro MPD diffractometer using Cu Kα radiation (λ = 1.5405 Å). Absorption spectra were recorded on a Perkin-Elmer (Lambda 2) UV-Visible spectrophotometer. Fluorescence spectra were recorded on a fluorolog-3 spectrofluorimeter F222 (Jobin Yvon) using a 450 W xenon lamp. The QY values were determined from the following equation:

QY(sample) = (F sample /F ref )(A ref /A sample )(n sample 2 /n ref 2 )QY (ref)
where F, A and n are the measured fluorescence (area under the emission peak), absorbance at the excitation wavelength and refractive index of the solvent respectively. PL spectra were spectrally corrected and quantum yields were determined relative to Rhodamine 6G (QY = 94%). 43 The fluorescence lifetimes were measured by a FluoroMax-4 spectrofluorometer (Jobin Yvon) using a NanoLED emitting at 372 nm as an excitation source with a nanoled controller module, Fluorohub from IBH, operating at 1 MHz. The detection was based on an R928P type photomultiplier from Hamamatsu.

Fourier transform infrared spectroscopy (FT-IR) was performed using a Brucker Vector 22 spectrometer.

Dynamic light scattering (DLS) was performed at room temperature using a Malvern zetasizer HsA instrument with a He-Ne laser (4 10 -3 W) at a wavelength of 633 nm. The QDs aqueous solutions were filtered through Millipore membranes (0.2 µm pore size). The data were analyzed by the CONTIN method to obtain the hydrodynamic diameter (d H ) and the size distribution in each aqueous dispersion of nanoparticles.

Results and discussion

A series of ZnSe:Mn QDs were synthesized in aqueous solution containing thioalkylacids (thioglycolic acid TGA, 3-mercaptopropionic acid MPA, 6-thiohexanoic acid MHA, mercaptosuccinic acid MSA or cysteine Cys). Thioalkyl acids contain two functional groups possessing lone pair electrons that might influence the nucleation, the growth and the stabilization of ZnSe:Mn QDs and hence affect their physical and optical properties. Since the synthesis of the d-dots is conducted at basic pH, the negative charges of thioalkyl acids will also impose a kinetic barrier for the growth due to Coulombic repulsions between the ligand shells of neighboring QDs and thus favour stable and small nanocrystals.

Preparation of Thioalkyl

Acid-capped ZnSe:Mn Nanocrystals. We first studied the synthesis and the optical properties of 3-mercaptopropionic acid-capped ZnSe nanocrystals.

Influence of the Zn/Se ratio and of the Zn 2+ concentration. Initial experiments performed with 4 % doping in Mn 2+ and a 2.7 mM ZnSO 4 solution demonstrated that good quality QDs (photoluminescence quantum yield = 3.5%) were obtained when using NaHSe in short supply (by ca. 0.9-fold) relative to the mole amount of the Zn 2+ precursor. A decrease (0.4 or 0.6) or an increase of the Zn/Se ratio (1.2 or 2.0) yielded ZnSe:Mn QDs with poor PL efficiency. It is well-known that Zn 2+ precursors have a low reactivity towards NaHSe. 44 Increasing the concentration of ZnSO 4 (5.0 or 10.0 mM) while maintaining a Zn/Se molar ratio of 1/0.9 did not effectively improve the PL properties of ZnSe:Mn QDs (PL QY < 2%).

Influence of the pH. To improve the PL properties of the Mn-doped ZnSe QDs, that is, improvement of the Mn 2+ -related 4 T 1 → 6 A 1 emission and suppression of the trap emission, we also investigated the influence of the pH of the Zn 2+ -MPA precursors solution using 4 % doping in Mn 2+ and a Zn 2+ /Se 2- /MPA molar ratio of 1/0.9/20. We found that when the pH of the Zn 2+ /MPA precursors solution was set at 10.3 (the pH increased to 10.4 after the injection of the NaHSe solution), the as-prepared nanocrystals obtained after 24 h heating showed a strong orange emission at ca 595 nm and a very weak trap emission. When the pH value of the Zn 2+ /MPA precursors was decreased below 10.0 (9.0 or 8.0), the trap emission increased and the PL QY decreased (<1%). We fixed the pH value of the Zn 2+ /MPA precursors solution at 10.3 in latter experiments.

Influence of the heating time. After the injection of NaHSe, small aliquots were removed from the reaction solution with a syringe to monitor the growth of nanocrystals via optical absorption. Using a Zn 2+ /Se 2-/MPA molar ratio of 1/0.9/20 and 4% Mn 2+ relative to Zn 2+ , a well-resolved peak at 343 nm was observed after 15 min heating at 100°C indicating the homogeneous nucleation of small crystallites (Figure 1). The absorption peak attenuated sharply for longer growth times becoming a shoulder due to the Ostwald ripening process (formation of larger nanocrystals at the expense of smaller ones) and gradually red-shifted to 378 nm, which was consistent with the growth of the ZnSe host. The subtle blue-shift in absorption spectra observed after 24 h reflux (375 nm) may be a consequence of the partial decomposition of the MPA ligand surrounding the ZnSe:Mn core, leading to a sulfur-enriched shell. with a Zn 2+ /Se 2-/MPA molar ratio of 1/0.9/20 and 4 % doping in Mn 2+ .

The photoluminescence of the clear pale yellow solution was also monitored in the course of the nanocrystals growth. All aliquots taken during the first 5 h of growth exhibited at first no PL. Upon resting at room temperature, a strong orange-red emission evolved over a matter of hours (in some cases, minutes) (Figure 2). PL QYs measured from nanocrystals obtained after 30 min or 2 h of growth at 100°C and respectively 6 and 21 h storage at room temperature were high (23 and 19%, respectively, without any post-preparative treatment and the absence of the broad trap-emission centered at ca. 480 nm) (Figure 2). The photoluminescence decreased however gradually and disappeared after 24 h for the dots synthesized during 30 min at 100°C (Figure 2a). It should also be pointed out that the stability of the d-dots obtained after short heating times (< 3 h) was found to be modest. A deterioration of these nanocrystals accompanied by precipitation was generally observed after 48 h storage at room temperature. Good quality powder X-ray diffraction (XRD) patterns could also not be obtained for MPA-capped ZnSe:Mn QDs prepared during short reaction times (< 3 h). As the reaction prolonged, widths of the diffraction peaks were found to decrease and were more resolved (vide infra). The stability of ZnSe:Mn d-dots was markedly increased by extending the heating time at 100°C. No sign of macroscopic aggregation or loss of fluorescence were observed for the dots obtained after 24 h heating upon storage at room temperature for up to 3 months.

One interesting feature during the early stages of growth (15 min to 2 h) is the significant blue-shift of the maximum wavelength of the Mn 2+ 4 T 1 → 6 A 1 emission from 609 to 595 nm after excitation at 350 nm. Negatively-charged MPA ligands strongly interact with surface zinc atoms of small doped-ZnSe nanocrystals formed at the beginning of the growth stage. As a consequence, the d-orbitals of the Mn center will experience a high electric field difference along different directions, which results in a higher splitting of the energy levels. 29 Once the thickness of the ZnSe shell increases around the doped-ZnSe cores, the lattice field around Mn 2+ ions becomes more symmetric leading to smaller splitting of the energy levels and a change of the PL color from orange-red to orange. Influence of the ligand. We also tested the impact of various thiols (TGA, MHA, MSA and Cys) on the PL properties of ZnSe:Mn QDs when the Zn 2+ /Se 2-/ligand molar ratio and the doping percentage were fixed at 1/0.9/20 and 4, respectively. All reactions were conducted at 100°C for 24 h. It was found that TGA-, MHA-, MSA-, and Cys-capped ZnSe:Mn QDs exhibited the Mn 2+ 4 T 1 → 6 A 1 transition at ca. 595 nm but the emissions were found to be broad (the full-width-at-half-maximum increased up to 94 nm while it was of ca. 60 nm for MPA-capped QDs). The PL QYs were also found to be weaker (< 1%) compared to nanocrystals prepared with MPA. The nanocrystal growth with TGA, MHA, MSA, and Cys ligands occurred finally with higher defect densities (increase of the ZnSe band-edge transition at ca. 460 nm). Compared to the other thiols used in this study, we suggest that MPA provides the better surface passivation of the QD crystalline lattice under our synthetic conditions.

Influence of the Mn 2+ concentration. The ratio of Mn 2+ relative to Zn 2+ was varied from 1 to 6% to evaluate the effects of Mn 2+ on the optical properties of the d-dots under the same experimental conditions (Zn 2+ /NaHSe/MPA = 1/0.9/20, heating time of 24 h). While the dopant PL peak position and spectral contour were found to be independent of the Mn 2+ concentration in the solution (PL emission at 596 nm and full-width-at-half-maximum of ca. 60 nm), an increase of the Mn 2+ concentration from 1 to 4% relative to Zn 2+ leads to an increase in the PL emission intensity, which also proves the successful incorporation of Mn 2+ into ZnSe QDs (Figure 3). Further increasing the concentration of Mn 2+ to 5 or 6% leads to a decrease of the PL QYs. The effect of Mn 2+ ion concentration on the PL properties of Mn 2+ -doped II-VI semiconductors QDs has been extensively studied in recent years and a so-called "concentration quenching effect" has been demonstrated. [45][46][47][48] When the concentration of Mn 2+ becomes higher than a certain threshold, the non radiative energy transfers between neighboring Mn 2+ dopant ions reduce and even annihilate the fluorescence. Under our synthetic conditions, the formation of such pairs of Mn 2+ dopant ions is probably observed when the amount of Mn 2+ ions in QDs is higher than 4 %.

Energy Dispersive X-ray (EDX) spectroscopy experiments (vide infra) conducted on ZnSe:Mn@MPA QDs after purification indicated that the Mn 2+ at. % was ca. 3%. This value is slightly higher than those used in organometallic syntheses (0.5-3% in growth solutions leading to nanocrystals doped with ca. 0.02-2.0 % Mn 2+ ions).

It was also surprising to observe that the 4 T 1 → 6 A 1 transition of Mn 2+ ions is red-shifted by ca. 25 nm compared to the MPA-capped ZnSe:Mn QDs prepared by Wang et al. 41 using an aqueous route. As demonstrated by several groups, [49][50][51] the position Mn 2+ emission peak is sensitive to the electric field imposed by the organic capping layer, to the Mn 2+ radial positions, to the heating time, and to the Mn 2+ concentration inside the nanocrystals. Wang and coworkers used a theoretical concentration of Mn 2+ twice lower than ours (4 versus 2) and heated the reaction mixture for 5 h, while the growth of our nanocrystals was performed during 24 h. A possible explanation for the red-shift observed with Mndoped ZnSe dots prepared under our synthetic conditions is that the crystal field for each dopant ion became more symmetric in long range due to the thicker ZnSe shell formed. As a result, the crystal field splitting of Mn 2+ d-orbitals became smaller, which resulted in a PL red-shift. 400 500 600 700 Mn 2+ ion (0.83 nm) is larger than that of the Zn 2+ ion (0.74 nm). 52 To confirm that the crystal lattice of the ZnSe host nanocrystals was intact, powder XRD pattern of the MPA-capped ZnSe:Mn QDs was recorded (Figure 4a). The XRD pattern of the nanocrystals exhibited characteristic peaks at ca. 27.5, 45.6, and 54.3° corresponding to the (111), ( 220) and (311) reflecting planes of cubic zinc blende ZnSe (JCPDS Card N° 80-0021). The peaks are broadened due to the finite crystalline size. These results prove that the ZnSe:Mn d-dots have a cubic zinc blende structure similar to undoped ZnSe QDs and confirm that the crystal lattice of the ZnSe host material is intact after doping with Mn 2+ ions. A crystallite size of 3.5 ± 1.0 nm was estimated from the diffraction line broadening, value which correlated well with the results obtained from TEM analysis (vide infra). Finally, it is worth noting that the XRD patterns of the ZnSe:Mn QDs are located on higher angle compared to that of the pure ZnSe material probably because the MPA ligand partially decomposed by heating leading to a sulfur-enriched ZnS shell with smaller lattice parameters compared to ZnSe. The formation of alloyed doped-ZnSe(S)

0,0 2,0x10 6 4,0x10 6 6,0x10 6 8,0x10 6 1,
QDs is likely after 24 h heating at 100°C. where E is the energy of the first exited state, Eg = 2.80 eV is the band gap energy of bulk ZnSe, m e and m h are the effective masses of the electron and hole in ZnSe (m e = 0.15m o , m h = 0.66 m o , where m o = 9.11 x 10 -28 g is the free-electron mass), respectively, ε = 9.2 is the semiconductor dielectric constant, ħ = 6.58 x 10 -16 eV is the reduced Planck constant, and e = 1.6 x 10 -9 C is the electron charge. Using the absorption peak for MPA-capped Mn-doped ZnSe QDs after 24 h of growth yields a nanocrystal diameter of 2.7 ± 0.4 nm. The deviation of 0.4 nm corresponds to ± 10 nm for the determination of the onset of absorption. The PLE spectrum corresponding to the 600 nm emission exhibits maximum at about 345 nm corresponding to the band gap excitation of host ZnSe nanocrystals.

E = Eg + π 2 h 2 2R 2 1 m e + 1 m h - 1.8e 2 εR
300 400 500 600 700

0,0 0,2 0,4 0,6 0,8 1,0 1,2 Absorbance (a.u.) Wavelength (nm)
PL Intensity (a.u.) Figure 7 shows a transmission electron microscopy (TEM) micrograph of 3.1 ± 0.3 nm diameter Mndoped ZnSe QDs, corresponding to the final spectrum of Figure 1 with 3.5% PL QY (labelled a), which is in good accordance with the value calculated from the quantum confinement equation. Dynamic lightscattering (DLS) measurements indicate that the QDs are well-dispersed in water and that their hydrodynamic diameters range from 5 to 11 nm (average hydrodynamic diameter of 7.6 nm). The hydrodynamic diameter is larger than the diameter of the inorganic core due to the solvatation layer around the QDs in aqueous solution. Energy dispersive spectroscopy (EDS) measurements revealed that the atomic composition of the d-dots had a value similar to that of the ratios of precursors (see details in the core/shell subsection). An average diameter of 3.1 nm represents ca 88 ZnSe structural units in each monocrystal (using a theoretical lattice parameter a of 0.56 nm, JCPDS Card N° 80-0021). Using 4% Mn 2+ dopant means ca 14 Mn 2+ ions per nanocrystal.

Epitaxial Growth of the ZnS Shell.

A synthesis of core/shell ZnSe:Mn/ZnS QDs was also developed with the goal to improve the chemical, PL stabilities and QY of the QDs but also with the hope that the higher band gap ZnS shell with a small lattice mismatch with ZnSe (ca. 5%) will provide a diffusion barrier by stabilizing the Mn 2+ doping ions. Although, there is no report on the preparation of ZnSe:Mn/ZnS d-dots, the preparation of ZnSe/ZnS d-dots in organic 31,[53][54][55][56] or aqueous 57,58 solution has already been investigated by a few research groups. The ZnS shell was successfully introduced at the periphery of ZnSe:Mn cores by thermal decomposition at basic pH of Zn 2+ /MPA complexes. The PL intensity of Mn 2+ that corresponds to the 4 T 1 → 6 A 1 transition within the 3d 5 configuration of Mn 2+ increased gradually upon heating at 100°C and reached its maximum value (PL QY = 9%) after 10 h (Figure 8a). No improvement was observed by further heating. For ZnSe:Mn@MPA nanocrystals under the same heating conditions but in the absence of the Zn 2+ /MPA complex, there is almost no change of PL QY. This suggests that the enhancement of PL QY is related to the growth of the ZnS shell. The enhanced QY of ZnS-capped ZnSe:Mn d-dots is a direct consequence of an effective surface passivation by which nonradiative recombination that compete with the regular luminescence are reduced significantly. A visual comparison of non-capped and of the ZnS-capped QDs is shown as inset in Figure 8b. Both samples are fluorescent under a handheld UV lamp providing 366 nm multiband irradiation, and it is obvious that the ZnS-passivated d-dots are much brighter. As shown in Figure 8b, the PL is shifted at 610 nm while it is located at 600 nm for the ZnSe:Mn core. Since Mn 2+ had a higher affinity for ZnS than for ZnSe, 30 we suppose that surface Mn 2+ ions of the ZnSe:Mn core diffuse into the ZnS lattice in which they are subjected to the electric field of the deprotonated MPA ligand. In the XRD pattern (Figure 4b), the peak positions shift to higher angles towards the position of ZnS, with peak widths and shapes being nearly unchanged, which proves the formation of core/shell ZnSe:Mn/ZnS QDs. This deviation could also result from the compression of the ZnS shell on the lattice planes of the doped-ZnSe core. Using TEM, the average diameter of core/shell nanocrystals obtained was found to be 4.3 ± 0.5 nm (average hydrodynamic diameter = 8.0 nm) (Figure 9). EDX measurements of the core ZnSe:Mn and of the core/shell ZnSe:Mn/ZnS QDs confirmed the nature of the materials (Figure 10). Zn, Se, and Mn peaks arise from the doped core while the detected S mainly comes from the capping MPA ligand. The presence of silica (peak at ca. 1.7 keV) may originate from a contamination of the samples by glassware during their storage. The growth of the ZnS shell around the ZnSe cores was clearly demonstrated by EDX experiments. With the deposition of the ZnS shell, the Zn/Se ratio increased from 1.02 in the ZnSe:Mn cores to 9.88 in the core/shell nanocrystals. In the meantime, the Zn/S atomic percentage decreased from 5.17 to 1.03. The overcoating of ZnSe:Mn d-dots with ZnS was further confirmed by measuring the slow PL decay of the nanocrystals before and after surface modification. Indeed, Mn 2+ ion emission lifetime is dependent on the diameter of the nanocrystals and is longer for QDs with thicker overcoating layer or thicker diffusion region. 59 The slow PL decay curves of Mn 2+ ions for the core ZnSe:Mn and the core/shell ZnSe:Mn/ZnS QDs are shown in Figure 11. The slow decay times, attributed to the emission of the single isolated Mn 2+ ions in a cubic site, [60][61][62] were 0.62 and 1.39 ms for the core and core/shell QDs, respectively. Interestingly, the shell dependence of the PL lifetime is related to the PL QY. After introduction of the shell, the PL QY increased by 257% (from 3.5 to 9%) while the PL lifetime was found to increase by 224%. at pH = 12). The PL QY at neutral pH is 8.9 %. A rapid decay of PL was observed when the pH was below 6.0 (PL QY = 2.5% at pH = 4). Under pH = 4.0, the QDs became unstable in the acidic aqueous environment, they precipitated and the PL disappeared.

Conclusions

We have developed a new synthetic method to prepare in aqueous solution highly luminescent and monodisperse Mn-doped ZnSe d-dots and demonstrated that molar ratios of precursors, Mn 2+ concentration, pH value, heating time and stabilizer play an important role in determining the stability and optical properties of these nanoparticles. Among the tested stabilizer, 3-mercaptopropionic acid is superior in terms of high PL QY and low trapping emission. The PL properties were markedly improved by overcoating the Mn-doped ZnSe core with a ZnS shell. The core/shell ZnSe:Mn/ZnS QDs had a PL QY of ca. 9% in water at neutral pH, a good crystallinity and a diameter of ca.

100°C

Highly luminescent manganese-doped zinc selenide (ZnSe:Mn) quantum dots were synthesized in aqueous solution at 100°C using the nucleation-doping strategy. The influence of various experimental variables, including Mn 2+ concentration, Zn 2+ /Se 2-/3-mercaptopropionic acid molar ratio, and reaction time have been systematically investigated to obtain quantum dots with high optical quality and pure dopant emission. Overcoating a high bandgap ZnS shell around the ZnSe:Mn cores allows significant improvement of the photoluminescence properties of the nanocrystals. A 2.5-fold increase in photoluminescence quantum yield (from 3.5 to 9%) and of Mn 2+ ion emission lifetime have been observed after surface passivation.

Figure 1 .

 1 Figure 1.Temporal evolution of the UV-vis spectra of ZnSe:Mn@MPA nanocrystals grown at 100°C

Figure 2 .

 2 Figure 2.Temporal evolution of the PL spectra of ZnSe:Mn@MPA nanocrystals prepared with a Zn 2+ /Se 2-/MPA molar ratio of 1/0.9/20 and 4 % doping in Mn 2+ and grown at 100°C (a) during 30 min,

3. 2 .

 2 Size, Shape, Crystal Structure and PL Properties of MPA-capped ZnSe:Mn d-dots. The dopants might also change the lattice structure of the ZnSe host material because the diameter of the

Figure 4 .

 4 Figure 4. XRD patterns of ZnSe:Mn@MPA QDs (a) before and (b) after introduction of the ZnS shell.

Figure 5

 5 Figure 5 presents UV-vis, PL, and PL excitation (PLE) spectra of ZnSe:Mn dots d-dots (~ 3 nm in

Figure 5 .Figure 6 .

 56 Figure 5. Optical absorption (solid line), PL excitation (doted line) and PL emission (dashed line) after

Figure 7 .

 7 Figure 7. (a) TEM micrograph of ZnSe:Mn@MPA d-dots, (b) the corresponding particle size

Figure 8 .Figure 9 .

 89 Figure 8. (a) Evolution of PL spectra during the overcoating of ZnSe:Mn QDs with the ZnS shell. (b)

Figure 10 .

 10 Figure 10. EDX spectra of ZnSe:Mn@MPA QDs (a) before and (b) after the introduction of the ZnS

Figure 11 . 3 . 4 .

 1134 Figure 11. PL decay curves recorded at room temperature of ZnSe:Mn and ZnSe:Mn/ZnS QDs at a

Figure 12 .

 12 Figure 12. (a) Fluorescence spectra of ZnSe:Mn/ZnS QDs in aqueous solution at different pH values.

2.1. Reagents and Materials.

  Zinc sulfate heptahydrate (ZnSO 4 . 7H 2 O, 99.99%), zinc acetate dihydrate (Zn(OAc) 2 . 2H 2 O, 98+%), manganese acetate tetrahydrate (Mn(OAc) 2 . 4H 2

  Influence of the concentration of Mn 2+ relative to Zn 2+ on PL spectra of ZnSe:Mn@MPA nanocrystals prepared with a Zn 2+ /Se 2-/MPA molar ratio of 1/0.9/20 and grown at 100°C for 24 h.

	1,6x10	7
	1,4x10	7	2%
	0x10 1,2x10 PL Intensity (a.u.)	7 7	3% 4% 5% 6%
		Wavelength (nm)
	Figure 3. Excitation wavelength is 350 nm.

CH 2 CO 2 H NaHSe pH = 10.3 100°C Zn(OAc) 2 HSCH 2 CH 2 CO 2 Na

  4.3 nm. It is hoped thatZnSe:Mn/ZnS core/shell QDs prepared by this convenient procedure will find diverse biological applications in areas where water-soluble QDs with high PL efficiency and long-term stability are needed.

			ZnS
	ZnSO 4
	+ Mn(OAc) 2	ZnSe:Mn	ZnSe:Mn
	+
	HSCH 2

Acknowledgment. This project is funded by the "EPST-Universités Lorraines". MG acknowledges partial funding of this research by the French Embassy in Poland. We would also like to thank Dr.

Ariane Boudier (EA 3452, Nancy-University) for FT-IR experiments.