Automatic Calibration of a Two-Axis Rotary Table for 3D Scanning Purposes

Rotary tables are often used to speed up the acquisition time during the 3D scanning of complex geometries. In order to avoid manual registration of the point clouds acquired with different orientations, automatic algorithms to compensate the rotation were developed. Alternatively, a proper calibration of the rotary axis with respect to the camera system is needed. Several methods are available in the literature, but they only consider a single-axis calibration. In this paper, a method for the simultaneous calibration of both axes of the table is proposed. A checkerboard is attached to the table, and several images with different poses are acquired. An optimization algorithm is then setup to determine the orientation and the locations of the two axes. A metric to assess the calibration quality was also defined by computing the average mean reprojection error. This metric is used to investigate the optimal number and distribution of the calibration poses, demonstrating that the optimum calibration results are achieved when a wider dispersion of the calibration poses is adopted.


Introduction
Image-based modeling is the process of reconstructing the three-dimensional geometry and the surface appearance of a real object from its two-dimensional images. In typical applications [1] images are acquired in a controlled environment. In order to achieve a 360°view of a complex object, the scanning must be repeated from different perspectives. Then, acquisition performed from different viewpoints needs to be registered, to reconstruct the full object geometry. The easiest solution to this issue is the manual point cloud registration, which is generally based on the 3-2-1 procedure: the first rough registration is manually provided by selecting three or more corresponding points on the clouds. Then, an iterative closest-point algorithm is run to refine the manual registration [2]. This procedure is robust and versatile, and it does not impose any constraints in terms of scanning perspective, but it requires long processing times. The procedure can be accelerated by gluing markers on the target surface, which can help the user to manually select the corresponding points, or can even allow for automatic algorithms [3]. Nevertheless, the marker placement may be cumbersome and also impossible to use where the surface cannot be altered. An alternative solution to this problem is represented by the use of turntables, which have proven to be effective in the reconstruction of complex shapes [4,5]. The object to be scanned can be placed on a turntable, and a stereo-camera system acquires images, while both orientation and placement of the item are changed within the work space by rotating the table by known angles. The use of a rotary table imposes a constraint to the possible scanning orientations, which can help the development of automatic algorithms to register the point clouds [6,7]. In addition, knowing the precise position and orientation of the turntable axis can further improve the registration step [8][9][10]. To this extent, a proper calibration of the acquisition system is necessary to obtain satisfying results.
The camera calibration problem has been extensively studied by researchers, and several approaches have been proposed. Linear methods [11,12] adopt a simple pinhole camera model and incorporate no distortion effects. Despite the algorithm being fast and non-interactive, neglecting camera distortion implies that lens distortion effects cannot be compensated for. Conversely, non-linear techniques [13][14][15] use general camera models, which also incorporate distortion parameters. At first, a relation between parameters is established, and then an iterative solution is found by minimizing properly defined error terms. Typically, the iterative procedure requires an initial guess, which often dictates the quality of the convergence. However, the most common approach is represented by two-step techniques [16][17][18], where in the first step, a direct solution of some parameters is given, and then an iterative procedure is followed in order to compute the other parameters as well as to reduce the error of the direct solution. The research on stereo camera systems calibration pushed forward the calibration of non-conventional optical setups based on mirrors (catadioptric systems), such as [19][20][21][22] which solve the problem of spherical mirror calibration to achieve an omnidirectional scanning.
Although technical literature well addresses the camera calibration problem, the extraction of the rotation axis of the turntable with reference to the fixed camera coordinate frame still deserves attention. Many different approaches have been found dealing with single-axis calibration: Mülayim et al. [23] proposed a multi-image approach capable of extracting the rotation axis of a single-axis turntable with respect to the camera coordinate frame. Li et al. [24] calibrated the rotation axis by means of a test sphere with a known diameter mounted on the turntable, similarly to [25]. Bi et al. [26] used a cuboid block and edge detection to assess the rotation axis. Langming et al. [27] developed a procedure based on the scanning of a cylindrical specimen. Ye et al. [28] used a planar specimen acquired in different orientations, basing the calibration on the intersection lines between subsequent planar surfaces. Chen et al. [29] formalized the problem as a constrained global optimization problem. However, these techniques do not generalize when a multi-axis turntable is considered, since each axis has to be calibrated separately, one after the other. Furthermore, no strict information about the placement of the calibration target as well as the number of images to be acquired is given in order to obtain acceptable results.
In this study, a novel approach for the calibration of a two-axis turntable is proposed. At first, a mathematical model for the rotation of a point attached to the turntable has been formalized and exploited in order to provide a minimum least-squares estimation of the two rotation axes as well as their distance with respect to the reference camera. Afterwards, a further analysis, based on optimal conditions to calibrate the rotary table, was conducted in order to optimize the estimates of the calibration process. Optimal conditions are considered to be the number and the orientation of placements of the calibration board with respect to the stereo-camera system. Qualitative criteria are finally provided based on experimental results.

Materials
With reference to Figure 1 the image acquisition system consisted of a fixed stereo-camera system, a rotary table, and a computer. The stereo-camera system, shown in Figure 1a, is characterized by two digital cameras with a resolution of 1600 pixel × 1200 pixel (The ImagingSource® DMK 51BU02), equipped with lenses having a 16 mm focal length, and a 1024 pixel × 768 pixel multimedia DLP projector (OPTOMA EX330e, resolution XGA 1024 pixel × 768 pixel). Figure 1b shows the two-axis rotary table developed in this work. The overall system schematic, reporting the main reference frames, is represented in Figure 2. The first axis is parallel to the camera image plane while, by construction, the second is perpendicular to the first; the two axes are also assumed to intersect each other. It is reasonable to make such an assumption, since the rotary table was designed as a CNC (Computer Numerical Control) machine with high dimensional accuracy. An asymmetrical checkerboard of 7 × 10 squares (each square having size of 12 mm × 12 mm), resulting in M = 54 inner corners, was attached Sensors 2020, 20, 7107 3 of 21 to the turntable during the image acquisition process. This pattern allows for extraction of the point sets as they are placed in corners on the squares (excluding those on the external side of the checkerboard). Therefore, information on the table configuration, through the arrangement of the points in space, can then be deduced.  with the addition of the Adafruit Motor Shield V2. The overall system was controlled via MATLAB ® on a PC platform where a Graphic User Interface (GUI) has been developed in order to simplify the entire workflow, providing an interface that allows to control the settings of cameras, image acquisition, table calibration, and analysis results. The software is available as a free package at [30], and a brief description of the GUI is given in Appendix A.

Calibration Approach
The geometrical construction of the system was exploited for the implementation of a suitable calibration procedure of the two-axis turntable. We assume that the centers of the tables are disposed such that their axes intersect perpendicularly at point q. That bears the advantage to consider a single point of application of rotation for both axes. Hence, just a single position vector from the cameras is enough. Since the orientation of the camera system, chosen as a global reference, was different from those of the tables, axis-angle notation was used to describe the rotation of each table.
The general rotation of a point around an axis outside of the origin, as shown in Figure 3, is represented as follows: where p is a generic point positioned on a system rotating around an axis with unit vector ω going through point q, p is the rotated point, and R ω (θ) is the axis-angle parametrization of the rotation matrix in which θ is the angle of rotation. Equation (1) can be expressed as a rigid-body transformation in homogeneous coordinates as follows: ibration Approach he geometrical construction of the system was exploited for the implementation of a suitable tion procedure of the two-axis turntable. We assume that the centers of the tables are disposed hat their axes intersect perpendicularly at point q. That bears the advantage to consider a single of application of rotation for both axes. Hence, just a single position vector from the cameras ugh. Since the orientation of the camera system, chosen as a global reference, is different from of the tables, axis-angle notation was used to describe the rotation of each table. he general rotation of a point around an axis outside of the origin, as shown in Figure 3, esented as follows: p is a generic point positioned on a system rotating around an axis with unit vector ω going gh point q, p is the rotated point and R ω (θ) is the axis-angle parametrization of the rotation in which θ is the angle of rotation. Equation (1) can be expressed as a rigid-body transformation ogeneous coordinates as follows: quation (2) can be applied to the studied system. As shown in Figure 2, the axis of the horizontal s defined as ω 1 , the axis of the vertical as ω 2 , θ 1 is the angle of rotation around ω 1 , and θ 2 the around ω 2 . Since ω 2 lies in the reference system of the horizontal table, the global transformation be defined in the following order: Equation (2) can be applied to the studied system. As shown in Figure 2, the axis of the horizontal table is defined as ω 1 , the axis of the vertical as ω 2 , θ 1 is the angle of rotation around ω 1 , and θ 2 the Sensors 2020, 20, 7107 5 of 21 angle around ω 2 . Since ω 2 lies in the reference system of the horizontal table, the global transformation must be defined in the following order: The axis-angle representation provides the advantage of being a direct function of the components of the axis unit vector for both tables. The two axis unit vectors ω 1 and ω 2 , together with the global coordinates of their intersection point q from the camera frame, denoted in Figure 2 as {S c }, are used as input parameters for an optimization procedure that finds the solutions through a nonlinear least-squares approach.
First, an initial reference pose X R is acquired, which is associated to angle configuration (θ 1 = where both X i and X R are expressed in homogeneous coordinates. However, since measured point clouds X i are considered, which are affected by various error sources, the following equations can be written: . . where δX i ∈ R 4×M (i = 1, . . . , N) are the residual errors between measured and ideal coordinates of corresponding point clouds originated by X R , and where M are the points (extracted via software) on the checkerboard. Finally, a nonlinear least-squares estimation of the parameters is established to minimize the mismatch between measured and ideal point cloud coordinates. The optimization function is therefore based on Equation (5), where the table rotation angles are known by counting the steps provided to the stepper motors (since the rotational speed is low, it is reasonable to neglect step losses). On the other hand, q, ω 1 , and ω 2 are unknown and should be estimated. A formalization of the nonlinear least-squares estimation, minimizing the sum of squared Euclidean error norms, is provided as follows: where i is the pose index examined (i = 1, . . . , N), j is the point index in the checkerboard (j = 1, . . . , M), and p i j is the corresponding rotated point (Cartesian coordinates) measured on the checkerboard. Point p i j (ω 1 , ω 2 , q; (θ 1 , θ 2 ) i ) is obtained from the rotation of the reference pose by θ 1,i around axis ω 1 and θ 2,i around axis ω 2 , both intersecting at the point q. The optimization must also satisfy the following constraints: where (7a) indicates the perpendicularity between the axes, whereas (7b) and (7c) impose the condition for ω 1 , ω 2 of having unit length. For instance, the proposed method was employed for the calibration of a two-axis rotary table, using the hardware previously described in Section 2; hence, M = 54 points per grid were considered. At first, 100 stereo-pair images were acquired in different poses (θ 1 , θ 2 ) i , i = 1, . . . , N of the rotary table. Then, using an estimation of (ω 1 , ω 2 , q) yielded by (6) applied to all the N = 100 stereo-pair images, and by taking the inverse of (4), the measured grids X i , i = 1, . . . , N were reprojected onto the reference grid X R . Finally, the N · M point-to-point errors committed were computed. Figure 4 shows the error cloud in the 3D space, expressed in the camera frame {S c }, whereas Figure 5 shows these errors in the three Cartesian planes. The absolute error was lower than 1 mm in every direction.

Influence of the Number of Calibration Images
Calibration quality depends both on the number K of poses used and on their spatial location, i.e., (θ 1 , θ 2 ) i , i = 1, . . . , K. Hence, it is useful to investigate the influence of the number K of images to be captured towards a quality and fast calibration. To this extent, a proper metric was introduced, as discussed further in this section.

Data Set Definition
Any i-th image, i = 1, . . . , K, can be captured in the rotary table joint configuration Space S, . This is the space any pair of angular positions (θ 1 , θ 2 ) i , i = 1, . . . , |S|, belongs to, i.e., S = [θ 1,min , θ 1,max ] × [θ 2,min , θ 2,max ]. In order to assess the influence of the number of images K, a test run was conducted. First of all, a set of images was captured in S (Acquisition set A, A ⊂ S). Then-excluding an image S R which is taken as Reference-it was partitioned into two subsets: the Calibration set C, C ⊂ A, and the Test set T, T ⊂ A, so that C ∪ T = A \ S R ∧ C ∩ T = ∅. C is the set from which to select the desired subset C H of K images to calibrate the turntable with, whereas T is the set subsequently used for quantifying the calibration error achieved with C H . With respect to C, for each K, K = 1, . . . , |C|, a series of b K subsets C H , C H ⊆ C could be extracted, where H is a set of indices i, K is the size of H, and b K is the following binomial coefficient: In the present study, a S = [−36°, 36°] × [−90°, 90°] was considered. Such S ensures that the checkerboard remains visible from the stereo-camera system during the entire acquisition process. To discretize S, a uniform sampling was performed, obtaining 101 images: the first one was associated to (0°, 0°) and was taken as reference, whereas the remaining 100 were equally spaced in steps of size (δθ 1 , δθ 2 ) = (8°, 20°) within S, as in Figure 6. Lastly, the obtained A was partitioned into two equally sized and spaced sets C and T. In more detail: . . , N C , that denotes the pairs of angular positions (θ 1 , θ 2 ) i , i = 2r + 1, of the subset of images to calibrate with; • K = |H|, K ∈ {1, . . . , N C }, is the size of H, i.e., the number of images considered; The i-th pair of angular positions is given by where      m = i + 8 10 n = m mod 2 · (i + 9 − 10m) + (m + 1) mod 2 · (10m

Error Metric Definition
The checkerboard, shown in Figure 7, consists of 7 × 10 squares, so the number of inner-corners M is (7 − 1) × (10 − 1) = 54. Figure 7 shows the turntable with the attached checkerboard in two different poses: Figure 7a is the reference pose at (θ 1 , θ 2 ) 1 = (0°, 0°), whereas Figure 7b represents a generic i-th configuration, which is taken as target at (θ 1 , θ 2 ) i . Given the i-th pose expressed in the coordinate frame {S c }, i = 1, . . . , N A , and recognizing the i-th Equation in (5), it is possible to rotate by (θ 1 , θ 2 ) i its measured target grid X i (blue dots in Figures 7b and 8a) onto the reference grid X R = X 1 (red dots in Figures 7a and 8), thus obtaining the reprojected grid X R = X 1 , as shown in Figure 8b. 1.5 Hence, given a subset C H within the Calibration set C of K images, and based on the results described in the previous section, it is possible to calibrate the rotary table, obtaining an estimate (ω 1 ,ω 2 ,q) H of (ω 1 , ω 2 , q). A series of metrics was arranged to quantify the reprojection accuracy obtained to reproject all N T grids of the Test set T onto the reference grid X 1 using the calibration results acquired with the subset C H of K images. The subscript H highlights the dependency of the two axes, ω 1 and ω 2 , and of their intersection point, q, on the particular subset of images used, thus (ω 1 ,ω 2 ,q) H . Considering only the grids belonging to T, for a given H, the point-to-point error between a j-th measured point p 1 j on the reference grid X 1 and the reprojection of the j-th point The arithmetic mean of j,2s,H over M points yields the 2s-th averaged point-to-point error committed using the (ω 1 ,ω 2 ,q) H that was estimated through the calibration previously made with the grids X h , Finally, averaging all of the N T averaged point-to-point errors over the entire Test set T gives the T-averaged point-to-point error, yielding the reprojection accuracy: Additionally, the root mean square (RMS) of these errors was defined as The number of possible subsets C H of size K that can be employed for calibration is given by the binomial coefficient b K in (8). First of all, in order to have a rough estimation of the effect of K on the calibration accuracy, and to reduce the computational time, for each K, K = 2, . . . , N C only N R K = 1000 combinations of C H are randomly extracted, obtaining N C − 1 sets R K , and N R K is the size of R K , i.e., N R K = |R K |. The case K = 1 is not considered from now on since it is not sufficient to calibrate a two-axis rotary table. Then, all the R K sets, K = 2, . . . , N C , were used to execute (N C − 1) · N R K = 49, 000 independent table calibrations. Subsequently, their respective T-averaged point-to-point errors H given by Equation (13), H | C H ∈ R K , were computed. Finally, averaging them over R K , where K = 2, . . . , N C , returns N C − 1 R K -averaged point-to-point errors: Additionally, the root mean square (RMS) of these errors was defined as An overview of the errors, and their RMS, is presented in Figure 9 for K = 2, . . . , N c . As can be seen, both E K and σ K decrease when K increases. While the effect of K on E K seems to saturate for K greater than 10, it can be noted that E K still decreases with K. Thus, if the poses are randomly chosen, greater values of K are more likely to produce high-quality calibrations. Nevertheless, even a small K can provide high-quality calibration results, if the poses are properly chosen as demonstrated in the following.

Methods Comparison
In order to verify its effectiveness, the developed algorithm was compared with the circles method proposed by Chen et al. [29]. Hence, several combinations of images were extracted from the Calibration Set, and the calibration of the two-axis rotary table was carried out using both methods. It is worth mentioning, being that the circles method is a single-axis procedure, the vertical axis ω 1 was firstly calibrated, and subsequently the horizontal axis ω 2 was considered. For each combination H of images, the T-averaged point-to-point error H and its standard deviation σ H were computed by considering all the images of the Test Set. The results are reported in Table 1, comparing the two methods for the tested calibration sets. It is possible to observe that the proposed method returned almost constant values of H and σ H (i.e., 0.33-0.36 mm and 0.16-0.18 mm, respectively), which slightly decreases with the increasing of K. On the contrary, a different behavior was noticed in the circles procedure case, where H decreased but not σ H , since this method is influenced appreciably by the specific calibration poses. In addition, the range of values of both H was considerably wider (i.e., 0.52-1.70 mm) than the results of the proposed method. Hence, the proposed method seems to be more robust and reliable with respect to the circles algorithm.

Optimal Calibration Set
On the basis of the data obtained in the previous section, a visual inspection was conducted to find out whether a certain pattern in the joint configuration Space S yielded preferred calibration results. No recurrent pattern (e.g., "X"-like, cross-like, circle-like, square-like, etc.) was evident among either of the studied subsets, which had led to quality calibrations of the tables or among the ones with poorer calibrations. Nevertheless, a direct correlation between dispersion of subsets and calibration quality has been found: for any given H, the more uniformly distributed the cloud of K points (θ 1 , θ 2 ) h ∈ C H ⊂ S, h ∈ H, the lower its calibration error H . For instance, in Figure 10, two opposite situations are reported. A wide dispersion can be noted in (a), corresponding to H = 0.3445 mm, while a lower dispersion can be noted in (b), which corresponds to H = 0.4239 mm.

Complete Graph Index
To quantify the scattering of any desired set of arbitrary dimension K, an index (I CG ) was devised based on the calculations of the arithmetic mean of all the weights of a weighted complete graph (CG). The subset of Figure 10a can be considered as an example. If its points were regarded as the vertices of a graph, it would be possible to think of their edges as being weighted by their Euclidean distances; hence, they can be summed and finally averaged by the number of links. Among all the weighted undirected graphs, the weighted complete graph (CG) was chosen, since it has the property to take into account all the pairs of points. In Figure 11 the CG constructed from the subset of Figure 10a is shown; it is possible to see its K points (nodes of CG, in red) and their K 2 links (edges of CG, in blue). Given K points v t = (θ 1 , θ 2 ) t , t = 1, . . . , K, it is advantageous to normalize them so that the index does not depend on the length of S ([θ 1,min , θ 1,max ] × [θ 2,min , θ 2,max ]). The normalized angular positionŝ v t = (θ 1 ,θ 2 ) t , t = 1, . . . , K, are obtained by applying the following normalizations: (θ 1 ,θ 2 ) t = (θ 1 , θ 2 ) t − (θ 1,min , θ 2,min ) (θ 1,max , θ 2,max ) − (θ 1,min , θ 2,min ) , where (θ 1 ,θ 2 ) t are the normalized coordinates ofv t along the ω 1 -axis and the ω 2 -axis, respectively. The weighted complete graph of order K considered is the pair (G, w), where . . , K} is the set of vertices (normalized angular positions); and the weight function is the Euclidean distance between any pair of vertices: The I CG index reported in (19) sums the weight function w evaluated in every pair of normalized points (v t ,v u ), t, u = 1, . . . , K, u > t and K > 1 (the case K = 1 must be avoided since at least two images are required to calibrate a two-axis rotary table), then divided by √ 2 times the number of pairs (|E|): Sensors 2020, 20, 7107

of 21
The dividing factor √ 2 |E| is useful for normalizing the index value with respect to the number of chosen images K, so that I CG ∈ [0, 1] ∀K ∈ N, K > 1. For instance, the subset in Figure 10a has a I CG = 0.4036, while the less scattered subset in Figure 10b gives a I CG of 0.2684.

Index-Based Analysis and Experimental Results
Let us now consider the joint configuration Space S, the Acquisition set A, the Calibration set C, the Test set T, and all the b K subsets C H ⊆ C presented in Section 4, as in Figure 6. Given a C H set of K images, the index is In order to examine the relationship between H and I CG,H , K = 7 has been selected, since it guarantees a relatively low binomial coefficient (Equation (8)), being b 7 = 50 7 = 99884400. Firstly, I CG,H is calculated for all the b 7 possible combinations of C H , obtaining the array I b 7 . Then, this is sorted in ascending order of I CG,H , obtaining the ordered array I b 7 . In order to reduce the computing time, I b 7 was downsampled by a factor of ∆ = 1000, obtaining N = b 7 ∆ = 99885 calibrations (ω 1 ,ω 2 ,q) H , equally spaced between min I b 7 and max I b 7 . The corresponding T-averaged point-to-point errors H were computed, as expressed in Equation (13), and shown in Figure 12a for all the N unsorted calibrations. Finally, these errors are graphed in function of I CG,H . Such plot is shown in Figure 12b: increasing the index, the errors tend to decrease. Moreover, it is clear that the greater the value of I CG,H , the lower the expected value of H and its variance.- The same analysis was conducted for K = 2, K = 4, and K = 46 to assess whether the parameter I CG,H was reliable regardless the value of K. At first, I CG,H was computed and used for ordering all the b 4 , b 46 , and b 2 combinations, respectively. Then, at every K, the N C H subsets were selected by uniformly downsampling by a factor-summed and rounded to the nearest integer index-of ∆ = b 4 N , b 46 N , and 1 among all the b 4 , b 46 , and b 2 combinations of C H , respectively. Finally, each combination of images was used for estimating (ω 1 , ω 2 , q), thus obtaining N (ω 1 ,ω 2 ,q) H . With reference to Figure 13, it is possible to observe that the trend of H , for K = 2, 4, 46, was the same as for K = 7 (Figure 12b), provided that the results were sorted in ascending order of I CG,H . It is worth noting that the higher the value of K, the more flattened the error range. In fact, while for K = 2 the errors ranged from around 0.3 mm to 3.7 mm, for K = 46 they were between 0.344 mm and 0.346 mm (as could be predicted considering that σ 2 σ 46 ). Additionally, Figure 13 confirms the effectiveness of I CG,H as a method to find the sets which yielded better calibration results, i.e., lower calibration errors. Moreover, considering as a reference the E K obtained by using all the N C images, i.e., E 50 = 0.3446 mm (Figure 9), and given a subset H, it is possible to define the percent relative error as It is worth noting that, if the H corresponding to the highest I CG,H was considered, η H was noticeably low for all the tested K, as shown in Table 2. As expected, the higher K, the lower η H . However, η H was lower than 10 % also for K = 2. In particular, for K = 2, there were four subsets corresponding to the highest I CG,H : for the worst η H = 8.72 %, whereas for the best η H = 1.28 %. Hence, high-quality calibration results can be achieved even with K = 2, provided that the proposed index I CG,H is used to choose the calibration set.  It is reasonable to use the index as a criterion for selecting the angular positions (θ 1 , θ 2 ) h within the joint configuration Space S, that is the K placements of the checkerboard. The main advantage is that it is able to predict whether a set is going to achieve an effective calibration before acquiring the images. An additional benefit is that it is less time-consuming for a computer, compared to an entire calibration. For instance, the code was tested using a computer equipped with Intel ® i7-6700HQ CPU (Dell Inc., Round Rock, Texas, USA, 4 physical cores at 2.60 GHz and 8 logical cores), and 32 GB of RAM, then running MATLAB ® R2020a (The MathWorks Inc., Natick, MA, USA) with parallel pool active on Microsoft Windows ® 10. It was observed that the time required to calculate b 7 times the index I CG,H was about 10 times shorter than the time needed for N calibrations. Thus, it is possible to assume that computing the index is 10 · ∆ = 10000 times faster than repeating the calibration for all the b 7 possible configurations.

Conclusions
This paper presents a novel approach for the automatic calibration of two-axis rotary tables through optical stereo systems. The proposed optimization algorithm allows for simultaneously calibrating both axes by imaging and processing a checkerboard attached to the rotary table in different poses. A test set of images was defined, which were not exploited during calibration, to assess the quality of the achieved calibration parameters by computing an average reprojection error and its RMS. Moreover, this paper investigated the effect of the number of calibration images on the results. It was shown that the higher the number of images, the lower the achieved mean error and RMS, if the calibration poses are randomly selected. Nevertheless, it was demonstrated that a proper choice of calibration poses can provide high-quality calibration results even with a limited number of images (thus, a shorter calibration time). An index based on the weighted complete graph, I CG , was defined, demonstrating that higher values of I CG (i.e., wider dispersion of the calibration poses) correspond to better calibration results. Thus, once the joint configuration Space S is defined, and given the number of calibration images to be used, it is possible to determine the ideal combination of poses by looking for the highest value of I CG . The experimental measurements demonstrated that, using this strategy, the worst calibration error achievable with only two calibration images is marginally higher, about 8 %, than the best calibration error achievable with 50 images, thus confirming the advantages of the proposed approach.  The Settings tab, shown in Figure A1a, allows the user to manage several hardware settings. In more detail, it is possible to • activate/deactivate the motors connected to the shield via Arduino and to specify its serial port; • activate/deactivate the camera(s) connected to the computer, adjusting several parameters such as exposure, image format, image resolution, and frame rate. Moreover, it is possible to test the stereo-camera system configuration performing an image preview or taking a one-shot acquisition.
The second tab, Turntable Control, is reported in Figure A1b. It allows the user to rotate the table within the working space in a target configuration determined by the pair of angular positions (θ 1 , θ 2 ). The directions of rotation are according to the rotary table model presented in Section 3, Figure 2. After reaching the desired target configuration, it is possible to set the current table pose as the reference pose X r for the image acquisition task by means of the Initialize button.
The Image Acquisition tab, shown in Figure A2a, allows the user to automatically capture the stereo-pair images within the joint configuration Space S needed for the calibration of the stereo-camera system or of the rotary table. As inputs, the user has to provide the following: 1.
the square dimensions of the checkerboard; 2. the desired image format; 3.
When the process is completed, the tab automatically creates a folder containing both the stereo-pair images and the associated input variables. It is worth mentioning that the first configuration to be acquired has to be the reference pose.
Subsequently, the Stereo-Camera Calibration tab enables the user to properly calibrate the stereo-camera system. With reference to Figure A2b, the user has to specify, as input, the path of a folder of images previously generated with the Image Acquisition tab. The software automatically recognizes the key variables and calibrates the stereo-camera system according to the Camera Calibration Toolbox for MATLAB ® . Thus, a stereo-camera calibration file is generated in the Camera Calibration folder (inside the Workspace folder) labeled with reference to the folder of images deployed in the task. Eventually, it is possible to obtain insights of a stereo-camera calibration file such as intrinsic and extrinsic parameters and geometrical displacement of both the N A grids and cameras in the 3D space by means of the Show Results button. Furthermore, the user is able to calibrate the two-axis rotary table by means of the Turntable Calibration tab, shown in Figure A3a. As inputs the user has to provide the following: 1.
a stereo-camera calibration file previously generated in the Stereo-Camera Calibration tab; 2.
a folder of images previously acquired with the aid of the Image Acquisition tab. The GUI fills the Calibration set panel with a number of checkboxes equal to the number of stereo-pair images contained in the specified folder, i.e., N A . By means of the checkboxes, the user is able to select the K stereo-pair images to be used in the calibration procedure, thus constituting the desired subset H; alternatively, the entire set of N A images is used by default. The rotary table calibration method is according to the one presented in Section 3. At the end of the process, a table calibration file is generated, and the estimated values of the unit vectors ω 1 and ω 2 , as well as the intersection point q, i.e., (ω 1 ,ω 2 ,q) H , are displayed in the tab.
The last tab offers an in-depth analysis of the errors committed during a turntable calibration. As input, only a folder of N A images-of which K are previously used for calibrating the rotary table-must be provided. At first, the computed T-averaged point-to-point error H , as well as its standard deviation σ H , are displayed. Then, as Figure A3b shows, the Error Analysis tab fills the Target grid panel with X i , i = 1, . . . , N A , Target grids. Finally, once a desired target grid has been selected by the user, a summarizing figure is displayed, as shown in Figure A4. Figure A4. A typical analysis figure displayed once a target grid (in this example, i = 2) is selected in the Error analysis tab of Figure A3b, showing the following: at the top, the reference grid (X 1 , in red) and the target grid (in blue) before (X 2 ) and after (the reprojected grid X 1 ) the overlapping, respectively; at the center, their point-to-point errors ( j,i,H , j = 1, . . . , 54 , i = 2, and for a certain H); at the bottom, the norm of these errors (| j,i,H |, blue asterisks), the T-averaged point-to-point error ( H , red line) and its standard deviation (±σ H , green lines). A (ω 1 ,ω 2 ,q) H previously obtained from the Turntable Calibration tab was used.
It is worthwhile to mention that, for each tab, various errors are properly treated by MATTA, as well as warnings, and that numerous informative messages stating the current status of the running process are displayed.