Djedjiga Benouioua
email: djedjiga.benouioua@femto-st.fr

Denis Candusso

Fabien Harel

P Picard

Xavier Francois

X François

On the issue of the pemfc operating fault identification: generic analysis tool based on voltage pointwise singularity strengths

The purpose of this article is to study the portability of a non-intrusive and free of any external / internal disturbance 16 diagnosis tool devoted to the monitoring of the State of Health (SoH) of PEM Fuel Cell (PEMFC) stack. The tool is based on a 17 thorough analysis of the stack voltage signal using a multifractal formalism and wavelet leaders. It offers well-suited signatures 18 indicators on the SoH of the Fuel Cell. Some relevant descriptors extracted from these patterns (singularity features) are used in 19 the frame of Machine Learning approaches to allow the PEMFC fault identification. The proposed diagnosis strategy is evaluated 20 with two different PEMFC stacks. The first one is designed for automotive applications and the second one is dedicated to stationary 21 use (micro combined heat and power -µCHP application).

22 23 24 The obtained classification results show that the proposed 1 PEMFC diagnosis tool allows identifying simple operating 2 faults and more complicated situations that contain several 3 fault types. The diagnosis are realized on two different test

Introduction

Fuel cells (FC) are considered as a promising alternative way for energy conversion [1]. To ensure their durability, reliability and safety, many fault diagnosis and fault tolerant control methods have been proposed. These methods can be classified into two groups: model-based methods [2-5] and data-driven methods [6][7][8]. The methods of the first group are very cumbersome and complex because they require an in-depth knowledge of the multi-physical mechanisms (thermal, electrical, electrochemical, and fluidic ones) which can occur in FC systems. They are based on numerous parameters governing its operation, and their values are difficult to estimate. Some model-based methods allows a deeper understanding of the FC physics but they might be not suitable to provide an accurate / quantitative description of the FC performances. Hence, the data-driven techniques attract more and more attention because of their simplicity regarding the implementation and the good performances obtained without profound system structure knowledge. They are supported by efficient signal processing methods as: Fourier transform [9], multi-resolution analysis [10,11], singularity analysis [12,13].

For fault identification and isolation tasks, some works use Electrochemical Impedance Spectra (EIS) as normal or faulty operation signatures to supply artificial intelligence algorithms

71

Then, the portability of the proposed tool is discussed. Main

Experimental

41

• an electric / electronic management sub-system,

42

• an electronic load.

44

The monitoring and the control of the FC testbenches

60

The set of the tests performed, with the two stacks and 61 according to the degrading operating parameters, is reported 62 in Table 3.

63 64 The Hölder exponent is defined as the largest such as:

ℎ( ) = sup : ∈ ( ) .
In practice, when ℎ( ) takes a value which tends to zero, the signal exhibits a strong singularity at (Fig. 2).

28

The singularity strength regarding the variability of the 29 regularity of ( ) vs. is usually described through the so- We assess the generalization and the usefulness of the proposed VSS by establishing the diagnosis on two databases issued from the two investigated PEMFC stacks.

For the stack PEMFCAuto, 10 scenarios are analyzed (normal and abnormal conditions including one fault or more complex situations with 2 or 3 faults occurred simultaneously).

For the stack PEMFCµCHP, 9 operating conditions (normal and abnormal conditions) are studied.

In this aim, we use 30 voltage profiles for each FC operating condition; each voltage profile covers 1000 voltage samples acquired at a frequency fa =11Hz (Fig. 4) using the monitoring data system of the FC experimental test bench.

To perform the VSS, the analyzing wavelet is selected as Daubechies 3 ('Db3') function (with 3 vanishing moments, Fig. 4).

49

For the second stack PEMFCµCHP, 9 classes are discriminated 50 with a success of 95.5 %.

52

The best results of the classification rates of the different 53 operating conditions studied are embedded in confusion 54 matrixes given by Tables 4 and5, for the stacks PEMFCAuto

(

  based for instance on fuzzy logic [14] or neural networks [15]) or conventional pattern recognition approaches (based on Support Vector Machines (SVM) [16], k-nearest neighbors 51 (knn) [13] methods). 52 This work aims at studying the portability of an innovative data 53 driven approach dedicated to PEMFC diagnosis, named 54 singularity analysis. This method consists in analyzing the 55 pointwise singularities stamped in the stack voltage signal for 56 various FC operating conditions. The singularity features are 57 then summarized in the form of concave arcs estimated thanks 58 to a set of mathematical equations, baptized multifractal 59 formalism [17-19]. The advanced analysis tool, named Voltage 60 Singularity Spectrum (VSS) is then obtained using a non-61 intrusive manner and without affecting in any way the FC 62 operation. Indeed, no external additional AC-solicitation has to 63 be superimposed to the existing DC load current as it is the case 64 in the usual EIS operation mode. 65 This paper is organized as follow. Section 2 deals with the 66 experimental work and environment conducted with two 67 PEMFC stacks. In section 3, a brief mathematical foundation 68 of the singularity measurement is given. In section 4, we show 69 how it is possible to make the singularity spectrum combined 70 with Machine Learning techniques as a PEMFC diagnosis tool.

  25

Fig. 1 .

 1 Fig. 1. Pictures of the two investigated PEMFC stacks:

45

  parameters are done through National Instruments materials 46 and a dedicated software. A Human-Machine Interface (HMI) 47 was also developed in-lab using Labview TM . 48 The purpose of the experimental protocol was to introduce 49 different controlled health states into the PEMFC stacks by 50 configuring various operating parameters. Two PEMFC health 51 states were set: no stack failure state (when the FC is operating 52 under normal parameters) and with stack failure state (when 53 the FC is operating under abnormal / severe parameters). 54 Different scenarios of FC system failure were associated to the 55 changing of different setting physical parameters through the 56 control-command interface of the testbenches, namely:57 cathode stoichiometry rate (FSC), anode stoichiometry rate 58 (FSA), gas pressure (P), cooling circuit temperature (T) and 59 relative humidity level (RH).

  Multifractals are the distribution of singularities, all lying on interwoven sets of varying fractal dimensions [20]. Multifractals arise in a variety of physical signals in the forms of self-similar sub-sets of samples, such as Multifractional Brownian Motion (MBM) signal [21], signals in turbulence [22], and temporal series in finance [23]. The singularity analysis, also named "multifractal analysis" allows the characterization of data by describing globally and geometrically the fluctuations of local regularity, usually measured by means of the Hölder exponent ℎ. It describes the time ( ) based fluctuations of the signal given by a function ( ). This is achieved by comparing the local variations of ( ) around fixed time position , against a local power law behavior: ( ) is said to belong to ( ) with ≥ 0 if there are a positive constant and a polynomial , satisfying deg ( ) < , such that ( ) -( ) ≤ | -| .

Fig. 2 .

 2 Fig. 2. Illustration of two pointwise regularity measures in the 26

30Fig. 3 .

 3 Fig. 3. Example of a singularity spectrum (ℎ).

Fig. 4 . 33 Fig. 5 .

 4335 Fig. 4. Stack voltage signal scanned by Daubechies wavelet (Db3) for VSS computing.

17 Table 4 .

 174 Confusion matrix of the good classification rates obtained with MRMR and kNN from the VSS computed with the stack voltage signals of the PEMFCAuto. The studied FC operating conditions are: Ref: normal conditions -DFSC: cathode flow fault (slight air starvation) -DFSA: anode flow fault (slight H 2 starvation) -DP: gas pressure fault (lower gas pressure) -DT: cooling circuit temperature fault (lower stack temperature) -DCO: carbon monoxide poisoning (Table 5. Confusion matrix of the good classification rates obtained with MRMR and kNN from the VSS computed with the stack voltage signals of the PEMFCµCHP. The studied FC operating conditions are: Ref : normal conditions -DFSC ↗ : cathode flow fault (air over-supply) -DFSC ↘ : cathode flow fault (slight air starvation) -DFSA ↗ : anode flow fault (H 2 over-supply) -DFSA ↘ : anode flow fault (H 2 starvation) -DT ↗ : cooling circuit temperature fault (higher temperature) -DT ↘ : cooling circuit temperature fault (lower stack temperature) -DRH ↗ : dew temperature fault (higher dew temperature) -DRH ↘ : dew temperature fault (: : ; ≡ #$% & , : < ≡ )* ↗ + , : = ≡ )* ↘ + , : > ≡ DFSA ↗ + , : ? ≡ DFSA↘ + , : @ ≡ DT ↗ + , : A ≡ DT↘ + , : B ≡ DRH ↗ + , : C ≡ DRH↘ + . 20 As we can see, the proposed diagnosis strategy identifies 21 successfully several complex operating faults (i.e. slight 22 deflections from the nominal operating conditions, and even 23 combination of faults) for both PEMFCAuto and PEMFCµCHP.

6

  i.e. more or less severe, deviations from the normal conditions) 29 with considering situations that combine 2 or 3 faults 30 simultaneously. Two stacks are investigated: a PEMFCAuto 31 stack that is designed for automotive applications and a 32 PEMFCµCHP stack, dedicated to stationary ones. 33 Singularity features extracted from the estimated VSS are 34 classified using Machine Learning approach, named k-Nearest stands, for different stack sizes, powers and technologies, with 5 different targeted environments of power application. In this paper, we demonstrate that the singularity analysis of 7 voltage signal offers a generic diagnosis tool for the PEMFC 8 SoH monitoring. One key-point of the method is that the VSS 9 can be estimated from the "free" evolution of the stack voltage 10 and without affecting in any way the FC operation: no external 11 additional solicitation is required to reveal the SoH patterns.

  was done thanks to the French ANR 15 project "DIAPASON 2"and in the framework of the 16 "Decentralized energy production" project, directed by 17 EFFICACITY, the French R&D Institute for urban energy 18 transition. 19 and Lottin O., One-dimensional Model of Oxygen Transport 28 Impedance Accounting for Convection Perpendicular to the 29 Electrode. Fuel Cells 12, 5, 848-861, 2012.

30[ 4 ]

 4 Chevalier, S., Auvity, B., Olivier, J.C., Josset, C., Trichet D. and 31 Machmoum M., Multiphysics dc and ac models of a pemfc for 32 the detection of degraded cell parameters. International Journal 33 of Hydrogen Energy, 38, 26, 11609-11618, 2013.

34[ 5 ]

 5 Petrone, R., Zheng, Z., Hissel, D., Péra, M.C., Pianese, C., 35 Sorrentino, M., Becherif, M., et al, A review on model-based 36 diagnosis methodologies for PEMFCs. International Journal of 37 Hydrogen Energy, 38, 7077-7091, 2013.

38[ 6 ]

 6 Zheng, Z., Petrone. R., Péra, M.C., Hissel, D., Becherif, M., 39 Pianese, C., et al. A review on non-model based diagnosis 40 methodologies for PEM fuel cell stacks and systems.

1 Synopsis of the investigated PEM Fuel Cells 2

	3	In our study, two PEMFC stacks are experimented to evaluate
	4	the portability of the proposed diagnosis tool. The first one is
	5	an 8 cells stack designed for automotive applications and
	6	manufactured by CEA (Alternative Energies and Atomic
	7	Energy Commission). The second one is a 12 cells stack
	8	dedicated to stationary application (micro combined heat and
	9	power -µCHP application). It is designed and marketed by
	10	Riesaer Brennstoffzellentechnik GmbH and Inhouse
	11	Engineering GmbH, Germany.
	12	The first stack (PEMFCAuto) is made of metallic gas distributor
	13	plates. The electrode active surface of a cell is equal to 220
	14	cm 2 . It is fed by air at cathode and pure hydrogen (H2) at anode.
	15	A summary of the FC nominal operating conditions is given in

Table 1 .

 1 The stack operates with a nominal current of 110 A. dioxide -CO2) simulating a reformat. It is made of graphite gas 20 distributor plates. The electrode's active surface is 196 cm 2 . The stack operates with a nominal current of 80 A. A picture 22 of the stack is given in Fig. 1(b). A summary of the FC

	21

nominal 23 operating parameters and other main characteristics are given 24 inTable 2 .

 2

Table 1 .

 1 PEMFCAuto nominal operating conditions.

	30
	Coolant flow: deionized water	2 l/min
	Anode stoichiometry rate (H2)	1.5
	Cathode stoichiometry rate (air)	2
	Absolute pressure for H2 and air inlets	150 kPa
	Max. anode -cathode pressure gap	30 kPa
	Temperature of the cooling circuit	80°C

Table 2 .

 2 The PEMFCµCHP nominal operating conditions.

	31
		Coolant flow: deionized water	3 l/min
		Anode stoichiometry rate (H2 & CO2 mix)	1.3
		Cathode stoichiometry rate (air)	2
		Absolute pressure for H2 and air inlets	150 kPa
		Max. anode -cathode pressure gap	20 kPa
		Temperature of the cooling circuit	75°C
		Anode and cathode relative humidity rate	50 %
		Nominal Current	80 A
	32
	33	Experimental process

The above described PEMFCs were experimented with 34 testbenches developed in the FC platform of Belfort, France. It 35 includes mainly: 36 • a complete gas conditioning sub-system including gas 37 humidifiers at anode and cathode, 38 • a test stand section dedicated to the control of the 39 temperature inside the stack and including the FC primary 40 water circuit.

Table 3 .

 3 The set of FC operating conditions applied during the experimentation of the two stacks. The underlying parameter values

	Parameters	Ref	DFSC	DFSA	DP	DT	DH	DCO
	value

conclusions are given in Section 5.

The second stack (PEMFCµCHP) is fed by air at cathode, and at

anode by a fuel mixture (75 % of H2 and 25 % of carbone

correspond to the introducing faults. The notation 'ND' means that the FC experiment is not done.

In this study, VSS are computed on voltage signals acquired

under different FC operating conditions (normal and abnormal,

Neighbors (kNN).

International Journal of Hydrogen Energy, 38, 21, 8914-8926.


Alexandr Klimchik
email: alexandr.klimchik@mines-nantes.fr

Stéphane Caro

Anatol Pashkevich

Optimal pose selection for calibration of planar anthropomorphic manipulators

The paper is devoted to the calibration experiment design for serial anthropomorphic manipulators with arbitrary number of links. It proposes simple rules for the selection of manipulator configurations that allow the user to essentially improve calibration accuracy and reduce identification errors. Although the main results have been obtained for the planar manipulators, they can be also useful for calibration of more complicated mechanisms. The efficiency of the proposed approach is illustrated with several examples that deal with typical planar manipulators and an anthropomorphic industrial robot.

Introduction

The standard engineering practice in industrial robotics assumes that the closed-loop control technique is applied only on the level of servo-drives, i.e. for actuating the manipulator joints. In such systems, the Cartesian space control is based on the open-loop method that incorporates numerous direct/inverse kinematic transformations derived from the manipulator geometric model. These transformations define correspondence between the manipulator joint coordinates and the Cartesian coordinates of the end-effector. Hence, to achieve desired accuracy, manipulator geometric model employed in the control algorithm should be carefully tuned (calibrated) to take into account manufacturing tolerances and parameter variations from manipulator to manipulator [START_REF] Roth | An overview of robot calibration[END_REF].

The problem of robot calibration has already been well studied and it has been in the focus of research community for many years [START_REF] Elatta | An overview of robot calibration[END_REF]. In general, the calibration process is divided into four sequential steps [START_REF] Mooring | Fundamentals of manipulator calibration[END_REF]: modeling, measurements, identification and compensation. The first step focuses on the design of the appropriate (complete but non-redundant [START_REF] Pashkevich | Computer-aided generation of complete irreducible models for robotic manipulators[END_REF]) mathematical model. At the second step, related measurements (calibration experiments) are carried out using commercially available or custom-made equipments [START_REF] Nubiola | Absolute robot calibration with a single telescoping ballbar[END_REF][START_REF] Santolaria | Self-alignment of on-board measurement sensors for robot kinematic calibration[END_REF]. The third step usually deals with the identification of the Denavit-Hartenberg parameters [START_REF] Hollerbach | A survey of kinematic calibration[END_REF], which may provoke numerical instabilities for the manipulators with collinear successive axes considered in this paper. For this particular but very common case, some authors (Hayati [8], Stone [START_REF] Stone | Kinematic modeling, identification, and control of robotic manipulators[END_REF], and Zhuang [START_REF] Zhuang | A complete and parametrically continuous kinematic model for robot manipulators[END_REF]) proposed some modifications, but here we will use a more straightforward approach that is more efficient for the planar manipulators. The last step is aimed at compensating identified parameter variations [START_REF] Chen | Spatial statistical analysis and compensation of machining errors for complex surfaces[END_REF][START_REF] Klimchik | Compliance error compensation technique for parallel robots composed of non-perfect serial chains[END_REF][START_REF] Shi | Kinematic modeling and calibration of a flexure based hexapod nanopositioner[END_REF].

Among numerous publications devoted to the robot calibration, there is a very limited number of works that directly addresses the issue of the identification accuracy and reduction of the calibration errors. In particular, Ikits and Hollerbach [START_REF] Ikits | Kinematic calibration using a plane constraint[END_REF] used noise amplification index to estimate the errors in the identified parameters of Puma 560 robot. Mirman and Gupta [START_REF] Mirman | Compensation of Robot joint variables using special Jacobian matrices[END_REF] proposed compensation algorithm using position-independent parameter error values. In [START_REF] Slamani | Modeling and assessment of the backlash error of an industrial robot[END_REF] the authors assessed backlash error for an ABB IRB 1600 6-dof serial industrial robot. Five different observability indexes were compared in [START_REF] Sun | Observability index selection for robot calibration[END_REF] and the authors detected that all of them are related to each other. In [START_REF] Sun | Active robot calibration algorithm[END_REF] the determinant-based observability index was used to evaluate the performance of active robot calibration algorithm applied to a 6-dof PUMA 560 robot. In further comparison study, Hollerbach et al. [START_REF] Hollerbach | Model identification[END_REF] proposed to treat all calibration methods as closed-loop ones and introduced the calibration index that categorizes all calibration methods in terms of number of equations per pose. Zhuang et al. [START_REF] Zhuang | Optimal planning of robot calibration experiments by genetic algorithms[END_REF] used the condition number of the identification Jacobian to compare the identification accuracy impact of different measurement configurations. In [START_REF] Daney | Choosing measurement poses for robot calibration with the local convergence method and Tabu search[END_REF] the authors used iterative one-by-one pose search algorithm in order to minimize the influence of measurement noise on the identification of geometrical parameters of a Gough-Stewart Platform.

It is clear that the calibration accuracy may be straightforwardly improved by increasing the number of experiments (with the factor 1/ √ m, where m is the number of experiments [START_REF] Atkinson | Optimum experimental designs[END_REF]). Besides, using diverse manipulator configurations for different experiments looks also intuitively promising and perfectly corresponds to some basic ideas of the classical theory of experiment design [START_REF] Hollerbach | Model identification[END_REF] that intends to use experiments that are as much distinct as possible. However, the classical results are mostly obtained for very specific models (such as the linear regression) and cannot be applied here directly due to the non-linearity of the relevant equations. There are few research works in the literature that deal with the optimal pose selection for robot calibration. For example, Borm and Menq [START_REF] Borm | Determination of optimal measurement configurations for robot calibration based on observability measure[END_REF] have investigated the implications of different observability measures in the robot position error and influence of the measurement configurations number on the final accuracy. They concluded that the number of measurements is less important than proper selection of measurement configurations. In [START_REF] Khalil | Identifiable parameters and optimum configurations for robots calibration[END_REF] the authors defined the set of optimal measurement configurations by minimizing the condition number of the observation matrix. Daney [START_REF] Daney | Optimal measurement configurations for Gough platform calibration[END_REF] used the constrained optimization algorithm based on the minimization of the singular values root-product for optimal measurement configurations selection for Gough platform calibration. The noise amplification index was used in [START_REF] Nahvi | The noise amplification index for optimal pose selection in robot calibration[END_REF] to quantify measurement configurations and to select the best one. In [START_REF] Klimchik | Design of experiments for calibration of planar anthropomorphic manipulators[END_REF] authors used D-optimality criteria to determine optimal measurement configurations for planar anthropomorphic manipulators. Zhuang et al. [START_REF] Zhuang | Optimal selection of measurement configurations for robot calibration using simulated annealing[END_REF] applied simulating annealing to obtain optimal or near optimal measurement configurations, which minimize at least one of two considered performance measures. Imoto et al. [START_REF] Imoto | Optimal kinematic calibration of robots based on maximum positioning-error estimation (theory and application to a parallel-mechanism pipe bender)[END_REF] proposed to use the end-effector position accuracy after calibration as a performance measure in order to generate measurement configurations. Similar idea was used in [START_REF] Klimchik | Optimal selection of measurement configurations for stiffness model calibration of anthropomorphic manipulators[END_REF] where the authors introduced test configurations related to the technological process, which allowed them to define the performance measure as the positioning accuracy after calibration that is also related to the weighted trace of the covariance matrix.

As follows from detailed analysis, all previous works in the area of calibration experiment design provide user with an iterative scheme that aims at minimizing an objective function that depends on the singular values of the identification Jacobian (condition number, for instance). However, this approach does not consider directly the identification accuracy and may lead to some unexpected results (where the condition number is perfect, but the parameter estimation errors are rather high). Besides, it requires very intensive and time consuming computations caused by a poor convergence and high dimension of the search space (number of calibration experiments multiplied by the manipulator joint number). Hence, to apply this technique in industry, strong mathematical background and good experience in the numerical optimization are required. It is obvious that practical engineers need some type of a "rule of thumb", which allows them to select measurement configurations without tedious computations.

In this paper, the problem of optimal design of the calibration experiments is studied for the case of a planar manipulator with arbitrary number of links. Such manipulators do not cover all architectures used in practice, but nevertheless this model allows us to derive some very useful analytical expressions and to propose some simple practical rules defining optimal configurations with respect to the calibration accuracy. In the following sections, particular attention will be given to planar manipulators with 2 and 3 d.o.f. that are essential components of all existing anthropomorphic robots. Practical significance of the obtained results will be illustrated by a case study that deals with the calibration experiment design for a 6-d.o.f. KUKA industrial robot, which is presented as a set of simple planar sub-manipulators.

The remainder of this paper is organized as follows. Section 2 defines the research problem and contains basic assumptions. Section 3 presents a motivation example that shows the importance of measurement pose selection in robot calibration. In Section 4, the identification algorithm is presented. Section 5 deals with the evaluation of the identification accuracy. Section 6 contains the main theoretical contributions that allow the user to generate desired measurement configurations without straightforward numerical optimization, using proposed rule of thumb. Sections 7 and 8 illustrate advantages of the developed approach and contain some simulation results. In Sections 9 and 10, the proposed technique is extended to the case of spatial manipulator and is applied to a 6-dof serial industrial robot. Section 11 contains discussion where weak and strong sides of the developed approach are considered. Finally, Section 12 summarizes the main contributions of the paper.

Problem statement

Let us consider a general planar serial manipulator consisting of n rigid links connected by the corresponding number of revolute joints. For this manipulator, the end-effector position (x, y) can be defined as follows:

x = l 1 cos q 1 + l 2 cos (q 1 + q 2 ) + • • • + l n cos (q 1 + q 2 + • • • + q n) y = l 1 sin q 1 + l 2 sin (q 1 + q 2 ) + • • • + l n sin (q 1 + q 2 + • • • + q n) (1)
where l 1 , l 2 , . . ., l n are the link lengths, q 1 , q 2 , . . ., q n are the actuated joint coordinates, n is the number of links. In practice, the actual values of the link length l i and the joint coordinates q i differ from the nominal ones l i 0 and q 0 i by some offsets l i and q i to be identified:

l i = l 0 i + l i ; q i = q 0 i + q i (2)
For further convenience, let us introduce the notations

 0 i = i k=1 q 0 k ;  i = i k=1 q k (3)
that allow us to rewrite (1) as

x = l 0 1 + l1 • cos  0 1 + Â1 + • • • + l 0 n + ln • cos  0 n + Ân y = l 0 1 + l1 • sin  0 1 + Â1 + • • • + l 0 n + ln • sin  0 n + Ân (4)
Below, the system (4) will be used to generate the set of calibration equations where the offset variables { l i , i = 1, n} and { Â i , i = 1, n} are treated as unknowns.

To find the desired offsets, a number of experiments are carried out providing a set of Cartesian coordinates x k , y k and corresponding joint angles q k 1 , q k 2 , . . ., q k n that theoretically satisfy the system of Eq. ( 4). However, due to measurement errors, the number of experiments should be excessive and the set of the calibration equations cannot be satisfied simultaneously. Hence, the identification procedure may be treated as the best fitting of the experimental data by the geometrical model (2), i.e. by minimizing the corresponding positional residuals.

To take into account the impact of the measurement noise, the calibration equations derived from (4) can be written in the following form:

x k = n i=1 l 0 i + l i • cos  0k i +  k i + ε k x y k = n i=1 l 0 i + l i • sin  0k i +  k i + ε k y ( 5
)
where k = 1, m, the subscript 'i' defines the manipulator joint/link number, the superscript 'k' indicates the experiment number, the superscript '0' denotes the nominal values of the corresponding variables, m is the total number of experiments. Besides, it is assumed that the measurement errors ε k x , ε k x are independent identically distributed random values (iid) with zero mean and the standard deviation . It should be noted that, for some technical reasons, the measurement errors in the joint variables q i are neglected because of their small impact [START_REF] Bernhardt | Robot calibration[END_REF].

In the frame of this model, the desired offsets { l i , i = 1, n} and { Â i , i = 1, n} computed from the system (5) can be also treated as random variables. Following the main ideas of the statistical estimation theory, it is natural to assume that the mean values of the obtained offsets are equal to their actual values (unbiased estimate) and their standard deviations are as small as possible, i.e.

E l i = l i ; Var l i → min; i = 1, n E q i = q i ; Var q i → min; i = 1, n (6)
where the hat symbol ' ' denotes the estimates of the corresponding parameters, E(.) and Var(.) are the expectation and variance, respectively. This generally leads to a multiobjective optimization problem whose solution highly depends on the selection of the manipulator joint coordinates q k 1 , q k 2 , . . ., q k n used in the calibration experiments. So, the main problem addressed in this paper is to find the set of manipulator configurations (defined by the joint coordinates) that insure the lowest dispersions of l i and q i or certain compromise between them. To make these results attractive for the practicing engineers, the main attention will be paid to the development of simple techniques ("rule of thumb") allowing us to generate the desired configurations without tedious computations, which are required while applying other existing methods (which involve time consuming numerical optimization in multidimensional space).

To develop the desired techniques, let us sequentially describe the identification algorithm, evaluate the related identification errors (via the covariance matrix) and develop some optimality conditions (imposed on the joint variables) that guarantee a diagonal structure of the covariance matrix with the smallest elements describing the parameter variances. The latter allows us to reformulate the conventional optimality criteria used in the literature and to present the final solution in a simple geometrical form.

Motivation example

To demonstrate importance of optimal measurement poses selection in robot calibration, let us present first a simple motivation example. For the purpose of simplicity, let us limit our study to a 2-dof planar manipulator with two revolute actuated joints. For this manipulator, the Cartesian coordinates x, y of the end-effector can be computed as:

x = l 1 • cos q 1 + l 2 • cos (q 1 + q 2 ) y = l 1 • sin q 1 + l 2 • sin (q 1 + q 2 ) (7)
where l 1 , l 2 are the link lengths, q 1 , q 2 are the actuated joint coordinates. It should be noted that usually real robot parameters l 1 , l 2 and q 1 , q 2 differ from their nominal values and this manipulator has four geometrical parameters to be identified: (i) two link length deviations l 1 , l 2 , and (ii) two joint encoder off-sets q 1 , q 2 . In order to show the significance of proper measurement configurations selection, let us examine two plans of calibration experiments: Plan (i)

q 1 = 30 o , q 2 -100 o ; q 1 = -150 o , q 2 -90 o Plan (ii) q 1 = 30 o , q 2 -90 o ; q 1 = 30 o , q 2 = 90 o
The first plan has been selected intuitively, while the second one has been generated using results obtained in this work. It should be mentioned that for the considered manipulator each measurement provides two Cartesian coordinates x i , y i . So, the minimum number of measurement configurations is equal to two, it insures the identifiability of all desired parameters l 1 , l 2 and q 1 , q 2 . Using these measurement configurations, simulation of calibration experiments has been done assuming that the measurement noise is Gaussian with zero mean and dispersion = 0.1 mm. Simulation results for l 1 = 0.6 m, l 2 = 0.4 m are summarized in Fig. 1, where the root-mean-square errors of the end-effector position after calibration throughout the robot workspace are presented. The figure clearly shows that, in the worst manipulator configuration, the intuitive Plan (i) provides a positioning error equal to 2.29 mm while the proposed Plan (ii) reduces the worst positioning error down to 0.14 mm, i.e. by a factor of 16. Hence, this simple example clearly shows that selection of measurement configurations is a very important issue in robot calibration. In fact, poorly chosen measurement configurations may have the negative effect and reduce the robot accuracy after calibration. This motivates careful planning of the calibration experiments in order to increase calibration efficiency, which is in the focus of the paper. In the frame of this paper, the calibration efficiency refers to the ability to reduce the impact of measurement errors on the robot accuracy without increasing the number of experiments.

Identification algorithm

To find the desired parameters using the noise corrupted measurements, the least square technique is usually applied. This approach aims at minimizing the square sum of the residuals in Eq. ( 5), which is the simplest way to minimize all residuals simultaneously. For the considered problem, the least square objective can be written as:

m k=1 n i=1 l 0 i + l i • cos  0k i +  i -x k 2 + m k=1 n i=1 l 0 i + l i • sin  0k i +  0 i -y k 2 → min (8)
where x k , y k and  0k 1 ,  0k 2 , . . .,  0k n denote, respectively, the end-effector Cartesian coordinates and the manipulator joint angles for the kth experiment, l 0 i is the nominal length of the ith link, while l i and  i are the offsets to be identified.

Since the unknown variables  i , i = 1, n are incorporated in the trigonometric functions, the optimization problem (8) cannot be solved analytically. On the other hand, in practice the offsets l i ,  i are relatively small. For this reason, the geometrical model ( 4) can be linearized with respect to the parameters to be identified and presented in a convenient matrix form as

P k = P 0k + J k • (9)
where the vector =  1 , . . .,  n , l 1 , . . ., l n T aggregates all unknown parameters, P k = (x k , y k ) T is the end-effector position obtained from the kth experiment, P 0k = (x 0k , y 0k ) T is the corresponding end-effector position computed using the nominal values of the link lengths l 0 i , i = 1, n and given coordinates  0k

1 , Â 0k 2 , . . ., Â 0k n defining the manipulator configuration for the kth experiment, i.e.

x 0k = n i=1 l 0 i cos  0k i ; y 0k = n i=1 l 0 i sin  0k i , k = 1, m, (10)
and J k is the calibration Jacobian, which can be computed by differentiating the system (4) with respect to . After relevant computations the latter can be presented in the form [START_REF] Chen | Spatial statistical analysis and compensation of machining errors for complex surfaces[END_REF] where

J k = J k x J k xl J k y J k yl 2×2n
J k x = -l 0 1 sin  0k 1 . . . l 0 n sin  0k n 1×n J k xl = cos  0k 1 . . . cos  0k n 1×n J k y = l 0 1 cos  0k 1 . . . l 0 n cos  0k n 1×n J k yl = sin  0k 1 . . . sin  0k n 1×n (12)
Using the above notation and defining P k = P k -P 0k , the leastsquare objective (8) can be rewritten as

m k=1 J k • -P k T J k • -P k → min (13)
which leads to an analytical solution

= m k=1 J k T J k -1 • m k=1 J k T P k (14)
that will be further used for the accuracy analysis. It should be noted that in the case of relatively large l i and  i , the linearization procedure can be applied several times with iterative modifications of the desired geometrical parameters and the configuration variables

l 0 i ← l 0 i + l i ; Â 0k i ← Â 0k i + Â i ; ∀i, k (15)
This allows us to find a numerical solution of the original optimization problem [START_REF] Hayati | Robot arm geometric link parameter estimation[END_REF], which cannot be solved analytically.

It should be mentioned that there exist different versions of the above algorithm in the literature, which differ in the type of original measurement data (Cartesian coordinates, distance to the reference point/line, end-effector orientation, etc.) and numerical optimization techniques (gradient search, simulated annealing, genetic algorithms, etc.). Detailed reviews on these issues can be found in [START_REF] Le | Joint calibration of multiple sensors[END_REF][START_REF] Beyer | Practical robot calibration with ROSY[END_REF][START_REF] Liu | Geometric parameter identification of a 6-DOF space robot using a laser-ranger[END_REF][START_REF] Park | Laser-based kinematic calibration of robot manipulator using differential kinematics[END_REF][START_REF] Liu | Improving the positioning accuracy of a neurosurgical robot system[END_REF][START_REF] Santolaria | Kinematic parameter estimation technique for calibration and repeatability improvement of articulated arm coordinate measuring machines[END_REF][START_REF] Lightcap | Improved positioning accuracy of the PA10-6CE robot with geometric and flexibility calibration[END_REF][START_REF] Santolaria | Laser tracker-based kinematic parameter calibration of industrial robots by improved CPA method and active retroreflector[END_REF][START_REF] Williams | Kinematic error calibration of the gantrytau parallel manipulator[END_REF]. In particular, Le et al. [START_REF] Le | Joint calibration of multiple sensors[END_REF] grouped sensors to produce 3D data and variety of geometric constraints that allowed them to calibrate all parameters simultaneously. In [START_REF] Beyer | Practical robot calibration with ROSY[END_REF], to increase the accuracy of serial and parallel robots, the authors developed the Robot Optimization System (ROSY) which uses a dedicated measuring tool with two digital CCD cameras. An alternative geometric parameter identification method that is based on a laser-ranger attached to the end-effector was used in [START_REF] Liu | Geometric parameter identification of a 6-DOF space robot using a laser-ranger[END_REF]. Kinematic calibration based on differential techniques that uses measurement data from structured laser module and stationary camera was used in [START_REF] Park | Laser-based kinematic calibration of robot manipulator using differential kinematics[END_REF] to estimate parameters of 7-DOF humanoid manipulator arm. In [START_REF] Liu | Improving the positioning accuracy of a neurosurgical robot system[END_REF], the authors applied a backpropagation neural network to compensate the joint errors of the neurosurgical robot system. Santolaria et al. [START_REF] Santolaria | Kinematic parameter estimation technique for calibration and repeatability improvement of articulated arm coordinate measuring machines[END_REF] utilized the ball bar gauge measurement device and identification technique whose objective function includes terms that are regarding repeatability and measurement accuracy. To improve the static positioning accuracy of the PA10-6CE robot, in [START_REF] Lightcap | Improved positioning accuracy of the PA10-6CE robot with geometric and flexibility calibration[END_REF] a 30parameter model incorporating elastostatic ones has been used. The screw measurement method that utilizes laser tracker with an active target and identification procedure that is based on circle point analysis was used in [START_REF] Santolaria | Laser tracker-based kinematic parameter calibration of industrial robots by improved CPA method and active retroreflector[END_REF] to identify the kinematic parameters of an industrial robot. In [START_REF] Williams | Kinematic error calibration of the gantrytau parallel manipulator[END_REF], the authors presented a 12-parameter error kinematic model for the three linear actuators of the Gantry-Tau robot and used three types of measurement equipments (laser interferometers, linear encoders and double-ball bars) to calibrate the linear actuators.

It is also worth mentioning that in practice, there exist non-negligible measurement noise, which is incorporated in the variables P k indicating the difference between the actual end-effector position (obtained from measurements) and those computed by using the nominal geometrical model. Hence, for the efficient calibration experiment planning, it is reasonable to evaluate the measurement noise impact on the desired manipulator parameters.

Identification accuracy

To estimate the measurement noise impact on the identification results, let us assume that the manipulator parameters differ from the nominal ones by and the measurements include the additive random errors ε k = ε k x , ε k y T (independent identically distributed with zero mean and standard deviation ). This allows us to present the measurement data as

P k = P 0k + J k • + k (16)
and, after substitution into [START_REF] Ikits | Kinematic calibration using a plane constraint[END_REF], to express the identified values as

= + m k=1 J k T J k -1 • m k=1 J k T k (17)
where the second term is the stochastic component in the model parameter estimation. As follows from this expression, the considered identification algorithm provides the unbiased estimates

E ( ) = (18)
where E (.) denotes the mathematical expectation. Besides, the corresponding covariance matrix defining the identification accuracy, can be expressed as

cov( ) = m k=1 J k T J k -1 E m k=1 J k T k • m k=1 k T J k m k=1 J k T J k -1 (19)
Further, taking into account the statistical properties of the measurement errors (see assumption at the beginning of this section)

E ε k ε jT = 0 2×2 , k / = j
the desired covariance matrix [START_REF] Hollerbach | Model identification[END_REF] can be presented as

cov ( ) = 2 m k=1 J k T J k -1 (21)
Therefore, the impact of the measurement errors (i.e. the "quality" of the selected measurement configurations q k 1 , q k 2 , . . ., q k n ) is completely defined by the matrix sum m k=1 J k T J k , which is constructed from the identification Jacobians.

Improvement of the identification accuracy via design of calibration experiments

To analyze the above derived matrix [START_REF] Daney | Choosing measurement poses for robot calibration with the local convergence method and Tabu search[END_REF] in detail, let us take into account that the manipulator model ( 4) includes two types of identified parameters ( l i and  i ). So, the sum m k=1 J k T J k can be presented in the block-wise form as

m k=1 J k T J k = L T C L L T S S T L C (22)
where L = diag(l 0 1 , l 0 2 , . . ., l 0 n ), while the matrices C, S depend on the measurement configurations in the following way:

C = m k=1 cos  k i - k j i = 1, n j = 1, n ; S = m k=1 sin  k i - k j i = 1, n j = 1, n (23)
where

 k i = i s=1 q k s .
It should be mentioned that the diagonal terms of the sub-matrices C, S can be computed straightforwardly

c ii = m; s ii = 0; ∀i = 1, n. (24)
and do not depend on the measurement configurations. However, the off-diagonal elements highly depend on the selection of the relevant joint angles q k 1 , q k 2 , . . ., q k n . It is clear that, using terminology from the classical design of experiments, the sum m k=1 J k T J k can be treated as the "information matrix," which may be optimized using several different criteria [START_REF] Atkinson | Optimum experimental designs[END_REF]. The most common ones are so-called A-and D-optimality criteria, which allow us to maximize the matrix trace and determinant, respectively. However, the first criterion cannot be applied to the considered problem because the trace of the information matrix [START_REF] Atkinson | Optimum experimental designs[END_REF] does not depend on the design variables q k 1 , q k 2 , . . ., q k n that define the plan of the calibration experiment. In fact, as follows from expressions ( 22)-( 24), here the trace is constant and is equal to

trace m k=1 J k T J k = m × n + m n i=1 l 0 i 2 (25)
Besides, this trace is composed of the elements of different units, which is not acceptable from an engineering view point.

In contrast, the second criterion leads to the optimization problem that produces reasonable and practically acceptable solutions. In this case, using the block-wise formulas of Schur [START_REF] Gantmakher | The theory of matrices[END_REF], the desired determinant can be expressed as

det m k=1 J k T J k = (det L) 2 • det C • C T -S • S T (26)
where det (L) = const. As follows from the relevant analysis (see Appendix A), the latter expression reaches its maximum value when all elements of the matrix S are equal to zero. This allows us to reformulate the considered problem in the following way

q k 1 , q k 2 , . . ., q k n |k = 1, m = argmax q det C (27)
subject to the equality constraint S = 0. It should be noted that det C > 0 here. Summarizing these results (see Appendix A for details), the optimality conditions of the desired calibration plan can be written as the set of n (n -1) trigonometric equations

m k=1 cos i s=1 q k s - j s=1 q k s = 0, ∀ i > j m k=1 sin i s=1 q k s - j s=1 q k s = 0, ∀ i > j . ( 28
)
in m • (n -1) variables q k 2 , . . ., q k n |k = 1, m , where i = 1, n and j = i + 1, n. It is worth mentioning that the joint variables q k 1 do not appear in Eq. ( 28), so the calibrations plans are invariant with respect to the coordinates of the first joint. If those optimality conditions are satisfied, the information matrix will be strictly diagonal since S = 0 and

C = m I m k=1 J k T J k = m L 2 0 0 I (29)
and, as follows from [START_REF] Daney | Choosing measurement poses for robot calibration with the local convergence method and Tabu search[END_REF], it ensures the minimum values of the variances Var ( l i ) and Var  i . Hence, to find the desired set of optimal measurement configurations, it is necessary to solve the underdetermined system of trigonometric Eq. [START_REF] Zhuang | Optimal selection of measurement configurations for robot calibration using simulated annealing[END_REF], where the number of variables m • (n -1) is higher than the number of equations n (n -1). As a matter fact, m is usually greater than n because the original identification problem (see Section 3) deals with estimation of 2n parameters l i ,  i using 2m Eq. ( 4). So, the minimum number of experiments m ensuring the parameters identifiability is equal to n. However, in practice m > n, because it allows reducing the measurement noise impact.

It can be easily proved that the underdetermined system of trigonometric Eq. ( 28) has an infinite number of solutions. For instance, even in the case of n = 2 and unreasonably small number of measurements m = 2, the complete set of solutions can be expressed as

q 1 1 = ˛1; q 1 2 = q2 1 = ˛2; q 2 2 = ˇ+ (30)
where ˛1, ˛2, ˇare arbitrary angles within the manipulator joint limits. In more realistic case, when n = 2 and m = 3, the desired set of solutions

q 1 1 = ˛1; q 1 2 = q2 1 = ˛2; q 2 2 = ˇ+ 2 3
q 3 1 = ˛3; q 3 2 = ˇ-2 3 ( 31
)
is defined via four arbitrary angles ˛1, ˛2, ˛3, ˇ. More details concerning the solutions of the considered system for different values of m and n are presented in Table 1.

It is worth mentioning that in practice, instead of generating the complete set of solutions of the system [START_REF] Zhuang | Optimal selection of measurement configurations for robot calibration using simulated annealing[END_REF], it is enough to find an appropriate partial one (since all of them are equivalent with respect to the identification accuracy and provide the same covariance matrix). To find this solution, it useful to take advantage of the additive structure of the equations defining the optimality conditions for the measurement plans (see Table 1). In fact, these equations can be easily separated into several independent groups

Table 1

Optimality conditions for measurement pose selection: calibration of typical planar manipulators.

Manipulator Equations defining optimal measurement poses 2-Links manipulator

m k=1 cos q k 2 = 0; m k=1 sin q k 2 = 0, 3-Links manipulator m k=1 cos q k 2 = 0; m k=1 sin q k 2 = 0; m k=1 cos q k 3 = 0; m k=1 sin q k 3 = 0; m k=1 cos q k 2 + q k 3 = 0; m k=1 sin q k 2 + q k 3 = 0 4-Links manipulator m k=1 cos q k 2 = 0; m k=1 sin q k 2 = 0; m k=1 cos q k 3 = 0; m k=1 sin q k 3 = 0; m k=1 cos q k 2 + q k 3 = 0; m k=1 sin q k 2 + q k 3 = 0 m k=1 cos q k 4 = 0; m k=1 sin q k 4 = 0; m k=1 cos q k 3 + q k 4 = 0; m k=1 sin q k 3 + q k 4 = 0; m k=1 cos q k 2 + q k 3 + q k 4 = 0; m k=1 sin q k 2 + q k 3 + q k 4 = 0
comprising non-identical variables. For instance, in the simplest case n = 2, the optimality conditions

m k=1 cos q k 2 = 0; m k=1 sin q k 2 = 0 (32)
can be rewritten as

m 1 k=1 cos q k 2 + m k=m 1 +1 cos q k 2 = 0 m 1 k=1 sin q k 2 + m k=m 1 +1 sin q k 2 = 0 (33)
and replaced by stronger conditions

m 1 k=1 cos q k 2 = 0; m k=m 1 +1 cos q k 2 = 0 m 1 k=1 sin q k 2 = 0; m k=m 1 +1 sin q k 2 = 0 ( 34
)
which is obviously sufficient but not necessary with respect to [START_REF] Le | Joint calibration of multiple sensors[END_REF]. Nevertheless, the desired information matrix ( 29) is achieved and the identification accuracy does not suffer from this intuitive simplification. It is clear that similar approach can be also applied in general case (for n≥2).

The above presented idea allows us to propose a heuristic approach for obtaining partial solutions of (28) for any given m using general solutions for the problems of lower dimensions. Below, this techniques will be referred to as the "superposition of low-dimensional plans". In practice, the most attractive and sufficient (for our problem) is the decomposition of the sums (28) into several sub-sums of size m 1 = 2 and m 2 = 3, which is equivalent to the presentation of experiment number m in the following form m = 2k 1 + 3k 2 [START_REF] Park | Laser-based kinematic calibration of robot manipulator using differential kinematics[END_REF] where k 1 , k 2 are some integers. It is obvious that the partitioning [START_REF] Park | Laser-based kinematic calibration of robot manipulator using differential kinematics[END_REF] is not unique (and exists if m≥2), but any of possible presentations is suitable here. For example, for n = 2 the desired solution can be easily composed of k 1 measurement configurations of the type (30) and k 2 configurations of the type [START_REF] Bernhardt | Robot calibration[END_REF]. Besides, in some cases it is possible to compose optimal plans using one type of configurations only: either type (30) (i.e. k 2 = 0 and k 1 = m/2) or type (31) (i.e. k 1 = 0 and k 2 = m/3). In the case when m is a multiple of 6, both mentioned plans of experiments can be used. For instance, for m = 6 the plans 2 × 3 and 3 × 2 are acceptable.

To evaluate visually the diversity of the obtained measurement configurations, Table 2 contains several case studies corresponding to 2-, 3-and 4-link manipulators. As follows from them, the optimal configurations differ essentially from each other, in order to ensure both parameters identifiability and low identification errors. This perfectly suits the basic ideas of the design of experiment theory.

In general, the proposed approach can be treated as a simple rule of thumb allowing to generate easily optimal poses for calibration of planar anthropomorphic manipulators. Applying this rule, it is possible to obtain the measurement configurations that provide the following covariance matrix for the parameters of interest l i , Â i |i = 1, n (as well as for the parameters

l i , q i |i = 1, n ) cov ( ) = 2 m • •L -2 0 0 I (36)
Table 2 Optimal measurement configurations for typical planar manipulators.

which allows to estimate the identification accuracy via the variances

Var ( q i ) = 2 m • l 2 i ; Var ( l i ) = 2 m ( 37
)
where is the s.t.d. of the measurement noise and i = 1, n. These results show that, if the optimal calibration plan is used, the identification errors for the linear parameters l i will depend on the number of experiments only, while the errors for the angular parameter q i will also depend on the link lengths. In more details, the results related to the identification accuracy are presented in Table 3. Hence, using the proposed rule of thumb one can easily generate the desired measurement configurations corresponding to the D-optimality criteria, which yields essential reduction of the identification errors. It is proved that, from the algebraic point of view, these configurations satisfy the optimality conditions, which require relevant sums of sin(.) and cos(.) of the joint coordinates to be equal to zero. Using geometrical approach, this rule can be represented on the unit circle, where the sum of corresponding unit vectors must be zero. The latter allows user to find the desired configurations easily. The advantages of this approach are illustrated in the following Sections.

Simulation study

To demonstrate the efficiency of the proposed technique, let us present some simulation results that deal with two-, three-and four-links manipulators and employ different calibration plans. The first plan (a conventional one) is based on the regular grid in the joint coordinate space. The second calibration plan implements the proposed strategy, which is described in detail in Section 6. The simulation study has been carried out for different number of experiments (m = 4. . .20); it was assumed that the s.t.d. of the measurement errors is equal to = 0.1 mm. To obtain reliable results, the simulations have been repeated 10,000 times and averaged.

The manipulator parameters used for the simulations are presented in Table 4. It was assumed that the link lengths are in the range of 100. . .260 mm, while the geometric errors in the link lengths may be from 0.4 to 1.5 mm and the joint off-sets are 0.3. . .0.7 • . This table also contains the summary of the simulation results for 4 and 20 calibration experiments, which are in good agreement with expressions [START_REF] Santolaria | Kinematic parameter estimation technique for calibration and repeatability improvement of articulated arm coordinate measuring machines[END_REF].

To show benefits of the proposed calibration plans, Fig. 1 contains more detailed results corresponding to the conventional and proposed approaches (for the case of 4-links planar manipulator). As follows from them, the optimal pose selection for the calibration experiments allows us to reduce the identification errors by a factor of 1.5-2. For instance, for the number of calibration experiments m = 10 the proposed technique allows us to obtain the s.t.d. of the identification accuracy 0.032 mm for the parameters L 1 , . . ., L 4 and from 0.007 • to 0.016 • for the parameters q 1 , . . ., q 4 . In contrast, for the conventional technique, corresponding values are 0.045. . .0.060 mm and 0.010 • . . .0.027 • , respectively. Thus, a simple and practically convenient "rule of thumb" developed in this paper

J k T J k Identification accuracy 2-Links manipulator diag (m • l1 2 , m • l2 2 , m, m), q 1 = √ m×l 1 ; q 2 = √ m×l 2 ; L 1 = √ m ; L 2 = √ m 3-Links manipulator diag(m • l1 2 , m• l2 2 , m • l3 2 , m, m, m) q 1 = √ m×l 1 ; q 2 = √ m×l 2 ; q 3 = √ m×l 3 ; L 1 = √ m ; L 2 = √ m ; L 3 = √ m 4-Links manipulator diag (m • l1 2 , m• l2 2 , m• l3 2 , m • l4 2 , m, m, m, m) q 1 = √ m × l1 ; q 2 = √ m × l2 ; q 3 = √ m × l3 ; q 4 = √ m × l4 ; L 1 = √ m ; L 2 = √ m ; L 3 = √ m ; L 4 = √ m Table 4
Identification accuracy for manipulator geometric parameters: simulation study ( = 0.1 mm).

Manipulator

L 1 = 0.050 mm, q 1 = 0.011 o L 2 = 0.050 mm, q 2 = 0.016 o L 1 = 0.022 mm, q 1 = 0.005 o L 2 = 0.022 mm, q 2 = 0.007 o 3-Links manipulator L1 = 260 mm, L1 = 1.5 mm, q1 = 0.5 o L2 = 180 mm, L2 = -0.6 mm, q2 = -0.5 o L3 = 120 mm, L3 = -0.4 mm, q3 = 0.7 o L 1 = 0.050 mm, q 1 = 0.011 o L 2 = 0.050 mm, q 2 = 0.016 o L 3 = 0.050 mm q 3 = 0.024 o L 1 = 0.022 mm, q 1 = 0.005 o L 2 = 0.022 mm, q 2 = 0.007 o L 3 = 0.022 mm q 3 = 0.011 o 4-Links manipulator L1 = 260 mm, L1 = 1.5 mm, q1 = 0.5 o L2 = 180 mm, L2 = -0.6 mm, q2 = -0.5 o L3 = 120 mm, L3 = -0.4 mm, q3 = 0.7 o L4 = 100 mm, L4 = 0.7 mm, q4 = -0.3 o L 1 = 0.050 mm, q 1 = 0.011 o L 2 = 0.050 mm, q 2 = 0.016 o L 3 = 0.050 mm q 3 = 0.024 o L 4 = 0.050 mm q 4 = 0.029 o L 1 = 0.022 mm, q 1 = 0.005 o L 2 = 0.022 mm, q 2 = 0.007 o L 3 = 0.022 mm, q 3 = 0.011 o L 4 = 0.022 mm, q 4 = 0.013 o
has obvious advantages and is very attractive from an engineering and practical view point.

Comparison analysis

To illustrate the advantages of the developed technique, let us consider an example of 3-link manipulator whose variations in joint coordinates are limited to the range of ±100 • . These limitations are quite common for industrial robots. The considered manipulator includes three linear and three angular geometrical parameters that should be identified by means of calibration. For comparison purposes three different plans of experiments are considered:

(i) regular plan within the joint limits; (ii) random plan within the joint limits; and (iii) plan of experiments that satisfies optimality conditions [START_REF] Zhuang | Optimal selection of measurement configurations for robot calibration using simulated annealing[END_REF] and takes into account the joint limits.

It should be noted that the approach for optimal measurement configurations selection proposed in [START_REF] Borm | Determination of optimal measurement configurations for robot calibration based on observability measure[END_REF], which also insures conditions [START_REF] Zhuang | Optimal selection of measurement configurations for robot calibration using simulated annealing[END_REF], is not applicable here, since joint limits do not allow us to generate regular plan of experiments within the whole joint space, i.e. ±180 • for each variable (Fig. 2).

For all plans of experiments, the same number of measurements m = 64 has been used (this value has been chosen to have four different angles q 1 , q 2 and q 3 for the case of regular plan). In the plans (i) and (ii), all measurement configurations are different, while in plan (iii) the measurements have been repeated 16 times for 4 optimal configurations. The link lengths have been assigned to 1.25 m, 1.10 m and 0.23 m (those values corresponds to the link lengths of Kuka industrial robot used in experimental study), the measurements noise parameter was = 0.1 mm. The identification accuracy of three examined plans of experiments is presented in Table 5. The results show that the proposed plan of calibration experiments allows us to increase essentially the parameters identification accuracy with respect to the random and regular plans (from 20% to 60%). Corresponding improvement of the robot positioning accuracy is presented in Fig. 3, which contains three plots showing scattering of the end-effector location after calibration due to random nature of the measurement errors in the following configuration q 1 = 0 • , q 2 = -90 o and q 3 = 90 o . The value in Fig. 3 denotes the dispersion of the end-effector errors computed using particular randomly generated data. An additional study has been performed for the density of points in the regular plan. The same calibration experiments has been performed with 2, 4, 8, 16 and 25 regular spaced points per joint. For each case, the improvement factor has been computed, which shows the benefits of increasing the number of measurement configurations on the parameters identification accuracy. For comparison purposes, the expected improvement factors have been also computed, which corresponds to repeating experiments in the same configurations several times. Simulation results are summarized in Table 6, they show that increasing the number of configurations up to 16,256 does not allow us to increase the identification accuracy more than by 2% comparing with the repeating experiments in 64 reference configurations. This confirms that points density augmentation in the regular plan does not leads to appreciable efficiency improvement, the main accuracy improvement is achieved here due to increasing the number of measurements.

Hence, this example shows that by using limited number of optimal measurement configurations and repeating experiments several times, it is possible to achieve better accuracy than using huge number of configurations distributed within the joint limits. So, the proposed approach speeds up the calibration process without loosing the identification accuracy, i.e. improves its efficiency.

Let us also show the benefits of the proposed approach for the measurement configuration selection compared to the techniques based on the conventional observability indices [START_REF] Sun | Observability index selection for robot calibration[END_REF]. It was shown before [START_REF] Joubair | Comparison of the efficiency of five observability indices for robot calibration[END_REF] that all these indices give similar results in terms of identification accuracy. Here, the comparison analysis is done using the same example of the 3-link planar manipulator with the joint limits ±100 • . For each performance measures O 1 -O 5 , four measurement configurations have been found, which are optimal in terms of the corresponding observability index. It should be mentioned that to generate optimal configurations using conventional techniques, each case requires solution of non-trivial optimization problem with 12 variables (4 measurement configurations). The obtained optimal configurations are presented schematically in Fig. 4. As follows from it, the obtained solutions are distributed in the robot workspace and differ from one case to another. It should be mentioned that all sets of measurement configurations ensure the same parameters identification accuracy that is equal to the accuracy of the developed technique (see optimal plan in Table 5). Moreover, the values of the observability indices are equal for all obtained plans of experiments (Table 7). Another important conclusion is that for the considered manipulator the positioning accuracy after calibration (using optimal measurement configurations) does not vary essentially within the robot workspace (see Table 7), the difference between maximal and root mean square errors is lower than 0.1%. Hence, the developed technique gives also good results comparing with the test-pose approach proposed in [START_REF] Klimchik | Design of calibration experiments for identification of manipulator elastostatic parameters[END_REF], but its obvious advantage is extra-simplicity (it does not require application any numerical optimization routine).

To make graphical presentation more clear, it is useful to take into account that the angle q 1 does not effect on the identification accuracy. So, in optimal plan of experiments, q 1 can be fixed to any arbitrary value. Hence, in order to make comparison analysis easier, it is reasonable to fix the angle q 1 to the same constant value in all obtained plans of experiments (q 1 = 0, for instance). After such modifications, the obtained optimal configurations become semisymmetrical (Fig. 5). Besides, in all cases (O 1 -O 5 ) the measurement configurations satisfy the proposed optimality conditions (28), i.e. all plans of experiments are optimal in terms of the developed rule. However, in contrast to the cases O 1 -O 5 , the developed technique allows us to generate optimal measurement configurations without solving optimization problem. This advantage becomes extremely important with increasing of the number of actuated joints and the number of measurement configurations. In particular, for a 6-dof manipulator and 20 measurement configurations, conventional techniques require solution of an optimization problem with 120 unknowns. This high-dimensional problem is very time consuming and no one can guaranty that the obtained solution the global optimum. In contrast, the developed technique allows user to get the optimum without solving the optimization problem.

Extension to the 3D case

In spite of numerous advantages, the proposed approach has essential limitations. They are related to the manipulator architecture, for which the optimality conditions are formulated and strictly proved. In fact, the manipulator topology described by Eq. ( 2) includes 2D serial planar mechanisms with revolute joints only. Although generalization for the case of prismatic joints is trivial,

Table 7

Performances measures for measurement configurations generated using different optimality criteria.

Performance measures extension to the 3D case is non obvious and requires additional efforts.

O1 : r √ 1 •...• r √ r → max O2 : r 1 → max O3 : r → max O4 : 2 r 1 → max O5 : 1 1 + • • • + 1 r -1 → max
The simplest way to extend the proposed "rule of thumb" to the 3D case is to apply the following procedure:

Step (a): decompose the spatial manipulator into a set of planar serial sub-chains;

Step (b): apply the developed rule to each sub-chain separately, without assigning certain values to the joint coordinates that can be selected arbitrary;

Step (c): aggregate the obtained sub-chain joint coordinates in order to find configurations of the entire manipulator (where some values are still arbitrary);

Step (d): apply the developed rule to the each set of the arbitrary coordinates, i.e. ensuring that the sums of sines and cosines are equal to zero for all of them.

To show the efficiency of this heuristic extension, let us consider several typical manipulator architectures used in industry:

(i) 3-Axis SCARA robot with 2 revolute and one prismatic joints (R z R z P z architecture). For this architecture, the rule of thumb can be applied straightforwardly since motions in the xy-plane and in z-direction are not coupled. This manipulator can be easily decomposed in two sub-chains corresponding to 2-link planar R z R z manipulator and 1-axis prismatic P z mechanism.

For the first one, the optimal solution is provided by Eq. [START_REF] Zhuang | Optimal selection of measurement configurations for robot calibration using simulated annealing[END_REF].

For the second one, the optimal solution is trivial and is defined by the classical design of experiment theory [START_REF] Atkinson | Optimum experimental designs[END_REF]. (ii) 4-Axis SCARA robot with 3 revolute and one prismatic joints (R z R z P z R z architecture). Taking into account that this architecture is nominally equivalent to R z R z R z P z one, here decomposition leads to R z R z R z and P z mechanisms. So, the optimal solutions are obtained in similar way: using results for 3-link planar manipulator and 1-axis prismatic mechanism. (iii) 2-Axis articulated robot with two revolute joints (R z R y architecture). In this case, decomposition yields two trivial sub-chains of R z and R y types. As follows from [START_REF] Zhuang | Optimal selection of measurement configurations for robot calibration using simulated annealing[END_REF], at the beginning of each separate sub-chain, the joint coordinates can be set arbitrary. However, applying Step (d) of the above presented procedure to the joint R y , one can get values of the corresponding joint coordinates ensuring zero sums of the sines and cosines. It can be easily verified that the latter operation yields a diagonal information matrix (without applying Step (d) to the first joint R z ). (iv) 3-Axis articulated robot with three revolute joints (R z R y R y architecture). Here, decomposition leads to a trivial sub-chain R z and a 2-link planar manipulator R y R y . They are treated in the same way as the above presented R z R y case, allowing the first joint R z to be set arbitrary. Similarly, obtained measurement configurations ensure diagonal structure of the information matrix.

Optimal measurement configurations for the above considers spatial manipulators are presented in Table 8, which also contains corresponding information matrices that are strictly diagonal. They give certain compromise between the identification accuracy of different manipulator parameters that is acceptable from a practical point of view, but obviously is not unique.

Hence, extension of the rule of thumb to the 3D case is rather simple, while it has been derived intuitively, using some heuristic ideas. Nevertheless, dedicated studies confirm that the results are quite good but perhaps not optimal with respect to some parameters. Therefore, the general (and strict) solution of the problem for the 3D-case is still an open question, which will be in the focus of further research. The virtue of these results in the case of 6-axis manipulator is shown in the next section.

Application example

Now let us apply the obtained results to the calibration of the 6 d.o.f. Kuka KR-270 industrial robot (Fig. 5), which obviously is not a planar manipulator studied in previous sections. However, it is attractive to try the proposed rule of thumb (see Section 9) heuristically, considering this spatial manipulator as set of planar kinematic sub-chains. In particular, it is natural to decompose the manipulator into planar sub-chains composed from the following joints and corresponding links: {joint no. 1}, {joints no.2, 3, 5} and {joints no.4, 6}. Using this decomposition and taking into account the joint limits of industrial robot KR-270, it is possible to generate four measurement configurations (Table 8) allowing us to identify nine desired parameters: the link length deviations l 1 , . . ., l 3 and the joint offsets q 1 , . . ., q 6 . The validity of these configurations for the entire 6-axes manipulator can be easily verified by evaluating the rank of the corresponding calibration Jacobian.

To obtain the desired identification accuracy (0.005 mm for l i and 0.01 deg for q i ), the calibration experiments for each measurement configuration should be repeated several times. Taking into account that conventional measurement systems used in industrial environment (laser trackers) provide the accuracy corresponding to the measurement noise parameter ≈ 0.03 mm, one can get the minimum number of calibration experiments (see Table 3) that in this case is m min = 40. It means that the calibration experiments should be repeated 10 times for all configurations presented in Table 8.

0 • 90 • 180 • q2 -30 • -60 • -90 • -120 • q3 90 • -90 • -90 • 90 • q4 0 • 180 • 90 • -90 • q5 120 • -120 • -60 • 60 • q6 -120 • 120 • 60 • -60 •
To show the efficiency of the proposed heuristic approach, Table 9 presents the identification accuracy for two sets of measurement configurations (calibration plans). The first one, Plan no. 1, corresponds to the regular distribution of the joint coordinates that is usually used in practice [START_REF] Dumas | Joint stiffness identification of six-revolute industrial serial robots[END_REF]. The second one, Plan no. 2, uses the measurement configurations (see Table 8) generated using the proposed rule of thumb. As follows from the obtained results, the proposed technique allowed us to improve the calibration accuracy by a factor of 1.7-4.2. The latter confirms usefulness of this heuristic idea: application of the results obtained for a planar case to spatial manipulators.

In more details, the advantages of the Plan no. 2 are confirmed by the following values. For the linear parameters, the regular plan ensures the identification accuracy 0.010. . .0.021 mm, while the proposed one allows us to reduce the identification errors to 0.005. . .0.006 mm. Similarly, for the angular parameters, application of the proposed technique leads to the reduction of the identification errors from 0.019. . .0.100 mrad to 0.008. . .0.050 mrad. Those results are also illustrated in Fig. 5, which presents stochastic simulation of the calibration experiments corresponding to two sets of measurement configurations. Corresponding graphs clearly show that the proposed technique provides unbiased estimates of the desired parameters and allows essentially reducing deviations of the identification errors. Hence, summarizing this Section, one may conclude that the developed approach, initially targeted to the planar manipulators, can be also efficiently applied to the spatial manipulators (using some heuristic decomposition of the corresponding kinematic chain).

Discussion

The main advantage of the results presented in this paper is related to essential reduction of computational efforts required for generation of optimal measurement configurations, which are used for geometrical calibration of robotic manipulators. In fact, instead of tedious numerical optimization, simple rules can be applied that yield the set of joint angles satisfying certain trigonometric equations (based on the sums of sines and cosines of the joint coordinates and their combinations, see Table 1). It is proved that, in the planar case, this solution ensures the minimization of the determinant of the covariance matrix (Dcriterion), which is often used for the calibration experiments planning. However, in the general case (for spatial manipulators), this statement is not valid anymore, while providing quite good solution completely acceptable in practice. The latter justifies the results utility for wider class of manipulator architectures (Table 10). Another important question, which is beyond the scope of this paper, is related to the selection of the optimization criteria used for the evaluation of the measurement configurations "quality". Here, the main idea came from the analysis of the information matrix that should be strictly diagonal (which corresponds to the maximum of its determinant). It is clear that other existing optimization criteria [START_REF] Zhuang | Optimal planning of robot calibration experiments by genetic algorithms[END_REF][START_REF] Daney | Choosing measurement poses for robot calibration with the local convergence method and Tabu search[END_REF][START_REF] Atkinson | Optimum experimental designs[END_REF][START_REF] Borm | Determination of optimal measurement configurations for robot calibration based on observability measure[END_REF][START_REF] Khalil | Identifiable parameters and optimum configurations for robots calibration[END_REF][START_REF] Daney | Optimal measurement configurations for Gough platform calibration[END_REF][START_REF] Nahvi | The noise amplification index for optimal pose selection in robot calibration[END_REF][START_REF] Klimchik | Design of experiments for calibration of planar anthropomorphic manipulators[END_REF] can be also used, but there is no guarantee that they can produce similar rule of thumb or corresponding solutions are better from engineering view point. Some details concerning comparison of different optimization criteria in calibration experiment planning can be found in our previous work [START_REF] Klimchik | Optimal selection of measurement configurations for stiffness model calibration of anthropomorphic manipulators[END_REF].

It should be noted that despite the fact that the original problem was formulated without paying attention to the joint limits, the obtained results can be also efficiently applied if the joint limits are wider than 180 o . Application of the proposed rule of thumb for this case has been presented in Section 8, where the achieved results are 20-60% better than those obtained with regular and random plans. It is clear that special attention should be paid to the problem where the joint ranges are lower than 180 o . In this case some of the optimality conditions (28) cannot be satisfied and a sub-optimal solution should insure minimum of corresponding residuals. However, the difficulty of applying the rule of thumb can only occur when at least two joints have narrow limits and it is required to minimize absolute values of several residuals (28) simultaneously. Nevertheless, this problem is not common in practice, since such joint ranges usually are typical for the first joint coordinate only. In fact, the joint angle q 1 is not included in the optimality conditions [START_REF] Zhuang | Optimal selection of measurement configurations for robot calibration using simulated annealing[END_REF], and it can be set arbitrary in the optimal plan developed above.

Besides, the developed rule of thumb can be also used to generate the measurement configurations if there are some constraints in the Cartesian workspace (work-cell limits, obstacles, etc.). Since any optimal plan defined by Eq. ( 28) has a number of redundant variables that can be set arbitrary, it is usually possible to shift the measurement configurations in the joint space to achieve an allowable zone in the Cartesian workspace.

Another question behind the main scope of this paper is related to the measurement noise impact in actuated joints. However, as follows from our experience, in practice this factor influence on the identification accuracy is much lower and almost negligible comparing with the measurement errors in Cartesian space. Moreover, it can be strictly proved that this type of errors does not change the optimal measurement configurations since their impact on the parameters identification accuracy (covariance matrix) is described by the constant diagonal term.

Nevertheless, in spite of the fact that obtained results have been efficiently applied to the 3D case and allowed us to improve the identification accuracy, there are several theoretical issues to be studied in future. In particular, special attention should be paid to the optimal pose selection for the calibration of serial kinematic chains with both revolute and translation joints mixed in arbitrary order (Figs. 6 and7).

Conclusion

The paper presents a new approach for the design of calibration experiments for robotic manipulators that essentially simplifies this procedure and eliminates tedious numerical optimization of measurement poses. The main result is expressed as a system of trigonometric equations describing the set of optimal configurations, which produces strictly diagonal structure of the information matrix. It corresponds to the D-optimality principle known in the design of experiment theory. Relevant technique is formulated in the form of the "rule of thumb" allowing practicing engineer to generate desired measurement configurations in a simple way. The validity of the obtained results and their practical significance were confirmed via simulation studies that deal with two-, threeand four-link planar manipulators. Compared to previous works, the obtained results can be treated as further development of the design-of-experiments theory that is adapted to the specific type of the non-linear models that arise in robot kinematics.

Although the main contributions have been obtained for the planar case, the developed rule has been heuristically generalized for the case of spatial manipulators and successfully applied to several articulated robots. An application example, which deals with the geometric calibration of a Kuka KR-270 industrial robot, also confirms the efficiency of the proposed technique. Nevertheless, a strict theoretical proof of this approach remains the subject of future work that will focus on the optimal pose selection for calibration of non-planar serial and parallel manipulators. for which it is required to find a maximum with respect to the angle coordinates q k 2 , q k 3 , q k 23 , k = 1, m. Expression (45) depends on the sum of cosines and sines of these arguments. For the convenience purposes these sums are substituted by single unknown parameters.

Using the above defined assumptions, solution of the optimization problem is reduced to finding the extremum of the inner function and to extracting solutions that correspond to the maximum. In order to find the extremum let us compute partial derivatives of expression (45) with respect to all unknown parameters. As a result, one can obtain the following system of equations

m × x = z × v + u × w m× y = z × w -u × v m × z = x × v + y × w m× u = x × w -y × v m × v = x × z -y × u m× w = y × z + x × u (46)
that should be solved with respect to x, y, z, u, v, w. This system of equations has the following solutions

⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ v = 0 w = 0 z = 0 u = 0 x = 0 y = 0 ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎭ ∪ ⎧ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎩ v 2 + w 2 = m 2 u = 0 z = × m x = × v y = × w ⎫ ⎪ ⎪ ⎪ ⎪ ⎬ ⎪ ⎪ ⎪ ⎪ ⎭ ∪ ⎧ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎩ z 2 + u 2 = m 2 w = 0 v = × m x = × z y = × u ⎫ ⎪ ⎪ ⎪ ⎪ ⎬ ⎪ ⎪ ⎪ ⎪ ⎭ ∪ ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ u = 0 v = 0 z = 1 × m w = 2 × m x = 0 y = 1 2 m ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎭ (47)
where = ±1, 1 = ±1, 2 = ±1.

Substituting the first solution of (47) in (45) one can obtain that m 3 is the maximum. The second and the third systems of (47) have two solutions that bring to zero expression (45). It is clear that such solutions do not ensure the maximum of the squared function. Similarly, the fourth solution does not bring maximize expression (45).

As a result, a necessary and sufficient condition to maximize (45) is That can be interpreted on the plane as a zero sum of unit vectors based on the optimal angles.

Fig. 1 .

 1 Fig. 1. Positioning accuracy of 2-dof manipulator after calibration using different plans of experiments (measurement noise = 0.1 mm): plan (i) q1 = 30 • , q2 = -100 • and q1 = -150 • , q2 = -90 • ; plan (ii) q1 = 30 • , q2 = -90 • and q1 = 30 • , q2 = 90 • .

  manipulator L1 = 260 mm, L1 = 1.5 mm, q1 = 0.5 o L2 = 180 mm, L2 = -0.6 mm, q2 = -0.5 o

Fig. 2 .

 2 Fig. 2.Comparison of the identification accuracy for different calibrations plans (case of 4-links planar manipulator): --äre analytical curves corresponding to the proposed optimal plan (expressions (36)); "x" are experimental results corresponding to the optimal calibration plan (simulation); "o" are experimental results corresponding to the calibration plan with regular grid (simulation).

Fig. 3 .

 3 Fig. 3. Geometrical error compensation efficiency for different calibration plans of experiments (3-dof planar manipulator, identification of six geometrical parameters, 64 measurement configurations, measurement noise parameter = 0.1 mm, manipulator configuration q1 = 0 • , q2 = -90 • , q3 = 90 • ).

Fig. 4 .

 4 Fig. 4. Schematic presentation of optimal measurement configurations for 3-link planar manipulator obtained using the developed technique and conventional observability indices.

Fig. 5 .

 5 Fig. 5. Schematic presentation of modified optimal measurement configurations for a 3-link planar manipulator obtained using the developed technique and conventional observability indices (with q1 = 0).

Fig. 6 .

 6 Fig. 6. Industrial robot Kuka KR-270 TM.

Fig. 7 .

 7 Fig. 7. Comparison of the identification accuracy for two types of calibration plans: (i) regular plan [44]; (ii) proposed plan (measurement noise parameter ≈ 0.03 mm, m min = 40).

  [START_REF] Joubair | Comparison of the efficiency of five observability indices for robot calibration[END_REF] defines optimality condition for calibration plan of 2-link manipulator.A.2. Optimality condition for a 3-link manipulatorFor the 3-link manipulator, let us introduce following notation x = to rewrite the optimization problem[START_REF] Lightcap | Improved positioning accuracy of the PA10-6CE robot with geometric and flexibility calibration[END_REF] asm 3 + 2yzw + 2xzv -2yuv + 2xuz -m x 2 + y 2 + z 2 + u 2 + v 2 + w 2 2

Table 3

 3 Identification accuracy for manipulator geometric parameters.

	Manipulator	m k=1

Table 5

 5 Comparison of parameters identification accuracy for different calibration plans.

	Parameters	Identification accuracy
		Random plan	Regular plan	Optimal plan
	Â1 (mrad)	0.012	0.012	0.010
	Â2 (mrad)	0.017	0.016	0.011
	Â3 (mrad)	0.075	0.067	0.054
	L1 (mm)	0.015	0.015	0.013
	L2 (mm)	0.019	0.018	0.013
	L3 (mm)	0.017	0.015	0.013

Table 6

 6 Improvement of parameter identification accuracy due to increasing the measurement configuration density for regular plan.

		Number of regular spaced points per joint
		2	4 *	8	16	25
	Total number of configurations	8	64	512	4096	15,625
	Parameters	Improvement factors
	Â1	0.33	1	2.86	8.12	15.9
	Â2	0.32	1	2.88	8.18	16.0
	Â3	0.33	1	2.86	8.12	15.9
	L1	0.33	1	2.86	8.12	15.9
	L2	0.32	1	2.88	8.18	16.0
	L3	0.33	1	2.86	8.12	15.9
	Expected	0.35	1	2.83	8.00	15.6

Table 8

 8 Optimal measurement configurations for typical spatial manipulators.

Table 9

 9 Optimal measurement configurations for calibration of industrial robot KUKA KR-270.

	Joint coordinate	Configuration no. 1	Configuration no. 2	Configuration no. 3	Configuration no. 4
	q1	-90 •

Table 10

 10 Identification accuracy for the regular (no. 1) and proposed (no. 2) plans of calibration experiments.

	Parameter	Identification accuracy	Improvement factor
		Plan no. 1	Plan no. 2
	l1 [mm]	0.015	0.005	3.0
	l2 [mm]	0.021	0.005	4.2
	l3 [mm]	0.010	0.006	1.7
	q2 [mrad]	0.021	0.008	2.6
	q3 [mrad]	0.019	0.009	2.1
	q5 [mrad]	0.100	0.050	2.0

I 2×2 , k = j (20)

* This case have been chosen as a reference one.

Acknowledgements

The work presented in this paper was partially funded by the Region "Pays de la Loire", France, by the project ANR COROUSSO, France and FEDER ROBOTEX, France (Grant no. ANR-2010-SEGI-003-02-COROUSSO).

Appendix A. Appendix A

Analytical solutions for the optimality conditions for 2-and 3link manipulators

The optimization function for the D-optimality criterion is based on the determinant of the information matrix and for the considered problem can be presented as:

where the matrices C, S are defined by Eq. ( 23). For illustrative purposes analytical solutions for the cases of 2-and 3-link manipulators are proposed in the following sub-sections.