Reduction of two-photon holographic speckle using shift-averaging

Holographic speckle is a major impediment for the emerging applications of multiphoton holographic projection in biomedical imaging, photo-stimulation and micromachining. Time averaging of multiple shifted versions of a single hologram (“shift-averaging”) is a computationallyefficient method that was recently shown to deterministically eliminate holographic speckle in single-photon applications. Here, we extend these results and show, computationally and experimentally, that in two-photon holographic excitation shift-averaging also reduces holographic speckle better than “random” averaging of multiple calculated holograms. ©2011 Optical Society of America OCIS codes: (030.6140) Speckle; (070.6120) Spatial light modulators; (090.2870) Holographic display; (190.4180) Multiphoton processes. References and links 1. N. J. Jenness, K. D. Wulff, M. S. Johannes, M. J. Padgett, D. G. Cole, and R. L. Clark, “Three-dimensional parallel holographic micropatterning using a spatial light modulator,” Opt. Express 16(20), 15942–15948 (2008). 2. 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(3), 1810–1816 (2001). 3. C. Lutz, T. S. Otis, V. DeSars, S. Charpak, D. A. DiGregorio, and V. Emiliani, “Holographic photolysis of caged neurotransmitters,” Nat. Methods 5(9), 821–827 (2008). 4. L. Golan, I. Reutsky, N. Farah, and S. Shoham, “Design and characteristics of holographic neural photostimulation systems,” J. Neural Eng. 6(6), 066004 (2009). 5. R. W. Gerchberg and W. O. Saxton, “A practical algorithm for the determination of phase from image and diffraction plane pictures,” Optik (Jena) 35, 237–246 (1972). 6. T. S. McKechnie, “Speckle Reduction,” in Laser Speckle and Related Phenomena, J. C. Dainty, ed. (SpringerVerlag, 1975). 7. E. Papagiakoumou, V. de Sars, D. Oron, and V. Emiliani, “Patterned two-photon illumination by spatiotemporal shaping of ultrashort pulses,” Opt. Express 16(26), 22039–22047 (2008). 8. J. Amako, H. Miura, and T. Sonehara, “Speckle-noise reduction on kinoform reconstruction using a phase-only spatial light modulator,” Appl. Opt. 34(17), 3165–3171 (1995). 9. J. W. Goodman, “Some fundamental properties of speckle,” J. Opt. Soc. Am. 66(11), 1145–1150 (1976). 10. L. Golan and S. Shoham, “Speckle elimination using shift-averaging in high-rate holographic projection,” Opt. Express 17(3), 1330–1339 (2009). 11. Y. Pan, X. Xu, S. Solanki, X. Liang, R. B. A. Tanjung, C. Tan, and T.-C. Chong, “Fast CGH computation using S-LUT on GPU,” Opt. Express 17(21), 18543–18555 (2009). 12. J. P. Parry, R. J. Beck, J. D. Shephard, and D. P. Hand, “Application of a liquid crystal spatial light modulator to laser marking,” Appl. Opt. 50(12), 1779–1785 (2011). 13. V. Nikolenko, B. O. Watson, R. Araya, A. Woodruff, D. S. Peterka, and R. Yuste, “SLM microscopy: scanless two-photon imaging and photostimulation using spatial light modulators,” Front. Neural Circuits 2, 5 (2008). 14. V. R. Daria, C. Stricker, R. Bowman, S. Redman, and H. A. Bachor, “Arbitrary multisite two-photon excitation in four dimensions,” Appl. Phys. Lett. 95(9), 093701 (2009). 15. G. Bautista, M. J. Romero, G. Tapang, and V. R. Daria, “Parallel two-photon photopolymerization of microgear patterns,” Opt. Commun. 282(18), 3746–3750 (2009). 16. E. Papagiakoumou, F. Anselmi, A. Bègue, V. de Sars, J. Glückstad, E. Y. Isacoff, and V. Emiliani, “Scanless two-photon excitation of channelrhodopsin-2,” Nat. Methods 7(10), 848–854 (2010). 17. S. Shoham, “Optogenetics meets optical wavefront shaping,” Nat. Methods 7(10), 798–799 (2010). 18. J. W. Goodman, Introduction to Fourier Optics, 3rd ed. (Roberts & Company, 2005). 19. J. Chen, S. M. Morris, T. D. Wilkinson, J. P. Freeman, and H. J. Coles, “High speed liquid crystal over silicon display based on the flexoelectro-optic effect,” Opt. Express 17(9), 7130–7137 (2009). #155182 $15.00 USD Received 22 Sep 2011; revised 9 Nov 2011; accepted 11 Nov 2011; published 5 Dec 2011 (C) 2011 OSA 19 December 2011 / Vol. 19, No. 27 / OPTICS EXPRESS 25891


Introduction
Computer Generated Holographic (CGH) projection using Liquid-Crystal Spatial Light Modulators (LC-SLMs) is a powerful strategy for generating intense multi-focal dynamical light distributions in two or three dimensions, which has recently found applications in a variety of technological fields, including micro-fabrication [1], optical trapping [2] and neural stimulation [3,4].SLMs have the added advantage of being dynamic; rapidly-changing sequences of CGHs may be displayed, with rates of up to several kHz.Holographic speckle is a major source of image distortion associated with the projection of patterns consisting of contiguous shapes (rather than sparse diffraction-limited spots), as required by certain applications, such as microfabrication and neural stimulation [3,4].In the reconstruction plane, a rectangular lattice of sinc-shaped Point Spread Functions (PSFs) is observed.Ideally, the complex amplitude of each such PSF should be user-defined.However, most phase retrieval algorithms, such as the Gerchberg-Saxton algorithm [5], only present restrictions on the magnitude of the field in the SLM and the reconstruction plane, with the phases remaining as free parameters.In practice, this results in random interferences between the contiguous PSFs leading to speckled light patterns.A common theme for many of the methods developed to solve or circumvent this problem is time-averaging [6].For example, using a rotating diffuser, it is possible to physically average many instances of the resultant pattern, each with a different pseudo-random phase function, thus reducing speckle [7].This, however, requires another optical element in the path, which could prove costly in terms of image quality.Another solution is to calculate N different holograms for the same desired pattern and display them in rapid succession [8], causing the speckle contrast to be reduced by a factor of N 1/2 [9].Of course, this method is marred by significantly heavier computational demands per pattern.
Recently, we introduced shift-averaging [10], a new method for speckle elimination in which a single hologram is calculated and undergoes several periodic two-dimensional cyclic shifts, which are sequentially displayed (as in regular time averaging).We found that cyclic shifts can be chosen that ensure the deterministic elimination of the mathematical cross-terms describing the interference between each pair of PSFs in the time-averaged intensity of the reconstruction plane pattern.In practice, this deterministic elimination (instead of random averaging in conventional time averaging) yields a major reduction of the intensity speckle contrast even when averaging a small number of holograms.Combined with the associated major reduction in computational burden, this has lead to the rapid adoption of this technique in applications ranging from holographic televisions [11], to one-photon neural stimulation [4] to laser marking [12] (which is likely to increase as faster SLMs become increasingly available).
One of the most exciting recent developments in holographic projection is the emergence of major applications where nonlinear optical excitation is used to achieve tight axial sectioning.These emerging applications include two-photon holographic fluorescence microscopy [13,14], photo-polymerization [1,15] and activation of light-gated ion channels for stimulating neurons in the brain [16].The problem of holographic speckle is even more critical in nonlinear projection than it is in one-photon projection: partially because speckle contrast is much stronger in <I 2 > than it is in <I> [7] and partially because the use of image patches is more crucial in nonlinear excitation where high numerical aperture objectives are used and individual diffraction limited spots have submicron dimensions that are too small to excite "large" targets like cells [16,17].In certain applications of two-photon CGH, holographic speckle has arguably emerged as the absolute limiting factor, and has led researchers to explore "smoother" alternatives like generalized phase contrast [16].A simple speckle-reduction solution like shift-averaging, could potentially be useful in tackling this challenge, however, shift-averaging deterministically eliminates speckle-forming cross-terms only in the time-averaged light intensity.In contrast, the expression for the time-averaged squared-intensity, which determines the integrated efficacy of two-photon phenomena, has a significantly more complex mathematical structure with a dramatically higher number of cross-terms, and there is no a priori reason to assume that it would become speckle-free under shift averaging.
In this paper, we study the performance of shift-averaging for reducing speckle noise in the squared intensity pattern.In section 2, the mathematical ground for one-photon shiftaveraging is first laid, followed by an analysis of the effect of shift-averaging on speckle noise in the two-photon case, which is then demonstrated by a simulation in section 3. Section 4 presents the experimental results obtained using a two-photon holographic projection system, which are discussed, along with the other findings, in section 5.

One-photon case
The following mathematical derivations directly follow [10].For simplicity, we assume a square SLM containing M × M pixels and define the phase-only hologram as a 2D matrix: exp( ); , 1, 2,.., , where mn φ are the phase values displayed on the SLM.The 2D discrete Fourier transform (DFT) of the SLM plane, a complex-valued 2D matrix, can be defined as A plane wave interacting with the SLM would have its amplitude transmitted according to , , where ⊗ indicates the convolution operation and where the definitions of the rect and delta are according to [18].If we also assume a wave with a unit-amplitude, the electrical field in the reconstruction plane may be written as ( , ) ( ( , )) ( , ), assuming f kl values are periodic in k and l, with periodicity M. The series of functions S kl are assumed to include all the mathematical information describing the wave's propagation through optical apertures, such as that of the SLM or that of any lens in the system.It also includes a slowly decaying sinc envelope brought about by the finite size of the SLM's pixel.For example, if the initial wave overfills the SLM's aperture, but is not truncated by any lenses that follow, S kl would become ( , ) sinc( , ) sinc( , ); if, on the other hand, the wave were to overfill a lens, S kl would include a Bessel term.However, the exact constituents of S kl are implicit in the rest of the analysis, provided that it maintains the form of a narrow Point Spread Function (PSF).It follows from Eqs. ( 4) and ( 5) that the reconstruction plane can be thought of as an infinite 2D lattice of such PSFs, each with its own phase introduced by the f kl values.In the case of Eq. ( 5), each such PSF is equal to zero in all lattice points save the one corresponding to it, rendering the reconstruction on the lattice points exact and interference-free.However, in every other point in the reconstruction plane, an interference occurs between the contributions of an infinite number of PSFs from all lattice points; assuming the desired reconstructed shape is that of a contiguous patch, the randomness of the phase values of f kl causes rapid spatial fluctuations in the electrical field, which is manifested as speckle noise.
Because the PSFs are spatially narrow functions, it is possible to approximate the electrical field's value at a specific point ( , ) u v , by taking into account the contributions of the lattice points in its immediate neighborhood only, as shown in Fig. 1(a).Assuming a square neighborhood consisting of c by c such PSFs, the electrical field in each point can be approximated by , , , where the indices k and l are as shown in Fig. 1(a).The square modulus of the approximation in Eq. ( 6) describes the intensity at a point ( , ) where we have defined two sets of indices: [ , ] k l and [ , ] r s .Assuming N holograms, theoretically resulting in the same reconstruction pattern but each accompanied with its own speckle patterns, and assuming they are displayed sequentially and periodically on the SLM at a fixed rate, the average intensity recorded over one period of display is given simply by averaging several instances of Eq. ( 7): where the index a refers to the a -th hologram.Using a one-photon system, a perfect timeaveraging-based speckle-elimination method would use N holograms that would eliminate all cross-terms, that is, the terms which include more than a single lattice points.Formulated mathematically, this requirement is akin to writing ( , ) .
In a previous paper [10], we have shown that by calculating only a single hologram and shifting its contents cyclically and periodically in a specific manner, the speckle is eliminated completely, in theory, and almost completely, in practice.Cyclic shifting, illustrated in Fig. 1(b), introduces a linear phase addition:  , , , turns the left side of Eq. ( 9) into where the last equation is true because the sums are of the roots of unity.The above result indicates that the specific two-dimensional cyclic, periodic shifting used for time-averaging purposes (shift-averaging, for short) theoretically eliminates all cross-terms from the average intensity of the electrical field in the reconstruction plane.

Two-photon case
The previous subsection's results are confined to phenomena which are proportional to the field's intensity.However, in two-photon phenomena, which are proportional to the square of the intensity, using a time-averaging approach will yield an effect proportional to the average of the squared intensities, which is different from a square of the averaged intensity (a square of the expression in ( 8)).Therefore, the specific choice of shifts defined by (12) does not guarantee a speckle-free reconstruction in the two-photon case, and its effect requires a different analysis that will be presented in this subsection.
The squared intensity of a single holographic reconstruction is given by ( ) ( ) where two new sets of indices were introduced: [ , ] g h and [ , ] p q .It is noteworthy that the number of terms for the squared intensity is now c 8 , indicating a much more complex speckle pattern.Again, upon time-averaging, Eq. ( 14) becomes ( )  The ideal averaging method would eliminate all c 8 terms except those for which all the first indices from each of the four pair presented above coincide, as well as all of the second indices.In other words, ideal averaging would eliminate all cross-terms, leaving us with c 2 terms, each corresponding to a single lattice point.Mathematically put, ideally

I u v S S S S f f f f N
However, trying to apply shift-averaging with the above-mentioned shifts in this case results in It has been shown, then, that using the shift-averaging method that worked well for the intensity case, performs differently for the case of the squared intensity: rather than being left with two four-dimensional Kronecker's deltas, as is ideal, we are instead left with two twodimensional ones.Effectively, this means that not all the cross-terms are averaged out: instead of c 8 terms, out of which one would ideally remain with c 2 , here we are left with all the crossterms satisfying ( ) ( ) The number of pairs ( ) . Thus, the number of quadruplets satisfying each equation is: , of which c terms are desirable diagonal terms.Hence, the overall proportion of unwanted terms which remain after shift averaging is ( ) a fraction on the order of c −2 .Figure 2 shows a graphic representation of this term as a function of c; even at a low c = 4 (i.e., using 16 holograms) the fraction of remaining crossterms is roughly 3%, indicating a major improvement potential in the expected SNR over the case of random averaging of c 2 samples.So, while shift-averaging cannot eliminate all existing cross-terms in the quantity of relevance for two-photon phenomena (the field's squared intensity), it does average out the vast majority of these cross-terms.

Simulation
The previous result provides an initial indication that shift averaging may considerably ameliorate two-photon holographic speckle phenomena.To further evaluate this effect we calculated and shift-averaged holograms producing two different patterns (a square patch and multiple cross-like patches, see Fig. 3).For the different projected patterns we evaluated the speckle contrast, defined as the ratio between the standard deviation and the mean value of the signal of interest G (for example, one-or two-photon fluorescence), across the pattern's area: As can be seen in Fig. 3(a), two-photon holographic speckle is much higher than in the one-photon case: without any averaging, the speckle contrast is about twice as large.However, as little as 4x4 averaging (16 holograms) significantly improves the two-photon holographic speckle, reducing it by a factor of 4 to 6.Moreover, compared to regular timeaveraging methods, not only is two-photon shift averaging beneficial in its significantly lower computational demands, but also in terms of the better speckle amelioration achieved using the same number of holograms, as shown in Fig. 3(b).This improvement appears particularly strong in the case of the pattern consisting of a number of patches -a difference that probably has to do with the smaller size of the contiguous patches.A point in a bigger patch is affected by more neighbors, requiring a larger c-by-c neighborhood to achieve the same speckle reduction effect obtained with a smaller patch.

Experimental results
The performance of two-photon shift-averaging was then experimentally verified using a multiphoton SLM-based holographic projection system portrayed in Fig. 4.An ultrafast Titanium-Sapphire laser (MaiTai WB, Spectra-Physics), emitting 180-fs 800nm laser pulses was expanded and its polarization matched to that of the SLM (XY Phase, Boulder Nonlinear Systems), before having its width modified by another telescope (L1 and L2) to match the objective lens's back aperture (60x, Apo NIR; Nikon).A CCD (GC1380H, Prosilica) was used to image the reconstruction plane.Figure 5 shows an example of shift averaging performed in the optical system portrayed above.A hologram corresponding to a pattern of several circular patches was calculated and displayed on the SLM.Upon projection onto a thin slide of fluorescein, the projected pattern gives rise to a very noisy two-photon fluorescence image.When applying 4x4 shift averaging, however, the high-frequency noise is mostly gone and the signal appears to be much smoother, accompanied by a speckle contrast reduction of about one order of magnitude.

Discussion
In this work, we have analyzed the performance of the shift-averaging method in the case of two-photon projection, where the observed excitation phenomenon is proportional to the squared intensity.We have shown analytically that shift-averaging eliminates many of the unwanted speckle terms in the average squared intensity pattern, although some unwanted terms remain.Our simulations show that shift-averaging equals or outperforms regular averaging in terms of holographic speckle contrast.This is achieved with the calculation of a single hologram per projected pattern, saving computation time compared to the multiplehologram averaging approach.The effectiveness of shift-averaging for ameliorating speckle was also experimentally demonstrated by imaging two-photon fluorescence excited by holographic patterns.
The applicability of two-photon shift averaging has some important limitations.One limitation of this approach is that it relies on the hologram being shift invariant except for the induced phase modulations in the output light pattern, and is therefore applicable only in light projection of 2D light patterns.This, however covers most of the recent applications of optically sectioned holographic nonlinear excitation [13,16], where the 3rd dimension is obtained by moving the objective.
Secondly, as the SLM used in this work could operate at a refresh rate of 100Hz at most, shift averaging using 16 holograms (4x4) reduces the maximal refresh rate to roughly 6Hz.Generally, time averaging methods are more suitable for fast SLMs, such as the binary, ferroelectric (rather than nematic) SLMs [10], which could operate at refresh rates of several kHz.Unfortunately, ferroelectric SLMs have a relatively low diffraction efficiency, and the practical implementation of rapid time-averaging in holographic two photon projection may await the development of projection tools that combine the benefits of high diffraction efficiency and fast refresh rates [19].

Fig. 1 .
Fig. 1.(a) Illustration of the c by c neighborhood of an arbitrary point in the (u,v) plane.(b) Illustration of two-dimensional cyclic shifts.Red ovals indicate a reference area in the hologram.

where 1 a d and 2 ad
indicate the amount of shift along each direction in the hologram plane.Choosing the following parameters:

Fig. 2 .
Fig. 2. Fraction of terms that remain out of c 8 initial terms after shift averaging.

Fig. 3 .
Fig. 3. Simulation of the performance of the shift-averaging method for one-photon and twophoton phenomena.(a) Performance for a single square patch (left panes) and a number of small patches (right panes).(b) Comparison of speckle contrast reduction for one-and twophoton shift-averaging and two-photon regular averaging as a function of the number of shifts in each direction, c.

Fig. 4 .
Fig. 4. Sketch of the optical system.BE -Beam Expander; λ/2 -half-wave plate; SLM -Spatial Light Modulator; M -mirror; Lenses L1 and L2 form a telescope, while L3 is used to image the plane onto the CCD; DM -dichroic mirror; OL -Objective Lens; B -zero-order point blocker.

Fig. 5 .
Fig. 5. (a) Demonstration of shift averaging in two-photon fluorescence.Using 16 holograms (4x4 shift averaging) the speckle contrast was reduced almost tenfold.Scale bar is 10µm.(b) Comparison of the cross-section along the diameters of the different patches before (left) and after (right) shift-averaging.