Issue 
Int. J. Metrol. Qual. Eng.
Volume 15, 2024



Article Number  21  
Number of page(s)  14  
DOI  https://doi.org/10.1051/ijmqe/2024016  
Published online  04 October 2024 
Research Article
A virtual CMM to estimate uncertainties
^{1}
CEA/DAM/ Gramat, 46500 Gramat, France
^{2}
University of Angers − LARIS, 62 avenue Notre Dame du Lac, 49000 Angers, France
^{*} Corresponding author: jeanfrancois.manlay@cea.fr
Received:
18
July
2023
Accepted:
14
August
2024
Coordinates measuring machines (CMMs) have become classical measuring instruments. Nevertheless, the uncertainties associated to measurements results obtained by CMMs, are often a global estimation. This work focuses on the development of a virtual CMM with Python, in order to estimate uncertainties of all measured dimensions, on basic parts. The originality of the approach consists in randomizing only the nine linear parameters of the CMM compensation matrix, which allows calculating the twelve other ones, and avoiding complex covariance calculations. The model takes into account CMM defects, probing uncertainties and some parameters of the part and the environment. Method uncertainty, difficult to be modeled, is treated as a set of recommendations. This model can be used as a preprocessing evaluation of uncertainty, or postprocessing of a real measurement.
Key words: Coordinate measuring machine / uncertainty / virtual machine
© J.F. Manlay et al., Published by EDP Sciences, 2024
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
1 Introduction
Each measurement result should be associated with its own uncertainty. In coordinate metrology, the uncertainty calculation appears as a complex domain, and it is rarely evaluated. Many users confuse measurement uncertainty with tolerance of the CMM acceptance test. Most advanced users can participate in intercomparisons or use specific software (UES: Uncertainty Evaluation Software) to estimate their uncertainties [1]. One of the most wellknown UES is VCMM (Virtual coordinate measuring machine) [2], but a few others also exist [3], or are being developed in labs [4]. Understanding the uncertainties on CMMs comes through having some knowledge about their functioning and their calibration [5]. The uncertainties of several parameters, such as part, environment [6] or methods are important to evaluate. The methods used are partially similar to classical evaluation of uncertainties in mechanics.
In this study, we will only focus on bridge CMMs, using contacting probes. We decided to follow another development approach, ensuring to be as close as possible to our CMM system and our measurements. The developed model allows calculating uncertainties on distances, angles, form... It is a new contribution to the evaluation of uncertainties on CMM results.
2 Methods
2.1 CMM defects description
A bridge CMM is a measuring instrument which allows taking hits with a probe that moves along three axes, nominally perpendicular in pairs. Each axis is considered separately from the other two.
During the move along an axis, the moving device is subjected to different defects, as shown in Figure 1.
Those defects are separated into two classes, linear defects and rotational defects. Linear defects are represented by the trueness of the scale and the path straightness of the air bearings [7]. Rotational defects are roll, pitch and yaw. There are six defects per axis, so eighteen defects, to which are added perpendicularity defects between axes two by two.
All those defects are corrected by a compensation matrix, usually called the “21 parameters compensation matrix”.
Fig. 1 Representation of CMM defects along an axis during its displacement. 
2.1.1 Compensation matrix
The compensation matrix is used to correct reading values on the scales, allowing giving an accurate position of the probe. Displacements along each axis are considered separately, and the different parts of the CMM are considered as rigid features, keeping their shape [8].
The compensation is done as a field of moments, described in Figure 2, where $\overrightarrow{{E}_{i}\left({O}_{j}\right)}$ represents defect of the i axis, named L_{ix}, L_{iy}, L_{iz} below, and $\overrightarrow{{W}_{i}}$ represents rotational defects of the i axis around j axes R_{ix}, R_{iy}, R_{iz.} P is the center ball of the probe.
The geometrical defects of a prismatic link can be represented by three rotations and two translations. Then, a displacement of P along Y axis can be described by the displacement of O_{1} and the effect of rotations on the lever arm between P and O_{1} (1).
$$\overrightarrow{{E}_{Y}\left(P\right)}}={\displaystyle \overrightarrow{{E}_{Y}\left({O}_{1}\right)}}+{\displaystyle \overrightarrow{{W}_{Y}}}\Lambda {\displaystyle \overrightarrow{{O}_{1}P}$$
with:
$$\overrightarrow{{E}_{Y}\left({O}_{1}\right)}}\left(\begin{array}{c}\hfill {L}_{YX}\hfill \\ \hfill Y+{L}_{YY}\hfill \\ \hfill {L}_{YZ}\hfill \end{array}\right);{\displaystyle \overrightarrow{{W}_{Y}}}\left(\begin{array}{c}\hfill {R}_{YX}\hfill \\ \hfill {R}_{YY}\hfill \\ \hfill {R}_{YZ}\hfill \end{array}\right);{\displaystyle \overrightarrow{{O}_{1}P}}\left(\begin{array}{c}\hfill {O}_{1x}+{T}_{x}\hfill \\ \hfill {T}_{Y}\hfill \\ \hfill {O}_{1z}+{T}_{z}\hfill \end{array}\right).$$(1)
The same operation on both other axes gives:
$$\overrightarrow{{E}_{X}\left(P\right)}}={\displaystyle \overrightarrow{{E}_{X}\left({O}_{2}\right)}}+{\displaystyle \overrightarrow{{W}_{X}}}{\displaystyle \overrightarrow{{O}_{2}P}$$(2)
$$\overrightarrow{{E}_{Z}\left(P\right)}}={\displaystyle \overrightarrow{{E}_{Z}\left({O}_{3}\right)}}+{\displaystyle \overrightarrow{{W}_{Z}}}\Lambda {\displaystyle \overrightarrow{{O}_{3}P}}.$$(3)
Compensation along each axis is calculated as shown below [8,9], by developing equations (1), (2) and (3):
$${C}_{X}={L}_{XX}+{L}_{YX}+{L}_{ZX}+{R}_{YY}.\left({P}_{Z}+{T}_{Z}\right){R}_{YZ}{T}_{Y}+{R}_{XY}.\left({P}_{Z}+{T}_{Z}\right){R}_{XZ}{T}_{Y}+{R}_{ZY}.{T}_{Z}{R}_{ZZ}{T}_{Y}$$(4)
$${C}_{Y}={L}_{XY}+{L}_{YY}+{L}_{ZY}+{R}_{YZ}.\left({P}_{X}+{T}_{X}\right){R}_{YX}\left({P}_{Z}+{T}_{Z}\right)+{R}_{XZ}{T}_{X}.{R}_{XX}\left({P}_{Z}+{T}_{Z}\right)+{R}_{ZZ}.{T}_{X}{R}_{ZX}{T}_{Z}$$(5)
$${C}_{Z}={L}_{XZ}+{L}_{YZ}+{L}_{ZZ}+{R}_{YX}.{T}_{Y}{R}_{YY}\left({P}_{X}+{T}_{X}\right)+{R}_{XX}.{T}_{Y}{R}_{XY}{T}_{X}+{R}_{ZX}.{T}_{Y}{R}_{ZY}{T}_{X.}$$(6)
where P_{i} is the position read on the scale along i axis, and T_{i} are the offsets of the center ball of the probe from the calibration origin, at the end of the Zaxis.
Fig. 2 Compensation of the position of P by applying a field of moment on each axis one by one. 
2.1.2 Rotational defects
In this study, we chose to calculate the rotational defects from the linear ones.
The air bearings move on paths (granite or aluminium) which are not perfectly flat. For one of the moving parts (Fig. 3), three air bearings build a plane, and other two construct a line.
During the displacement of the moving part of the CMM along its guide ways, the plane and the line vectors change direction, because of the linear defects. The projections of those vectors define different rotations, roll, pitch and yaw (Fig. 4).
Extending this principle to the whole CMM is explained in Figures 5 and 6.
On the CMM presented, the bridge has the triangular shape. The plane vector built on the upper face is not parallel to CMM axes. The rotational defects are calculated from this vector projection on the two axes perpendicular to the move.
Fig. 3 Scheme of the bridge air bearings structure representing a sliding connection. 
Fig. 4 Representation of the link between the flatness of the guide ways and the rotation defects of the bridge 
Fig. 5 Relationship between flatness defects on guide ways, air bearing positions and rotational defects along X and Y axes. 
Fig. 6 Relationship between flatness defects on guide ways, airbearing positions and rotational defects along Z axis. 
2.1.3 Perpendicular defects
Different manufacturers of CMMs assign perpendicular defects as a constant value added to a rotation or as a slope of straightness. In both cases, those defects depend on linear defects, so they are calculable from them.
In the example given in Figure 7, the manufacturer uses a rule where the straightness must be zero at the start and at the end of the stroke axis. The measurement of the straightness of an axis along a second one during the calibration usually starts at zero, and provides deviations from it all along the axis. The slope of the line between extreme points defines the perpendicularity defect between both axes.
In this study, the perpendicular defects are added as a constant value to rotational defects. The perpendicularity defect between Z and X is added to R_{XY}, between Y and Z is added to R_{ZX}, and between X and Y is added to R_{XZ}.
Fig. 7 Relationship between straightness of Y axis measured along Z axis and perpendicularity defect between Y and Z. The slope defines the compensation value for the perpendicularity defect. 
2.1.4 Trueness defects
The trueness defects of each scale are corrected in the compensation matrix by an interferometer comparison between actual readings on the CMM and calibrated distances measured by the interferometer (Fig. 8). Those measurements give the values of L_{XX}, L_{YY} and L_{ZZ} in (1), (2) and (3).
Fig. 8 : Trueness calibration using an interferometer by comparison to the scale readings. 
2.1.5 Straightness and rotational calculations
In the compensation matrix, straightness defects are represented as a constant value for a given position along an axis. This value has to be calculated from air bearings positions on the guide ways.
In the application of a field of moments, there is a point M where the rotation vector is zero. Knowing the air bearings vertical position dP_{iz} allows calculating the coordinates of this point, and therefore the value of the straightness defects.
The first step involves the three air bearings constructing the plane (example on the bridge displacement on the granite):
$$\{\begin{array}{c}\hfill {\displaystyle \overrightarrow{dM}}.{\displaystyle \overrightarrow{z}}=\left({\displaystyle \overrightarrow{d{P}_{1}}}+{\displaystyle \overrightarrow{R}}\Lambda {\displaystyle \overrightarrow{{P}_{1}M}}\right).{\displaystyle \overrightarrow{z}}\hfill \\ \hfill {\displaystyle \overrightarrow{dM}}.{\displaystyle \overrightarrow{z}}=\left({\displaystyle \overrightarrow{d{P}_{2}}}+{\displaystyle \overrightarrow{R}}\Lambda {\displaystyle \overrightarrow{{P}_{2}M}}\right).{\displaystyle \overrightarrow{z}}\hfill \\ \hfill {\displaystyle \overrightarrow{dM}}.{\displaystyle \overrightarrow{z}}=\left({\displaystyle \overrightarrow{d{P}_{3}}}+{\displaystyle \overrightarrow{R}}\Lambda {\displaystyle \overrightarrow{{P}_{3}M}}\right).{\displaystyle \overrightarrow{z}}\hfill \end{array}$$(7)
$$with\left({\displaystyle \overrightarrow{R}}\Lambda {\displaystyle \overrightarrow{{P}_{i}M}}\right)=\left({\displaystyle \overrightarrow{R}}\Lambda {\displaystyle \overrightarrow{{P}_{i}O}}\right)+\left({\displaystyle \overrightarrow{R}}\Lambda {\displaystyle \overrightarrow{OM}}\right).$$(8)
By definition, $\left({\displaystyle \overrightarrow{R}}\Lambda {\displaystyle \overrightarrow{OM}}\right)={\displaystyle \overrightarrow{0}}$ so $\left({\displaystyle \overrightarrow{R}}\Lambda {\displaystyle \overrightarrow{OM}}\right).{\displaystyle \overrightarrow{z}}=0$
It gives:
$${R}_{yx}.{y}_{M}{R}_{yy}.{x}_{M}=0.$$(9)
Reporting in (9) in (7) gives:
$$\left[\begin{array}{ccc}\hfill {y}_{{P}_{1}}\hfill & \hfill {x}_{{P}_{1}}\hfill & \hfill 1\hfill \\ \hfill {y}_{{P}_{2}}\hfill & \hfill {x}_{{P}_{2}}\hfill & \hfill 1\hfill \\ \hfill {y}_{{P}_{3}}\hfill & \hfill {x}_{{P}_{3}}\hfill & \hfill 1\hfill \end{array}\right]\left[\begin{array}{c}\hfill {R}_{xy}\hfill \\ \hfill {R}_{yy}\hfill \\ \hfill {L}_{yz}\hfill \end{array}\right]=\left[\begin{array}{c}\hfill {z}_{{P}_{1}}\hfill \\ \hfill {z}_{{P}_{2}}\hfill \\ \hfill {z}_{{P}_{3}}\hfill \end{array}\right].$$(10)
Solving this set of equations (10) gives the values of the pitch of X, the roll of Y and straightness of Y measured along Z.
Then, the second step can be done, about two air bearings driving the straightness of Y measured along X:
$$\{\begin{array}{c}\hfill {\displaystyle \overrightarrow{dM}}.{\displaystyle \overrightarrow{x}}=\left({\displaystyle \overrightarrow{d{P}_{4}}}+{\displaystyle \overrightarrow{R}}\Lambda {\displaystyle \overrightarrow{{P}_{4}M}}\right).{\displaystyle \overrightarrow{x}}\hfill \\ \hfill {\displaystyle \overrightarrow{dM}}.{\displaystyle \overrightarrow{x}}=\left({\displaystyle \overrightarrow{d{P}_{5}}}+{\displaystyle \overrightarrow{R}}\Lambda {\displaystyle \overrightarrow{{P}_{5}M}}\right).{\displaystyle \overrightarrow{x}}\hfill \end{array}$$(11)
So:
$$\left[\begin{array}{cc}\hfill {y}_{{P}_{4}}\hfill & \hfill 1\hfill \\ \hfill {y}_{{P}_{5}}\hfill & \hfill 1\hfill \end{array}\right]\left[\begin{array}{c}\hfill {R}_{yz}\hfill \\ \hfill {L}_{yx}\hfill \end{array}\right]=\left[\begin{array}{c}\hfill {x}_{{P}_{4}}+{R}_{yy}.{z}_{{P}_{4}}\hfill \\ \hfill {x}_{{P}_{5}}+{R}_{yy}.{z}_{{P}_{5}}\hfill \end{array}\right].$$(12)
This calculation (12) defines the yaw of Y and the straightness of Y measured along X.
The same method is applied on the two other axes.
2.1.6 Temperature effect
The temperature can modify the CMM geometry because of difference between the temperature gradient within the table at the moment of the CMM calibration and the actual gradient. This can generate rotational defects such as pitch for a difference of temperature between upper and lower faces of the granite (Fig. 9).
Considering the curvature as a local portion of a circle, the pitch defect can be calculated:
$$\{\begin{array}{c}\hfill {L}_{1}=L\left(1+\frac{\alpha .\delta T}{2}\right)=\beta .\left(R1+\frac{H}{2}\right)\hfill \\ \hfill {L}_{2}=L\left(1\frac{\alpha .\delta T}{2}\right)=\beta .\left(R1\frac{H}{2}\right)\hfill \end{array}.$$(13)
So:
$$\{\begin{array}{c}\hfill \beta =\frac{L.\alpha .\delta T}{H}\hfill \\ \hfill R1=\frac{L}{\beta}\hfill \end{array}$$(14)
The versed sine V is:
$$V=\left(R2\right)\left(1\text{cos}\beta \right)withR2=R1+\frac{H}{2}.$$(15)
Therefore, the defect on the pitch is:
$${E}_{xz}=VR2\left(1\text{cos}\left(Asin\left(\frac{x}{R2}\right)\right)\right).$$(16)
Fig. 9 Representation of a temperature gradient between both faces of the granite. The expansion difference allows calculating the versed sine. 
2.2 Probing defects
The probing system creates different types of defects, which give another part of uncertainty.
2.2.1 Probe calibration
The probe calibration is done on a calibrated sphere, whose diameter and shape are in relationship with the meter definition.
The uncertainties on those values are often of a second order in front of probe calibration standard deviation, but should be taken into account in the result.
The probe calibration defines different parameters such as the real ball radius, the offsets from the master probe position (Fig. 10), and the standard deviation of sphere radii measurements.
This last parameter depends on the type of probe, and on the direction of probing. It directly defines a part of the uncertainty coming from the probe. In our model, it is a randomized value on the probe diameter, and on offsets (T_{x}, T_{y} and T_{z} in the compensation equations). A touch trigger probe calibration needs to be performed with enough points to accurately determine the standard deviation.
Fig. 10 Representation of a probe calibration. The results are a calculated ball radius and an offset deviation from theoretical position and a standard deviation of radii. 
2.2.2 Rotational head case
This type of head cannot give an accurate repositioning between rotations. It is a really useful tool, but it needs to be used with some concerns.
The head manufacturer gives a standard deviation of repositioning at a given distance of rotation axes (Fig. 11), which can be used as an uncertainty, depending on the length of the actual probe.
Those defects can be cancelled out by a probe calibration after each rotation, which decrease the interest of the head, but increase the accuracy.
Our model purposes both solutions to consider user's practices.
This case is also available by using a probe changer, which does not keep the accuracy of the previous calibration.
Fig. 11 Repositioning defects on a rotation head without calibration. The uncertainty depends on the length of the probing system. 
2.2.3 Probing force
This parameter is at the boundary between probe and part uncertainties.
It depends on the type of probe and on the material of the part.
Using Hertz formulas allows calculating the defect d depending on probe diameter, material parameters, and probing force:
$$d=\sqrt[3]{\frac{1}{R}{\left[\frac{3\pi F\left({K}_{1}+{K}_{2}\right)}{4}\right]}^{2}}$$(17)
$${K}_{i}=\frac{1{\vartheta}_{i}}{\pi {E}_{i}}$$
d : defect (mm)
R : probe radius (mm)
F : probing force (N)
ϑ_{i} : Poisson coefficient (part and probe − no dimension)
E_{i} : Young modulus (part and probe − Mpa).
Table 1 presents defect for a probing force of 0.03 N, for a contact between a sphere and a plane.
For now, this parameter is not taken into account by our model, because of our usual practices (ball diameter and materials). In the future, it will be implemented as a radius variation, depending on used materials.
Effect of the force probing on different materials.
2.3 Applying the model to the CMM
Assuming that the calibration of the CMM is perfect, so the compensation matrix compensates correctly geometrical imperfections, the residual defects come from the interferometer uncertainties and temperature variations.
Then, randomized defects are added on flatness guide ways and on trueness of scales, in relation to the interferometer uncertainties. We also randomized values of temperature.
This gives a set of linear values randomized along each axis, which allows calculating rotational and perpendicular values, implemented in (4), (5) and (6). Taking all these parameters into account, this allows calculating the position uncertainty of any measured coordinates on the CMM.
2.4 Part parameters
The most important part parameters in our model are the part temperature and roughness.
2.4.1 Temperature of the part
Many CMM's types of software propose a temperature compensation, which corrects the dimensions to provide results at the standard temperature (often 20 °C).
This correction assumes that the temperature and the coefficient of thermal expansion (CTE) of the material are known. It also depends on the homogeneity of the temperature inside the part and around the CMM. If the airconditioning is well set, the part must still stay a few hours in the room to get the right temperature (it depends on its mass).
The uncertainty of the correction depends on those parameters, and should be written:
$${U}_{\Delta L}^{2}={\left({u}_{L}.\alpha .\Delta T\right)}^{2}+{\left(L.{u}_{\alpha}.\Delta T\right)}^{2}+{\left(L.\alpha .{u}_{\Delta T}\right)}^{2}.$$
Considering the terms of ISO/TS 23165:2006 [10], in case of calibration for temperature compensated CMMs, the temperature measurement should not be considered as an uncertainty contributor (§A.3.2.4.3). On our CMM, the temperature is compensated, the most important uncertainty concerns the knowledge of the CTE of the measured part. Our model takes into account an uncertainty on temperature measurement as a randomized value, and an uncertainty on the CTE of ±1.10^{−6} K^{−1}.
2.4.2 Roughness of the part
The roughness of the part can be an important part of uncertainty, particularly when the manufacturing gives low surface quality (Fig. 12). The rougher the surface is, the larger the ball radius should be.
Table 2 shows the difference of flatness measured on different surfaces with different probes (real measurements).
In this study, roughness is defined by Ra, which is an usual parameter.
Using real profiles of roughness, we calculated the measured profile depending on the ball radius and a step of measure (Fig. 13).
Table 3 presents the dispersion of proof values calculated depending on the probe diameter, with a step of 10 μm along a standard rough Ra 3.2.
For this kind of roughness (Ra 3.2), using a ball diameter of 3 mm generates an uncertainty of ± 0.5 μm (k=2). If representative, this uncertainty is added to the probe radius deviation as a randomized value.
Fig. 12 Representation of the same target point on a rough surface measured with different ball radii. 
Influence of roughness and ball radius with flatness result.
Fig. 13 Influence of ball diameter on a standard profile Ra 3.2 (simulation). The measured height depends on the ball diameter and on the position of the probe in front of the real profile. 
dispersion of height measurement depending on the probe radius.
2.5 Method parameters
The method of measurement is a large part of uncertainty, unfortunately, it is hard to simulate.
We separate it into two different classes: one which is imposed and the second which needs to follow a rule.
2.5.1 Imposed parameters
Those parameters usually come from standards (ISO or ASME) and concern the algorithms to be used to fit the measured points with perfect features. It is particularly the case with associating features that represent datums. The associated features must be tangent to the set of measured points, on the free side of the material. ISO and ASME purpose two different algorithms to obtain the right result, it is the definition drawing that imposes which to be used.
Our model uses ISO association (Chebychev or minimax algorithm) to fit datums to measured points.
Unfortunately, standards do not fix the number of points to create a feature.
2.5.2 Rules
There are different rules to follow in order to minimize the uncertainty, such as the number of points per feature, their location, etc…
Those parameters should be described in the measurement result.
2.5.2.1 Number of points
The number of points depends on the size of the surface, the type of machining and the roughness.
A classical rule is to measure a minimum of three times the number of points that defines the feature. Since a plane is defined by three points, it should be measured with at least nine points.
Another rule follows in Table 4 [11].
Those rules define the location and the orientation of a feature. Defining its form needs more points. Using an odd number of points allows avoiding symmetrical issues. Figure 14 shows the influence of the selection of measured points on the flatness result. The reference flatness is obtained on a plane measuring 100 mm by 70 mm, with about 6000 points. For a same number of points, depending on their location on the surface, the flatness can vary from one to three. This figure also shows the influence of an outlier filtering (±3 standard deviations).
Recommended number of hits per feature.
Fig. 14 Influence of number of points on the flatness result for the same measured surface. 
2.5.2.2 Location of points
The location of points on the measured surface can give different results, as shown in Figure 15.
The location of points is also very important when measuring a circular part, such as a circle, a cylinder, a sphere… In those cases, measuring a feature on an angle lower than 180° increases the uncertainty a lot. Figure 16 represents the center coordinates and the radius deviations of two different circles. The first one is a perfect threelobed circle, the second one is a perfect circle with a 50 UPR defect. Both circle defects amplitude are ±2.5 μm.
The metrologist is free to select the number of points, their location, the probe diameter, the type of machining, the tolerance… Nevertheless, he must choose a set of parameters relationship to the considering accuracy.
Fig. 15 Influence of location of measured points depending on the distance between them. 
Fig. 16 center coordinates and radius deviations depending on portion angle of probing around the circle. 
2.6 Methodology of use
The proposed methodology retained to calculate uncertainties with our model consists in an application of previous parameters to the points coordinates measured on the part, in the CMM alignment. Each coordinate is modified by the compensation matrix, the probing parameters, the temperature…
The part is virtually moved in the CMM space (depending on axes strokes), 300 times: ten times along X, ten times along Y and three times along Z. This allows taking into account the different residual defects of the compensation matrix along each axis.
In addition, the compensation matrix is changed 100 times, so that the part is “measured” 30,000 times (Fig. 17).
For each measurement, the model calculates the datums, and then calculates the deviation along the surface vector of each point of a feature.
The global uncertainty is estimated by the sum of the absolute value of the average of deviations d_{i} with the standard deviation of deviations, multiplied by two, as described in [12].
$${d}_{i}=\left({\displaystyle \overrightarrow{O{M}_{{i}_{meas}}}}{\displaystyle \overrightarrow{O{M}_{{i}_{th}}}}\right).{\displaystyle \overrightarrow{{n}_{{i}_{th}}}}$$(15)
$${U}_{measj}=\pm {\left(\left{\displaystyle \stackrel{\u203e}{{d}_{i}}}\right+2{\sigma}_{d}\right)}_{measj}.$$(16)
Fig. 17 logic diagram of model implemented using Python. 
3 Results
First, we applied our methodology to elementary features, such as planes or circles, by entering perfect cloud of points. Then we extended the application to real cases.
In addition to the points' coordinates, the text file contains the vector associated with each point, the rotation head angles and a digit defining the type of feature.
3.1 Case of a plane
The cloud contains 900 points, spaced by 5 mm (Fig. 18). We also considered a case by filtering the dataset to keep only 25 points. The points define a perfect plane.
As it concerns a form defect, we use a refolded normal distribution, flatness cannot be negative (Fig. 19).
In this perfect plane case, the uncertainty of flatness could not be less than ±2.2 μm (k = 2).
The same approach on a measured plane (real coordinates) gives 56 μm ± 4 μm (k = 2). This result is obtained on a plane measured with 35 points, equally spaced on the surface, with a ball diameter of 5mm, using ISO 1101 standard [13].
Fig. 18 Examples of flatness results calculated from coordinates of a perfect plane. 
Fig. 19 statistical representation of 30,000 calculations of flatness on a perfect plane with the virtual CMM. 
3.2 Case of a circle
In addition to the form, the circle gets its own dimension, so the model can calculate another kind of uncertainty.
A perfect circle of Ø 150 mm measured horizontally on the CMM by 36 points gives those uncertainties:
2 μm ±1.2 μm for the roundness (k = 2)
150.001 mm ±0.8 μm for the least squares diameter (k = 2).
A real circle measured on an angle of 35° gives:
6.8 μm ±1.5 μm for the roundness (k = 2)
129.08 ±10 μm for the least squares diameter (k = 2).
As described in the rules, measuring a circle on a short portion gives a large uncertainty.
3.3 Case of a calibration sphere
In this application, we focus on the uncertainties coming from the calibration of different probes on the calibration sphere. In this case, the sphere coordinates do not change, the variation is applied to the rotation head angles nine times.
The result is a standard deviation of the measured sphere radii for each angle, and the defects of the sphere localization, which correspond to an offset variation (Fig. 20).
The results for a probe of 50 mm of length and 5 mm of diameter are:
±0.4 μm for the roundness
±0.8 μm for the location.
This example shows the importance of recalibrating the probe after a rotation when looking for maximal accuracy.
Fig. 20 Representation of points measured on the calibration sphere with different positions of the rotation head. 
3.4 Simulating an acceptance test
We also applied the model to the measurement of a step gauge in ISO 10360 conditions [14].
The cloud of points contains coordinates of five lengths of the gauge, in seven different positions (Fig. 21).
This application demonstrates the residual defects of the CMM using our model. As an example, the perpendicularity between axes during a simulated test, calculated from the four volumes measurements gives a defect of 6 μrad on one of them, while all the values are inside the tolerance of the acceptance test.
Fig. 21 example of simulated acceptance test. 
3.5 Comparison to other values
We also applied the model to a part which had been used by the past to different approach of uncertainty estimation [15,16].
This part is an ellipsoid, as shown in Figure 22. The uncertainty concerns the profile measurement of part on both sides, internal and external.
This part has been measured on more than 20 different CMMs, using different softwares, during an intercomparison. It defined a first approach of the uncertainty on profile measurements [15].
The uncertainty is calculated on the two profile dimensions, inner and outer. The interest of differencing both sides is to take into account the different probe offsets (Figure 23). The internal surface is measured with one vertical probe, the external surface is measured by three different horizontal probes. This induces another source of variability.
Figure 23 shows the deviations of the profiles in relation to two datums, a plane as primary reference and a circle as secondary reference. Those surfaces are processed by a model in accordance with ISO 5459 [17]. That explains why the points at the top of the part show a zero deviation, corresponding to two perfect features associated to the measured points.
In addition to previous values, the part data has also been tested with the software VCMM2 [18].
Table 5 shows the different results.
Our model seems to be close to VCMM2 and “practical method” for the internal measurement, and seems to reduce the uncertainty on the external part.
Fig. 22 intercomparison part used to estimate uncertainties [15;16]. 
Fig. 23 example of profile measurement treated by our virtual CMM on the intercomparison part. 
Comparison of different values of uncertainty coming from different evaluations.
4 Conclusion
We built a model of virtual CMM using Python. This model is based on an analysis of linear displacements of the CMM, which allows calculating rotational and perpendicular defects, so that covariance can be ignored.
This model is easily transposable to another CMM, because it depends on CMM strokes and distances between air bearings, assuming that the CMM compensation is perfect, except on the calibration uncertainty.
The results obtained seem to show an acceptable correlation to other results on a same part. The model has to be tested on more complicated parts.
Funding
This research received no external funding.
Conflicts of interest
JF. Manlay, A. Charki and A. Delamarre certifie that they have no financial conflicts of interest (eg., consultancies, stock ownership, equity interest, patent/licensing arrangements, etc.) in connection with this article.
Data availability statement
Data associated with this article cannot be disclosed due to legal reasons.
Author contribution statement
Conceptualization, JF. Manlay; Methodology, JF. Manlay; Software, JF. Manlay; Validation, JF. Manlay, A. Charki and A. Delamarre; Formal Analysis,JF. Manlay, A. Charki and A. Delamarre; Investigation, JF. Manlay, A. Charki and A. Delamarre; Data Curation, JF. Manlay; Writing – Original Draft Preparation, JF. Manlay; Writing – Review & Editing, JF. Manlay, A. Charki and A. Delamarre.
References
 ISO/TS 15 5304: 2008 Geometrical Product Specifications (GPS) − Coordinate measuring machines (CMM) − Technique for determining the uncertainty of measurement − Part 4: Evaluating taskspecific measurement uncertainty using simulation. [Google Scholar]
 D. Heisselmann, M. Franke, K. Rost, K. Wendt, T. Kistner, C. Schwehn, Determination of measurement uncertainty by MonteCarlo simulation, WSPC Proceedings 361–371 (2019) [Google Scholar]
 N. Van Gestel, Determining measurement uncertainties of feature measurements on CMMs, Thesis, Leuven (2011) [Google Scholar]
 P. Salacheep, Error compensation and uncertainty evaluation of CMMs based on kinematic error models and Gaussian processes, Thesis, London (2016) [Google Scholar]
 S. Hammad Mian, A. AlAhmari, New developments in coordinate measuring machines for manufacturing industries, Int. J. Metrol. Qual. Eng. 5, 101 (2014) [CrossRef] [EDP Sciences] [Google Scholar]
 Gaska, W. Harmatys, P. Gąska, M. Gruza, K. Gromczak, K. Ostrowska, Virtual CMMbased model for uncertainty estimation of coordinate measurements performed in industrial conditions, Measurement 98, 361–371 (2017) [CrossRef] [Google Scholar]
 Charki, K. Diop, S. Champmartin, A. Ambari, Numerical simulation and experimental study of thrust air bearings with multiple orifices, Int. J. Mech. Sci. 72, 28–38 (2013) [CrossRef] [Google Scholar]
 T. Coorevits, Contribution au développement de techniques d'autocalibrage appliquées aux machines à mesurer tridimensionnelles, Thesis, Lille (1990) [Google Scholar]
 J.F. Manlay, Fast checking of CMM geometry with a patented tool, in 17th International congress of metrology, Paris, Proceedings (2015) [Google Scholar]
 ISO/ TS 2316 5: 2006, Geometrical Product Specifications (GPS) − Coordinate measuring machines (CMM) test uncertainty [Google Scholar]
 J. Berthold, D. Imkamp, Looking at the future of manufacturing metrology: roadmap document of the German VDI/VDE Society for measurement and automatic control, J. Sens. Sensor Syst. (2013) [Google Scholar]
 E4136110N, Moyens de production − Agrément capabilité des moyens de mesure − CNOMO [Google Scholar]
 NF EN ISO 1101 April 2017, Geometrical Product Specifications (GPS) − Geometrical tolerancing − Tolerances of form, orientation, location and runout [Google Scholar]
 NF EN ISO 103602, January 2010, Geometrical Product Specifications (GPS) − Acceptance and reverification tests for coordinate measuring machines (CMM) − Part 2: CMMs used for measuring linear dimensions [Google Scholar]
 J.F. Manlay, CMM: why doing a new intercomparison? in 13th International Congress of Metrology, Lille, Proceedings (2007) [Google Scholar]
 J.F. Manlay, Une approche pratique des incertitudes de mesure en 3D, AFNOR éditions, BIVI I8022 (2011) [Google Scholar]
 NF EN ISO 5459, November 2011, Geometrical Product Specifications (GPS) − Geometrical tolerancing − Datums and datums systems [Google Scholar]
 G. Theveneau, M. Cadignan, G. Lanier, Incertitudes de mesure sur une machine à mesurer tridimensionnelle, revue Chocs Focus 8, 22–23 (2023) [Google Scholar]
Cite this article as: JeanFrançois Manlay, Abdérafi Charki, Anthony Delamarre, A virtual CMM to estimate uncertainties, Int. J. Metrol. Qual. Eng. 15, 21 (2024)
All Tables
All Figures
Fig. 1 Representation of CMM defects along an axis during its displacement. 

In the text 
Fig. 2 Compensation of the position of P by applying a field of moment on each axis one by one. 

In the text 
Fig. 3 Scheme of the bridge air bearings structure representing a sliding connection. 

In the text 
Fig. 4 Representation of the link between the flatness of the guide ways and the rotation defects of the bridge 

In the text 
Fig. 5 Relationship between flatness defects on guide ways, air bearing positions and rotational defects along X and Y axes. 

In the text 
Fig. 6 Relationship between flatness defects on guide ways, airbearing positions and rotational defects along Z axis. 

In the text 
Fig. 7 Relationship between straightness of Y axis measured along Z axis and perpendicularity defect between Y and Z. The slope defines the compensation value for the perpendicularity defect. 

In the text 
Fig. 8 : Trueness calibration using an interferometer by comparison to the scale readings. 

In the text 
Fig. 9 Representation of a temperature gradient between both faces of the granite. The expansion difference allows calculating the versed sine. 

In the text 
Fig. 10 Representation of a probe calibration. The results are a calculated ball radius and an offset deviation from theoretical position and a standard deviation of radii. 

In the text 
Fig. 11 Repositioning defects on a rotation head without calibration. The uncertainty depends on the length of the probing system. 

In the text 
Fig. 12 Representation of the same target point on a rough surface measured with different ball radii. 

In the text 
Fig. 13 Influence of ball diameter on a standard profile Ra 3.2 (simulation). The measured height depends on the ball diameter and on the position of the probe in front of the real profile. 

In the text 
Fig. 14 Influence of number of points on the flatness result for the same measured surface. 

In the text 
Fig. 15 Influence of location of measured points depending on the distance between them. 

In the text 
Fig. 16 center coordinates and radius deviations depending on portion angle of probing around the circle. 

In the text 
Fig. 17 logic diagram of model implemented using Python. 

In the text 
Fig. 18 Examples of flatness results calculated from coordinates of a perfect plane. 

In the text 
Fig. 19 statistical representation of 30,000 calculations of flatness on a perfect plane with the virtual CMM. 

In the text 
Fig. 20 Representation of points measured on the calibration sphere with different positions of the rotation head. 

In the text 
Fig. 21 example of simulated acceptance test. 

In the text 
Fig. 22 intercomparison part used to estimate uncertainties [15;16]. 

In the text 
Fig. 23 example of profile measurement treated by our virtual CMM on the intercomparison part. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.