1. Introduction

The design of phase-only diffractive optical elements (DOEs) which generate a desired intensity distribution in a single specified diffraction plane has received a great deal of attention [1]. Early work on computer generated holography has been used to generate specific three-dimensional (3-D) intensity distributions [2, 3, 4] for applications including holographic displays.

This study was inspired by recent interest in interference lithography for fabricating volume structures, such as photonic crystals, where the 3-D intensity distribution behind a diffraction grating is recorded in a photosensitive material [5]. Better control over the propagating wave-field via the design of the diffraction gratings provides equal control over the properties of the synthesized volume structure. In contrast to computer holograms for display purposes the 3-D intensity structure which is desired does not necessarily comply with the wave equation. This means, the design of suitable diffraction screens typically involves nonlinear optimization to determine the best compromise between several competing properties of the resulting intensity profile. The design of DOEs for generating 3-D intensity distributions has been investigated before, for instance in the context of non-diffracting beams [6, 7, 8, 9] and self-imaging [10, 11].

Our contribution describes an alternative approach to design DOE’s for generating periodic intensity distributions. Based on the fractional Talbot effect [12] and the concept of Talbot array illuminators (TAIs) [13], we employ the matrix formalism of the fractional Talbot effect [14]. This formalism was shown to provide a convenient framework for implementing numerical design algorithms if the output is specified in a single diffraction plane [15].

We show that it is straightforward to implement a similar algorithm for the design of volume intensity distributions. In particular the simulated annealing algorithm described in [15 ] can be used with little modifications to design binary optics DOEs with a small number of discrete phase levels. However, the analysis of the degrees of freedom which can be exploited defines a non-trivial problem, since the propagating wavefield is governed by the wave equation. We discuss how constraints imposed by the wave equation translate into inherent symmetries of periodic wavefields. The performance of the design algorithm then critically depends on the specified intensity to be compatible with these symmetries. Features as well as limitations of the design problem in 3-D are highlighted based on results obtained for a particular example, i.e. the design of TAIs with extended focal depth.

 

Fig. 1. Basic geometry of the optical setup.

Download Full Size | PDF

2. Degrees of freedom

For this study we restrict the discussion to a two-dimensional configuration, i.e. to optical signals which are functions of one transverse coordinate and the longitudinal position. For separable problems this automatically yields the solution to the 3D problem. For non-separable configurations the proposed procedure can be generalized accordingly. In addition, we constrain our discussion to the domain of paraxial optical signals. Figure 1 illustrates the geometry we consider. A periodic binary optics DOE with Q pixels of constant phase per period d is illuminated with a coherent plane wave. Transverse periodicity results in a longitudinal periodicity of the complex amplitude distribution with period z_T = 2d ²/λ, where l is the wavelength of the incident wave. This is the essence of the Talbot effect [16]. For our discussion we consider so-called fractional Talbot planes z_N,M = z_T M/N specified by integer numbers M and N, and we limit our attention to the case of even numbers Q = N/2. It can be shown that the complex amplitude u(x,z_N,M ) can be expressed as a superposition of Q copies of the complex amplitude g(x) in the grating plane z = 0, where each copy is shifted and modulated [14], i.e.

The c_q are the so-called Talbot coefficients. Given the discrete nature of the grating profile it is further possible to represent the complex amplitude at the grating by Q samples g_q = g(q_d /Q). The diffraction amplitude at the fractional Talbot planes can equally be expressed by Q samples u_n = u(nd/Q,z_N,M ), which are related to the grating samples via a linear transformation [14]

where elements D_n,q are the Talbot coefficients c_k , with k = (n - q) mod Q. For all planes (M,N) the matrix D is unitary and symmetric. In other words, the diffraction planes associated with these linear transforms are those where the number of equidistant samples necessary for uniquely describing the complex amplitude corresponds to the degrees of freedom of the optical signal and the diffraction grating. For M = 1 the Talbot coefficients can be calculated as [14]

For M > 1 the Talbot coefficients can conveniently be determined with the help of Eqs. (1) and (2), by computing (D)^M. Now we are prepared to estimate the degrees of freedom of a volume distribution of the intensity which is controlled with the DOE. From Eq. (1) we find that specifying the intensity distribution in any plane z_N,M provides us with Q equations,

which can be used to determine the set of Q ² unknowns g_q$g_{q\mathit{\prime }}^{*}$. This might imply that the intensity can be specified for Q out of TV planes before exhausting the available degrees of freedom. However, the g_q$g_{q\mathit{\prime }}^{*}$ are not independent variables, which becomes obvious if we re-write Eq. (4) by substituting g_q$g_{q\prime }^{*}$ = 2 cos(Δϕ _q,q′). Here, Δϕ _q,q′ is the phase difference between pixels g_q′ . and g_q . Since cases (q, q′) and (q′,q) result in equal intensities and terms q = q′ do not correspond to unknowns, we can specify at the most Q/2 - 1 diffraction planes. We note that Eq. (4) constitutes a nonlinear relationship between the phase only distribution of the DOE and the output signal. This makes the calculation of the number of independent parameters difficult. We therefore emphasize that our discussion only provides an estimate which should not be regarded as a strict condition. Additional constraints are imposed by the physics of wave propagation which limit the variety of intensity distributions which can be implemented. The matrix formalism of the fractional Talbot effect again is the convenient tool to investigate these constraints. In particular it is possible to analyze inherent symmetries of the periodic diffraction amplitude. For diffraction planes where N and M share a common denominator only a subset of matrix elements D_n,q is non-zero, i.e. the number of shifted and modulate copies of g(x) in Eq. (1) is smaller than Q. For instance, the complex amplitude at the plane M = N/2 is described by a single copy of the grating profile shifted by half the period d. This has two consequences. Firstly, as illustrated in Fig. 1, the intensity distributions within the first and the second half of the Talbot length are identical, except for a transverse shift of d/2, and cannot be specified independently. Secondly, at z = z_T /2 the complex amplitude again is phase-only regardless of the structure of the DOE. Obviously, no DOE can be designed to generate the specified 3D intensity distribution with high accuracy if the desired intensity is not consistent with these fundamental properties of the propagating wavefield. Similarly, we find that the wavefield at z = z_T /4 is described by two copies of the grating profile laterally shifted by d/2 and d/4, respectively. This implies that the intensity in this plane cannot exceed twice that of the incident field, no matter how the DOE is designed. This in turn means that we cannot hope to confine all of the intensity to a cross section smaller than d/2. Similar relationships can be deduced for 1/8,1/16,.... of the Talbot distance as well as for multiples of these fractional Talbot planes for which N and M are coprime.

3. Nonlinear optimization

Based on the assumptions in the previous section we implemented a simulated annealing algorithm for synthesizing intensity distributions in transverse and longitudinal direction. The optimization scheme was adopted from previous work on the design of Talbot array illuminators [15]. Simulated annealing has the advantage of being compatible with binary optics, where each pixels can take one of L discrete phase values. The set of discrete phase levels can conveniently be interpreted as a state of the thermodynamic system which is annealed. We use a standard implementation of simulated annealing [17, 18] with only small modifications. The cost function was computed as the mean squared deviation between the desired and the actual intensity distribution. Key for achieving satisfactory performance was the use of the relative change of the cost function as the energy function to compute the Boltzmann distribution for accepting an increase of the cost function. Convergence to the optimum configuration was ensured by introducing only small changes to the phase-only function of the DOE at each iteration.

 

Fig. 2. TAIs with extended focal depth. The intensity is specified in (a) one, (b) two, and (c) three diffraction planes. For the optimum DOE one period of the intensity in all specified planes (upper part) and the intensity distribution in x and z (lower part) is shown.

Download Full Size | PDF

Figure 2 shows the results for designing a TAI with enhanced depth of focus. Each color in the upper row corresponds to the intensity distribution in one diffraction plane. The intensity is normalized with respect to the intensity of the incident plane wave. The design goal was to confine the total intensity to one quarter of the transverse period d. In accordance with the results of the previous section we selected (N = 16,M = 2) for the location of the spot array, i.e. the associated diffraction plane is located at z_T /8. If the intensity is specified for this diffraction plane only, the solution of the design problem can be obtained analytically [19]. This also provides a convenient test case for fine-tuning the convergence properties of the optimization algorithm. Assuming L = 8 discrete phase levels the DOE can compress 100% of the intensity into one quarter of each period (see upper part of Fig. 2(a)). The intensity distribution along the z axis (lower part of Fig. 2(a)) reveals a sharp set of periodic spots at the specified diffraction plane marked with white triangles. The intensity of the propagating wavefront was obtained with a paraxial angular spectrum propagation method. The input vector contained four periods of the DOE where each pixel of the DOE was sampled with four points. If the same intensity distribution is used to specify the intensity in two planes M = (1,2), and in three planes M = (1,2,3) the increase of the longitudinal extension of the spots is clearly recognizable. Figs. 2(b) and (c) show the respective intensity distributions in all specified diffraction planes as well as the intensity distributions of the propagating wavefront.

A quantitative evaluation of the confined intensity as a function of the propagation distance is shown in Fig. 3. The average intensity within the designated quarter of the transverse period is plotted in Fig. 3(a). While specifying the the desired intensity in more than one plane reduces the maximum intensity by little more than 6%, the half width of the intensity distribution is increased about five times. The ideal plateau-like shape is not reproduced and part of the total intensity is confined with the designated volume. In part, the deviations from the specified distribution are caused by the phase quantization and the spatial quantization of the DOE. Similar to DOEs which are design for the Fraunhofer domain [20] either the number of pixels per period Q or the number of phase levels L can be increased to improve primarily the confined portion of intensity and to some extent the shape of the intensity distribution. Figure 3 also contains the result for three diffraction planes and L = 16 phase levels which results in a small, but distinct improvement for diffraction planes between M = 2 and M = 4. For comparison we also evaluated the part of the longitudinal intensity distribution which is confined to one half of the period (Fig. 3(b)). The behavior is similar to that in Fig. 3(a), however a slightly more distinct plateau shape can be observed. This result can be generalized in that the optimization can be performed for a high compression ratio in order to obtain acceptable results for the design goal with smaller compression ratio.

 

Fig. 3. Average intensity confined to (a) 1/4, (b) 1/2 of the transverse period. Solid lines refer to DOE with L = 8 discrete phase levels, the dashed green line corresponds to L = 16.

Download Full Size | PDF

4. Discussion and Summary

The example used in Figs. (2) and (3) primarily serves as an illustration of the behavior we observed with a larger set of simulations. The algorithm typically converges to an acceptable solution if for optimization in a single diffraction plane a perfect or near perfect design can be obtained. Consistent with our expectations the algorithm invariably provides poor results if the specified intensities do not match the inherent properties of periodic propagating waves. The nonlinear optimization procedure cannot be applied blindly, but requires a substantial amount of guidance, to define the optimization problem in accordance with the wave equation.

We also point to the behavior of the intensity in Fig. 2 for diffraction planes z > z_T /4. From Fig. 2(c) we find that the intensity is also confined to some extent. The intensity can be specified for additional planes in this range as well. The algorithm then determines the optimum configuration which is a reduced intensity confinement for planes z < z_T /4 and a better confinement for planes z > z_T /4. However, the resulting design is governed by the diffraction pattern at plane z_T /4 where the maximum compression ratio cannot exceed a factor of two.

Finally, we would like to emphasize that the assumptions we have used for this study do not limit the generality of this approach. This includes the strict periodicity of the intensity distribution. For sufficiently large Q it is possible to investigate the case, where the intensity distributions of adjacent periods are essentially isolated from one another. This in effect allows treatment of non-periodic problems within the framework of the fractional Talbot effect. In particular, this allows us to recast the design problem of non-diffracting beams [7, 8, 9] into the framework of the fractional Talbot effect. However, while we believe that a thorough investigation of this relationship would be beneficial, it is beyond the scope of this contribution.