## Abstract

We present a beam shaping technique in controlling the complex amplitude of an optical beam. The constraint on the amplitude of the output beam in the Gerchberg-Saxton algorithm is replaced with constraints both on the amplitude and phase of the output beam in the proposed method. The total areas of the constrained regions and free regions on the complex amplitude of the output beam in the proposed method are maintained. An output beam with arbitrary complex amplitude can be realized with the proposed method. The computing result from the proposed method is a phase-only distribution, which can be fabricated as diffractive optical element for higher diffraction efficiency. Both simulations and experiments are present and the effectiveness of the proposed method is verified.

© 2015 Optical Society of America

## 1. Introduction

Laser beams with desired intensity and phase distributions are pursuit in many areas such as optical trapping, optical sensing, laser machining, and so on. However, a conventional laser emits a beam with a fixed distribution such as Gaussian profile only. To obtain an optical beam with desired intensity or phase profile we need use a diffractive optical element (DOE) or computer-generated hologram to modulate the wavefront of a laser beam. The DOEs are generally designed as phase-only to achieve higher diffraction efficiency. Phase distribution of a phase-only DOE can be obtained by running a beam shaping algorithm such as iterative [1], genetic [2], simulated annealing [3], and etc. Commercial softwares [4] for solving the beam shaping problems are also available. Among the algorithms, the Gerchberg-Saxton (GS) algorithm [5] is the most widely used for solving the phase-retrieval problems. Since the GS algorithm is simple and effective for beam shaping, many modified algorithms [6–8] are also proposed to improve diffraction efficiency of designed DOEs. However, the GS-based algorithms can control only amplitude or phase of the complex amplitude of an output beam. In the fields such as optical trapping, we need not only the desired intensity distributions of optical beams to trap particles in the desired positions but also specified phase distributions to rotate, transport, and/or sort particles. A vortex beam can rotate the trapped particles automatically as the beam possesses a helical phase [9]. Roichman et al. proved that an optical line trap with uniform intensity and gradient phase can transfer particles along the optical line automatically owing to the optical forces arising from the phase gradient of the beam [10]. Since control of both the amplitude and phase of an optical beam is important for many applications, researchers proposed various algorithms to beam shape the complex amplitude of an arbitrary optical beam [11–17]. For example, non-iterative algorithms were proposed to generate a beam whose intensity and phase were prescribed along arbitrary 3D curves [12,13]. Rigorous approaches were also proposed to directly calculate phase-only distributions of DOEs for beam shaping of complex-amplitudes of output beams [14–16]. In [17] a modified GS algorithm which put constraints on the optical fields of non-zero intensity in the output domain was proposed to generate an arbitrary phase-gradient optical beam. However, in [17] the locations and profiles of the constrained regions have to be manually chosen, and an adjustable constant, which modifies the amplitude in the iteration, is obtained with tedious calculation. In this letter a method is proposed to control both the amplitude and phase of an arbitrary output beam in a beam shaping system. The proposed method is based on the GS algorithm, but the constraint on the amplitude of the output beam in the GS is replaced with constraints both on the amplitude and phase distributions of the output beam instead. Thus, the computed result from the proposed method remains a phase-only distribution, which can be fabricated as relief structure with current micro-fabrication technology. The proposed method is simple in operation, efficient in beam shaping, and fast in computing. Moreover, arbitrary complex amplitudes for an optical beam can be realized conveniently.

## 2. Method

A beam shaping system mainly comprises an input beam in the input plane, an output beam in the output plane, and a beam propagation method for the beam propagating between the planes. For the GS-based system, amplitudes of the input and output beams on the input plane and output plane are given as constraints, and phases of the input and output beams are free for value assignment. The beam propagation method can be near field or far field diffraction, depending on the optical layout of the exact system. In the GS-based beam shaping system, the beam propagates forward and backward between the input and output planes iteratively. When the given condition is satisfied, the iteration terminates, and the phase distribution of the input plane will be saved as the computing result for the fabrication of a phase-only DOE. In the GS-based system, both the amplitudes in the input beam and output beam are constrained, and both the phases in the input and output beams are free for calculation. Hence, the constrained areas in the input and output planes are equal, and the free areas in the input and output planes are also equal. In the proposed method, the constraint on the amplitude of the output beam is replaced with a combination of two constraints on the amplitude and phase of the output beam. The sum of the constrained domains of the amplitude and phase of the output beam is kept as not greater than the sum of the constrained domain of the amplitude plane in the GS algorithm, so the freedom in the proposed method is not less than that in the GS algorithm. In Fig. 1 a schematic diagram shows the freedoms and constraints of the proposed beam shaping system. The blank region stands for the freedom, and the gridding region stands for the constraint.

In Fig. 1, all the four domains ‘a’, ‘p’, ‘A’, and ‘P’ have equal size of $m{}_{x}\cdot {m}_{y}$ pixels, where $m{}_{x}$ and $m{}_{y}$ are the pixel numbers in the $x$ and $y$ directions, respectively. In the conventional system, there are $2m{}_{x}\cdot {m}_{y}$constrained pixels, which are equally allocated to the amplitude domains of the input and output planes, and $2m{}_{x}\cdot {m}_{y}$free pixels, which are equally allocated to the phase domains of the input and output planes. In the proposed system, the input plane still has $m{}_{x}\cdot {m}_{y}$ constrained pixels in the amplitude domain and $m{}_{x}\cdot {m}_{y}$ free pixels in the phase domain. But in the output plane, the constrained $m{}_{x}\cdot {m}_{y}$ pixels are split into two parts, one is allocated to the amplitude domain and another is allocated to the phase domain. The sum of the constrained pixels in the amplitude and phase domains in the output plane should not be greater than $m{}_{x}\cdot {m}_{y}$, but sums of the constrained pixels in the amplitude and phase domains of the output plane can be unequal. The remained pixels in the amplitude and phase domains of the output plane are set as freedoms. Thus, the total freedoms in the beam shaping system of the proposed method are maintained. Since both the amplitude and phase in the output plane are partly restrained, diffraction efficiency of the computed DOE by the proposed method could be affected.

In the proposed algorithm, we use modified complex amplitude to replace the calculated complex amplitude of the output beam in the iteration. The modified complex amplitude of the output beam $CAmp$ is expressed in Eq. (1).

The flow chart of the proposed method is similar to that of the GS algorithm and is shown in Fig. 2, where $n$ and $N$ are the number of the iteration run and the total iteration number, respectively. The amplitude of the input beam $am{p}_{t}$ is constrained and usually set as a uniform intensity, which can be realized by collimating a laser beam. The initial phase of the input beam $ph$ is generated randomly, and $n$ is set starting from 1.

In the iteration, $T$ and ${T}^{-1}$represent the transform functions of the forward propagation and backward propagation of the beam, respectively. Generally, the free-space propagation can be classified as near field or far-field diffraction, depending on the parameters such as pitch size of the DOE, wavelength of the beam, propagation distance between the input plane and the output plane, and so on. Hence, the transforms can be expressed by Eq. (2), where we use Fourier transforms and the angular plane wave spectrum (APWS) method [18] to simulate the far-field and near-field diffractions, respectively.

We use diffraction efficiency $\eta $ and non-uniformity $\epsilon $ both expressed in Eq. (3) to evaluate beam shaping effects of the reconstructed intensity and phase of a computed DOE, respectively.

## 3. Simulations and experiments

In the simulations, the parameters are preset as follows, the sampling points $Nx=Ny=512$, wavelength $\lambda $ = 532 nm, pixel size $pitch=15\mu m$, and length (equal to width) of the DOE $L=512\times pitch$. Firstly, we design an optical line beam with uniform intensity and phase-gradient along a line as the target output beam. We run the proposed algorithm with the Fourier transforms and the APWS method, respectively. Intensity and phase distributions of the target output beam are shown in Figs. 3(a) and 3(b), respectively. The intensity and phase distributions of the beam reconstructed by the Fourier transforms are shown in Figs. 3(c) and 3(d), respectively. The intensity and phase distributions of the beam reconstructed by the APWS method are shown in Figs. 3(e) and 3(f), respectively. The phases in Figs. 3(d) and 3(f) are normalized to the range from 0 to 2$\pi $. It can be seen that both the intensity and phase distributions of the reconstructed beams are compatible with those of the target beam.

For a better view we also show the plots of the reconstructed intensity and phase profiles of the output beams in Fig. 4. The reconstructed intensity and phase distributions by the Fourier transforms are shown in Figs. 4(a) and 4(b), respectively, and the reconstructed intensity and phase distributions by the APWS method are shown in Figs. 4(c) and 4(d), respectively. The curve in each plot of Fig. 4 cuts through the center of the reconstructed profile of the output beam and along the direction of the phase gradient. It is worth mentioning that, for a convenient view, we only show phases of the reconstructed beams in the region corresponding to the non-zero intensity region of the target beam, as phases in the zero-intensity region of a beam have no physical values. It can be seen that the intensity of the reconstructed beam by the Fourier transforms fluctuates remarkably whereas the phase is in good agreement with the designed phase gradient. On the contrary, the intensity of the reconstructed beam by the APWS method has better uniformity. Although the phase has a spike, the value is very close to 2$\pi $ and can be normalized to 0. So the phase is still in good agreement with the designed. Phase distributions of the DOEs computed by the Fourier transforms and the APWS method are shown in Figs. 5(a) and 5(b), respectively. The phase values range from 0 to 2$\pi $, corresponding to the scaled grey values from 0 to 256. We use equations of the diffraction efficiency and non-uniformity in Eq. (3) to evaluate the performance of the computed DOEs. Diffraction efficiencies for the DOEs computed by the Fourier transforms and the APWS method are less than 5%. However, if we calculate diffraction efficiency of the reconstructed light in the constraint region of the output plane only, i.e., the ratio between the light energy of the reconstructed output beam in the constrained region and the light energy of the target output beam in the constrained region, the diffraction efficiencies for the Fourier transforms and the APWS method are about 68% and 97%, respectively. Actually, the diffraction efficiency varies with complex amplitude distribution of a target output beam. We choose a window function corresponding to the non-zero intensity regions of the target beam to calculate the non-uniformities of phases. The non-uniformities of the phases for the beams reconstructed by the Fourier transforms and the APWS method are about 0.5% and 2%, respectively, which exhibit high fidelity between the reconstructed phases and the designed ones. The computed non-uniformities also vary with actual phase distribution of a target output beam.

In the experimental beam shaping configuration, a continuous wave laser beam is emitted from a diode-pumped laser (Coherent, Genesis MX532-1000 STM), expanded and collimated by a spatial filter and a convex lens. The laser beam has a wavelength of 532 nm and a maximum optical power of 1 W. Then the light is modulated and reflected by a spatial light modulator (SLM, Phase series from BNS), which is loaded with the phase distribution of a computed DOE through a personal computer interface. The SLM has 512 × 512 pixels and the pixel area is $15\mu m\times 15\mu m$. Finally, the light is received by the charge-coupled device (CCD) of a beam profiler (Newport, LBP-2-USB). The experimental setup is schematically shown in Fig. 6.

CCD captured intensity profiles of Fourier transform reconstructed line beams with phase gradients from 0 to 2$\pi $, 0 to 20$\pi $, 2$\pi $ to 0, and 20$\pi $ to 0 are shown in Figs. 7(a)-7(d), respectively. From Fig. 7 we can observe that although the desired maximum phases and directions of phase gradients vary greatly, the intensity distributions of the reconstructed beams are almost unchanged. However, the intensity distributions have a rectangular frame with lower intensity in the center region rather than the uniform intensity in the desired region. This could be resultant from the intrinsic problem of the GS-based algorithms in beam shaping top-hat intensity distributions.

For the verification of the reconstructed phase distribution of an output beam, we choose a vortex beam as the target beam since the vortex beams are widely used in optical rotation of particles induced by the helical phases of the beams. The target intensity and phase of the vortex beam is obtained by a Fourier transform or the APWS propagation of an input beam which has the uniform intensity and the helical phase $\mathrm{exp}(i\cdot l\cdot \theta )$, where $l$ is the topological charge and $\theta $ is the azimuthal angle, respectively. The target intensity and phase of the helical beam diffracted by the Fourier transform are shown in Figs. 8(a) and 8(b), respectively. The topological charge of the vortex beam is 10. The reconstructed intensity profile of the output beam by the computed DOE is shown in Fig. 8(c), and the phase distribution of the computed DOE is shown in Fig. 8(d). Note that, for a better view we display only the central regions of the intensity distributions in Figs. 8(a) and 8(c). In a similar way, we can design the target intensity and phase of the vortex beam obtained by the APWS method. The target intensity and phase are shown in Figs. 9(a) and 9(b), respectively. The topological charge of the target vortex beam is also 10. The reconstructed intensity profile of the computed DOE is shown in Fig. 9(c), and the computed DOE is shown in Fig. 9(d). We can observe from Fig. 9(c) that the phase distribution of the DOE is the same as an analytic solution to generation of a vortex beam, and the effective region has been expanded to the whole DOE. Hence, the computed DOE can also achieve high diffraction efficiency. Using the proposed algorithm, we have generated the vortex beam based on both the amplitude and phase distributions deduced from the propagation of an input beam. On the contrary, the common way to generate a vortex beam is using an analytic expression for the phase distribution of DOE. With the proposed method, one can conveniently design various forms of special vortex beams if the desired amplitudes and phases of the vortex beams are given as the targets.

We can observe from Figs. 8 and 9 that the intensity distributions of the target and the reconstructed are in good agreement. We also find from simulations that the reconstructed phases and the target phases are matching. It is worth mentioning that, for the phase distributions of DOEs generated by the APWS based algorithm, we can make use of the effective region, whose size is equal to that of the constrained regions in the amplitudes, as a truncated DOE for the fabrication, i.e., if a DOE is 512 × 512 pixels and the corresponding effective region is 362 × 362, we can use the region of 362 × 362 only for the beam shaping. On the other hand, if we want to design a 512 × 512 DOE, we can design the DOE with a size of 724 × 724 by the APWS based method firstly and then truncate the DOE into a size of 512 × 512. The simulations find that the truncated DOE can beam shape the designed complex amplitude of the output beam with high diffraction efficiency.

To verify the existence of helical phase in the reconstructed vortex beams, we incorporate the beams into optical tweezers to manipulate particles. The optical tweezers mainly comprise the above mentioned beam shaping system and an inverted optical microscope (Olympus, IX71). The SLM-modulated light is re-directed and focused into the sample stage by a series of geometric optical components and a 100 × oil-immersed objective lens. We load the phase distribution of the computed DOE obtained by the Fourier transforms onto the SLM. The reconstructed image in the sample cell of the microscope is shown in Fig. 10(a). The particles to be trapped are copper beads with diameter from 1 to 5$\mu m$ and immersed in de-ionized water. When the power of laser is increased to 650 mW, we observe that two particles are trapped and rotated, and if we change the sign of the topological charge of the phase in the computed DOE, the rotation direction of the trapped particles changes accordingly. Four captured video frames showing the trapped particles in rotation are displayed in Fig. 10(b), where the dashed rectangles show the position of the focused beam. Thus, the reconstructed vortex beams by the DOEs are proved to possess not only the doughnut shape but also the helical phase, which are in accordance with designed vortex beams.

## 4. Conclusion

To summarize, we have proposed a beam shaping algorithm based on the GS algorithm. The proposed method puts constraints both on the amplitude and phase of an output beam in the iteration and modifies the constraint conditions both in the amplitude and phase of the output beam, so the complex amplitude of the output beam can be tailored simultaneously. The simulations and experiments have verified the effectiveness of the proposed method. The reconstructed intensity and phase distributions of the computed DOEs are in good agreement with those of the designed ones. The proposed method is simple and convenient for implementation, and the running time is almost the same short as that of the GS algorithm. The proposed method can find extensive applications in optical trapping, optical sensing, laser machining, and so on.

## Acknowledgments

The research was financially supported by the National Natural Science Foundation of China (Grant No. 61178017). The authors thank Shubo Cheng for his assistance in the acquisition of experimental images.

## References and links

**1. **J. R. Fienup, “Phase retrieval algorithms: a personal tour [Invited],” Appl. Opt. **52**(1), 45–56 (2013). [CrossRef] [PubMed]

**2. **G. Zhou, Y. Chen, Z. Wang, and H. Song, “Genetic local search algorithm for optimization design of diffractive optical elements,” Appl. Opt. **38**(20), 4281–4290 (1999). [CrossRef] [PubMed]

**3. **S. Kirkpatrick, C. D. Gelatt Jr, and M. P. Vecchi, “Optimization by simulated annealing,” Science **220**(4598), 671–680 (1983). [CrossRef] [PubMed]

**5. **R. W. Gerchberg and W. O. Saxton, “A practical algorithm for the determination of phase from image and diffraction plane pictures,” Optik (Stuttg.) **35**(2), 237–246 (1972).

**6. **J. S. Liu and M. R. Taghizadeh, “Iterative algorithm for the design of diffractive phase elements for laser beam shaping,” Opt. Lett. **27**(16), 1463–1465 (2002). [CrossRef] [PubMed]

**7. **H. Zhai, F. Liu, X. Yang, G. Mu, and P. Chavel, “Improving binary images reconstructed from kinoforms by amplitude adjustment,” Opt. Commun. **219**(1–6), 81–85 (2003). [CrossRef]

**8. **S. H. Tao and X. C. Yuan, “Practical implementation of the phase-quantization technique in an iterative Fourier-transform algorithm,” Appl. Opt. **43**(10), 2089–2092 (2004). [CrossRef] [PubMed]

**9. **K. T. Gahagan and G. A. Swartzlander Jr., “Optical vortex trapping of particles,” Opt. Lett. **21**(11), 827–829 (1996). [CrossRef] [PubMed]

**10. **Y. Roichman, B. Sun, Y. Roichman, J. Amato-Grill, and D. G. Grier, “Optical forces arising from phase gradients,” Phys. Rev. Lett. **100**(1), 013602 (2008). [CrossRef] [PubMed]

**11. **A. Jesacher, C. Maurer, A. Schwaighofer, S. Bernet, and M. Ritsch-Marte, “Full phase and amplitude control of holographic optical tweezers with high efficiency,” Opt. Express **16**(7), 4479–4486 (2008). [CrossRef] [PubMed]

**12. **S.-H. Lee, Y. Roichman, and D. G. Grier, “Optical solenoid beams,” Opt. Express **18**(7), 6988–6993 (2010). [CrossRef] [PubMed]

**13. **J. A. Rodrigo, T. Alieva, E. Abramochkin, and I. Castro, “Shaping of light beams along curves in three dimensions,” Opt. Express **21**(18), 20544–20555 (2013). [CrossRef] [PubMed]

**14. **J. A. Davis, D. M. Cottrell, J. Campos, M. J. Yzuel, and I. Moreno, “Encoding amplitude information onto phase-only filters,” Appl. Opt. **38**(23), 5004–5013 (1999). [CrossRef] [PubMed]

**15. **V. Arrizón, U. Ruiz, R. Carrada, and L. A. González, “Pixelated phase computer holograms for the accurate encoding of scalar complex fields,” J. Opt. Soc. Am. A **24**(11), 3500–3507 (2007). [CrossRef] [PubMed]

**16. **L. Z. Liu, K. O’Keeffe, D. T. Lloyd, and S. M. Hooker, “General analytic solution for far-field phase and amplitude control, with a phase-only spatial light modulator,” Opt. Lett. **39**(7), 2137–2140 (2014). [CrossRef] [PubMed]

**17. **Z. Z. Yuan and S. H. Tao, “Generation of phase-gradient optical beams with an iterative algorithm,” J. Opt. **16**(10), 105701 (2014). [CrossRef]

**18. **J. W. Goodman, *Introduction to Fourier Optics* (McGraw-Hill, 1996).