Since our work involves 3D forms, we focus primarily on findings related to samples with relief, including both 2.5D and 3D methods. Also called relief printing, 2.5D printing uses a similar layer-by-layer approach to that of 3D printing. However, unlike 3D printing, it cannot build true 3D objects and is only used to print texture variations or slightly raised surfaces [1]. As 3D printing is also relevant to such applications, a primary distinction between 3D and 2.5D printing is the use of support material: support is necessary for complex geometries in 3D printing but generally not required in 2.5D printing due to its limited depth [20].

But although 2.5D and 3D technologies are considered different, they do share some similarities. Indeed, both 2.5D and 3D technologies can use the same techniques and materials during printing [20]. This is true for color 2.5D printing and PolyJet 3D printing for example, where both techniques use inkjet UV-curable CMYKW ink in their printing process. This means that findings from studies using one technology or the other can be generalized interchangeably, as the printed samples share the same manufacturing process and the features/biases related to it. Accordingly, we include findings from studies conducted on naturalness of 2.5D printed samples as well. In the present investigation, we used a 3D printer, and so we refer to our samples as 3D printed.

Naturalness assessment for objects with relief is not a straightforward process for observers, and it depends on a multitude of factors [14]. An investigation into naturalness perception of 2.5D samples depicting different materials such as wood, stone, and so on was conducted by Kadyrova et al. [15]. They found that surface texture, characterized by surface elevation and roughness, is a prominent factor in judging a 2.5D object’s naturalness. The aforementioned work was continued in another study [16] only for wood samples, where reliefs of lower elevations were printed. Together, the studies found that observers had a preference for midlevel elevations at 0.5 mm of maximum relief height to be the most natural looking. Van Hoey et al. [32] found that visual distortions caused by surface texture modifications on relief samples decreased the perceived quality among observers. But aside from the visual aspect, surface texture is also defined by its tactile attribute, roughness. Tymms et al. [30] demonstrated that the shape and size of the facets making up the surface of a 3D printed object widely varies the tactile perception and understanding of the said object. They also showed that surface texture depends not only on elevation but also on the positioning of surface facets.

Since surface texture depends on the positioning of the facets, added relief should be meaningful to improve tactile perception as shown by another study by Tymms et al. [29]. When based on 2D images, relief was added using simple grayscale images of these colored 2D images [15]. In contrast, Wang et al. [35] investigated better methods to add relief to printed samples. They introduced a methodology that combines texture analysis and semantic understanding to enable the reproduction of textures in 2.5D prints. This involves identifying texture features using filters in the frequency domain [22], such as a high-pass filter or Gabor filter. Wang et al. [35] concluded that incorporating semantic understanding into the printing process produces more realistic and perceptually rich 2.5D prints that closely resemble the original textures. Furthermore, Texture-Aware Error Diffusion Halftoning (TAED), like the algorithm developed by Li et al. [17], may also produce semantically relevant relief placement. In their investigation [17], Li et al. demonstrated that TAED outperformed traditional error-diffusion methods (Floyd–Steinberg [8] and Shiau [24]) both objectively (using 2D image quality metrics, such as the the Sparse Feature Fidelity, Peak Signal-to-Noise Ratio, Mean-Structural Similarity, and computation time) and subjectively (via psychometric experiments). Therefore, TAED offers an innovative approach to digital halftoning that aims to enhance image quality, particularly in preserving texture details.

Displacement maps, hereafter referred to as DM(s) in the paper, were generated using RGB reference images (images obtained from Pixabay [21], an open-source image database). Considering familiarity and everyday use, these reference images were selected from four different material categories: Wood, Stone, Fruit, and Cultural Heritage Object (two from each category), and are 591 × 591 pixels wide due to the technical limitations of the 3D printer (refer to Section 3.2 for more details). The cultural heritage images were chosen for their complex surface ornamentation, which presents a challenging scenario for texture reproduction, allowing us to assess the capabilities and limitations of the various algorithms under demanding conditions. The intricate details of the Alhambra ornamentation, for example, characterized by fine lines, sharp edges, and subtle variations in depth, provide a rigorous test for the algorithms’ ability to capture and reproduce fine details that are crucial to naturalness perception. However, we acknowledge a potential limitation associated with the use of the Alhambra image. The naturalness judgments in our study are based on the participants’ perception of how well the 3D printed samples match their internal representation of the depicted objects or materials. Some participants may be unfamiliar with the Alhambra and its specific ornamentation. This unfamiliarity could potentially influence their judgments, as they may not have a clear mental image of the original structure for comparison. Although this could introduce some variability into the naturalness ratings for the Alhambra samples, the complexity of the image still provides valuable insights into the performance of the DM algorithms.

DMs are grayscale images used in 3D design software to physically alter a 3D object’s surface [27]. Physical alteration of a 3D object’s surface generates more realistic micro-surfaces, enhancing its appearance through improved light–surface interaction, resulting in better shading and light reflection in the 3D model.

DMs are generated from the eight reference images shown in Figure 1. First, a pixel-wise transformation is applied to calculate a grayscale DM (Gray DM) from the 2D RGB image [26]. Thus, for an image with N*M pixels, the grayscale value at the pixel (i, j) is calculated using Eq. (1). In Eq. (1), Y represents the grayscale level and R, G, and B are the RGB values of the corresponding pixel in the 2D image. This is a standard RGB-to-grayscale transformation that represents grayscale in a manner aligned with human vision. Consequently, it assigns a higher coefficient to green, followed by red, and then blue based on the eye’s sensitivity to each color [26]. The Gray DM, therefore, serves as a starting point for texture extraction and DM generation. Six DMs were generated in MATLAB (version 2023b) using the Gray DM where different appearance attributes like color, roughness, and gloss were considered. Table I presents the description of the Gray DM and the six DMs generated. The DMs were generated using a specific set of algorithms. These DM-generating algorithms are not the only ones we could have used but were chosen because they allow us to have DMs that are related to the content of the image and therefore reproduce texture in a meaningful way. The choice of algorithms is part of the novelty of this work because as far as we know, they have not been used in 3D texture printing before. Therefore, the choice of algorithms is not based on a ranking of performance, it is rather based on the surface texture features we want to control: roughness, foreground highlighting, and edge highlighting. Consequently, we employ algorithms that enable us to control these features either individually or in combination. Furthermore, these algorithms are based on standard image processing functions that are mostly in-built in widely used image processing software like MATLAB and Python. For instance, we used binarization algorithms like imbinarize for image forefront highlighting. We also employed halftoning algorithms like Minimized Average Error Diffusion (MAE) and TAED (see Table I) for pixel-wise roughness addition. However, we had to limit the number of algorithms due to practical constraints: a large selection would have resulted in longer printing times and more complex and time-consuming subjective experiments for observers. However, although these algorithms are useful for feature extraction, they also have disadvantages. These range from information loss due to halftoning to a reliance on grayscale images, making the output highly sensitive to shading and gloss variations in the reference images. Hence, this is by no means a definitive list of algorithms to always use, and future research might uncover novel and better ways to extract and control surface texture from reference 2D images. Thus, the selection of algorithms for this study was primarily based on the variables we wanted to control that are specific to this study.

As an example, Figure 2 illustrates the DMs generated for the Alhambra image.

The 2D images were converted to 3D models by applying the DMs using the Displace modifier in the rendering software Blender (version 3.0) [2]. The midlevel was set to zero, ensuring that all micro-facets were added onto the surface rather than being etched into the model. The strength was adjusted to ensure that the 3D models had three different elevation levels: 0.75 mm, 1 mm, and 1.5 mm. The models were exported as .obj files and loaded as such into the 3D printer software. Figure 3 illustrates the transformation of the DM into different elevations depending on the grayscale values according to Eq. (2). In Eq. (2), d is the elevation height, Y is the grayscale value of the DM at a given pixel (i,j), and E is the maximum elevation to be printed (for Y = 1) in the 3D sample. Figure 4 illustrates the simulation of a model generated using a DM with three different elevation levels.

To limit the variables affecting naturalness reproduction to elevation, roughness, and surface texture, the color information (in the RGB color space) from the 2D images was UV-mapped (see Fig. 1) onto the 3D model’s surface before 3D printing. As we are working with cubes that are rather simple 3D shapes, no complex UV-unwrapping was needed. We separated the UV-facets of the cube manually. Then, we identified the UV-coordinates corresponding to each face of the cube and added the color information to the corresponding cube face by UV-mapping the reference 2D image onto the upper face of the cube.

With seven DMs, eight reference images, and three elevation levels, a total of 168 (7 × 8 × 3) samples were 3D-printed. For printing, we used the Stratasys J55 multimaterial PolyJet 3D printer [25], hereafter referred to as J55 in the paper, to print the 3D models. The ink used is CMYKW (cyan–magenta–yellow–key (black)–white) ink belonging to the Vero family of UV-curable ink [33]. This ink family is characterized by its translucent, specular, and scattering optical properties. The J55 can print with a resolution of 300 dpi (dots per inch) in the X and Y directions, and a resolution of 1354 dpi in the Z direction. This corresponds to a voxel of 0.08 × 0.08 × 0.018 mm in size. We used the GrabCAD Print 3D printing software for 3D model orientation, slicing, and printing. We printed the samples in high-quality mode, with a print time of 4.5 h per print project, with each print project containing 20 samples. To balance material use, printing time, and the need for sufficient detail clarity, the 3D samples were designed with dimensions of 50 mm × 50 mm × 5 mm. This size ensured that the printed features were discernible while minimizing material waste and avoiding excessively long printing times. Furthermore, we used translucent materials on the J55 printer, making it necessary to print thicker samples than usual to obtain an opaque object. The 3D printed sample dimensions and the DM pixel dimensions were kept equal to avoid any automatic stretching or compression of the DM pixels. This was achieved by maintaining a fixed pixel ratio of the 2D image before generating its corresponding DM. That is, for a maximum resolution of 300 dpi for an object of 50 mm × 50 mm in dimension, the corresponding number of pixels would be 591 × 591. After printing, the samples were verified against the expected outcome through visual inspection. We also verified that the printed elevations match the theoretical values by 3D-scanning some of the samples and comparing the sizes of the 3D models and the scanned models. Figure 5 shows the resulting printed samples of the Alhambra image.

A categorical judgment psychometric experiment was conducted to evaluate the perceived naturalness of the 3D printed samples. Fifteen observers (7 female, 8 male, average age = 29 years, standard deviation = 8 years, age range [22, 52]) participated in the experiment. None of the participants is an expert in 3D printing although all but two had prior knowledge of color science. The experiment adhered to the ethical guidelines for the “Protection of Research Subjects” set by the Norwegian National Committee for Research Ethics in Science and Technology [9]. Participants provided signed consent before commencing the experiment, and their anonymity was ensured. Although the number of participants in this experiment is in the acceptable range [11], we acknowledge that we are at the lower end of that range. This should be considered a potential source of bias in our study, and future work needs to heed this limitation.

For this experiment naturalness was defined as a matching between the memory of an object in an observer’s mind and the scene portrayed in each sample in this experiment [15]. Participants were shown all of the 168 3D printed samples divided into 8 sets of 21 samples each. Each set corresponds to the samples derived from each reference image shown in Fig. 1. The categorization took place in a Gretag Macbeth Judge II viewing booth under a CIE D50 illuminant. Ambient lighting conditions were set to darkness. Five categories, as listed below, were used to categorize the 3D printed samples.

Participants were not provided with any reference images to aid them in categorizing the samples. Participants relied on visual assessment of the samples and could also touch the samples (without lifting the sample from the booth floor) to feel the surface texture. Each observer spent approximately 75 minutes completing the experiment. Figure 6 shows the experimental setup.

To understand how participants categorized the naturalness of samples, it would be important to understand their thought process when judging a sample’s appearance both in terms of the naturalness definition used and the attributes considered. Therefore, each session was audio-recorded in its entirety with the participant’s consent. Participants were encouraged to describe all the steps they were following and to describe their thought process when ranking the naturalness of 3D samples along with a few statements when they were done categorizing each set of samples. The audio recordings were then transcribed, and a word frequency analysis was conducted to determine which attributes the observers relied on during their assessment.

Raw individual scores were recorded during the experiment and analyzed. With data collected being categorical and the difference between the categories being variable, we consider our data as ordinal [37]. This type of data should not be misinterpreted as interval data, as this can lead to conceptual and statistical inaccuracies. For instance, it is not necessarily true that a rating of 4 (Somewhat Natural) is twice as good as a rating of 2 (Un-natural) or that the improvement from 1 (Very Un-natural) to 2 (Un-natural) is equivalent to the improvement from 4 (Somewhat Natural) to 5 (Natural). Therefore, it is crucial to handle ordinal data appropriately to avoid misconceptions and incorrect interpretations. Cumulative Link Mixed Models (CLMMs) are a fitting method for analyzing ordinal data. They operate under the assumption that the observed ordinal variable is a categorization of an underlying continuous variable [5]. This enables CLMMs to provide the probability for each rating rather than a single prediction. The Bayesian approach, which can be implemented using the brms package in R, is used to fit these models [4]. This approach offers several benefits over traditional frequentist statistics, including the ability to incorporate prior knowledge, quantify uncertainty, and estimate complex models. This makes it a powerful tool for analyzing ordinal data. In this paper, the Bayesian approach was deemed more appropriate as naturalness perception and sample categorization can be highly subjective, that is, incurring random effects [28]. Therefore, we use it to estimate the likelihood of sample categorization.

Using Observer (referred to as ID), Algorithm, Elevation Level, and Original Image as variables, seven statistical models were generated using different variable combinations. We include these variables because we know from the literature that elevation [16] and algorithm [35] affect the naturalness and reproduction quality of 3D printed samples. We also include the participants and images as variables because, as shown in the literature [28], these variables incur random effects that arise from the participants’ understanding of naturalness and image content shown. We begin with simpler statistical models that include only a few variables. We then formulate the said models and test their predictive performance. After testing, model complexity is increased by adding more variables or by incorporating interactions between existing variables. Through this iterative process and with our selection of variables already determined, we construct and evaluate various statistical models based on their predictive power. Once we build a model with strong predictive power, that is, model prediction values match the observed data, we use the said model for result analysis. All variables in this model are expected to influence the naturalness rating. The effect of each variable on the model can then be examined to test our initial hypotheses.

We start with models that assume that variables are independent at first, but we do not neglect the probability of interaction between the different variables in the naturalness rating given during the subjective experiment. Therefore, we also include models that take into consideration this variable interaction between a select number of variables. Table II gives a detailed explanation of the statistical models used.