Laser-wavefront measurement is essential in the quality analysis of laser beams [1]. According to the actual measurement, the laser wavefront in general is highly influenced by tilt, defocus, and astigmatism. The phase-diversity phase retrieval (PDPR) [2, 3] algorithm is one of the critical methods of wavefront detection, using a nonlinear optimization algorithm to try to match the collected laser-spot image to a calculated image based on scalar diffraction propagation.

Nonlinear optimization requires reasonable starting guesses. An excellent initial estimate improves the accuracy of the PDPR algorithm. Currently, the starting guesses are usually set to zero or random values [4, 5]. However, this leads to an issue known as the “Capture-range problem”: The nonlinear optimization algorithm falls into a local minimum during iteration when the starting guesses for the wavefront are too far from the ground-truth wavefront [4, 6]. Turman [7] proposed a geometrical-optics method to obtain wavefront tilt data for the phase retrieval (PR) algorithm, and it was modified by Jurling [4, 8, 9] and Carlisle [10]. The method improves the capture range when segment tilt errors are large. Paine [6, 11], Cao et al. [12], and Zhao et al. [13] proposed methods to overcome the capture-range problem by calculating a piston for segmented systems. However, they all ignored the influence of other low-order aberrations. The application of the currently proposed particle-swarm optimization (PSO) [14] may result in an incorrect solution, due to the stagnation problem. Machine-learning methods require large representative datasets for cumbersome model training [15–19]. Therefore, the “Capture-range problem” with the PDPR algorithm remains to be solved with high accuracy and less time.

In this paper, we propose a geometrical-optics initial-estimate method to improve the capture-range problem in the normal PDPR algorithm for laser-wavefront detection. Two spot images in the laser’s optical path are collected to build geometrical-optics theoretical models. We calculate the low-order aberration based on the models and set it as the initial estimate for the PDPR algorithm. The proposed method provides a good idea for selecting the starting guesses in the laser-wavefront restoration using the PDPR method, and has strong practicability.

We use the Zernike polynomials of Fringe [20] to represent the laser wavefront. The low-order Zernike terms are extremely important in wavefront aberration. In this section, we present two geometrical-optics theoretical models to calculate the low-order aberrations of the laser spot. These aberrations include x-tilt, y-tilt, defocus, 0° astigmatism, and 45° astigmatism. We analyze the propagation path of the laser in the model and carefully derive the theoretical formula for aberration calculations. It should be noted that we only consider a single aberration and ignore all others in each derivation process.

We take x-tilt aberration as an example for analysis, since the derivation process for x-tilt and y-tilt aberrations are the same. As shown in Fig. 1, the light spot propagates from Fig. 1(a) to Fig. 1(b). Due to wavefront distortion, the position of the light spot in Fig. 1(b) is shifted from its position in Fig. 1(a). According to geometrical optics and the geometric relationships in Fig. 1,

Figure 1. Schematic of the x-tilt term. (a) and (b) laser spot images. l, image width; T_x, x-tilt aberration; WF, the line formed by the intersection of the wavefront and the xz-plane; ∆z, propagation distance; ∆d, the deviation of the centroids of the two light spots on the x-axis.

(1)$$\tan \theta =\left\{\begin{array}{c}\Delta d/\Delta z\\ T_x/\left(l/2\right)\end{array}\right.,$$

where θ is the angle between the ray and the z-axis and the angle between the WF (the line formed by the intersection of the wavefront and the xz-plane in Fig. 1) and the x axis, ∆d is the deviation of the centroids of the two light spots on the x axis, ∆z is the propagation distance, T_x is the x-tilt aberration, and l is the image width. We can obtain

(2)$$T_x=\frac{\Delta d\times l}{2\Delta z},$$

and T_x is also the x-tilt term’s coefficient in the Zernike polynomial.

The defocus and astigmatism terms share a geometrical-optics model for derivation and calculation. The astigmatism term is divided into separate terms for 0° astigmatism and 45° astigmatism. The term of highest degree in the Zernike polynomial for defocus and astigmatism is quadratic. So, we use the quadratic model.

We take the derivation of the defocus term as an example to explain the model shown in Fig. 2. The light spot propagates from Fig. 2(a) to Fig. 2(b). The spot radii in Fig. 2(a) and Fig. 2(b) are different due to wavefront distortion. The light is emitted from point A to point B on the xz plane. We set point O of Fig. 2(a) as the origin of the optical system. The term of highest degree in the Zernike polynomial for defocus is quadratic, and there is no first-order term, so WF of Fig. 2 is symmetrical about the z-axis. We set WF of Fig. 2 as

Figure 2. Schematic of the defocus and astigmatism terms. WF, the curve formed by the intersection of the wavefront and xz-plane.

(3)$$w=a{x}^2+c,$$

where a and c are coefficients. We set F (w, x) = ax² + c − w and calculate the partial derivative

(4)$$\left\{\begin{array}{c}\frac{\partial F}{\partial x}=2ax,\\ \frac{\partial F}{\partial w}=-1.\end{array}\right.$$

According to Eq. 4, the normal vector at point (r_x1, $a{r}_{x_1}^2$ + c) is n = (2ar_x1, − 1), where r_x1 is the radius of the light spot on the x axis in Fig. 2(a). Point (r_x1, $a{r}_{x_1}^2$ + c) is not necessarily point A, but the propagation direction of the light AB is the same as the direction of n. The slope of line AB is expressed by k,

(5)$$k=-\frac{1}{2a{r}_{x_1}}.$$

According to geometrical optics and the geometric relationships in Fig. 2,

(6)$$\left\{\begin{array}{c}k=\frac{\Delta z}{\Delta r_x},\\ \Delta r_x=r_{x_2}-r_{x_1},\end{array}\right.$$

where ∆z in Fig. 2 is the propagation distance and r_x2 is the radius of the light spot on the x axis in Fig. 2(b). According to Eqs. 5 and 6, we get

(7)$$\left\{\begin{array}{l}a=-\frac{\Delta r_x}{2r_{x_1}\Delta z},\\ c=-a{r}_{x_1}{}^2.\end{array}\right.$$

So, we obtain w_max and w_min. The defocus aberration on the xz-plane is D_x = (w_max − w_min) / 2. Similarly, the defocus aberration of the yz-plane can be calculated, and is expressed by D_y. We denote the defocus term coefficient as D, and

(8)$$D=\frac{D_x+D_y}{2}.$$

We denote the 0°-astigmatism term as A₀ and calculate A_0x and A_0y according to the calculation method of D_x and D_y. Astigmatism causes the difference in points of convergence between meridional and sagittal beams.

(9)$$A_0=A_{0_x}-A_{0_y}.$$

The calculation method for the 45°-astigmatism term has one more step than the calculation method for the 0°-astigmatism term. We need to rotate the two light-spot images counter-clockwise by 45 degrees in advance when calculating the 45°-astigmatism term, and the subsequent calculation steps are the same.

In this section, we experiment with initial estimates for the simulation. The experiment aims to verify the accuracy of the theory presented in the second section.

The principle of the simulation experiment is shown in Fig. 3. We simulate a Gaussian laser spot’s intensity image I₁ with resolution of 256 × 256 pixels. The intensity image I₂ can be obtained with W₁ and I₁ by angular spectrum propagation (ASP). We apply the proposed algorithm to I₁ and I₂ to obtain the predicted wavefront map W₂.

Figure 3. Schematic of the aberration model’s validation. W₁, wavefront map of simulation; I₁, intensity image of simulation; ASP, angular spectrum propagation; I₂, intensity image after ASP; W₂, predicted wavefront map after the proposed method; λ, wavelength; ∆z, propagation distance; ps, pixel size.

The radius of the I₁ laser spot is a random value in the range 30–80 pixels. Six-term Zernike polynomials are used for reconstructing the wavefront W₁, with piston set to zero. The second through sixth terms are randomly distributed in the range −0.5 λ~0.5 λ that is detected from the GY-10 He-Ne laser. The purpose of this experiment is to verify the theoretical correctness of the algorithm, so we do not add noise to I₁ and I₂. In addition, λ = 632.8 nm and ∆z = 5 mm. The pixel size is set to 6.45 μm × 6.45 μm. We take 100 runs of the same simulation experiment and record the difference between each group’s Zernike coefficients W₁ and W₂. Figure 4 presents the mean absolute value of the difference for each term.

Figure 4. The errors between the low-order aberrations calculated by the proposed method and the true values.

As shown in Fig. 4, the maximum error did not exceed 0.03 λ. The calculated low-order aberrations were closed to the true value, verifying the proposed method’s correctness. In addition, the reason for the large defocus aberration may be that the shape of the laser wavefront is approximated as spherical when calculating the defocus aberration.

After that, we apply the initial value predicted by the proposed method to the PDPR algorithm. Additionally, the accuracy of the wavefront sensors is no less than 0.01, so it is inappropriate to use the value detected by the sensors as the true value. Therefore, simulation experiments are carried out to verify whether the proposed method can improve the PDPR algorithm.

As shown in Fig. 5, W₃ could be restored with I₁ and I₂ by the PDPR algorithm. We calculated the root mean square (RMS) error of W₁ and W₃ to show the accuracy of PDPR wavefront recovery. We used the limited-memory Broyden-Fletcher-Goldfarb-Shanno (LBFGS) algorithm as the nonlinear optimization algorithm of the PDPR algorithm.

Figure 5. Schematic of the proposed method applied to phase-diversity phase retrieval (PDPR). I₁ and I₂, intensity map before and after angular spectrum propagation (ASP); W₁, wavefront map before ASP; W₃, wavefront map recovered by PDPR; IP, initial iteration value of the limited-memory Broyden-Fletcher-Goldfarb-Shanno (LBFGS) algorithm [21] is predicted by the proposed method; IZ, initial iteration value of the LBFGS algorithm is zero; IR, initial iteration value of the LBFGS algorithm is random.

Traditionally, zeros or random numbers are generally used as the initial values for the PDPR algorithm. In this article, we use the values calculated by the method as the initial values for PDPR.

The values for λ, ∆z, and ps are consistent with the initial-estimate experiment. Fifteen-term Zernike polynomials are used for reconstructing the wavefront W₁, with piston set to zero. Subsequently, we detect the laser wavefront of the GY-10 laser at the different positions using Beam Wave 500 (PhaseView, Verrieres Le Buisson, France) to obtain the range of the various aberrations. According to the test results in Beam Wave 500 (Phaseview), the second through sixth terms are randomly distributed in the range −0.5 λ~0.5 λ, and the others are in −0.05 λ~0.05 λ. The background noise of a camera [Beam Wave 500 (Phaseview)] is added to I₁ and I₂, to be closer to reality. The simulation experiment is conducted 100 times. All of the experiments are carried out on a computer with Intel i7-6700 CPU (3.40 GHz) and 16 GB of RAM. We record the difference between the Zernike coefficients W₁ and W₂; Figure 6 presents the mean absolute value of the difference.

Figure 6. Comparison of Zernike-coefficient differences. y-axis: The difference between the true Zernike coefficients and the restored Zernike coefficients from the limited-memory Broyden-Fletcher Goldfarb-Shanno (LBFGS) algorithm with different initial values.

As shown in Fig. 6, the error for initial iteration value of the LBFGS algorithm predicted by the proposed method (IP) is closer to 0 in the noisy environment. We compare the RMS error of the restored wavefront obtained from different initial values and the ground-truth wavefront to verify the influence of different initial values on the wavefront. The Zernike coefficients from 2^nd to 15^th orders are combined to yield the RMS value. Table 1 shows the mean RMS error for 100 experiments. For a more intuitive comparison, Fig. 7 shows the wavefront map and residuals map in one of the experiments under the background noise of the Phaseview Beam Wave 500’s camera.

TABLE 1. System performance for different initial values in a noisy environment.

Initial Value	RMS Residual Error (λ)			Total Time(s)
	Max	Min	Mean
IP	0.0588	0.0065	0.0304	19.657
IZ	0.1319	0.0129	0.0516	24.143
IR	0.1802	0.0057	0.0470	26.557

Figure 7. Maps of wavefronts and residuals. (a) Ground-truth wavefront map. (b)–(d) are the restored wavefront maps for IP, IZ, and IR. (e)–(g) are the maps of the residuals between a and (b)–(d). IP, initial iteration value of the limited-memory Broyden-Fletcher-Goldfarb-Shanno (LBFGS) algorithm predicted by the proposed method; IZ, initial iteration value of the LBFGS algorithm is zero; IR, initial iteration value of the LBFGS algorithm is random.

From Table 1, it can be seen that the RMS error for IP was smaller than those for initial iteration value of the LBFGS algorithm is zero (IZ) and initial iteration value of the LBFGS algorithm is random (IR). The iteration time for IP is reduced by about 28.96% (5.69 seconds) compared to the other initial values. As shown in Fig. 7, the wavefront maps generated by different initial values in Figs. 7(b)–7(d) are similar to the ground-truth wavefront map in Fig. 7(a). However, we can find that our method achieves smaller residuals than others, by comparing to Figs. 7(e)–7(g). Therefore, the proposed method can obviously improve the accuracy of wavefront reconstruction.

In this paper, two geometric models were constructed to calculate the low-order aberration of a laser beam’s wavefront, which provided a high-accuracy initial value for the PDPR method. We conducted simulation experiments. In noise-free simulations, the maximum errors between the calculated low-order aberrations and the true value were less than 0.03 λ. In a noisy environment, the wavefront with the initial value calculated by our method was closer to the ground-truth wavefront, compared to initial values of zero or a random number. Also, the iteration time with our method was reduced by 28.96%. The proposed method was proven to solve the capture-range problem.

The authors declare no conflict of interest.

All of the data have been listed in this work; no extra data were generated in the presented research.

111 project of China (D21009); Ministry of Science and Technology of the People’s Republic of China (2018YFB 1107600); Jilin Scientific and Technological Development Program (20160204009GX, 20170204014GX).