## Abstract

Benefitting from the development of commercial spatial light modulator (SLM), holographic optical tweezers (HOT) have emerged as a powerful tool for life science, material science and particle physics. The calculation of computer-generated holograms (CGH) for generating multi-focus arrays plays a key role in HOT for trapping of a bunch of particles in parallel. To realize dynamic 3D manipulation, we propose a new tilted-plane GS algorithm for fast generation of multiple foci. The multi-focal spots with a uniformity of 99% can be generated in a tilted plane. The computation time for a CGH with 512×512 pixels is less than 0.1 second. We demonstrated the power of the algorithm by simultaneously trapping and rotating silica beads with a 7×7 spots array in three dimensions. The presented algorithm is expected as a powerful kernel of HOT.

© 2020 Optical Society of America under the terms of the OSA Open Access Publishing Agreement

## 1. Introduction

Optical tweezers has emerged as a powerful and flexible tool for manipulating biological and colloidal particles in a non-contact manner [1,2]. Although the standard Gaussian single-beam trap remains a useful tool, simultaneous manipulation of a large number of microscopic objects in three-dimensional (3D) space is appealing for parallel optical tweezers [3–5]. Originally, multiple optical traps were implemented by splitting or rapidly scanning a single beam [6,7]. The number of the trapped objects in both of the cases is limited to less than 10 for the system complexity and cycle time. On the other hand, interference of two or more coherent beams can produce a large number of local intensity maximum, but the positions of the local intensity maximum are difficult to change [8]. Holographic optical tweezers (HOT) provide another convenient way to produce optical trap arrays [9–11], which usually utilize liquid crystal spatial light modulators (SLM) to shape the input field into multiple foci. More importantly, the number and the position of foci can be easily varied by changing the computer-generate hologram (CGH) loaded on the SLM without any mechanical motions. The tailor-made optical trap arrays have exciting and surprising applications, such as continuously sorting fluid-borne particles, controlling microscopic interactions of colloidal particles, and rotationally melting the substrate-induced crystal [12–14].

To date, many algorithms have been developed to calculate the CGH, including non-iterative algorithms [15] and iterative algorithms [16]. Compared to the non-iterative algorithms, the iterative algorithms are more suitable for producing complex trap arrays. One of the well-known iterative algorithms is the Fourier-transform-based Gerchberg-Saxton (GS) algorithm [17], in which the CGH that yields the desired fields in the focal region is obtained by searching iteratively between the entrance plane and the focal plane. By utilizing phase induction technique [18], GS algorithm has been utilized to realize dynamic manipulation of multiple optical traps at video or higher frame-rates, but these traps are confined to the focal plane. To extend the algorithm to 3D space, one should introduce defocus and refocus factors in the iterative processing and superpose phase patterns of light fields from multiple discrete transverse planes [19,20]. In this way, the multi-plane 3D GS algorithm can be used to create multiple foci spreading over the whole focal volume. In addition, the 3D Fourier-transform (FT)-based algorithm was also developed to generate 3D multiple foci with a high uniformity (nearly 99%) [21,22]. However, by repeating the complex 3D FT operation in the iteration, the algorithm suffers from slow processing speed and high computation cost. Therefore, the CGHs need to be pre-calculated and then displayed in sequence, which prevents the HOT from dynamic 3D applications.

To realize dynamic 3D manipulation, we propose a tilted-plane GS algorithm for fast generation of multiple foci in a tilted-plane. The tilted plane method has appeared in different contexts like digital holography and optical security [23–26]. Our work focuses mainly on beam shaping and provides a GS algorithm suitable for generation of multi-focus arrays in an arbitrary tilted-plane. By using this algorithm, the time required to compute a CGH with 512×512 pixels for 3D manipulation of multiple traps can be reduced to 0.1 s. We demonstrated the effectiveness of the algorithm by manipulating silica beads in three dimensions and realizing dynamic manipulation of high-quality multiple optical traps. The study offers new insights into applications such as tomographic phase microscopy [27] that entails dynamic manipulation of multiple living cells in suspension.

## 2. Methods

#### 2.1 Principle

In HOT, the manipulation of trapped particles is realized by refreshing the holograms loaded on the SLM. So, efficient computation of the hologram is essential for applications involving instantaneous dynamical process. Here, we propose a new tilted-plane GS algorithm that can calculate a hologram in an efficient way to generate a large number of optical traps in a tilted plane. The novel optical traps can be applied to manipulate particles in three dimensions.

The principle of the tilted-plane GS algorithm is similar to that of the conventional GS algorithm, in which the desired hologram is obtained by iterating back and forth between the input field in the SLM panel and the target field in the focal volume as shown in Fig. 1(a). However, the conventional GS algorithm confines the target field to the focal plane of the objective [28,29], which limits its applications in general scenarios. Here, we extend the plane of target field to a general plane. As illustrated in Fig. 1(b), we denote such an arbitrary plane by the normal vector **n**(*θ*_{0}, *ϕ*_{0}) with (*θ*_{0}, *ϕ*_{0}) referring to the spherical polar angles. For instance, the vector **n**(0, 0) corresponds the focal plane. To understand the basic principle of the tilted-plane GS algorithm, we begin with a brief review of the GS algorithm for beam shaping. As shown in Fig. 1(a), a SLM-modulated input field is focused by an objective lens to obtain the focused field. Denote the input field on the SLM and the target field near the focus by *E _{H}*(

*ξ*,

*η*) and

*E*(

_{T}*x*,

*y, z*), respectively. Here, (

*ξ*,

*η*) are the coordinates on the SLM and (

*x*,

*y*,

*z*) are the coordinates in the focal region. The focal field

*E*(

_{T}*x*,

*y, z*) is connected to the input field

*E*(

_{H}*ξ*,

*η*) by the diffraction integral [30,31]

*C*is a constant and (

*k*,

_{x}*k*) =

_{y}*k*(

*ξ*/

*f*,

*η*/

*f*) are the transverse wavenumbers with

*k*and

*f*respectively referring to the wavenumber and the focal length of the objective; the longitudinal wavenumber

*k*= (

_{z }*k*

^{2}-

*k*

_{x}^{2}-

*k*

_{y}^{2})

^{1/2}and the domain of integration is

*θ*denoting the maximum acceptance angle for the objective lens. According to Eq. (1), one can iterate between the illumination field

_{max}*E*(

_{H}*ξ*,

*η*) and the focal field

*E*(

_{T}*x*,

*y, z*= 0) by using the fast Fourier transform (FFT) and its inverse to design a hologram associated with the target field in the focal plane efficiently.

To design foci arrays in an arbitrary tilted plane, an appropriate coordinate transform is required to make FFT applicable for the integral in Eq. (1). The new coordinate system denoted by (*x’*, *y’*, *z’*) is obtained from the old one by a rotation matrix R_{T}_{ }= R* _{z}*(

*ϕ*

_{0})R

*(*

_{y}*θ*

_{0}), as shown in Fig. 1(b). The matrix R

*consists of a rotation that rotating the old coordinate system around the*

_{T}*y*-axis through an angle of

*θ*

_{0}followed by another rotation around the

*z*-axis through an angle of

*ϕ*

_{0}. Then, the new coordinates (

*x’*,

*y’*,

*z’*) are related to the original coordinates (

*x*,

*y*,

*z*) as

*z'*-axis is normal to the

*x'y'*-plane (

*z'*= 0). After applying the transformations in Eqs. (3) and (4), the integral in Eq. (1) now has new integration variables (

*k*,

_{x’}*k*) while the corresponding ‘longitudinal’ wavenumber

_{y’}*k*can take two values: +(

_{z’}*k*

^{2}-

*k*

_{x'}^{2}-

*k*

_{y'}^{2})

^{1/2}or -(

*k*

^{2}-

*k*

_{x'}^{2}-

*k*

_{y'}^{2})

^{1/2}. As a result, the integral domain in the coordinates (

*k*,

_{x’}*k*) has two components

_{y’}*D*

_{+}and

*D*

_{-}as

*k*|=(

_{z'}*k*

^{2}-

*k*

_{x'}^{2}-

*k*

_{y'}^{2})

^{1/2}and (

*k*,

_{x}*k*,

_{y}*k*) expressed in (

_{z}*k*,

_{x’}*k*, |

_{y’}*k*|) via the inverse of Eq. (4). With Eqs. (3)–(5), the integral in Eq. (1) may be rewritten as

_{z'}**k**

_{⊥}=(

*k*,

_{x′}*k*) and

_{y′}**x**

_{⊥}=(

*x′*,

*y′*). The integral domain $\mathbb{D}$={

**k**

_{⊥}, |

*k*|<

_{x'}*k*, |

*k*|<

_{y'}*k*}; χ

_{+}(

**k**

_{⊥}) is the characteristic function of

*D*

_{+}(

*θ*

_{0},

*ϕ*

_{0}), that is, χ

_{+}(

**k**

_{⊥}) = 1 for

**k**

_{⊥}∈

*D*

_{+}(

*θ*

_{0},

*ϕ*

_{0}) and zero for

**k**

_{⊥}∉

*D*

_{+}(

*θ*

_{0},

*ϕ*

_{0}), and likewise for χ

_{-}(

**k**

_{⊥}).

*E*

_{H}_{+}and

*E*

_{H}_{-}are two components of the input field corresponding to the cases

*k*>0 and

_{z'}*k*<0, respectively. The details of derivation from Eq. (1) to Eq. (5) can be found in our previous work [32]. Note that in Eqs. (1) and (5), we have ignored the vector nature of the fields and have dropped the (

_{z'}*k*/

_{z}*k*)

^{1/2}factor associated with energy conservation in tight focusing. In fact, as demonstrated in the following, such a treatment has negligible effect on the quality of generated multiple traps over the entire focal volume. By Eq. (6), the field

*E*in the

_{T}*x'y'*-plane is now expressed as an angular spectrum integral over

*k*and

_{x′}*k*, which can be efficiently calculated via inverse FFT. However, the angular spectrum function is a sum of

_{y′}*E*

_{H}_{+}and

*E*

_{H}_{-}, which are not determined only from the field

*E*(

_{T}*x′*,

*y′*, 0). For definiteness, we assume that

*E*is symmetry in

_{H}*k*, that is,

_{z′}*E*

_{H}_{+}(

*k*,

_{x′}*k*,

_{y′}*k*) =

_{z′}*E*

_{H}_{-}(

*k*,

_{x′}*k*,

_{y′}*-k*), leaving the integral for

_{z′}*E*in the arbitrary plane as

_{T}*E*as well as the hologram to be loaded onto the SLM can be calculated as done in the conventional GS algorithm.

_{H}#### 2.2 Tilted-plane GS algorithm

Based on the integral in Eq. (7), a tilted-plane GS algorithm of finding the input optical field *E _{H}* corresponding to a predesigned tilted-plane focal field

*E*can be implemented. The flow chart of the tilted-plane GS algorithm is presented in Fig. 2. The target (prescribed) field consists of multiple focal spots in the

_{T}*x'y'*-plane: ${I_T} = \sum\nolimits_{i = 1}^M {{I_i}\delta ({x^{\prime} - x{^{\prime}_i},y^{\prime} - y{^{\prime}_i}} )} $, where

*M*denotes the number of focal spots and

*I*denotes the desired peak intensity of the

_{i}*i*-th focal spot. Meanwhile, the hologram field and the focal field are expressed as

*E*

_{H}_{,j}

*=*

_{ }*A*

_{H}_{,j}exp(

*iφ*

_{H}_{,j}) and

*E*

_{T}_{,j}=

*A*

_{T}_{,j}exp(

*iφ*

_{T}_{,j}). In the algorithm, the magnitude of the input field

*A*

_{H}_{,j}is assumed to be constant and the amplitude of the target focal field ${A_{T,0}} = \sqrt {{I_T}} $ while the phase

*φ*of the input field is the unknown to be determined. The

_{H}*x'y'*-plane is defined by the angles (

*θ*

_{0},

*ϕ*

_{0}). The algorithm is performed by propagating the beam forward and backward: the two processes numerically calculated by the inverse FFT and FFT as shown in Fig. 2. The details of the algorithm are listed as follows:

- (1) The iterative process starts with a constant phase
*φ*_{H}_{,0}. - (2) Set the input field outside of the domains
*D*_{+}(*θ*_{0},*ϕ*_{0}) and*D*_{-}(*θ*_{0},*ϕ*_{0}) be zero, equivalent to the constraints by χ_{+}and χ_{-}. - (3) The focal ﬁeld in the
*x'y'*-plane is computed by the inverse FFT of the constrained input field, after which the field’s magnitude is replaced by the target magnitude proﬁle adjusted by a weighting factor while the phase distribution is left unchanged. The weighting factor (WF) will be discussed in the next section. - (4) An improved phase distribution is obtained by letting the reshaped focal ﬁeld propagate backward to the hologram region through FFT. Then the input field is renewed by setting the magnitude to unity.

As shown in Fig. 2, steps (2)-(4) form a closed loop. The iteration stops when the focal field converges to the intensity distribution *I _{T}* and the phase

*φ*(

_{H}*k*,

_{x′}*k*) is obtained. In the experiment, the phase to be loaded onto the SLM is expressed in terms of (

_{y′}*k*,

_{x}*k*). Thus the projection after step 4 is needed for obtaining the actual phase hologram. By loading the hologram onto the SLM, the desired field in the focal region is generated.

_{y}#### 2.3 Weighting factor

As mentioned in Section 2.2, a weighting factor *w* can be introduced to modify the target focal field [16,33], which is defined as

*j*is the iteration number. Then, we can balance the difference among the multiple foci during the iteration to generate a high-quality field distribution. To evaluate the intensity uniformity of the focal field, a uniformity factor [21] is defined as where

*I*and

_{max}*I*are the maximum and minimum peak intensities of the foci, respectively. As an illustration, we consider a 7×7 focus array uniformly distributed in a region of 60×60 µm

_{min}^{2}in the

*xz*-plane. Figures 3(a) and 3(b) show the vectorial calculation [34] of the fields resulting from the holograms with and without considering the weighting factor, respectively. In Fig. 3(c), we plot the normalized intensity distribution of the foci scanned along the blue and red lines in Figs. 3(a) and 3(b). It is obvious that the intensity uniformity of the multiple foci with weighting factor is significantly superior to that without introducing the weighting factor. The uniformity factors for the two cases are measured to be 99% and 60%, respectively, that is, the intensity uniformity of the point arrays is improved by 65% with WF. This result demonstrates the necessity of introducing a balance factor in the iterative process.

## 3. Results

To demonstrate the performance of the presented algorithm, we compare it with the multi-plane 3D GS algorithm, a 3D beam shaping algorithm. While 3D beam shaping algorithms such as 3D FT-based algorithm can produce high-quality optical fields in a tilted plane [35], they are commonly time-consuming when dealing with a large scale of 3D data unless GPU acceleration technique is used [36]. Figure 4 shows the intensity distributions of several on-axis trap arrays (*M* = 3, 5, 7) computed by the vectorial diffraction integration using multi-plane 3D GS algorithm (Fig. 4(a)) and the tilted-plane GS algorithm (Fig. 4(b)). The quality of the focal fields is accessed by the uniformity factor *u*. In the simulation, the numerical aperture (NA) of the objective is set to be 1.27 immersed into water (index of refraction 1.33). The input beam is assumed to have a wavelength of *λ* = 1064 nm and be in circular polarization. From Fig. 4(a), it is clear that the intensity uniformity of the fields given by the multi-plane 3D GS algorithm tends to drop with increasing the number of foci. As a comparison, the tilted-plane GS algorithm gives rise to an unaffected intensity uniformity factor of 99%, owing to the fact that it treats the multi-foci pattern in different *z* positions as an entity. Such an operation can reduce the crosstalk between different foci. However, when the separation distance of two foci is too small, the crosstalk appears, degrading the intensity uniformity. We find the focal field quality can maintain above 95% even if the separation distance is reduced to 2µm with the proposed algorithm.

We next discuss the time efficiencies of the multi-plane 3D GS algorithm and the tilted-plane GS algorithm when generating an array of *M* on-axis foci. All the CGHs are calculated with 512×512 pixels on a personal computer (Intel Core i5-4570 K CPU @3.2 GHz and 8GB RAM) with MATLAB 2017b running on Windows 7 x64. We use the “big-*O*” notation to describe computational time [37]. Assuming the tilted-plane GS algorithm has time complexity of *O*(1), the multi-plane 3D GS algorithm then has time complexity of *O*(*M*). The resulting total computation time of the multi-plane GS algorithm almost increases linearly with the value of *M*. It should note that the multi-plane GS algorithm is a little faster than the tilted-plane GS algorithm when the focus number *M* = 1. This is attributed to the extra time consumed in coordinate transforms in the tilted-plane GS algorithm (∼0.05 s). When the number of foci *M* is larger than 2, the tilted-plane GS algorithm will become much faster than the multi-plane 3D GS algorithm with an almost unchanged computing time of ∼0.1 seconds.

To test the practical performance of the tilted-plane GS algorithm, we built a HOT system as illustrated in Fig. 5. A linearly polarized CW laser beam with power of 5 W and wavelength of *λ* = 1064 nm is expanded and collimated by a telescope consisting of Lens 1 and Lens 2, and is then reflected onto a phase-only spatial light modulator (SLM, Pluto ΙΙ, HoloEye Photonics AG, Germany) by a specially designed triangle reflector. The required shaping beam is accomplished by loading the CGH onto the SLM. The modulated beam then passes through a 4*f* system consisting of Lens 3 and Lens 4. A quarter-wave plate (QWP) is placed behind the Lens 4 to convert the linearly polarized beam into a circularly polarized beam. Then the beam is reflected by a beam splitter into the back aperture of an infinity-corrected objective (60×, NA1.27, water-immersion, E Plan, Nikon Inc., Japan) to form the desired 3D optical traps in the focal volume. The water-immersion objective is adopted here to minimize the spherical aberration introduced by the glass-water interface [38]. The objective also serves to collect the image of the sample. A CCD camera (DCC3240M, Thorlabs Inc., USA) with resolution of 1280×1024 pixels and maximal frame rate of 60 fps is used to record the manipulation process of the trapped particles.

The sample used for optical trapping in the experiment is 4 µm-diameter silica sphere suspended in water. A 7×7 focal spots array is generated in an arbitrary tilted plane (defined by the angles *θ*_{0} and *ϕ*_{0}), which can simultaneously trap multiple beads in the array and manipulate them in 3D space. The tilted-plane GS algorithm is used for the real-time refresh of the CGHs producing 7×7 spots array in different tilted (*θ*_{0}, *ϕ*_{0}) planes. To reduce the aberration, the generated phase pattern ${\varphi _H}$ is overlaid with a premeasured aberration correction pattern ${\varphi _{aberr}}$ [39], giving an effective phase $\varphi $ displayed on the SLM by

To access the actually realized uniformity, the 3D intensity profile of the trap arrays is measured by axially scanning a front-surface mirror by which the light is sent back into the objective and then onto the CCD camera. The 3D intensity profile is reconstructed from a stack of scanned two-dimensional images. We define the intensity of the *i*-th focus as the average of those of points inside an 880nm×880nm×1800nm cube centered at the intensity peak location [16]. By this definition, the experimental results are in good agreement with the theoretical estimations. For example, in a configuration of traps shown in Fig. 5(b), the uniformity is ∼94%. The deviation lies in several facts: (a) the objective cannot capture all the reflected light from the mirror. (b) non-uniform spatial profile of the input beam.

In the experiment, the total refreshing time for the SLM to update a CGH is ∼0.125s, including 0.1s calculation of tilted-plane GS algorithm and 0.025s of SLM displaying. This gives a maximum refresh rate of ∼8 frame per second. In the Visualization 1, we show three dynamic actions of the trap array: (i) the 7×7 trap array initially in the focal plane is rotated around *y*-axis; (ii) the 7×7 trap array in the (*θ*_{0}, *ϕ*_{0}) = (10°, 0°) plane is rotated around *z*-axis; (iii) the 7×7 trap array in the (*θ*_{0}, *ϕ*_{0}) = (10°, 45°) plane is rotated around its normal vector. In the whole process (∼ 12.5s) presented in Visualization 1, the tilted-plane GS algorithm generated 100 CGHs.

Figure 6 shows a sequence of the video frames extracted from the Visualization 1. For each frame, the left part shows the image of the trapped beads and the right part illustrates the foci array in the focal volume. Figure 6(a) shows the initial state of the trapped beads in the dynamic process (i), in which the 7×7 beads are trapped in the focal plane. Then the trap array is rotated around the *y*-axis, and all the beads are seen to stay in the array in the rotated plane as shown in Figs. 6(b) [(*θ*_{0}=17.5°, *ϕ*_{0}=0°)] and 6(c) [(*θ*_{0}=-17.5°, *ϕ*_{0}=0°)]. We see that the images of the beads take on 7 different levels of blurriness, implying that the beads are distributed in 7 different transverse planes off the focal plane. Figures 6(d) and 6(e) give two snapshots in the dynamic process (ii), where the same particle marked with black dotted circle show the same blurring level as the beads are rotated around *z*-aixs. In Figs. 6(f)–6(h), three snapshots show the dynamic process (iii), where the trap array is rotated around its normal (*θ*_{0}, *ϕ*_{0}) = (10°, 45°). It is seen that the beads are trapped in the same plane and rotated stably. The experimental results demonstrate that it is realizable to trap multiple particles in an arbitrary plane and dynamically rotate them in 3D space with our tilted-plane GS algorithm.

## 4. Conclusion

In summary, we have presented a tilted-plane GS algorithm targeted to 3D optical manipulations. With this algorithm, a large number of point traps can be generated in a tilted plane with a substantially high uniformity of 94%. Besides, the algorithm can effectively reduce the problem to a 2D one, thus enabling a dynamic 3D rotation of multiple trapped particles. We also realized the dynamic manipulation and rotation of 7×7 silica beads in 3D space. As an efficient method for arbitrary tilted-plane beam shaping, the algorithm can be further combined with other techniques like multifocal confocal microscopy, tilted light sheet microscopy, and parallel laser beam micro-fabrication to extend its applications.

## Funding

National Natural Science Foundation of China (11704405, 11974417, 61522511, 81427802); National Key Research and Development Program of China (2017YFC0110100).

## Disclosures

The authors have no relevant financial interests in this article and no potential conflicts of interest to disclose.

## References

**1. **M. Padgett and R. Di Leonardo, “Holographic optical tweezers and their relevance to lab on chip devices,” Lab Chip **11**(7), 1196–1205 (2011). [CrossRef]

**2. **Y. Liang, S. Yan, Z. Wang, R. Li, Y. Cai, M. He, B. Yao, and M. Lei, “Simultaneous optical trapping and imaging in the axial plane: a review of current progress,” Rep. Prog. Phys. **83**(3), 032401 (2020). [CrossRef]

**3. **D. B. Conkey, R. P. Trivedi, S. R. P. Pavani, I. I. Smalyukh, and R. Piestun, “Three-dimensional parallel particle manipulation and tracking by integrating holographic optical tweezers and engineered point spread functions,” Opt. Express **19**(5), 3835–3842 (2011). [CrossRef]

**4. **M. Righini, A. S. Zelenina, C. Girard, and R. Quidant, “Parallel and selective trapping in a patterned plasmonic landscape,” Nat. Phys. **3**(7), 477–480 (2007). [CrossRef]

**5. **R. L. Eriksen, V. R. Daria, and J. Gluckstad, “Fully dynamic multiple-beam optical tweezers,” Opt. Express **10**(14), 597–602 (2002). [CrossRef]

**6. **K. Sasaki, M. Koshioka, H. Misawa, N. Kitamura, and H. Masuharah, “Pattern-formation and ﬂow-control of ﬁne particles by laser-scanning micromanipulation,” Opt. Lett. **16**(19), 1463–1465 (1991). [CrossRef]

**7. **J. E. Molloy, J. E. Burns, J. C. Sparrow, R. T. Tregear, J. Kendrickjones, and D. C. S. White, “Single-moleculemechanics of heavy-meromyosin and s1 interacting with rabbit or drosophila actins using optical tweezers,” Biophys. J. **68**(4), 298S–303S (1995).

**8. **M. P. MacDonald, L. Paterson, K. Volke-Sepulveda, J. Arlt, W. Sibbett, and K. Dholakia, “Creation and manipulation of three-dimensional optically trapped structures,” Science **296**(5570), 1101–1103 (2002). [CrossRef]

**9. **M. Reicherter, T. Haist, E. U. Wagemann, and H. J. Tiziani, “Optical particle trapping with computer-generated holograms written on a liquid–crystal display,” Opt. Lett. **24**(9), 608–610 (1999). [CrossRef]

**10. **J. E. Curtis, B. A. Koss, and D. G. Grier, “Dynamic holographic optical tweezers,” Opt. Commun. **207**(1-6), 169–175 (2002). [CrossRef]

**11. **B. Jia, H. Lin, and M. Gu, “Dynamic generation of Debye diffraction-limited multifocal arrays for direct laser printing nanofabrication,” Opt. Lett. **36**(3), 406–408 (2011). [CrossRef]

**12. **P. T. Korda, M. B. Taylor, and D. G. Grier, “Kinetically locked-in colloidal transport in an array of optical tweezer,” Phys. Rev. Lett. **89**(12), 128301 (2002). [CrossRef]

**13. **K. Mangold, P. Leiderer, and C. Bechinger, “Phase transitions of colloidal monolayers in periodic pinning arrays,” Phys. Rev. Lett. **90**(15), 158302 (2003). [CrossRef]

**14. **M. Brunner and C. Bechinger, “Phase behavior of colloidal molecular crystals on triangular light lattices,” Phys. Rev. Lett. **88**(24), 248302 (2002). [CrossRef]

**15. **M. Montes-Usategui, E. Pleguezuelos, J. Andilla, and E. Martín-Badosa, “Fast generation of holographic optical tweezers by random mask encoding of Fourier components,” Opt. Express **14**(6), 2101–2107 (2006). [CrossRef]

**16. **R. Di Leonardo, F. Ianni, and G. Ruocco, “Computer generation of optimal holograms for optical trap arrays,” Opt. Express **15**(4), 1913–1922 (2007). [CrossRef]

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

**18. **H. Kim, M. Kim, W. Lee, and J. Ahn, “Gerchberg-Saxton algorithm for fast and efficient atom rearrangement in optical tweezer traps,” Opt. Express **27**(3), 2184–2196 (2019). [CrossRef]

**19. **G. Sinclair, J. Leach, P. Jordan, G. Gibson, and E. Yao, “Interactive application in holographic optical tweezers of a multi-plane Gerchberg-Saxton algorithm for three-dimensional light shaping,” Opt. Express **12**(8), 1665–1670 (2004). [CrossRef]

**20. **E. H. Waller and G. V. Freymann, “Multi foci with diffraction limited resolution,” Opt. Express **21**(18), 21708–21713 (2013). [CrossRef]

**21. **H. Ren, H. Lin, and M. Gu, “Three-dimensional parallel recording with a Debye diffraction-limited and aberration-free volumetric multifocal array,” Opt. Lett. **39**(6), 1621–1624 (2014). [CrossRef]

**22. **G. Whyte and J. Courtial, “Experimental demonstration of holographic three-dimensional light shaping using a Gerchberg-Saxton algorithm,” New J. Phys. **7**, 117 (2005). [CrossRef]

**23. **H. T. Chang, C. Y. Chen, J. S. Lin, and W. J. Li, “Image reconstruction by applying Fresnel transform on phase-only computer generated hologram at tilted planes,” Digital Holography and Three-Dimensional Imaging, Optical Society of America, DW2A, 11 (2015).

**24. **H. T. Chang, Y. T. Wang, and C. Y. Chen, “Angle multiplexing optical image encryption in the Fresnel transform domain using phase-only computer-generated hologram,” Photonics **7**(1), 1 (2019). [CrossRef]

**25. **K. Matsushima, H. Schimmel, and F. Wyrowski, “Fast calculation method for optical diffraction on tilted planes by use of the angular spectrum of plane waves,” J. Opt. Soc. Am. A **20**(9), 1755–1762 (2003). [CrossRef]

**26. **C. Chang, J. Xia, and Y. Jiang, “Holographic image projection on tilted planes by phase-only computer generated hologram using fractional Fourier transformation,” J. Disp. Technol. **10**(2), 107–113 (2014). [CrossRef]

**27. **M. Habaza, B. Gilboa, Y. Roichman, and N. T. Shaked, “Tomographic phase microscopy with 180° rotation of live cells in suspension by holographic optical tweezers,” Opt. Lett. **40**(8), 1881–1884 (2015). [CrossRef]

**28. **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]

**29. **M. Johansson and J. Bengtsson, “Robust design method for highly efficient beam-shaping diffractive optical elements using an iterative-Fourier-transform algorithm with soft operations,” J. Mod. Opt. **47**(8), 1385–1398 (2000). [CrossRef]

**30. **M. Leutenegger, R. Rao, R. A. Leitgeb, and T. Lasser, “Fast focus field calculations,” Opt. Express **14**(23), 11277–11291 (2006). [CrossRef]

**31. **B. R. Boruah and M. A. A. Neil, “Focal field computation of an arbitrarily polarized beam using fast Fourier transforms,” Opt. Commun. **282**(24), 4660–4667 (2009). [CrossRef]

**32. **Y. Cai, Z. Wang, Y. Liang, F. Ren, B. Yao, M. Lei, and S. Yan, “Direct calculation of tightly focused field in an arbitrary plane,” Opt. Commun. **450**, 329–334 (2019). [CrossRef]

**33. **S. P. Poland, N. Krstajić, R. D. Knight, R. K. Henderson, and S. M. Ameer-Beg, “Development of a doubly weighted Gerchberg-Saxton algorithm for use in multibeam imaging applications,” Opt. Lett. **39**(8), 2431–2434 (2014). [CrossRef]

**34. **X. Hao, C. Kuang, T. Wang, and X. Liu, “Effects of polarization on the de-excitation dark focal spot in STED microscopy,” J. Opt. **12**(11), 115707 (2010). [CrossRef]

**35. **H. Chen, Y. Guo, Z. Chen, J. Hao, J. Xu, H. Wang, and J. Ding, “Holographic optical tweezers obtained by using the three-dimensional Gerchberg-Saxton algorithm,” J. Opt. **15**(3), 035401 (2013). [CrossRef]

**36. **S. Bianchi and R. Di Leonardo, “Real-time optical micro-manipulation using optimized holograms generated on the GPU,” Comput. Phys. Commun. **181**(8), 1444–1448 (2010). [CrossRef]

**37. **J. Schmidhuber, “A fixed size storage O (n3) time complexity learning algorithm for fully recurrent continually running networks,” Neural Comput. **4**(2), 243–248 (1992). [CrossRef]

**38. **S. N. S. Reihani, M. A. Charsooghi, and H. R. Khalesifard, “Efficient in-depth trapping with an oil-immersion objective lens,” Opt. Lett. **31**(6), 766–768 (2006). [CrossRef]

**39. **Y. Liang, Y. Cai, Z. Wang, M. Lei, Z. Cao, Y. Wang, M. Li, S. Yan, P. R. Bianco, and Y. Bao, “Aberration correction in holographic optical tweezers using a high-order optical vortex,” Appl. Opt. **57**(13), 3618–3623 (2018). [CrossRef]