Optimized holographic optical traps

Marco Polin

Dept. of Physics and Center for Soft Matter Research, New York University, New York, NY 10003

Kosta Ladavac

Department of Physics, James Franck Institute and Institute for Biophysical Dynamics, The University of Chicago, Chicago, IL 60637

Sang-Hyuk Lee

Yael Roichman

David G. Grier

Dept. of Physics and Center for Soft Matter Research, New York University, New York, NY 10003

Abstract.

Holographic optical traps use the forces exerted by computer-generated holograms to trap, move and otherwise transform mesoscopically textured materials. This article introduces methods for optimizing holographic optical traps' efficiency and accuracy, and an optimal statistical approach for characterizing their performance. This combination makes possible real-time adaptive optimization.

A single laser beam brought to a focus with a strongly converging lens forms a type of optical trap widely known as an optical tweezer (1). Multiple beams of light passing simultaneously through the lens' input pupil focus to multiple optical tweezers, each at a location determined by the associated beam's angle of incidence and degree of collimation as it enters the lens. Their intersection at the input pupil yields an interference pattern whose amplitude and phase corrugations characterize the downstream trapping pattern. Imposing the same modulations on a single incident beam at the input pupil would yield the same pattern of traps. Such wavefront modification can be performed by a computer-designed diffractive optical element (DOE), or hologram.

Holographic optical trapping (HOT) uses computer-generated holograms (CGHs) to project arbitrary configurations of optical traps (2); (3); (4); (5); (6), and so provides exceptional control over microscopic materials dispersed in fluid media. Holographic micromanipulation provides the basis for a rapidly growing field of applications in the physical and biological sciences as well as in industry (7).

This article describes refinements to the HOT technique that help to optimize the traps' performance. It also introduces self-consistent and statistically optimal methods for characterizing their performance. Section I describes modifications to the basic HOT optical train that compensate for practical limitations of dynamic holography. Section II discusses a direct search algorithm for HOT CGH computation that is both faster and more accurate than commonly used iterative refinement algorithms. Together, these modifications yield marked improvements in the holographic traps' performance that can be quantified rapidly using techniques introduced in Section III. These techniques are based on optimal statistical analysis of trapped colloidal spheres' thermally-driven motions, and lend themselves to simultaneous real-time characterization and optimization of entire arrays of traps through digital video microscopy. Such adaptive optimization is demonstrated experimentally in Section IV.

§ I. Improved optical train

Figure 1. Simplified schematic of a holographic optical tweezer optical train before and after modification. (a) A collimated beam is split into multiple beams by the DOE, each of which is shown here as being collimated. The diffracted beams pass through the input pupil of an objective lens and are focused into optical traps in the objective's focal plane. The undiffracted portion of the beam, shown here with the darkest shading, also focuses into the focal plane. (b) The input beam is converging as it passes through the DOE. The DOE collimates the diffracted beams, so that they focus into the focal plane, as in (a). The undiffracted beam comes to a focus within the coverslip bounding the sample. (c) A beam block can eliminate the undiffracted beam without substantially degrading the optical traps.

Figure 1(a) depicts a conventional HOT implementation in which a collimated laser beam is modified by a computer-designed DOE, and thereafter propagates as a superposition of independent beams, each with individually specified wavefront characteristics (4); (6). These beams are relayed to the input pupil of a high-numerical-aperture lens, typically a microscope objective, which focuses them into optical traps. Although a transmissive DOE is shown in Fig. 1, comparable results are obtained with reflective DOEs. The same objective lens used to form the optical traps also can be used to create images of trapped objects. The associated illumination and image-forming optics are omitted from Fig. 1 for clarity.

Practical DOEs only diffract a portion of the incident light into the intended modes and directions. Some of the incident beam may not be diffracted at all, and the undiffracted portion typically forms an unwanted trap in the middle of the field of view (8). This “central spot” has been removed in previous implementations by spatially filtering the diffracted beam (8); (9). Practical DOEs also tend to project spurious “ghost” traps into symmetry-dictated positions within the sample. Spatially filtering a large number of ghost traps generally is not practical, particularly in the case of dynamic holographic optical tweezers whose traps move freely in three dimensions. Projecting holographic traps in the off-axis Fresnel geometry automatically eliminates the central spot (10), but limits the number of traps that can be projected, and also does not address the formation of ghost traps.

Figure 1(b) shows a basic improvement that minimizes the central spot's influence and effectively eliminates ghost traps. Rather than illuminating the DOE with a collimated laser beam, a converging beam is used. This moves the undiffracted central spot upstream of the objective's normal focal plane. The intended traps can be moved back to the focal plane by incorporating wavefront-shaping phase functions into the hologram's design (6). Deliberately decollimating the input beam allows the central spot to be projected outside of the sample volume, thereby ensuring that the undiffracted beam lacks both the intensity and the intensity gradients needed to influence the sample's dynamics.

An additional consequence of the traps' displacement relative to the modified optical train's focal plane is that most ghost traps also are projected out of the sample volume. This is a substantial improvement for processes such as optical fractionation (11); (12) and optical ratchets (13); (14), which require defect-free intensity distributions.

Even though the undiffracted beam may not create an actual trap in this modified optical train, it still can exert radiation pressure on parts of the sample near the center of the field of view. This is a particular problem for large arrays of optical traps in that the central spot, which typically receives a fixed proportion of the input beam, can be brighter than the intended traps. Illuminating the DOE with a diverging beam (15) reduces the undiffracted beam's influence by projecting some of its light out of the optical train. In a thick sample, however, this has the deleterious effect of projecting both the weakened central spot and the undiminished ghost traps into the sample.

These problems all can be mitigated by placing a beam block as shown in Fig. 1(c) in the intermediate focal plane within the relay optics to spatially filter the undiffracted portion of the beam. The trap-forming beams focus downstream of the beam block and therefore are only partially occluded, even if they pass directly along the optical axis. This has little effect on the performance of conventional optical tweezers and can be compensated by increasing the occluded traps' relative brightness.

§ II. Algorithms for HOT CGH Calculation

Holographic optical tweezers' efficacy is determined by the quality of the trap-forming DOE, which in turn reflects the performance of the algorithms used in their computation. Previous studies have applied holograms calculated by simple linear superposition of the input fields (3), with best results being obtained with random relative phases (4); (6), or with variations (4); (5); (6) on the classic Gerchberg-Saxton and adaptive-additive algorithms (16). Despite their generality, these algorithms yield traps whose relative intensities can differ greatly from their design values, and often project an unacceptably large fraction of the input power into ghost traps. These problems can become acute for complicated three-dimensional trapping patterns, particularly when the same hologram also is used as a mode converter to project multifunctional arrays of optical traps (4); (6). This section describes faster and more effective algorithms for HOT DOE calculation based on direct search optimization.

The holograms used for holographic optical trapping typically operate only on the phase of the incident beam, and not its amplitude. Such phase-only holograms, also known as kinoforms, are far more efficient than amplitude-modulating holograms, which necessarily divert light away from the traps. Quite general trapping patterns can be achieved with kinoforms because optical tweezers rely for their operation on intensity gradients and not on phase variations. The challenge is to find a phase pattern in the input plane of the objective lens that encodes the desired intensity pattern in the focal volume.

Most approaches to designing phase-only holograms are based on scalar diffraction theory, in which the complex field $E(\boldsymbol{r})$ , in the focal plane of a lens of focal length $f$ is related to the field, $u(\boldsymbol{\rho})\,\exp(i\varphi(\boldsymbol{\rho}))$ , in its input plane by a Fraunhofer transform,

$E(\boldsymbol{r})=\int u(\boldsymbol{\rho})\,\exp(i\varphi(\boldsymbol{\rho}))\,\exp\left(-i\frac{k\boldsymbol{r}\cdot\boldsymbol{\rho}}{f}\right)\, d^{2}\rho.$

(1)

Here, $u(\boldsymbol{\rho})$ and $\varphi(\boldsymbol{\rho})$ are the real-valued amplitude and phase, respectively, at position $\boldsymbol{\rho}$ in the input pupil, and $k=2\pi/\lambda$ is the wavenumber of light of wavelength $\lambda$ . If $u(\boldsymbol{\rho})$ is taken to be the amplitude profile of the input laser beam, then $\varphi(\boldsymbol{\rho})$ is the kinoform encoding the intensity distribution $I(\boldsymbol{r})=\left|E(\boldsymbol{r})\right|^{2}$ . Finding the kinoform $\varphi(\boldsymbol{\rho})$ to project a particular pattern, $I(\boldsymbol{r})$ , is nontrivial because the inherent nonlinearity of Eq. (1) defies straightforward inversion. Even so, kinoforms may be estimated through indirect search algorithms. The particular requirements of holographic trapping lend themselves to especially fast and effective computation.

Simply computing the superposition of input beams required to create a desired trapping pattern, disregarding the resulting amplitude variations, and retaining the phase as an estimate for $\varphi(\boldsymbol{\rho})$ turns out to be remarkably effective (3); (4); (6). Fidelity to design is particularly good if the input beams are assumed to have random relative phases. Not surprisingly, such randomly phased superpositions yield traps with widely varying intensities as well as a great many ghost traps.

The process of refining such an initial estimate begins by noting that most practical DOEs, including those projected with SLMs, consist of an array of discrete phase pixels, $\varphi _{j}=\varphi(\boldsymbol{\rho}_{j})$ , centered at locations $\boldsymbol{\rho}_{j}$ , each of which can impose any of $P$ possible discrete phase shifts on the incident beam. The field in the focal plane due to such an $N$ -pixel DOE is, therefore,

$E(\boldsymbol{r})=\sum _{{j=1}}^{N}u_{j}\,\exp(i\varphi _{j})\, T_{j}(\boldsymbol{r}),$

(2)

where the transfer function describing the light's propagation from input plane to output plane is

$T_{j}(\boldsymbol{r})=\exp\left(-i\frac{k\boldsymbol{r}\cdot\boldsymbol{\rho}_{j}}{f}\right).$

(3)

Unlike more general holograms, the desired field in the output plane of a holographic optical trapping system consists of $M$ discrete bright spots located at $\boldsymbol{r}_{m}$ :

	$\displaystyle E(\boldsymbol{r})$	$\displaystyle=\sum _{{m=1}}^{M}E_{m}\,\delta(\boldsymbol{r}-\boldsymbol{r}_{m}),\quad\text{with}$		(4)
	$\displaystyle E_{m}$	$\displaystyle=\alpha _{m}\,\exp(i\xi _{m}),$		(5)

where $\alpha _{m}$ is the relative amplitude of the $m$ -th trap, normalized by $\sum _{{m=1}}^{M}\left|\alpha _{m}\right|^{2}=1$ , and $\xi _{m}$ is its (arbitrary) phase. Here $\delta(\boldsymbol{r})$ represents the amplitude profile of the focused beam of light, and may be approximated by a two-dimensional Dirac delta function. For simplicity, we may also approximate the input beam's amplitude profile by a top-hat function with $u_{j}=1$ within the input pupil's aperture, and $u_{j}=0$ elsewhere. In these approximations, the field at the $m$ -th trap is (6)

$E_{m}=\sum _{{j=1}}^{N}{\mathbf{K}}_{{j,m}}^{{-1}}\,{\mathbf{T}}_{{j,m}}\,\exp(i\varphi _{j}),$

(6)

with ${\mathbf{T}}_{{j,m}}=T_{j}(\boldsymbol{r}_{m})$ . We introduce the inverse operator ${\mathbf{K}}_{{j,m}}^{{-1}}$ because the hologram $\varphi _{j}$ may modify the wavefronts of each of the diffracted beams it creates in addition to establishing its direction of propagation. Such wavefront distortions are useful for creating three-dimensional arrays of multifunctional traps. However, they also distort the traps' otherwise sharply peaked profiles in the focal plane, which were assumed in Eq. (4). The inverse operators correct for these distortions so that even generalized traps can be treated discretely.

For example, a trap can be displaced a distance $z$ away from the focal plane by curving the input beam's wavefronts into a parabolic profile

$\varphi _{z}(\boldsymbol{\rho},z)=\frac{k\rho^{2}z}{2f^{2}}.$

(7)

The operator that displaces the $m$ -th trap to $z_{m}$ is (4); (6)

${\mathbf{K}}^{z}_{{j,m}}=\exp(i\varphi _{z}(\boldsymbol{\rho}_{j},z_{m})).$

(8)

Its inverse, ${{\mathbf{K}}^{z}_{{j,m}}}^{{-1}}={{\mathbf{K}}^{z}_{{j,m}}}^{\ast}$ returns the $m$ -th trap to best focus in the focal plane.

Similarly, a conventional TEM ${}_{{00}}$ beam can be converted into a helical mode with the phase profile

$\varphi _{\ell}(\boldsymbol{\rho},\ell)=\ell\theta,$

(9)

where $\theta$ is the azimuthal angle around the optical axis and $\ell$ is a winding number known as the topological charge. Such corkscrew-like beams focus to ring-like optical traps known as optical vortices, which can exert torques as well as forces (17); (18); (19); (20). The topology-transforming kernel (6)

${\mathbf{K}}_{{j,m}}=\exp(i\varphi _{\ell}(\boldsymbol{\rho}_{j},\ell _{m}))$

(10)

can be composed with ${\mathbf{T}}_{{j,m}}$ in the same manner as the displacement-inducing ${\mathbf{K}}^{z}_{{j,m}}$ to convert the $m$ -th trap into an optical vortex. A variety of analogous phase-based mode transformations have been described, each with applications to single-beam optical trapping (7), all of which can be applied to each trap independently in this manner.

Calculating the fields only at the traps' positions greatly reduces the computational burden of HOT CGH refinement. It also eliminates the need to account for the beams' propagation through intermediate trapping planes when designing three-dimensional patterns (4). Unlike more general FFT-based algorithms (5), this restricted approach does not directly optimize the field between the traps. If the converged amplitudes match the design values, however, no light is left over to create ghost traps.

Applying Eq. (6) directly in an iterative refinement algorithm (6) also has drawbacks. In particular, only the $M-1$ relative phases $\xi _{m}$ in Eq. (5) can be adjusted when inverting Eq. (6) to solve for $\varphi _{j}$ . Having so few free parameters severely limits the improvement over simple superposition that can be obtained. Equation (6) suggests an alternative approach that not only is far more effective, but also is substantially more efficient.

The operator ${\mathbf{K}}_{{j,m}}^{{-1}}\,{\mathbf{T}}_{{j,m}}$ describes how light in the mode of the $m$ -th trap propagates from the $j$ -th phase pixel on the DOE to the trap's projected position $\boldsymbol{r}_{m}$ . Changing the pixel's value by $\Delta\varphi _{j}$ therefore changes each trap's field by

$\Delta E_{m}={\mathbf{K}}_{{j,m}}^{{-1}}\,{\mathbf{T}}_{{j,m}}\,\exp(i\varphi _{j})\,\left[\exp(i\,\Delta\varphi _{j})-1\right].$

(11)

If such a change were to improve the overall pattern, we would be inclined to retain it, and to seek other such improvements. This is the basis for direct search algorithms. The simplest involves selecting a pixel at random from a trial phase pattern, changing its value to any of the $P-1$ alternatives, and computing the effect on the projected pattern. Quite clearly, there is a considerable computational advantage in calculating changes only at the $M$ traps' positions, rather than over the entire focal plane. The updated trial amplitudes then are compared with their design values and the proposed change is accepted if the overall error is reduced. The process is repeated until the result converges to the design or the acceptance rate for proposed changes dwindles.

Effective and efficient refinement by the direct search algorithm depends on the choice of metric for quantifying convergence. The standard cost function, $\chi^{2}=\sum _{{m=1}}^{M}(I_{m}-\epsilon I^{{(D)}}_{m})^{2}$ , assesses the mean-squared deviations of the $m$ -th trap's projected intensity $I_{m}=\left|\alpha _{m}\right|^{2}$ from its design value $I^{{(D)}}_{m}$ , assuming an overall diffraction efficiency of $\epsilon$ . It requires an accurate estimate for $\epsilon$ and places no emphasis on uniformity in the projected traps' intensities. An alternative proposed by Meister and Winfield (21),

$C=-\left<I\right>+f\sigma,$

(12)

avoids both shortcomings. Here, $\left<I\right>$ is the mean intensity at the traps and

$\sigma=\sqrt{\frac{1}{M}\,\sum _{{m=1}}^{M}(I_{m}-\gamma I^{{(D)}}_{m})^{2}}$

(13)

measures the deviation from uniform convergence to the design intensities. Selecting

$\gamma=\frac{\sum _{{m=1}}^{M}I_{m}\, I^{{(D)}}_{m}}{\sum _{{m=1}}^{M}\left(I^{{(D)}}_{m}\right)^{2}}$

(14)

minimizes the total error and optimally accounts for non-ideal diffraction efficiency (21). The weighting fraction $f$ sets the relative importance attached to diffraction efficiency versus uniformity, with $f=0.5$ providing a generally useful balance.

A direct binary search proceeds with any candidate change that reduces $C$ being accepted, and all others being rejected. In a worst-case implementation, the number of trials required for convergence should scale as $NP$ , the product of the number of phase pixels and the number of possible phase values. In practice, this estimate is accurate if $P$ and $N$ are comparatively small and if the starting phase function is either uniform or purely random. Much faster convergence can be obtained by starting from the a randomly phased superposition of input beams. In this case, convergence typically is obtained within $N$ trials, even for fairly complex trapping patterns, and thus requires a computational effort comparable to the initial superposition.

Figure 2. (a) Design for 119 identical optical traps in a two-dimensional quasiperiodic array. (b) Trapping pattern projected without optimizations using the adaptive-additive algorithm. (c) Trapping pattern projected with optimized optics and adaptively corrected direct search algorithm. (d) Bright-field image of colloidal silica spheres 1.53 $\mathrm{\upmu}$ $\mathrm{m}$ in diameter dispersed in water and organized in the optical trap array. The scale bar indicates 10 $\mathrm{\upmu}$ $\mathrm{m}$

As a practical demonstration, we have implemented a quasiperiodic array of optical traps, which is challenging because it has no translational symmetries. The traps are focused with a $100\times$ NA 1.4 S-Plan Apo oil immersion objective lens mounted in a Nikon TE-2000U inverted optical microscope. The traps are powered by a Coherent Verdi frequency-doubled diode-pumped solid state laser operating at a wavelength of 532 $\mathrm{n}$ $\mathrm{m}$ . Computer-generated phase holograms are imprinted on the beam with a Hamamatsu X8267-16 parallel-aligned nematic liquid crystal spatial light modulator (SLM). This SLM can impose phase shifts up to $2\pi~\mathrm{rad}$ at each pixel in a $768\times 768$ array. The face of the SLM is imaged onto the objective's 5 $\mathrm{m}$ $\mathrm{m}$ diameter input pupil using relay optics designed to minimize aberrations. The beam is directed into the objective with a dichroic beamsplitter, which allows images to pass through to a low-noise charge-coupled device (CCD) camera (NEC TI-324AII). The video stream is recorded as uncompressed digital video with a Pioneer 520H digital video recorder (DVR) for processing.

Figure 2(a) shows the intended planar arrangement of 119 holographic optical traps designed by the dual generalized method for generating quasiperiodic lattices (22). Even after adaptive-additive refinement, the hologram resulting from simple superposition with random phases fares poorly for this aperiodic pattern. Figure 2(b) shows the intensity of light reflected by a front-surface mirror placed in the sample plane. This image reveals extraneous ghost traps, an exceptionally bright central spot, and large variability in the intended traps' intensities. Imaging photometry on this and equivalent images produced with different random relative phases for the beams yields a typical root-mean-square (RMS) variation of more than 50 percent in the projected traps' brightness. The image in Fig. 2(c) was produced using the modified optical train described in Sec. I and the direct search algorithm described in Sec. II, and suffers from none of these defects. Both the ghost traps and the central spot are suppressed, and the apparent relative brightness variations are smaller than 5 percent, a factor of ten improvement. Figure 2(d) shows 119 colloidal silica spheres, $2a=1.5\pm 0.3~\mathrm{\upmu}\mathrm{m}$ in diameter (Bangs Labs, lot 5238), dispersed in water at $T=27^{^{{\circ}}}\mathrm{C}$ and trapped in the quasiperiodic array.

Figure 3. A three-dimensional multifunctional holographic optical trap array created with the direct search algorithm. (a) Refined DOE phase pattern. (b), (c) and (d) The projected optical trap array at $z=-10~\mathrm{\upmu}\mathrm{m}$ , $0~\mathrm{\upmu}\mathrm{m}$ and $+10~\mathrm{\upmu}\mathrm{m}$ . Traps are spaced by $1.2~\mathrm{\upmu}\mathrm{m}$ in the plane, and the 12 traps in the middle plane consist of $\ell=8$ optical vortices. (e) Performance metrics for the hologram in (a) as a function of the number of accepted single-pixel changes. Data include the DOE's overall diffraction efficiency as defined by Eq. (15), the projected pattern's RMS error from Eq. (16), and its uniformity, $1-u$ , where $u$ is defined in Eq. (17).

To place the benefits of the direct search algorithm on a more quantitative basis, we augment standard figures of merit with those introduced in Ref. (21). In particular, the DOE's theoretical diffraction efficiency is commonly defined as

${\cal Q}=\frac{1}{M}\,\sum _{{m=1}}^{M}\frac{I_{m}}{I^{{(D)}}_{m}},$

(15)

and its root-mean-square (RMS) error as

$e_{\text{rms}}=\frac{\sigma}{\max(I_{m})}.$

(16)

The resulting pattern's departure from uniformity is usefully gauged as (21)

$u=\frac{\max(I_{m}/I^{{(D)}}_{m})-\min(I_{m}/I^{{(D)}}_{m})}{\max(I_{m}/I^{{(D)}}_{m})+\min(I_{m}/I^{{(D)}}_{m})}.$

(17)

Figure 3 shows results for a HOT DOE encoding 51 traps, including 12 optical vortices of topological charge $\ell=8$ , arrayed in three planes relative to the focal plane. The excellent results in Fig. 3 were obtained with a single pass of direct-search refinement. The resulting traps, shown in the bottom three images, again vary from their planned relative intensities by less than 5 percent. In this case, the spatially extended vortices were made as bright as the point-like optical tweezers by increasing their requested relative brightness by a factor of 15. This single hologram, therefore, demonstrates independent control over three-dimensional position, wavefront topology, and brightness of all the traps. Performance metrics for the calculation are plotted in Fig. 3(b) as a function of the number of accepted single-pixel changes, with an overall acceptance rate of 16 percent. Direct search refinement achieves greatly improved fidelity to design over randomly phase superposition at the cost of a small fraction of the diffraction efficiency and roughly doubled computation time. The entire calculation can be completed in the refresh interval of a typical liquid crystal spatial light modulator.

§ III. Optimal characterization

Gauging a HOT system's performance numerically and by characterizing the projected intensity pattern does not provide a complete picture. The real test is in the projected traps' ability to localize particles. A variety of approaches have been developed for measuring the forces exerted by optical traps. The earliest involved estimating the hydrodynamic drag required to dislodge a trapped particle (23). This has several disadvantages, most notably that it identifies only the marginal escape force in a given direction and not the trap's actual three-dimensional potential. Complementary information can be obtained by measuring a particle's thermally driven motions in the trap's potential well (24); (25); (26). For instance, the measured probability density $P(\boldsymbol{r})$ for displacements $\boldsymbol{r}$ is related to the trap's potential $V(\boldsymbol{r})$ through the Boltzmann distribution

$P(\boldsymbol{r})\propto\exp(-\beta V(\boldsymbol{r})),$

(18)

where $\beta^{{-1}}=k_{B}T$ is the thermal energy scale at temperature $T$ . Similarly, the power spectrum of $\boldsymbol{r}(t)$ for a harmonically bound particle is a Lorentzian whose width is the viscous relaxation time of the particle in the well (24); (27).

Both of these approaches require amassing enough data to characterize the trapped particle's least probable displacements, and therefore oversample the trajectories. Oversampling is acceptable when data from a single optical trap can be collected rapidly, for example with a quadrant photodiode (24); (25); (26); (28). Tracking multiple particles in holographic optical traps, however, requires the area detection capabilities of digital video microscopy (29), which yields data much more slowly. Analyzing video data with optimal statistics (30) offers the benefits of thermal calibration by avoiding oversampling.

An optical trap is accurately modeled as a harmonic potential energy well (25); (26); (27); (28),

$V(\boldsymbol{r})=\frac{1}{2}\,\sum _{{i=1}}^{3}\kappa _{i}r_{i}^{2},$

(19)

with a different characteristic curvature $\kappa _{i}$ along each axis. This form also is convenient because it is separable into one-dimensional contributions. The trajectory of a colloidal particle localized in a viscous fluid by a harmonic well is described by the one-dimensional Langevin equation (31)

$\dot{x}(t)=-\frac{x(t)}{\tau}+\xi(t),$

(20)

where the autocorrelation time $\tau=\gamma/\kappa$ , is set by the particle's viscous drag coefficient $\gamma$ and by the curvature of the well, $\kappa$ . The Gaussian random thermal force, $\xi(t)$ , has zero mean, $\left<\xi(t)\right>=0$ , and variance

$\left<\xi(t)\,\xi(s)\right>=\frac{2k_{B}T}{\gamma}\,\delta(t-s).$

(21)

If the particle is at position $x_{0}$ at time $t=0$ , its trajectory at later times is given by

$x(t)=x_{0}\,\exp\left(-\frac{t}{\tau}\right)+\int _{0}^{t}\xi(s)\,\exp\left(-\frac{t-s}{\tau}\right)\, ds.$

(22)

Sampling such a trajectory at discrete times $t_{j}=j\,\Delta t$ , yields

$x_{{j+1}}=\phi\, x_{j}+a_{{j+1}},$

(23)

where $x_{j}=x(t_{j})$ ,

$\phi=\exp\left(-\frac{\Delta t}{\tau}\right),$

(24)

and where $a_{{j+1}}$ is a Gaussian random variable with zero mean and variance

${\sigma _{a}^{2}}=\frac{k_{B}T}{\kappa}\,\left[1-\exp\left(-\frac{2\Delta t}{\tau}\right)\right].$

(25)

Because $\phi<1$ , Eq. (23) is an example of an autoregressive process (30), which is readily invertible.

In principle, the particle's trajectory $\{ x_{j}\}$ can be analyzed to extract $\phi$ and ${\sigma _{a}^{2}}$ , and thus the trap's stiffness, $\kappa$ , and the particle's viscous drag coefficient $\gamma$ . In practice, however, the experimentally measured particle positions $y_{j}$ differ from the actual positions $x_{j}$ by random errors $b_{j}$ , which we assume to be taken from a Gaussian distribution with zero mean and variance $\sigma _{b}^{2}$ . The measurement then is described by the coupled equations

$x_{j}=\phi\, x_{{j-1}}+a_{j}\text{\quad and \quad}y_{j}=x_{j}+b_{j},$

(26)

where $b_{j}$ is independent of $a_{j}$ . We still can estimate $\phi$ and ${\sigma _{a}^{2}}$ from a set of measurements $\{ y_{j}\}$ by first constructing the joint probability

$p(\{ x_{i}\},\{ y_{i}\}|\phi,{\sigma _{a}^{2}},{\sigma _{b}^{2}})=\\ \prod _{{j=2}}^{N}\left[\frac{\exp\left(-\frac{a_{j}^{2}}{2{\sigma _{a}^{2}}}\right)}{\sqrt{2\pi{\sigma _{a}^{2}}}}\right]\prod _{{j=1}}^{N}\left[\frac{\exp\left(-\frac{b_{j}^{2}}{2{\sigma _{b}^{2}}}\right)}{\sqrt{2\pi{\sigma _{b}^{2}}}}\right].$

(27)

The probability density for measuring the trajectory $\{ y_{j}\}$ , is then the marginal (30)

	$\displaystyle p(\{ y_{j}\}\|\phi,{\sigma _{a}^{2}},{\sigma _{b}^{2}})=\int p(\{ x_{j}\},\{ y_{j}\}\|\phi,{\sigma _{a}^{2}},{\sigma _{b}^{2}})\, dx_{1}\cdots dx_{N}$		(28)
	$\displaystyle=\frac{(2\pi{\sigma _{a}^{2}}{\sigma _{b}^{2}})^{{-\frac{N-1}{2}}}}{\sqrt{{\sigma _{b}^{2}}\,\det({\mathbf{A}}_{\phi})}}\exp\left(-\frac{1}{2{\sigma _{b}^{2}}}(\boldsymbol{y})^{T}\left[{\mathbf{I}}-\frac{{\mathbf{A}}_{\phi}^{{-1}}}{{\sigma _{b}^{2}}}\right]\,\boldsymbol{y}\right),$		(29)

where $\boldsymbol{y}=(y_{1},\dots,y_{N})$ with transpose $(\boldsymbol{y})^{T}$ , ${\mathbf{I}}$ is the $N\times N$ identity matrix, and

${\mathbf{A}}_{\phi}=\frac{{\mathbf{I}}}{{\sigma _{b}^{2}}}+\frac{{\mathbf{M}}_{\phi}}{{\sigma _{a}^{2}}},$

(30)

with the tridiagonal memory tensor

${\mathbf{M}}_{\phi}=\begin{pmatrix}\phi^{2}&-\phi&0&0&\cdots&0\\ -\phi&1+\phi^{2}&-\phi&0&\cdots&\vdots\\ 0&-\phi&1+\phi^{2}&-\phi&\cdots&\vdots\\ 0&0&-\phi&\ddots&\cdots&\vdots\\ \vdots&\vdots&\cdots&-\phi&1+\phi^{2}&-\phi\\ 0&0&\cdots&0&-\phi&1\end{pmatrix}.$

(31)

Calculating the determinant, $\det({\mathbf{A}}_{\phi})$ , and inverse, ${\mathbf{A}}_{\phi}^{{-1}}$ , of ${\mathbf{A}}_{\phi}$ is greatly facilitated if we artificially impose time translation invariance by replacing $M_{\phi}$ with the $(N+1)\times(N+1)$ matrix that identifies time step $N+1$ with time step 1. Physically, this involves imparting an impulse, $a_{{N+1}}$ , that translates the particle from its last position, $x_{N}$ , to its first, $x_{1}$ . Because diffusion in a potential well is a stationary process, the effect of this change is inversely proportional to the number of measurements, $N$ . With this approximation,

$\displaystyle\det({\mathbf{A}}_{\phi})$	$\displaystyle=\prod _{{n=1}}^{{N}}\left\{\frac{1}{{\sigma _{b}^{2}}}+\frac{1}{{\sigma _{a}^{2}}}\,\left[1+\phi^{2}-2\phi\,\cos\left(\frac{2\pi n}{N}\right)\right]\right\}$	(32)
and
$\displaystyle({\mathbf{A}}_{\phi}^{{-1}})_{{\alpha\beta}}$	$\displaystyle=\frac{1}{N}\,\sum _{{n=1}}^{N}\frac{{\sigma _{a}^{2}}{\sigma _{b}^{2}}\,\exp\left(i\frac{2\pi}{N}n(\alpha-\beta)\right)}{{\sigma _{a}^{2}}+{\sigma _{b}^{2}}\,\left[1+\phi^{2}-2\phi\,\cos\left(\frac{2\pi n}{N}\right)\right]},$	(33)

so that the conditional probability for the measured trajectory, $\{ y_{j}\}$ , is

$p(\{ y_{j}\}|\phi,{\sigma _{a}^{2}},{\sigma _{b}^{2}})=(2\pi)^{{-\frac{N}{2}}}\\ \times\prod _{{n=1}}^{{N}}\Big\{{\sigma _{a}^{2}}+{\sigma _{b}^{2}}\,\left[1+\phi^{2}-2\phi\,\cos\left(\frac{2\pi n}{N}\right)\right]\Big\}^{{-\frac{1}{2}}}\\ \times\exp\left(-\frac{1}{2{\sigma _{b}^{2}}}\,\sum _{{n=1}}^{N}y_{n}^{2}\right)\\ \times\exp\left(\frac{1}{2{\sigma _{b}^{2}}}\,\frac{1}{N}\,\sum _{{m=1}}^{N}\frac{\tilde{y}_{m}^{2}\,{\sigma _{a}^{2}}}{{\sigma _{a}^{2}}+{\sigma _{b}^{2}}\,\left[1+\phi^{2}-2\phi\,\cos\left(\frac{2\pi m}{N}\right)\right]}\right),$

(34)

where $\tilde{y}_{m}$ is the $m$ -th component of the discrete Fourier transform of $\{ y_{n}\}$ . This can be inverted to obtain the likelihood function for $\phi$ , ${\sigma _{a}^{2}}$ , and ${\sigma _{b}^{2}}$ :

$L(\phi,{\sigma _{a}^{2}},{\sigma _{b}^{2}}|\{ y_{i}\})=-\frac{N}{2}\,\ln 2\pi\\ -\frac{1}{2{\sigma _{b}^{2}}}\,\sum _{{n=1}}^{N}y_{n}^{2}+\frac{{\sigma _{a}^{2}}}{2{\sigma _{b}^{2}}}\,\frac{1}{N}\,\sum _{{n=1}}^{{N}}\frac{\tilde{y}_{n}^{2}\,{\sigma _{a}^{2}}}{{\sigma _{a}^{2}}+{\sigma _{b}^{2}}\,\left[1+\phi^{2}-2\phi\,\cos\left(\frac{2\pi n}{N}\right)\right]}\\ -\frac{1}{2}\,\sum _{{n=1}}^{N}\ln\left({\sigma _{a}^{2}}+{\sigma _{b}^{2}}\left[1+\phi^{2}-2\phi\,\cos\left(\frac{2\pi n}{N}\right)\right]\right).$

(35)

Best estimates $(\hat{\phi},\hat{{\sigma _{a}^{2}}},\hat{{\sigma _{b}^{2}}})$ for the parameters $(\phi,{\sigma _{a}^{2}},{\sigma _{b}^{2}})$ are solutions of the coupled equations

$\frac{\partial L}{\partial\phi}=\frac{\partial L}{\partial{\sigma _{a}^{2}}}=\frac{\partial L}{\partial{\sigma _{b}^{2}}}=0.$

(36)

§ III.1. Case 1: No measurement errors ( ${\sigma _{b}^{2}}=0$ )

Equations (36) can be solved in closed form if ${\sigma _{b}^{2}}=0$ . In this case,

$\hat{\phi}_{0}=\frac{c_{1}}{c_{0}}\text{\quad and \quad}\hat{{\sigma _{a}^{2}}}_{0}=c_{0}\,\left[1-\left(\frac{c_{1}}{c_{0}}\right)^{2}\right],$

(37)

where

$c_{m}=\frac{1}{N}\sum _{{j=1}}^{{N}}y_{j}\, y_{{(j+m)\bmod N}}$

(38)

is the barrel autocorrelation of $\{ y_{j}\}$ at lag $m$ . The associated statistical uncertainties are

$\Delta\hat{\phi}_{0}=\sqrt{\frac{\hat{{\sigma _{a}^{2}}}_{0}}{Nc_{0}}}\text{, \quad and \quad}\Delta\hat{{\sigma _{a}^{2}}}_{0}=\hat{{\sigma _{a}^{2}}}_{0}\,\sqrt{\frac{2}{N}}.$

(39)

In the absence of measurement errors, $c_{0}$ and $c_{1}$ constitute sufficient statistics for the time series (30) and thus embody all of the information that can be extracted.

§ III.2. Case 2: Small measurement errors ( ${\sigma _{b}^{2}}\ll{\sigma _{a}^{2}}$ )

The analysis is less straightforward when ${\sigma _{b}^{2}}\ne 0$ because Eqs. (36) no longer are simply separable. The system of equations can be solved approximately if ${\sigma _{b}^{2}}\ll{\sigma _{a}^{2}}$ . In this case, the best estimates for the parameters can be expressed in terms of the error-free estimates as

$\displaystyle\hat{\phi}$	$\displaystyle\approx\hat{\phi}_{0}\,\left\{ 1+\frac{{\sigma _{b}^{2}}}{\hat{{\sigma _{a}^{2}}}_{0}}\,\left[1-\hat{\phi}_{0}^{2}+\frac{c_{2}}{c_{0}}\right]\right\}$
and
$\displaystyle\hat{{\sigma _{a}^{2}}}$	$\displaystyle\approx\hat{{\sigma _{a}^{2}}}_{0}-\frac{{\sigma _{b}^{2}}}{\hat{{\sigma _{a}^{2}}}_{0}}\, c_{0}\,\left[1-5\,\hat{\phi}_{0}^{4}+4\,\hat{\phi}_{0}^{2}\;\frac{c_{2}}{c_{0}}\right],$	(40)

to first order in ${\sigma _{b}^{2}}/{\sigma _{a}^{2}}$ , with statistical uncertainties propagated in the conventional manner. Expansions to higher order in ${\sigma _{b}^{2}}/{\sigma _{a}^{2}}$ involve additional correlations, and the exact solution involves correlations at all lags $m$ . If measurement errors are small enough for Eq. (40) to apply, the computational savings relative to other approaches can be substantial, and the amount of data required to achieve a desired level of accuracy in the physically relevant quantities, $\kappa$ and $\gamma$ , can be reduced dramatically.

Figure 4. Power dependence of (a) the trap stiffness, (b) the viscous drag coefficient and (c) the viscous relaxation time for a 1.53 $\mathrm{\upmu}$ $\mathrm{m}$ diameter silica sphere trapped by an optical tweezer in water.

The errors in locating colloidal particles' centroids can be calculated from knowledge of the images' signal to noise ratio and the optical train's magnification (29). Centroid resolutions of 10 $\mathrm{n}$ $\mathrm{m}$ or better can be attained routinely for micrometer-scale colloidal particles in water using conventional bright-field imaging. In practice, however, mechanical vibrations, video jitter and other processes may increase the measurement error. Quite often, the overall measurement error is most easily assessed by increasing the laser power to the optical traps to minimize the particles' thermally driven motions. In this case, $y_{j}\approx b_{j}$ , and ${\sigma _{b}^{2}}$ can be estimated directly.

§ III.3. Trap characterization

The stiffness and viscous drag coefficient can be estimated simultaneously as

$\frac{\kappa}{k_{B}T}=\frac{1-\hat{\phi}^{2}}{\hat{{\sigma _{a}^{2}}}}\text{, \quad and \quad}\frac{\gamma}{k_{B}T\,\Delta t}=-\frac{1-\hat{\phi}^{2}}{\hat{{\sigma _{a}^{2}}}\,\ln\hat{\phi}},$

(41)

with error estimates, $\Delta\kappa$ and $\Delta\gamma$ , given by

	$\displaystyle\left(\frac{\Delta\kappa}{\kappa}\right)^{2}$	$\displaystyle=\left(\frac{\Delta\hat{{\sigma _{a}^{2}}}}{\hat{{\sigma _{a}^{2}}}}\right)^{2}+\left(\frac{2\hat{\phi}^{2}}{1-\hat{\phi}^{2}}\right)^{2}\,\left(\frac{\Delta\hat{\phi}}{\hat{\phi}}\right)^{2}\text{ \quad and}$		(42)
	$\displaystyle\left(\frac{\Delta\gamma}{\gamma}\right)^{2}$	$\displaystyle=\left(\frac{\Delta\hat{{\sigma _{a}^{2}}}}{\hat{{\sigma _{a}^{2}}}}\right)^{2}+\left(\frac{2\hat{\phi}^{2}}{1-\hat{\phi}^{2}}+\frac{1}{\ln\hat{\phi}}\right)^{2}\,\left(\frac{\Delta\hat{\phi}}{\hat{\phi}}\right)^{2}.$		(43)

If the measurement interval, $\Delta t$ , is much longer than the viscous relaxation time $\tau=\gamma/\kappa$ , then $\phi$ vanishes and the error in the drag coefficient diverges. Conversely, if $\Delta t$ is much smaller than $\tau$ , then $\phi$ approaches unity and both errors diverge. Consequently, this approach does not benefit from excessively fast sampling. Rather, it relies on accurate particle tracking to minimize $\Delta\hat{\phi}$ and $\Delta\hat{{\sigma _{a}^{2}}}$ . For trap-particle combinations with viscous relaxation times exceeding a few milliseconds and typical particle excursions of at least 10 $\mathrm{n}$ $\mathrm{m}$ , digital video microscopy provides the resolution needed to simultaneously characterize multiple optical traps (29).

In the event that measurement errors can be ignored ( ${\sigma _{b}^{2}}\ll{\sigma _{a}^{2}}$ ),

$\displaystyle\frac{\kappa _{0}}{k_{B}T}$	$\displaystyle=\frac{1}{c_{0}}\,\left[1\pm\sqrt{\frac{2}{N}\left(1+\frac{2c_{1}^{2}}{c_{0}^{2}-c_{1}^{2}}\right)}\right]$
and
$\displaystyle\frac{\gamma _{0}}{k_{B}T\,\Delta t}$	$\displaystyle=\frac{1}{c_{0}\,\ln\left(\frac{c_{0}}{c_{1}}\right)}\,\left(1\pm\frac{\Delta\gamma _{0}}{\gamma _{0}}\right)$	(44)

where

$N\,\left(\frac{\Delta\gamma _{0}}{\gamma _{0}}\right)^{2}=2+\frac{1}{c_{0}^{2}-c_{1}^{2}}\,\left[\frac{c_{0}^{2}+2c_{1}^{2}\,\ln\left(\frac{c_{1}}{c_{0}}\right)-c_{1}^{2}}{c_{1}\ln\left(\frac{c_{1}}{c_{0}}\right)}\right]^{2}.$

(45)

These results are not reliable if $c_{1}\lesssim{\sigma _{b}^{2}}$ , which arises when the sampling interval $\Delta t$ is much longer or much shorter than the viscous relaxation time, $\tau$ . Accurate estimates for $\kappa$ and $\gamma$ still may be obtained in this case by applying Eq. (40).

As a practical demonstration, we analyzed the thermally driven motions of a single silica sphere of diameter $1.53~\mathrm{\upmu}\mathrm{m}$ (Bangs Labs lot number 5328) dispersed in water and trapped in a conventional optical tweezer. With the trajectory resolved to within $6~\mathrm{n}\mathrm{m}$ at 1/30 s intervals, 1 minute of data suffices to measure both $\kappa$ and $\gamma$ to within 1 percent error using Eqs. (41). The results plotted in Fig. 4(a) indicate trapping efficiencies of $\kappa _{x}/P=\kappa _{y}/P=142\pm 3~\mathrm{p}\mathrm{N}\mathrm{/}{\mathrm{\upmu}\mathrm{m}\;\mathrm{W}}$ . Unlike $\kappa$ , which depends principally on $c_{0}$ , $\gamma$ also depends on $c_{1}$ , which is accurately measured only for $\tau\gtrsim 1$ . Over the range of laser powers for which this condition holds, we obtain the expected $\gamma _{x}/\gamma _{0}=\gamma _{y}/\gamma _{0}=1.0\pm 0.1$ , as shown in Fig. 4(b). The viscous relaxation time becomes substantially shorter than our sampling time for higher powers, so that estimates for $\gamma$ and its error both become unreliable, as expected.

§ IV. Adaptive Optimization

Optimal statistical analysis offers insights not only into the traps' properties, but also into the properties of the trapped particles and the surrounding medium. For example, if a spherical probe particle is immersed in a medium of viscosity $\eta$ far from any bounding surfaces, its hydrodynamic radius $a$ can be assessed from the measured drag coefficient using the Stokes result $\gamma=6\pi\eta a$ . The viscous drag coefficients, moreover, provide insights into the particles' coupling to each other and to their environment. The independently assessed values of the traps' stiffnesses then can serve as a self-calibration in microrheological measurements and in measurements of colloidal many-body hydrodynamic coupling (32). In cases where the traps themselves must be calibrated accurately, knowledge of the probe particles' differing properties gauged from measurements of $\gamma$ can be used to distinguish variations in the traps' intrinsic properties from variations due to differences among the probe particles.

These measurements, moreover, can be performed rapidly enough, even at conventional video sampling rates, to permit real-time adaptive optimization of the traps' properties. Each trap's stiffness is roughly proportional to its brightness. So, if the $m$ -th trap in an array is intended to receive a fraction $\left|\alpha _{m}\right|^{2}$ of the projected light, then instrumental deviations can be corrected by recalculating the CGH with modified amplitudes:

$\alpha _{m}\rightarrow\alpha _{m}\,\sqrt{\frac{\sum _{{i=1}}^{N}\kappa _{i}}{\kappa _{m}}}.$

(46)

Analogous results can be derived for optimization on the basis of other performance metrics. A quasiperiodic pattern similar to that in Fig. 2(c) was adaptively optimized for uniform brightness, with a single optimization cycle yielding better than 12 percent variance from the mean. Applying Eqs. (41) to data from images such as Fig. 2(d) allows us to correlate the traps' appearance to their actual performance.

With each trap powered by 3.4 mW, the mean viscous relaxation time is found to be $\tau/\Delta t=1.14\pm 0.11$ . We expect reliable estimates for the viscous drag coefficient under these conditions, and the result $\gamma/\gamma _{0}=0.95\pm 0.10$ with an overall measurement error of 0.01, is consistent with the manufacturer's rated 10 percent polydispersity in particle radius. Variations in the measured stiffnesses, $\left<\kappa _{x}\right>=0.38\pm 0.06~\mathrm{p}\mathrm{N}\mathrm{/}{\mathrm{\upmu}\mathrm{m}}$ and $\left<\kappa _{y}\right>=0.35\pm 0.10~\mathrm{p}\mathrm{N}\mathrm{/}{\mathrm{\upmu}\mathrm{m}}$ , can be ascribed to a combination of the particles' polydispersity and the traps' inherent brightness variations. This demonstrates that adaptive optimization based on the traps' measured intensities also optimizes their performance in trapping particles.

§ V. Summary

The quality and uniformity of the holographic optical traps projected with the methods described in the previous sections represent a substantial advance over previously reported results. We have demonstrated that optimized and adaptively optimized HOT arrays can be used to craft highly structured potential energy landscapes with excellent fidelity to design. These optimized landscapes have potentially wide-ranging applications in sorting mesoscopic fluid-borne objects through optical fractionation (11); (12), in fundamental studies of transport (33); (9), dynamics (13); (14) and phase transitions in macromolecular systems, and also in precision holographic manufacturing.

We have benefitted from extensive discussion with Alan Sokal. This work was supported by the National Science Foundation under Grant Number DBI-0233971 with additional support from Grant Number DMR-0451589. S.L. acknowledges support from the Kessler Family Foundation.

References

(1)
A. Ashkin, J. M. Dziedzic, J. E. Bjorkholm and S. Chu. “Observation of a single-beam gradient force optical trap for dielectric particles.” Opt. Lett. 11, 288–290 (1986).
(2)
E. R. Dufresne and D. G. Grier. “Optical tweezer arrays and optical substrates created with diffractive optical elements.” Rev. Sci. Instrum. 69, 1974–1977 (1998).
(3)
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, 608–610 (1999).
(4)
J. Liesener, M. Reicherter, T. Haist and H. J. Tiziani. “Multi-functional optical tweezers using computer-generated holograms.” Opt. Commun. 185, 77–82 (2000).
(5)
E. R. Dufresne, G. C. Spalding, M. T. Dearing, S. A. Sheets and D. G. Grier. “Computer-generated holographic optical tweezer arrays.” Rev. Sci. Instrum. 72, 1810–1816 (2001).
(6)
J. E. Curtis, B. A. Koss and D. G. Grier. “Dynamic holographic optical tweezers.” Opt. Commun. 207, 169–175 (2002).
(7)
D. G. Grier. “A revolution in optical manipulation.” Nature 424, 810–816 (2003).
(8)
P. T. Korda, G. C. Spalding, E. R. Dufresne and D. G. Grier. “Nanofabrication with holographic optical tweezers.” Rev. Sci. Instrum. 73, 1956–1957 (2002).
(9)
P. T. Korda, M. B. Taylor and D. G. Grier. “Kinetically locked-in colloidal transport in an array of optical tweezers.” Phys. Rev. Lett. 89, 128301 (2002).
(10)
A. Jesacher, S. Furhpater, S. Bernet and M. Ritsch-Marte. “Size selective trapping with optical ”cogwheel” tweezers.” Opt. Express 12, 4129–4135 (2004).
(11)
K. Ladavac, K. Kasza and D. G. Grier. “Sorting by periodic potential energy landscapes: Optical fractionation.” Phys. Rev. E 70, 010901(R) (2004).
(12)
M. Pelton, K. Ladavac and D. G. Grier. “Transport and fractionation in periodic potential-energy landscapes.” Phys. Rev. E 70, 031108 (2004).
(13)
S.-H. Lee, K. Ladavac, M. Polin and D. G. Grier. “Observation of flux reversal in a symmetric optical thermal ratchet.” Phys. Rev. Lett. 94, 110601 (2005).
(14)
S.-H. Lee and D. G. Grier. “Flux reversal in a two-state symmetric optical thermal ratchet.” Phys. Rev. E 71, 060102(R) (2005).
(15)
A. Jesacher, S. Fürhapter, S. Bernet and M. Ritsch-Marte. “Size selective trapping with optical ”cogwheel” tweezers.” Opt. Express 12, 4129–4135 (2004).
(16)
V. Soifer, V. Kotlyar and L. Doskolovich. Iterative Methods for Diffractive Optical Elements Computation (Taylor & Francis, Bristol, PA, 1997).
(17)
H. He, N. R. Heckenberg and H. Rubinsztein-Dunlop. “Optical particle trapping with higher-order doughnut beams produced using high efficiency computer generated holograms.” J. Mod. Opt. 42, 217–223 (1995).
(18)
H. He, M. E. J. Friese, N. R. Heckenberg and H. Rubinsztein-Dunlop. “Direct observation of transfer of angular momentum to absorptive particles from a laser beam with a phase singularity.” Phys. Rev. Lett. 75, 826–829 (1995).
(19)
K. T. Gahagan and G. A. Swartzlander. “Optical vortex trapping of particles.” Opt. Lett. 21, 827–829 (1996).
(20)
N. B. Simpson, L. Allen and M. J. Padgett. “Optical tweezers and optical spanners with Laguerre-Gaussian modes.” J. Mod. Opt. 43, 2485–2491 (1996).
(21)
M. Meister and R. J. Winfield. “Novel approaches to direct search algorithms for the design of diffractive optical elements.” Opt. Commun. 203, 39–49 (2002).
(22)
J. L. Aragón, G. G. Naumis and M. Torres. “A multigrid approach to the average lattices of quasicrystals.” Acta Cryst. A58, 352–360 (2002).
(23)
K. Svoboda, C. F. Schmidt, B. J. Schnapp and S. M. Block. “Direct observation of kinesin stepping by optical trapping interferometry.” Nature 365, 721–727 (1993).
(24)
L. P. Ghislain, N. A. Switz and W. W. Webb. “Measurement of small forces using an optical trap.” Rev. Sci. Instrum. 65, 2762–2768 (1994).
(25)
F. Gittes, B. Schnurr, P. D. Olmsted, F. C. MacKintosh and C. F. Schmidt. “Microscopic viscoelasticity: Shear moduli of soft materials determined from thermal fluctuations.” Phys. Rev. Lett. 79, 3286–3289 (1997).
(26)
E. L. Florin, A. Pralle, E. H. K. Stelzer and J. K. H. Horber. “Photonic force microscope calibration by thermal noise analysis.” Appl. Phys. A 66, S75–S78 (1998).
(27)
K. Berg-Sørensen and H. Flyvbjerg. “Power spectrum analysis for optical tweezers.” Rev. Sci. Instrum. 75, 594–612 (2004).
(28)
F. Gittes and C. F. Schmidt. “Interference model for back-focal-plane displacement detection in optical tweezers.” Opt. Lett. 23, 7–9 (1998).
(29)
J. C. Crocker and D. G. Grier. “Methods of digital video microscopy for colloidal studies.” J. Colloid Interface Sci. 179, 298–310 (1996).
(30)
G. E. P. Box and G. M. Jenkins. Time Series Analysis: Forecasting and Control (Holden-Day, San Francisco, 1976).
(31)
H. Risken. The Fokker-Planck Equation. Springer series in synergetics (Springer-Verlag, Berlin, 1989), 2nd ed.
(32)
M. Polin, D. G. Grier and S. R. Quake. “Anomalous vibrational dispersion in holographically trapped colloidal arrays.” Phys. Rev. Lett. 96, 088101 (2006).
(33)
P. T. Korda, G. C. Spalding and D. G. Grier. “Evolution of a colloidal critical state in an optical pinning potential.” Phys. Rev. B 66, 024504 (2002).