Evaluation of interpolation methods for generating maps in cultural heritage chemical applications

Journal of Cultural Heritage 62 (2023) 293–303 Contents lists available at ScienceDirect Journal of Cultural Heritage journal homepage: www.elsevier.com/locate/culher Original article Evaluation of interpolation methods for generating maps in cultural heritage chemical applications Domingo Martín a , ∗, Germán Arroyo a , Juan Ruiz de Miras a , Luis López a , María Rosario Blanc b , Philippe Sarrazin c , Juan Carlos Torres a a Department of Software Engineering, University of Granada, Spain b Department of Analytical Chemistry, University of Granada, Spain c eXaminart, CA, USA a r t i c l e i n f o Article history: Received 14 June 2022 Accepted 5 June 2023 Keywords: XRF mapping Interpolation Data analysis Painting study a b s t r a c t The different non-invasive techniques that have been developed for the study of works of art in Cultural Heritage have become an indispensable tool for researchers and practitioners. In particular, the creation of images showing the spatial distribution of chemical elements and pigments, called maps, helps give a better understanding of the artwork. While high-cost devices can measure the artwork at many positions with high resolution, the cheapest and most common devices are often used manually producing a small number of measurements. The solution is to use interpolation methods. In this article we present a statistical study of the feasibility of using interpolation, we discuss the accuracy of the results and propose the best solutions and a scheme of work. Additionally, we provide all the data developed and programs for future use. ©2023 The Author(s). Published by Elsevier Masson SAS on behalf of Consiglio Nazionale delle Ricerche (CNR). This is an open access article under the CC BY-NC-ND license ( http://creativecommons.org/licenses/by-nc-nd/4.0/ ) 1. Introduction The study of materials is a very important research area. Particularly in Cultural Heritage, it is a key element of the study and restoration of artwork. This information allows the expert, mainly of material sciences and painting conservation and restoration areas, to know what pigments were used, to identify possible deterioration or previous restorations, to date the artwork, to perform more accurate restorations, to prescribe conservation methods, etc. The evolution of technology has permitted great advancements in the study of materials, from invasive procedures that needed to take a sample of the artwork to non-invasive methods that allow study without affecting the works of art. Non-invasive techniques use illumination with radiation (varying the energy) and capture the interaction with the object in different ways. One of the most common techniques in Cultural Heritage is XRay Fluorescence (XRF), which allows an elemental and chemical analysis of materials. The data of the chemical elements obtained by the device are reviewed by the expert. This data in conjunction with the corresponding spatial positions where the measurements ∗Corresponding author. E-mail address: [email protected] (D. Martín) . were made, allow the creation of maps, which can permit the expert to identify the pigments. The creation of maps and their accuracy depends on the specifications of the scanning device: some devices are automatically controlled and placed and can provide thousands of measurements at a very high spatial resolution, the technique being often called Macro-XRF, while others are manually placed and have a lower spatial resolution, usually generating a few dozen measurements. The capabilities and quality of the devices usually depend on the price: for example, chinese XRF handheld devices can cost around 15,0 0 0 to 25,0 0 0$, a well-known brand around 50,0 0 0$, and a more complex XRF plus XRD device costs around 175,0 0 0$ (all of them without a motorized placement). 1 So, while some museums and institutions can afford high cost devices, researchers and professionals with a low budget use portable devices operated by hand or bulky devices that are more difficult to place and move. Computerized devices can produce as many measurements as the user wants, only limited by time. In the case of manually placed devices, the artwork is measured at only a few representative positions selected by the expert providing information about 1 Information obtained from Internet. https://doi.org/10.1016/j.culher.2023.06.004 1296-2074/© 2023 The Author(s). Published by Elsevier Masson SAS on behalf of Consiglio Nazionale delle Ricerche (CNR). This is an open access articleundertheCC BY-NC-ND license ( http://creativecommons.org/licenses/by-nc-nd/4.0/ ) D. Martín, G. Arroyo, J. Ruiz de Miras et al. Journal of Cultural Heritage 62 (2023) 293–303 chemical compositions only at the selected positions, having no information about the spatial distribution. One solution to this problem is the use of interpolation methods: the device produces values at some positions of a surface and the interpolation method provides the values for the rest of positions by applying a color table, what is called a map (See some examples in Figs. 7 and 8 ). The measurements can be taken at regular intervals in a grid but, as we have commented, it is more common to take them at scattered positions. Therefore, we are interested in interpolation methods that can use those disperse positions. This idea is not new and it has been used in many scientific areas. Particularly interesting and related to our goal are the methods developed in the geophysical fields [1–3] . There are also some solutions in the chemical field for Cultural Heritage, but most cases they are proprietary programs that usually accompany the scanning device, or resolve a particular problem, or the programs are not publicly available. Based on all these considerations, we want to evaluate the reliability of interpolation methods with XRF data and identify those that produce the best results. In order to achieve this we need: • To select a representative set of interpolation methods for scattered positions. We have chosen a set of methods that cover most of the families of interpolation procedures that exist, focusing on methods that have been published, that are easily available, or have enough information to be implemented. • To apply the selected methods to reference datasets. Two datasets will be used to produce the interpolation maps. • To carry out a statistical analysis of the results to check their reliability and accuracy. We will use the Mean Squared Error, MSE (see Section 3 ), to compute how good the results are and establish the best methods depending on the applications. This study will not only produce important results about the use of interpolation with XRF data, but we will also provide a procedure for working with scattered data of any source to produce reliable results. Finally, we provide all the programs and data that will allow the comparison of different solutions for future work in the area of Cultural Heritage. 2. Previous work The development of non-invasive methods for the study of materials has been increasingly important in science. One example are “common” radiography, and the computed axial tomography, where X-rays are used to produce images of the body or other objects. More generally, we can discuss different techniques depending on the type of radiation used: ultraviolet radiation (UV) and infrared radiation (IR) [4] , multi-spectral [5–8] , Raman [9] , X-ray based techniques, such as the previously commented radiography, X-Ray Fluorescence (XRF) [10] , X-Ray Diffraction (XRD) [11] , electron emission [12,13] , etc. New techniques include, for example, multi-porous polycapillary systems [14–16] or pinhole cameras [17] . Each method may provide a different type of information [18] . The apparition of all these types of techniques has brought about a radical change in the way art pieces are studied [19–32] . While interpolation, as part of numerical analysis, is not new [33–35] , we are interested in its use, particularly for the creation of maps. Some of the methods that we are going to evaluate were developed in order to solve a particular problem. For example, the Kriging method for obtaining information about the presence of gold from the observations at a few positions (the mines) is based on the work of Krige [36] . Other important examples are those related with geoscience [37] . For example, Barnes [38] developed a method for producing weather maps based on atmospheric pressure. One thing in common with all these methods is that they can work using a set of measures obtained from scattered positions. This contrasts with other interpolation techniques that require the measurements to be arranged in a grid, for instance bilinear or bicubic interpolation, which are commonly used to scale images [39–41] . Some papers that proposed a similar approach to ours are [1–3] . Willmott et al. [1] presents a generic discussion of the methods for evaluating geophysical models. They use the Mean Square Error, MSE, and different variations of the Root of the Mean Square Error, RMSE, to assess the goodness of the method. The bootstrap process is used to estimate the confidence intervals, CI. Weng [2] explores the results of using six interpolation methods when elevation data are used. Bhowmik and Cabral [3] compares three interpolation methods, Spline, Inverse Distance Weighting, IDW, and Kriging, to create maps that describe temperature. While there is an increasing use of maps in chemical analysis and particularly in Cultural Heritage, e.g. [42–44] , we are not interested in when it is used as a means but rather when it is an end in itself. We particularly want to know if maps produced by interpolation procedures are valuable and accurate. Using this criteria, the search of previous studies related with Cultural Heritage is reduced to the paper of Martín-Ramos and Chiari [45] , who presented an interesting solution based on the Nearest-Neighbor type (see Section 3.2.1 ) for the creation of maps from general chemical information, not only XRF data. This method has been used in several recent studies [46–49] . The key feature of this simple approach is that the color of the artwork in the sampled positions is also used to produce the interpolation. That is, the method not only uses the coordinates of the positions where the measurements are taken, but also their values in color space. This is relevant because it is a different approach to other methods where the measured values used to compute the interpolation are related directly to the output (e.g. atmospheric pressure is used to create atmospheric pressure maps). This approach makes sense because it assumes that color is related to pigments, and this is usually true, but can produce unrealiable results in the case where there are several layers and the measured element is not in the visible one because the XRF technique analyses all the layers of the stratigraphy. 3. Methods Our main goal is to evaluate whether it is possible to interpolate images with a low MSE that accurately represent the quantitative spatial distribution of the pigments used in an artwork from a few scattered sampled positions. The MSE is defined in this way: Given m trusted/measured values V i and m computed values ˆ V i , then MSE = 1 m  m i =1 (V i −ˆ V i ) 2 . To achieve this goal we need to perform a statistical analysis whose statistical variable is the MSE that is committed by using a reduced set of trusted values to produce the interpolated values in the other of positions. The study pipeline is shown in Fig. 1 . Given the full set of measurments, the idea is to obtain a random subset that represents the measures taken by the expert. Then an interpolation function is created using those measurements, which is in turn used to obtain the interpolated values in the remaining positions of the dataset (where we have the known information). In this way, we can compute the error commited and perform the statistical study. The different parts of this study are discussed in more detail below. 3.1. Dataset We need to establish the dataset considered as our reference for checking the interpolation methods: the trusted dataset. Ideally, we would like to obtain datasets with many measurements 294 D. Martín, G. Arroyo, J. Ruiz de Miras et al. Journal of Cultural Heritage 62 (2023) 293–303 Fig. 1. Study pipeline. performed in a regular grid at the maximum resolution of the scan device. Such a dataset would have three advantages: the more measurements the better for a statistical study, it makes all zones equally important, and it allows us to obtain configurations that simulate real cases where the measurements are not on a grid. Although it is possible to create a custom piece of art and produce a good dataset, we have considered that it is much more useful to work on data from real ancient paintings. We have decided to use two datasets: “The Man” from The Miraculous Interventions of the Jiz ¯ o Bosatsu [50] , a 13th century painted Japanese handscroll, studied in several papers [43,44] . The dimensions of the full handscrol are 1431.9 cm ×30 . 5 cm ( W ×H ). The scanning was achieved with a Tracer 5 g (Bruker) handheld XRF with a graphene window and 1.2 mm collimator connected to an MPS-400E Mobile Art Scanner (Dewitt Systems) (data taken Fig. 2. The original and measurement positions of both datasets. from Clarke et al. [43] ). It has 1314 measurements, and the chemical elements for each one are: As, Ba, Ca, Cd, Co, Cr, Cu, Fe, Hg, K, Mn, Ni, Pb, Sb, Se, Sn, Ti and Zn. The original image and the all sample positions are shown in Fig. 2 . This dataset contains a lot of measurements at a high resolution because the data was captured with a mechanized device. It is also an example of how a large number of measurements with a high level of detail makes the interpolation process practically unnecessary. “The Transfiguration”, a restored painting that belongs to a private collection. It is a copy of the painting that is located in the Vatican. It is currently being studied and allegedly belongs to Raphael’s school [47] . Its dimensions are 63.3 cm ×93 . 2 cm ( W ×H ). In the case of this painting, the number of measurements is quite large for a scanning process performed manually, a grid of 11 ×15 positions, 165 in total. A Thermo Scientific Niton XL3t XRF Analizer was used. The original image and all of the sample positions are shown in Fig. 2 . For each position we have the values of these chemical elements: As, Ba, Ca, Cd, Co, Cr, Cu, Fe, Hg, K, Mn, Ni, Pb, Sb, Se, Sn, Ti and Zn. One advantage of this dataset is that it has a lot of different colors, which is the main parameter for MHD method (see Section 3.2.1 ). This dataset has far fewer measurements than the previous one, but it is more similar to the normal case we are interested in where only a few distributed positions are sampled. 295 D. Martín, G. Arroyo, J. Ruiz de Miras et al. Journal of Cultural Heritage 62 (2023) 293–303 3.2. Interpolation methods Among the numerous methods available for interpolating spatial data, we are interested in those that are able to handle scattered positions rather than in a grid that can also produce values outside of the convex hull (the smallest convex polygon that includes all the measured points [51] ), and that are available. The code or the implementation are public. In order to reduce the size of the comparison, we have focused on some families of methods, assuming that the selected ones are similar to others (e.g. inverse distance weighting is similar to radial basis functions and polyharmonic splines are a special case of radial basis functions), and those that are commonly used due to their good results. The selected algorithms are the following ones. 3.2.1. Nearest-neighbor interpolation This is a simple method that works in n dimension spaces and with multiple variables. Given a set of m 2D measured points, where each point is a vector of 3 values (the 2 position coordinates and the measured value, P i (x i , y i , V i ) , with i between 0 and m −1 ), we want to obtain the unknown value at any position, P s (x, y, ?) . This goal is achieved by finding the nearest point P i , and assigning its value to the studied position. The nearest point is computed using some type of metric, for example, the Euclidean distance. This method does not have parameters. A particular solution is the method of Martín-Ramos and Chiari [45] . This method uses both the normalized RGB color, (R, G, B ) with R,G,B values between 0 and 1, and their normalized position coordinates, (x, y ) with the x , y values between 0 and 1, to compute the distance. Each non measured position has five known parameters, the RGB and the position coordinates, and one unknown parameter, the value, (R, G, B, x, y, ?) . The method computes the distance of the unknown position to all the known ones using this formula: d =  (R −R i ) 2 + (G −G i ) 2 + (B −B i ) 2 + (x −x i ) 2 + (y −y i ) 2 . The value of the known measurement that produces the smallest distance is assigned to the non measured position. It is also possible to apply some modulation to the value by computing a factor that depends on the ratio between the minimum computed distance and the maximum distance of the hypercube: f = Max _ distance −d Max _ distance . As the Martín-Ramos and Chiari [45] paper does not assign any name to the method, and after checking the way it works, we label it in this paper as Minimum Hypercube Distance, MHD. Two important variants would be to use only the color information or the position information. 3.2.2. Barnes method The Barnes [38] method generates interpolated maps from pressure data gathered by monitoring stations at dispersed positions. Its approach uses the inverse distance to control a Gaussian function. This solution requires two steps: in the first step an initial result of the interpolation is produced and the second step refines the result of the initial interpolation in order to improve the final data. One disadvantage is that it requires manual adjustment of the input parameters, the radius of the selection being the more relevant. That is, given an unknown point, the method selects measured points that are at a distance less or equal to the radius. Given that points for the test are selected randomly (the process is explained in Section 4 ), we need a radius value that always guarantees that at least one point is selected. We have used two weight functions: w = e −distance κγ for Barnes, and w = rad ius −d istance rad ius + d istance for Cressman, a variation of Barnes. 3.2.3. Kriging method The Kriging method [36] is based on a weighted linear combination of the measured data. It assumes that the data collected from a certain population are correlated in space, which matches with the paintings’ realization. We have used the ordinary Kriging type. This method needs several parameters depending on the different types of functions, called variograms. If they are not provided, they are calculated automatically using an L1 norm minimization. The main problem with this method is that it requires solving an equation system, which in some cases is not possible or produces a great MSE. The following functions have been tested: Linear, Power, Gaussian, Spherical, Exponential and Hole-effect (check the formulas in Section 4.1 ). 3.2.4. Radial basis function The interpolation using radial basis functions (RBF) consists of constructing an interpolant from the weighted sum of RBF. A radial basis function [52] is a function that depends on the distance to one fixed point. For example, a Gaussian function complies with this definition. So the idea is to place a radial function at each point and adjust the weights of the functions and then combine them. The parameter, for gaussian or multiquadric functions defaults to average distance between known points. There is also a smooth parameter that adjust the smoothness of the approximation. By default it interpolates at the known points. We have tested these radial functions: Multiquadric, Inverse, Gaussian, Linear, Cubic, Quintic and Thin plate (check the formulas in Section 4.1 ). 4. Statistical analysis For our study we have followed a similar approach to Willmott et al. [1] . We are interested in methods that produce the lowest level of MSE. In such cases, the maps will have enough accuracy to be considered by the experts as a valid option to interpret the results. For each test, the process is as follows: 1. A set of m values at m different positions is obtained by scanning the surface of the painting. These are considered our measured or known values. In our case there are 1314 and 165 values for each detected element. We obtain the following information for each position and for each element: • Position coordinates, (x, y ) , normalized between 0 and 1: the x is divided by width and y is divided by height. It must be taken into account that this process produces some deformation in non-square shape artworks. • The value, V , which is normalized by dividing by the maximum of each element. • The RGB color components for each position, (R, G, B ) , only for MHD method, obtained from any good quality image that has been acquired with a photographic camera or other devices. Their values are normalized between 0 and 1 by dividing each value by the maximum. The most common format is to use a byte for each color component, so, the normalization is achieved by converting the integer value to float and then dividing by 255. The use of normalized values is not just for simplicity but also for compatibility with all of the interpolation methods: MHD needs normalized data before computing the distance. This normalization does not affect the results because the possible deformation is produced in all cases. 2. From known values we extract two disjoint sets: a test subset and an interpolation subset . The test subset is conformed by m values, obtained randomly, that are used by the interpolation method to obtain the results at the remaining positions, the interpolation subset. For example, with “The Transfiguration” reference set, if our test subset is composed of 50 values at their respective positions, 165 −50 = 115 remaining positions which 296 D. Martín, G. Arroyo, J. Ruiz de Miras et al. Journal of Cultural Heritage 62 (2023) 293–303 Fig. 3. Results for the “The Man”. The MSE of Fe and Hg with different interpolation methods is shown. From top to bottom: MHD, Barnes, RBF linear and Kriging exponential. The MSE bars show boot-strapped 95% confidence intervals (CIs). A bar has a green color if the test of Shapiro–Wilk is passed (Normal distribution) or pink if it is not passed. Each tick of the x axis shows the size of the percentage of points used, the test subset. The y axis shows the MSE plust the 95% of CI. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.) conform the interpolation subset. We have defined 10 cases to create test sets, with m representing the 5%, 10%, 20%, 30%, 40%, 50%, 60%, 70%, 80% and 90% of the total, 1314 and 165. Our hypothesis is that the fewer measurements are used the greater the MSE will be and vice versa. The random selection is made only on the set of available positions without taking color into account. 3. We use the MSE to compute the error. 4. For each value of m , it is necessary to repeat the process of computing the MSE a required number of times, the sample size. It is important to note that in general pigment distributions are not continuous and smooth, unlike other kinds of data, such as atmospheric pressure. The pigments we find in paintings are distributed unevenly following the strokes and shapes painted by the artists. This can produce areas of discontinuity, where two nearby positions can have very different values. So we do not have information about the distribution function. To determine the sample size given these conditions we have performed a previous study regarding the optimal value, testing different cases from 20 to 10 0 0. We determined that 100 times is a good compromise because it produces results similar to 10 0 0 while making the computing process much faster. After applying the interpolation method, we obtain a set of 100 measures of MSE for each m . We employed the non-parametric approach of bootstrapping [53] to construct the probability function and to obtain the Confidence Intervals (CI). From a statistical point of view, this is a procedure that follows APA recommendations [54] and addresses the recent critique of Null Hypothesis Significance Testing (NHST) within statistics and application domains [55–57] . The use of parametric tests for this type of study has been at least questioned in previous papers 297 D. Martín, G. Arroyo, J. Ruiz de Miras et al. Journal of Cultural Heritage 62 (2023) 293–303 Fig. 4. Results for “The Transfiguration”. The MSE of Fe and Hg with different interpolation methods is shown. From top to bottom: MHD, Barnes, RBF linear and Kriging exponential. The MSE bars show boot-strapped 95% confidence intervals (CIs). A bar has a green color if the test of Shapiro–Wilk is passed (Normal distribution) or pink if it is not passed. Each tick of the x axis shows the size of the percentage of points used, the test subset. The y axis shows the MSE. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.) [58–60] . In the case that the NHST-based statistics are preferred, the guidelines by Krzywinski and Altman [61] can be used to infer p-values from our results. To obtain additional information about the distribution function we also have applied the Shapiro–Wilk test [62] to determine whether the data follow or do not follow a Normal distribution. This process was applied to each of the selected interpolation methods. 4.1. Implementation For the implementation of the study we have used Python due to the availability of the implementation of most methods and procedures for interpolation. For each interpolation method we have used the following programs/libraries (also check the references for the parameters and functions): • Nearest-neighbor interpolation. We have implemented the MHD method • Barnes method [63] . • Kriging method [64] . • Radial Basis Function [65] . 5. Results Many results have been obtained and it is impossible to show them all in the main article. Therefore, we have selected a few chemical elements in order to comment on the values obtained 298 D. Martín, G. Arroyo, J. Ruiz de Miras et al. Journal of Cultural Heritage 62 (2023) 293–303 Fig. 5. Results for the “The Man”. The MSE of Cu, Fe, Hg and Pb with different interpolation methods is shown. Each tick of the x axis shows the size of the percentage of points used, the test subset. The y axis shows the MSE. Fig. 6. Results for “The Transfiguration”. The MSE of Cu, Fe, Hg and Pb with different interpolation methods is shown. Each tick of the x axis shows the size of the percentage of points used, the test subset. The y axis shows the MSE. and to establish consequences, leaving the whole set of results as additional material. We have selected 4 chemical elements of the 18 possibilities, Cu, Fe, Hg and Pb because they appear in common pigments. The individual results for “The Man” are shown in Fig. 3 and for “The Transfiguration” in Fig. 4 . Each graph shows the MSE with a bar for the different percentage of selected points. The bar marks the CI of 95% and the black sticks the standard deviation. The green color indicates that the data can be associated to a Normal distribution (after passing a Shapiro–Wilk test). Otherwise the color is pink. Figs. 5 and 6 show the results for the same elements but comparing some methods. In search of greater clarity we have only shown 4 methods, the four best ones. If we look at the graphs with all of the 18 methods 2 , there some methods that produce a level of MSE so high that it looks like noise. From the results, we can see that: 2 These graphs can be produced with the additional material • All the interpolation methods in both datasets produce results that are better than using a random process. This implies, as expected, that there is some correlation in the data. • The results for “The Man” are, in general, better than those for “The Transfiguration”. In the case of “The Transfiguration”, we performed some tests with a low number of points, e.g. 8 and 17, but even if we leave these out, statistically the results show a different picture. If we compare the resulst of Fig. 3 with Fig. 4 , we can see that: • Fe produces a larger MSE than Hg. This result is clearer for “The Transfiguration” but is also visible with “The Man”. –The results are very good even using a very low number of points. This is particularly clear in the graphs of “The Man”. In the case of “The Transfiguration” it is shown that in some cases the results are less precise (care must be taken with the different scales). We think that these results are due to the characteristics of each painting, specifically the number of colors, the pigments, and their variation: paintings with a low chromatic variation 299 D. Martín, G. Arroyo, J. Ruiz de Miras et al. Journal of Cultural Heritage 62 (2023) 293–303 Fig. 7. Example of interpolation maps for Pb element using an input image of 204 pixels ×357 pixels ( W ×H ) and all the measured values, 1314, of “The Man” dataset. will produce a low variation in measurements, which can be captured with a low number of measurements, and vice versa. – Some elements are more sensitive to the interpolation method than others. We can see that, for both datasets, Cu and Hg can produce lower MSE levels than Fe and Pb. – Some methods produce more consistent results than others. Related with the previous point, we can see that a method that produces very good results with one element, can produce not so good results with others (e.g. Barnes, see Fe in Fig. 5 ). There are other cases where the method produces good results for some configuration of numbers, test subset vs. interpolation subset, but not for others (e.g. RBF linear, see Hg in Fig. 6 ). We think that this can be due to two causes: some methods capture better the distribution of elements in the paint, and some methods can have stability problems as they need to solve equation systems (this can be observed in some results not shown here). –As expected (in most cases), the use of more information provides better results. If we take all these considerations into account, we can affirm that the use of interpolation is fully justified, allowing us to obtain, in most cases, maps with a very low level of MSE with the best methods. It is important to note that each painting or work of art will have its own characteristics and that it will be necessary to study individually which interpolation methods are best suited to the data. However, before indicating which methods have been shown to be the best and most consistent, other characteristics need to be studied, for example the need of parameters, the speed or its robustness. Some methods, e.g. Barnes, Kringing and RBF, have the possibility of manually adjusting parameters. This could be useful because it would allow us to produce better results, but it has the disad300 D. Martín, G. Arroyo, J. Ruiz de Miras et al. Journal of Cultural Heritage 62 (2023) 293–303 Fig. 8. Example of interpolation maps for Pb element using an input image of 756 pixels ×1112 pixels ( W ×H ) and all the measured values, 165, of “The Transfiguration” dataset. vantage of requiring interaction with the expert. For our study, we have preferred to use the same conditions, allowing each method to select the default values. In the case of Barnes and Cressman, the radius, which controls what points are selected, is defined by the user. We could have tried to create an algorithm to limit this distance depending on the number of points used, which would maximize the use of the functions as well as reduce the MSE conditions and the time computation but instead we have used the default distance 1, which will work for all cases. Another important characteristic is speed. We have checked the time needed for each computation 3 and the results are shown in Table 1 . The results show that even in the worst case, the times obtained confirm that we can work in real time or interactive time. 3 AMD Ryzen 7 2800X, 16 GB of 3200 MHz DDR4 memory. The fastest method is the linear RBF, while the slowest is MHD, although part of its inefficiency is probably due to the fact that it is programmed entirely in Python. Finally, some methods are very simple to program and even to parallelize, e.g. MHD, which is a great advantage. While our recommendation for obtaining a low MSE value would be Kriging (exponent) or especially RBF (linear) because of its speed and no need of parameters, we must take into account the problem with the ill-defined matrices. We should also consider the MHD method as a valuable alternative because it has not issues and it does not require manual adjustment of its parameters. To obtain a more complete picture, we have tried generating maps for a large set of unknown positions, a photo of the painting, in order to compare the results. We have used 3 different interpolation methods to produce maps computing the values for all the pixels of “The Man” image, and an image of 756 pixels ×1112 pix301