Fast Texture Synthesis for Discrete Example-Based Elements

Considering the problem of discrete texture synthesis and the time for texturing, this paper proposes a novel framework for synthesizing texture images based on discrete example-based elements. We start with extracting texture feature distribution from exemplars and then produce discrete elements based on the cluster algorithm. After initializing a texture image, we propose a texture optimization algorithm based on heuristic searching to improve the quality of the texture image. Final, we use a texture transfer method based on Convolutional Neural Network (CNN) to stylize the optimized texture image. Our results show that the proposed texture synthesis method can significantly improve the quality of discrete texture synthesis and effectively shorten the time for texture generation.


I. INTRODUCTION
The extraction and analysis of texture is the key issue in the technical system of computer vision [1]. However, the technologies for further applying extracted texture need to be explored. Such applications include texture mapping, rendering and synthesis [2]. Nowadays, many studies and methods for texture mapping and rendering have been discussed (see [17] to [31]), but it is still worthy to discuss how to realize effective texture synthesis based on different texture properties that can be divided into two classes: discrete texture and continuous texture. For this purpose, it is necessary to propose a method that can integrate the process of feature extraction, optimization, synthesis and even stylization [3]. To our best knowledge, there are many pixel-based or patchbased methods for texturing based on mathematical statistical models [4], such as Markov Random Field (MRF) model. However, these pixel-based or patch-based methods are more suitable for synthesizing continuous texture [2], as shown in Fig. 1, and most of them are quite time-consuming because such methods need to grow texture pixel by pixel or patch by patch. Indeed, methods for synthesizing discrete texture are still limited.
The associate editor coordinating the review of this manuscript and approving it for publication was Gianluigi Ciocca .
In this paper, to achieve fast texture synthesis for discrete texture, we propose an integrative method to achieve this goal, including five subsections. First, we obtain the feature distribution based on images gradient information. Second, discrete example-based elements can be extracted based on image feature distribution. Here, we adopt the mean shift algorithm to determine the number of discrete elements and extract them from the exemplar discretely. Next, we initialize texture images in a stochastic or a regular way and store their properties (such as respective position, angle and structure) with a special method to speed up the synthesis process in following subsections. Then, we propose an optimization algorithm to improve the quality of texture images based on heuristic searching. Final, we use the CNN-based texture transfer technology to stylize the optimized texture images. Compared with some existing methods (such as [18], [20] and [28]) for texturing that may produce some image seams and just optimize the local texture distribution (see Fig. 12), our method starts with preserving the integrity of elements in the exemplar, then optimizes the elements distribution globally, and can save more time for texturing process. Moreover, the overall synthesis process is automatic (users only need to prepare their own materials), and it is quite easy to transfer our model to a 3D surface according to user needs. The final stylized texture images (see Fig. 13) that produce some strong visual effects are also interesting. The patch-based method for texture synthesis [19]. (b) is the output image of continuous texture from (a). (d) is the output image of discrete texture from (c) that may present some image seams.
The rest of this paper is structured as follows. In Section II, we present an introduction of related work. Section III introduces the proposed process of texture synthesis. The results are presented in Section IV, and the conclusions are drawn in Section V.

II. RELATED WORK A. FEATURE EXTRACTION
Common feature extraction methods mainly include: statistical methods, transform-based methods and model-based methods.
For statistical methods, the feature distribution of texture images is extracted by some statistical models, such as a gray-level model, a gradient histogram model and a mapped pattern-based model. In [5], the gray level co-occurrence matrix (GLCM) is used to extract second or higher order feature from RGB images. Following [5], GLCM is further applied to other scenarios, such as extracting structural information from seismic attributes [6] and monitoring changes in baguettes [7]. Furthermore, the authors in [8] propose a texture descriptor for image segmentation based on the histogram of image gradient, and the research in [9] designs a method to recognize buildings according to the histogram of oriented gradient (HOG). In this paper, we also extract image feature based on image gradient changes. Further, users can adopt some methods, such as Log-Polar Transformation (LPT), to reduce the rotation effects.
For transform-based methods, each texture image can be transformed into a frequency space or a scale space to interpret its texture. In [10], a texture classification method is proposed by applying 1D or 2D Fourier filter to extract local frequency information. Following [10], the authors propose an effective method for texture representation based on wedge filters and local descriptors [11]. To extract more image information, the author in [12] uses a copula model to decompose texture information produced by Circularly Symmetric Gabor Wavelet (CSGW).
For model-based methods, some mathematical models are adopted. Among these models, the complex network model has better effect than others, such as a gravitational model [13] or an autoregressive model [14]. In [15], a vocabulary of words based on Complex Network (CN) is built to analyze texture. Following [15], a local spatial pattern mapping (LSPM) is proposed to classify texture images [16].

B. TEXTURE SYNTHESIS
Common methods for texture synthesis mainly include: pixel-based methods, patch-based methods and examplebased methods.
For pixel-based methods, the basic idea stems from [17] that designs a non-parameter synthesis method by extending a noise image pixel by pixel. Following [17], the authors in [18] propose a synthesis method based on MRF and adopt a tree-structured vector to speed up such process. In [19], a method called ''image quilting'' is proposed by minimizing the cost path through the overlapped areas horizontally and vertically to stitch input texture images. Such method is also effective for resizing images by reduction and expansion operation [20].
For patch-based methods, the main idea is to speed up the synthesis process by overlapping patches instead of growing pixels. In [21], the authors first extract the texture patches from exemplars and then copy them into a 3D surface seamlessly by a liner model. The research in [22] presents an analytical method for texture synthesis by decomposing an exemplar into patches and recomposing them on a 3D surface directly. Following [22], a graph-cut method is developed to compute the optimal patch size. The extracted patches can be copied into the output image dynamically [23], and the method in [23] is also suitable for video synthesis. Moreover, a search-based algorithm for texturing is proposed by constructing a image pyramid and using a cache reuse technology [24]. Compared with pixel-based methods, patch-based methods are more efficient.
For example-based methods, an input image is usually regarded as a complete and small element, and the synthesis process is similar to assemble these elements. In [25], a datadriven method for discrete texture synthesis is proposed, but such method requires users to determine their materials structures in advance. Following [25], the authors propose a global optimization algorithm to generate large building models based on small ones [26], and a tile-based method is proposed in [27] to generate facade building images. With the advent of neural networks, a method based on the long-short term memory (LSTM) that stems from Recurrent Neural Network (RNN) is proposed to produce regular texture [28], a feedforward network (FFN) is used to synthesize diverse texture by interpolating texture images generated by different FNN layers and maximizing stylized texture quality [29], [30], and the authors in [31] propose a conditional generative CNN-based (CGCNN) algorithm for synthesizing examplebased texture with non-local structures. Moreover, the research in [1] designs a novel method for growing texture over a 3D surface, and a camera-aided texturing method is proposed by integrating artistic tools into the texture synthesis system [32]. In this paper, the main idea for texture synthesis also stems from example-based methods.

III. METHODOLOGY
Our method mainly consists of five parts, as shown in Fig. 2. A feature extraction method is used to obtain discrete example-based elements in (b). In (c), we initialize a texture image in a stochastic way. The regular initialization method can be also used, as shown in Fig. 8. Then, we propose a heuristic algorithm to optimize the texture image (d). Final, a stylized image is obtained by texture transfer (e).  Patch-based method [23] for texture synthesis based on Fig. 3(a).

A. EXTRACT FEATURE FROM INPUT EXEMPLARS
Given an exemplar I (x, y), such as Fig. 3(a), if we directly use some patch-based algorithms (such as [20]) to synthesize texture, the final texture image may present some image seams, as shown in Fig. 4, and the texturing process is time-consuming. To preserve the integrity of elements in I (x, y) and speed up the synthesis process, we need to extract texture elements from I (x, y) discretely.
First, we obtain the blurred image g (x, y) from I (x, y) by Gaussian filter [33] to reduce pattern noise: where * means the convolution operation and σ 2 is the variance. Next, the gradient value G (x, y) and the gradient direction θ (x, y) of g (x, y) can be obtained by Sobel operator [8]: where g x = ∂g(x,y) ∂x is the horizontal gradient component and g y = ∂g(x,y) ∂y is the vertical gradient component. Here, we remove some pixels that are near an element edge by threshold T S to reduce the edge width: Then, two thresholds, T L and T H , are used to recognize the feature points based on the result of (4): is an uncertain point. If an uncertain point is adjacent to known feature points, its value is equal to 255; otherwise, its value is equal to 0. G (x, y) = 0 means that (x, y) is a non-feature point. Final, we can extract a pixel coordinate set of feature points P = {(x i , y i )} , i = 1, . . . , p based on gradient information by a border following algorithm, where p is the number of feature points, as shown in Fig. 5. Moreover, if necessary, users can add feature points to the set P manually to improve the accuracy of edge extraction when they are facing with some complex patterns, such as flower patterns. VOLUME 8, 2020

B. GENERATE DISCRETE EXAMPLE-BASED ELEMENTS
Based on the result of set P, we need to cluster these feature points (see Fig. 5) to determine the approximate number of elements. Therefore, some clustering algorithms, such as K-means clustering or mean shift clustering [34], can be used to determine the number of centroids (the pixel coordinates need to be integerized) in an exemplar that are equal to the number of discrete elements. If users want to use cluster algorithms, such as K-means clustering, the number of centroids should be determined in advance. Here, we adopt the mean shift algorithm to cluster feature points because such algorithm can search centroids based on the distribution of feature points automatically, and the main steps are listed as follows: Step1: Select one random point as a center point c i from the unclassified points (feature points).
Step2: Add all points around c i within r i to the set M i , where r i is the radius.
Step3: Compute each vector c i , m j , where m j ∈ M i , and add all c i , m j together to obtain s.
Step4: Let c i move along the direction of s, and the moving distance is s .
Step5: Repeat Step2-Step4 until s ≤ ε, where ε is the threshold to stop the loop, and record the current position of c i .
Step7: Compute an access frequency f c i for each point based on different center points c i , and let the current point belong to c i by arg max f c i .
Then, we divide the set P into single discrete elements p 1 , . . . , p j , as shown in Fig. 6, where j is the number of centroids. Each p j includes its pixel coordinate, its relative angle and distance between the centroid j and corresponding feature points, and the RGB intensities, which can speed up the process for texture synthesis in following subsections.

C. INITIALIZE TEXTURE IMAGES
There are many methods that can be used to initialize texture images, such as patch-based methods. Here, we adopt a stochastic pattern for initialization, as shown in Fig. 7. However, a regular pattern can also be used to initialize texture images, as shown in Fig. 8, and we will present the final texture images based on both two initialization methods.
In this paper, the Poisson cluster process (PCP) is chosen to model the stochastic position distribution of elements that are randomly selected from p 1 , . . . , p j , and we always copy selected elements into a output image completely [35]: where is the parent point process, n is the number of the parent point process, and i , i ∈ n is a finite set of the child point process. The modeling process is listed as follows and shown in Algorithm 1: Step1: Generate n ∼ Poisson πr 2 λ , where λ is the density function, r is the radius, and u 1 ∼ U (0, 1) , . . . , u n ∼ U (0, 1).
In addition, to speed up the initialization process, we align the centroid of selected element with the node generated by PCP or the regular pattern and reconstruct the element based on the relative angle and distance and the RGB intensities restored in p j .

D. OPTIMIZE TEXTURE BASED ON HEURISTIC SEARCHING
From Fig. 7 or Fig. 8, we find that some elements overlap over others, which cannot preserve the integrity of all elements. for j = 1, 2, . . . , n do 7: Generate u i,j ∼ U (0, 1) 8: Compute θ i,j = 2πu i,j and generate R i , θ i,j 9: end for 10: end for 11: return texture image To improve the quality of a texture image, we propose a heuristic algorithm to optimize a texture image in this subsection.
First, we formulate our optimization problem via the following energy function: where M is the number of parent nodes, N is the number of child nodes in each parent cluster, O S j , S k is the average overlap rate (AOR) between element j and element k, and the expression of O S j , S k can be expressed as: where S j and S k are the relative area of element j and element k. Second, we define the relative area S j or S k in (8) as the bounding box of element j or element k. Here, Graham algorithm is used to determine such bounding box, and its program interface can be found in CV2, as shown in Fig. 9.
Then, we can define our optimization problem as:  Final, it is hard to solve (9) directly because its gradient cannot be obtained. Therefore, we propose a heuristic algorithm to optimize (9) based on Coordinate Descent Optimization (CDO) that mainly consists of five steps, as shown in Algorithm 2, with a time complexity of O (MNN ).
Step1: Set each centroid j as an origin, divide the searching area for S j into four quadrants, and set initial parameters for the searching radius, the searching length in pixel coordinate frame and the AOR threshold to stop the searching process (set by users).
Step2: Before each iteration, we need to determine whether element j and element k are overlapped based on (8). If the value of (8) is greater than the AOR threshold, we will add the number pair (j, k) to the set {(j, k)}.
Step3: Before each searching process, we select (j, k) from the set {(j, k)}, then fix the element j position and choose one random direction to begin our searching process. Here, we adopt three transform methods to optimize the global texture distribution, including translation, scale and rotation. For example, we first fix the value of y k and search the possible value x k for the position (x k , y k ). Then, based on the new position of x k , we will search the possible value for y k .
Step4: After the searching process, we add the element k position to a set {(x k , y k )}, where (x k , y k ) is the possible position for element k, and determine the optimal position for element k by finding x * k , y * k that can minimize (8) in the set {(x k , y k )}. If so, (j, k) is removed from the set {(j, k)}.
Step5: Repeat Step2-Step4 until the set {(j, k)} is empty or the number of loop reaches its threshold.

E. TEXTURE TRANSFER
In this subsection, our aim is to stylize optimized texture images, and we adopt a neural network method to stylize texture images. Such method can provide the feature space   The experiment results of texture synthesis. For each discrete texture (left), the output images based on the pixel-based method in [18] and the patch-based method in [20] are on the middle, and the output image based on our method is on the right. produced by CNN layers that are very useful for feature representation.
Given an input image A (x, y) ( Fig. 10(g)) and an artistic image B (x, y) ( Fig. 11(b)), we use 16 convolutional layers and 5 average pooling layers from VGG-Network [36] to extract feature representations from A (x, y) and B (x, y) in different layers and define the squared-error loss between different feature representations in layer l as: where N l is the number of filters in layer l, M l is the size of feature map, F l i,j is the content loss produced by the filter i at position j in A (x, y), and P l i,j is the style loss produced by the filter i at position j in B (x, y). The derivative of (10) can be expressed as [37]: Therefore, a gradient descent algorithm (GDA) can be an effective method to optimize A (x, y) until it has the same feature representation as B (x, y), i.e., minimize (10) based on the convergence of (11), and the final texture image is shown in Fig. 11.

IV. SIMULATION RESULTS
For our simulations, the parent nodes are distributed in an area 500 × 500 (the size of output image), and the distribution of parent nodes satisfies PCP, where the distribution density λ is 1, the number of subprocess points n is 50, and the distribution radius r is 250 (usually half of the size of an output image). Table 1 lists the simulation parameters of texture synthesis and Table 2 lists the system parameters of texture synthesis.
All results in this section are independent experiments through Python (including NUMPY, CV2, SKLEARN and KERAS). The results of texture synthesis are presented in Fig. 12 and the corresponding results of texture stylization are presented in Fig. 13. Moreover, the results of time for texture synthesis are summarized in Table 3. As illustrated in Table 3, our proposed method can save nearly 1 hour and 2 hours compared with the pixel-based method in [18] and the patch-based method in [20], respectively. Here, the main steps in [18] are summarized as follows, with a time complexity that is proportional to the size of a sample image and an initial noise image: Step1: Input a sample image I s and generate a random noise image I n .
Step2: Set the searching radius r = 2 and the searching window Step3: Select a patch in I s that is similar to the local region in I n based on the convolution operation.
Step4: Repeat Step3 until all pixels in I n are updated.
The main steps in [20], with a time complexity that is proportional to the size of an output image, a block and an overlap region, include: Step1: Input a sample image I s , initialize an image I 0 and set the size for an overlap region S and a block B.
Step2: Copy one block B 1 from I s into I 0 randomly to begin.
Step3: Compute the error e = (B current − B next ) 2 at the overlap region S e to minimize E i,j = e + min E i−1,j−1 , E i−1,j , E i−1,j+1 , where i and j are the size of the overlap region S e and E is the cumulative error.
Step4: Repeat Step3 until I 0 are filled with blocks. Compared with our method, the methods in [18] and [20] are more time-consuming obviously.

V. CONCLUSION
With a focus on fast texture synthesis for discrete texture, this paper proposes a novel framework for texturing that mainly consists of five steps. We first extract texture feature from input exemplars and generate discrete elements based on the feature distribution and the mean shift algorithm. After initializing the texture image in a stochastic or regular way, we propose a texture optimization algorithm base on heuristic searching to improve the quality of optimized texture images. Final, a CNN-based method is used to stylize our texture images. The results indicate that our method can effectively optimize and speed up the process of discrete texture synthesis.
In the future, we will explore methods for extracting more complex patterns, like some feature matching methods, and take into account factors that may improve the extraction accuracy, such as illumination variation, viewpoint angles, color constancy, and contrast balancing. Moreover, it is important to design more robust optimization algorithms and apply our texturing method to other scenarios, such as remote sensing images fusion and artistic or industrial design.
XIAOHONG LIN received the B.S. degree from the School of Design, Hunan University, Hunan, China, in 1989. She is currently working with the Fuzhou Technology and Business College, Fuzhou, China. Her major research interests include industrial design, product design and 3D modeling, machine learning, and computer vision.
CHONGCHENG CHEN (Member, IEEE) received the M.A. degree in hydro and engineering geology from the Chengdu University of Technology, Chengdu, China, in 1992, and the Ph.D. degree in cartography and GIS from the Institute of Geographical Science and Resource, Chinese Academy of Science, Beijing, China, in 2000. He initiated the geographical information system development and application at Fuzhou University. He was the Co-Founder of the Spatial Information Research Center, Fuzhou, China, in 2001. He is currently dedicated to bridging the gap between innovative IT and earth science research under the umbrella of Digital Earth, and promote geospatial technology to inoculate with applications in broad fields of nature resource, environment, ecology, and government management since then. VOLUME 8, 2020