#### Abstract

SFA (Surface Fitting Algorithm) for continuous displacement is an important method for digital image correlation with antinoise ability and computational efficiency advantages in practical applications. In order to improve the algorithm accuracy and expand its application range, this paper tries to improve the SFA and studies the modified cubic surface fitting algorithm CTSFA (Corrected Three Surface Fitting Algorithm), which is suitable for solving the initial value of continuous displacement. Bilinear interpolation and adjacent interpolation are used to analyze the gray level at any integer-pixel position in the displacement matrix and the weight coefficient is given. The distance-weighted method is used to approximate the true initial displacement value of the continuum, and the algorithm suitable for digital image processing is extended to the continuum displacement solution. The cubic surface expression of the CTSFA programmatic application is solved by the least squares method, and the correlation coefficient of the power basis function is calculated. In the computer simulation of speckle test, the comparison between CTSFA and SFA on the calculation results of linear and nonlinear displacement fields shows that the calculated amount of CTSFA is basically the same as that of SFA, but the calculation accuracy is doubled. The study of analysing the Brazilian splitting test using CTSFA and SFA reveals that CTSFA is better than SFA in observing the development of cracks.

#### 1. Introduction

As a simple and direct method of calculating displacement, Digital Image Correlation (DIC) has the advantages of noncontact, full-field measurement, simple experimental equipment, and low environmental requirements, which can use an object’s own surface’s natural texture to obtain speckle images of deformation [1–5]. After the DIC method was first proposed [6, 7], numerous scholars conducted research on the accuracy and computational efficiency of different algorithms in-depth, and the advantages and disadvantages of different methods had been pointed out [8–12].

SFA is one of the most commonly used methods for DIC. It performs surface fitting on the correlation coefficient matrix, obtaining the subpixel displacement of the point on the surface by fitting the correlation coefficient. In theory, the method has been researched from different levels in depth [13–15], and the displacement positioning accuracy of the influence of various influencing factors has been evaluated through experiments [16–19]. In practical applications, matching the center of the subarea window by calculating the correlation coefficient through the SFA is limited to the integer pixel precision, without locating at the actual displacement point, which resulted in the calculation of the correlation coefficient is not accurate enough. In response to this problem, the initial displacement of SFA has been tried to solve by using linear interpolation and gradient method to obtained high precision continuum displacement [20–23].

This paper introduces the concept of distance weighting in digital image processing to modify SFA of initial displacement determination by CTSFA. CTSFA adopts the gray level interpolation algorithm which combines the nearest neighbor interpolation with bilinear interpolation. By moving the fitting center and calculating the weighted weights, the points participating in the surface fitting calculation become closer to the actual initial point. For the cubic surface fitting function of CTSFA, the least square method is used to solve the extreme coordinates of the correlation coefficient. The accuracy of the algorithm is verified by computer-generated simulated speckle images and cracked Brazilian disc tests before application.

#### 2. Algorithm Principle

##### 2.1. Basic Theory of SFA

Based on the actual data, the SFA obtains the analytical expression between the function and the variables *x* and *y*, and it can pass or approximate all data points, which can be distributed in the function on the spatial surface represented. Taking the two polynomials as an example, the fitting function is [24]

The function should satisfy the equations at the extreme points of the quadric:

Then according to equation (2) [24], the position of the extreme point of the fitting surface is obtained [24]:

After finding the location of the extreme point, the difference between the coordinates of the initial value point and the initial value point is a continuous displacement. Therefore, obtaining high-precision initial displacement and calculating the coefficients of polynomial fitting function are the two important of CTSFA.

##### 2.2. Determination of Initial Value of CTSFA Displacement

At first, CTSFA calculates the interpolation results of bilinear, nearest neighbor points and the gray variance of four neighbor pixels around the interpolation point which based on the traditional linear interpolation method when determining the initial displacement value. Then, it constructs the weighting coefficients through the gray variance. Finally, the two interpolation results are weighted together to get the final interpolation results. CTSFA not only considers the distance between interpolation points and neighboring points, but also takes the gray distribution characteristics of neighboring points into account, which can improve the subpixel accuracy of initial displacement effectively.

The specific steps for determining initial values of displacement in CTSFA are as follows:(1)Calculating bilinear interpolation and adjacent interpolation of the interpolation points . As shown in Figure 1, bilinear interpolation uses the gray-scale weighted interpolation of four neighboring points as the gray-scale value of the point, which can be decomposed into two one-dimensional linear interpolations expressed by Among the four adjacent pixels of the pixel which will be sought, the nearest neighbor pixel gray is assigned. Define *i* *+* *u* and (*i* and *j* are positive integers, *u* and are decimal greaters whose value more than 0 as well as less than 1) are the noninteger pixel coordinates, the value of the pixel whose gray is solved and shown in Figure 2. If falls in area A that means *u* < 0.5 and , the gray value of the pixel in the upper left corner is assigned to the pixel to be sought. Similarly, the falling in the B area gives the upper right corner, and so on.(2)Calculating the standard deviation *H* of four adjacent points around the interpolation point and *H* is expressed by The is expressed by(3)Weighting fusion of and is expressed by equation (8) [25, 26]:

In the equation (8), *W*_{d} and *W*_{n} are weighted coefficients and calculated by nonlinear weighting which can be shown by

The interpolation result is the initial value of the actual displacement required.

In Figure 3, linear interpolation in *x* direction is carried out at four points *a*_{1}, *a*_{2}, *a*_{3}, and *a*_{4} to get *A*_{1} and *A*_{2} and carried out at two points *A*_{1} and *A*_{2} to get *A* in *y* direction, and the gray value of *A* is the result of bilinear interpolation *.* After the result of four-point interpolation was calculated, is got though the four points of *b*_{1}, *b*_{2}, *b*_{3}, and *b*_{4} interpolate adjacent using the above method and the is the point *B* in Figure 3. Then, the *W*_{d} and *W*_{n} are calculated by using point *A* and point *B*. Finally, (point *C* in Figure 3) is obtained through equations (8)-(9), and the point *D* is the true initial displacement. The flow chart of CTSFA calculation is shown in Figure 4.

##### 2.3. Calculation of CTSFA Fitting

In order to improve the fitting effect of polynomial function, this paper adopts three times surface fitting function. The correlation coefficient of the target matching point (*x*, *y*) and the points in surrounding neighborhood, which are commonly used by any cubic surface function, can be expressed by power basis function by

In order to obtain the maximum value of the correlation coefficient and the corresponding coordinates, it is necessary to calculate the cubic surface fitting polynomial. The least squares method is an approximation theory that the surface is generally not obtained by a known point, but by minimizing the sum of the squares of the difference between the value of the sampled surface and the actual value of the fitted surface. Its main idea is that the sum of squared deviations of prediction data and actual values is minimal. Thus, the cubic surface equation of equation (11) is expressed as

This paper gets equation (12) through the given set of data points * m* times (*m* = *N*) polynomials:

After obtaining equation (12), the square deviation *Q* of the predicted data and the real value are calculated, and *Q* is shown by

Equation (14) is obtained by transforming the construction problem of the fitting polynomial of equation (10) into the extremum problem of the multivariate function, and satisfies :

Equation (15) is simplified by equation (14):

Equation (16) is further sorted out through equation (15):

The power basis function of CTSFA is solved after is determined, and the maximum of the correlation coefficients and the position coordinates at the maximum can be obtained. Then the actual displacement can be obtained by taking the difference between the maximum and the initial displacement.

#### 3. Experimental Verification

##### 3.1. Analysis of Digital Speckle Computer Simulation

In order to verify the performance of CTSFA, the accuracy is tested by using the computer simulated speckle experiment. A digital speckle image is selected as the reference frame, with the preset theoretical displacement field given. CTSFA and SFA are used for analyzing the experimental results through the reconstructed speckle field which reconstructed as the calculation frame.

The displacement of the speckle in the image can be precisely controlled by using the numerical method to simulate the speckle image of the uniform and nonuniform theoretical deformation field. The simulated speckle image size in Figure 5 is 128 × 128 pixels, and the intensity of the speckle is Gaussian. The speckle size is 3 pixels, and the number of speckles is 300. The size of calculating window in SFA and CTSFA is 41 × 41.

**(a)**

**(b)**

**(c)**

###### 3.1.1. Uniform Deformation Field Test

The effectiveness of CTSFA and SFA is tested by using the uniform deformation field speckle image of Figure 5(b). The displacement of the speckle field in Figure 5(b) is shown by

Corresponding calculations are performed using CTSFA and SFA, respectively. The image analysis area is *x* = 1∼100 pixels and *y* = 1∼100 pixels. Figure 6 shows the lateral displacement field obtained by CTSFA and SFA.

**(a)**

**(b)**

**(c)**

It is found from Figure 6 that the displacements in CTSFA and SFA are stratified, and the stratification values are consistent with the applied displacement function, indicating that both CTSFA and SFA can calculate the uniform displacement field. The displacement boundary in the *u*-field cloud map of CTSFA is smoother than that of SFA, and the displacement boundary of SFA is irregularly jagged, indicating that CTSFA is more accurate than SFA on calculation of uniform displacement field.

In Figure 7(a), the lateral displacement calculated by CTSFA is closer to the theoretical value, and the calculation result of SFA has a certain deviation from the theoretical value, indicating that the calculation accuracy and stability of CTSFA are better than SFA. In Figure 7(b), the relative error of CTSFA in homogeneous deformation field is about 1.2%, which is half of that of SFA.

**(a)**

**(b)**

###### 3.1.2. Nonuniform Deformation Field Test

The speckle image of nonuniform deformation fields is simulated numerically by CTSFA and SFA which is shown in Figure 5(c). The displacement of the speckle field in Figure 5(c) is shown by

CTSFA and SFA are used, respectively, to calculate nonuniform deformation fields and the area analyzed is *x* = 44∼56 pixel and *y* = 1∼100 pixel. Figure 8 shows the lateral displacement field calculated by CTSFA and SFA. Figure 9(a) shows the results of theoretical, CTSFA, and SFA. Figure 9(b) shows the error analysis of CTSFA and SFA.

**(a)**

**(b)**

**(c)**

**(a)**

**(b)**

In the *u*-field nephogram, the displacement boundary of CTSFA is smoother than that of SFA, which indicates that CTSFA is better than SFA in dealing with nonuniform displacement field. Especially in the case of large displacement, the value of CTSFA is closer to the theoretical value, while the value of SFA is lower, the main reason being that the initial value of displacement selected by CTSFA is more accurate.

The data on the *y* = 40 pixel interface are selected for error analysis. The abscissa coordinates in Figure 9 represent the pixel coordinates of *x*, and the ordinate coordinates represent the value of *u*.

In Figure 9(a), at the position *x* = 50 pixel where the gradient of deformation is large, the calculated result of SFA deviates from the theoretical value greatly, and CTSFA is closer to the theoretical value, indicating that the spatial resolution and displacement resolution of CTSFA are higher than that of SFA and CTSFA is more suitable for analyzing in the region of the large gradient of displacement. Figure 9(b) shows that the relative error of CTSFA is stable at about 3%, and the relative error of SFA is 3 times that of CTSFA.

##### 3.2. Test of Cracked Brazil Disk

In order to evaluate the accuracy and efficiency of CTSFA and SFA further, the test of cracked Brazilian disk taken by industrial CCD camera under static loading is carried out using CTSFA and SFA, and the evolution of lateral displacement field is analyzed. The experimental specimens choose sandstone specimens with diameter of 50 mm and thickness of 25 mm. A crack of 10 mm in length is prefabricated at the center of the specimen, and four sets of experiments are carried out. The speed of collecting photos by industrial CCD camera is 10 frames per second, and the testing machine is RLJW-2000. The loading rate of the test is 0.05 mm/min. The picture of the state of rock failure in the rock mechanical test is shown in Figure 10, and the curve of load-time obtained by the testing machine is shown in Figure 11.

**(a)**

**(b)**

Based on the characteristics of the loading curve, 7 typical moments in the whole loading process are selected for identification. The digital image correlation calculation based on CTSFA and SFA is used to calculate the above time points and compare the results. Among them, the marking point *a* is the reference point analyzed by digital image correlation, the marking points *b* and *c* are located in the linear elastic phase of the test, the marking points *d*, *e,* and *f* are located in the plastic phase of the test, and the marking point is located at the peak point of the test.

The calculation area and time of CTSFA and SFA are same, and the results are shown in Table 1. In the linear elastic stage, CTSFA observes the precursor of crack formation and SFA only observes the vertical compression of the specimen without showing the crack location. In the plastic stage, both CTSFA and SFA observe basically the extension of cracks, but the stratification of the lateral displacement of CTSFA is more obvious. At the peak load point, the crack propagation pattern observed by CTSFA conforms to the actual damage of the specimen. The crack observed by SFA appears jagged in the upper and lower sides of the precrack, indicating that the observation effect of CTSFA is better than SFA.

#### 4. Conclusion

This paper comes up with CTSFA based on the idea of distance weight, the validity of CTSFA is verified by digital speckle simulation and Brazilian disk tests, and the main conclusions are as follows:(1)With the characteristics of high computation efficiency and maneuverability, it is easy for CTSFA to implement programmed application. Computer simulation of digital speckle shows that CTSFA has good adaptability and stability, and it is suitable for both uniform and inhomogeneous deformation fields. The relative error of CTSFA has a low fluctuation, which is kept at 1.2% in the uniform deformation field and at 3.0% in the nonuniform deformation field. Additionally, the measurement accuracy of CTSFA is 3 times as high as that of SFA.(2)For the region with large displacement gradient in the continuum deformation field, the calculation result of SFA is lower than that of the theoretical value, and the relative error is larger. The relative error of CTSFA, however, is closer to the theoretical value. The spatial resolution of CTSFA is higher than that of SFA.(3)Using the programmed CTSFA and SFA to analyze the same Brazilian disk test, CTSFA can observe the crack propagation at the early stage of continuum failure, and the observed crack morphology is in good agreement with the actual situation. In other words, the performance of CTSFA in programmed application is better than that of SFA.

#### Data Availability

No data were used to support this study.

#### Conflicts of Interest

The authors declare that they have no competing financial interests.

#### Acknowledgments

The research described in this paper was financially supported by the Shandong Provincial Natural Science Foundation (no. ZR2019QEE026), National Natural Science Foundation of China (no. 51674160), Tai’shan Scholar Engineering Construction Fund of Shandong Province of China (no. ts201511026), and Taishan Scholar Talent Team Support Plan for Advantaged & Unique Discipline Areas.