We explore implementations of an experimental scheme, proposed by Gori and Santarsiero [1,2], for the generation of a wide range of spatially partially coherent optical beams. The scheme is based on a globally incoherent field with a specific spatially varying intensity distribution p(v) at the input plane of an optical system that performs a linear transformation described by an arbitrary spatial kernel function H(ρ,v). Partially coherent fields generated in this way are guaranteed to have genuine (non-negative definite) cross-spectral density functions (CSDs) W(ρ1,ρ2).

In our implementation the required spatial variation of the intensity distribution is created by first reducing the spatial coherence of a laser beam with a random rotating diffuser [3–5] and then imaging the scattered field onto an intensity-modulating metasurface. The kernel function of the optical system is taken to be of Fourier form, which is sufficient here since the fields to be demonstrated follow the Schell-model [6]. The generation scheme and its practical limitations are discussed in Section 2.

Among numerous possible options, we concentrate on a recently introduced class of azimuthally periodic and radially quasi-periodic partially coherent beams [7]. In a special case these beams have annular far-zone distributions, and in a certain limit they reduce to fundamental Laguerre–Gaussian correlated beams [8,9]. To make the present paper self-contained, we cover the main properties of these fields in Section 3.

The metasurfaces are realized as binary diffractive structures based on Bragg diffraction by surface-relief gratings as information carriers [10]. The spatial intensity information for coherence control is included by spatial-frequency modulation with a continuous variation of the local fill factor of the carrier grating [11,12]. The operating principle, design, and fabrication of the metasurfaces are the subject of Section 4.

The experimental results are presented in Section 5, along with simulations based on the theory in Section 3. Results of beam propagation over several distances Δz in the half-space z>z0 are presented, and the complex degrees of coherence at the plane z=0, which give rise to the propagation properties, are measured with a mirror-based wavefront folding interferometer (WFI) [13].

While only particular examples of partially coherent beams are demonstrated here, the generation scheme is applicable to any beam covered by the model in Refs. [1,2]. It can also be extended to vector beams. These issues are discussed in Section 6 before conclusions.

Gori and Santarsiero [1] introduced a class of genuine spatially partially coherent fields, for which the CSD has a non-negative definite form W(ρ1,ρ2)=∫−∞∞p(v)H∗(ρ1,v)H(ρ2,v)d2v.(1)Here, ρ1 and ρ2 are transverse position coordinates in a plane z=constant, p(v) is a non-negative function of a 2D variable v, and H(ρ,v) is a kernel of an arbitrary linear transformation. The frequency dependence of the CSD is implicit. A wide range of physically relevant CSDs is obtained if the kernel function is taken to be of the form H(ρ,v)=H(ρ)exp(−i2πρ·v),(2)which can be realized with a 2D Fourier-transforming optical system. Upon inserting Eq. (2) into Eq. (1) we get W(ρ1,ρ2)=H∗(ρ1)H(ρ2)∫−∞∞p(v)exp(−i2πΔρ·v)d2v,(3)where Δρ=ρ2−ρ1. This implies that the complex degree of spatial coherence depends only on the coordinate difference, resulting in fields which are of the Schell-model form. Fields that depart from this model require more sophisticated kernels, as discussed in Section 6.

The schematic setup for generating spatially partially coherent beams with genuine CSDs is illustrated in Fig. 1(a). One starts from a globally incoherent field at the input plane A of the linear system T, where the intensity distribution is proportional to the kernel function H(ρ,v). The weight function p(v) is embedded in the metasurface at A to modulate the intensity profile; then the CSD given by Eq. (1) is produced at the output plane B. If desired, the output field can be modulated by adding a complex-amplitude transmission filter with transmission function t(ρ) at plane B, which transforms the CSD as W′(ρ1,ρ2)=t∗(ρ1)t(ρ2)W(ρ1,ρ2)(4)and therefore adds freedom to the choice of realizable fields. In the simplest case the filter is just a hard or smoothly apodized aperture.

The system we use to generate the incoherent field at plane A is illustrated schematically in Fig. 1(b). In our implementation a linearly polarized HeNe laser beam (wavelength λ0=632.8nm) is focused to a width wD on a standard rotating diffuser to reduce its spatial coherence area. The scattered field is collimated by lens L1 and passes through a Gaussian apodization filter G to produce a Gaussian Schell-model (GSM) beam [5]. The transmitted beam is then demagnified by an afocal system consisting of lenses L2 and L3 to illuminate the metasurface with a GSM beam of radius wA.

Admittedly, the use of a rotating diffuser is a bulky and somewhat cumbersome approach to reduce spatial coherence. For future needs we therefore put forward an alternative technical solution: moving the diffuser randomly in the xy plane with a piezoelectric device. Random movement over a sufficient area is possible since piezoelectric devices allow fast 2D scanning over distances larger than the required coherence area. Strictly speaking, the field on A is not completely incoherent, but it is quasi-homogeneous. Therefore the conditions for the validity of Eq. (1) do not hold exactly, but non-negative definiteness is nevertheless ensured since the field is physically realizable.

We demonstrate the scheme in Fig. 1 by applying it to a recently introduced class of radially quasi-periodic and azimuthally periodic (RAP) beams [7]. These are generally Bessel-correlated beams with a Gaussian intensity profile at the source plane z=z0, which, in a special case, can have an annular distribution of radiant intensity. The main characteristics of RAP beams are briefly covered to make the paper self-contained and to summarize the formulas that enable numerical simulations of experimental results.

We introduce RAP beams by assuming that the function H(ρ) in Eq. (2) has a Gaussian form H(ρ)=exp(−ρ2/w02). The weight function p(v) is written in polar coordinates v=(v,ϕ) as p(v,ϕ)=pmexp(−2π2σ02v2)[1+CmJm(2πρ0v)cos(mϕ)](5)with pm=2πσ02{{1+Cmexp[−ρ02/(2σ02)]}−1whenm=0,1whenm≠0.(6)Here, σ0 and ρ0 are positive constants, Jm(z) denotes the Bessel function of the first kind and order m, C0=±1, and Cm=1 for m>1. Equations (1) and (2) result in CSDs of the form W(ρ1,ρ2)=S(ρ1)S(ρ2)μ(Δρ), where the spectral density at z=z0, using Eqs. (5) and (6), has a Gaussian profile S(ρ)=exp(−2ρ2/w02).(7)

The complex degree of coherence (DOC) for solutions m≥1 is μ(Δρ)=exp(−Δρ2/2σ02)×[1+Cm(−i)mcos(mΔφ)exp(−ρs22)Im(ρsΔρσ0)],(8)where Δρ=(Δρ,Δφ), ρs=ρ0/σ0, and Im(z) denotes the modified Bessel function of the first kind and order m.

A closed-form expression for the propagated CSD in the paraxial domain can be found with the standard Fresnel formula [7]. In particular, the spectral density, which we will measure in Section 5, evolves as $$\begin{gathered} S(ρ,z)=w02w2(Δz)exp[−2ρ2w2(Δz)]{1+Cm(−1)mcos(mφ) \\ ×exp[−R22q2c(Δz)]⁢Jm[2Rq1+q2Δz/zR1+(Δz/zR)2ρw0]}, \end{gathered}$$(9)where ρ=(ρ,φ), Δz=z−z0 is the propagation distance from plane B, q=σ0/w0, and R=ρ0/w0. Further, w(Δz)=w0[1+(Δz/zR)2]1/2(10)is the beam width at z>0, and the quantity c(Δz) is c(Δz)=1+(Δz/zG)21+(Δz/zR)2,(11)where zG=kw02/2 is the Rayleigh range of a fully coherent Gaussian beam, and zR=zG(1+q−2)−1/2(12)is the Rayleigh range of a partially coherent GSM beam [14] characterized by coherence parameter q.

Let us consider RAP beams with m=0 and C0=±1 in Eq. (5), which gives a weight function p(v,ϕ)=p0exp(−2π2σ02v2)[1±J0(2πρ0v)],(13)where p0=2πσ02B0 and B0=[1±exp(−ρs2/2)]−1.(14)

The spectral density at the source plane is still given by Eq. (7), but the complex DOC takes the form μ(Δρ)=B0exp(−Δρ2/2σ02)×[1±exp(−ρs22)I0(ρsΔρσ0)].(15)

The spectral density at z>0, S(ρ,z)=B0w02w2(Δz)exp[−2ρ2w2(Δz)]×{1±exp[−R22q2c(Δz)]⁢J0[2Rq1+q2Δz/zR1+(Δz/zR)2ρw0]},(16)takes an annular (donut-shaped) form upon propagation to the far field. It therefore shares some properties of other ‘hollow’ beams, including fully coherent higher-order Bessel beams [15,16], orbital angular momentum, azimuthally polarized beams [17,18], generalized Bessel–Gauss beams [19], partially coherent dark propagation-invariant beams [20,21], and various Bessel-correlated beams that may have donut shapes either initially or in the far-zone [22]. Donut-shaped beams have proven important in particle trapping [23] and stimulated emission depletion (STED) microscopy [24].

In the limiting case ρs→0, and when C0=−1, the complex DOC given in Eq. (15) reduces to μ(Δρ)=exp(−Δρ22σ02)(1−Δρ22σ02),(17)which is the DOC of a Laguerre–Gauss Schell-model (LGSM) source of order n=1 [8]. The spectral density of the LGSM beams takes the form S(ρ,z)=w02w2(Δz)exp[−2ρ2w2(Δz)]×[c(Δz)+2Δz2zR2ρ2(w02+σ02)(Δz2+zR2)2].(18)

If C0=1 and we again let ρs→0, Eq. (15) becomes the DOC of a GSM beam. Figure 2 illustrates some properties of the m=0 solutions not discussed in Ref. [7]. Cross sections of the DOC plotted in Fig. 2(a) show central peaks and side lobes of expanding widths when the value of ρs in the I0 term in Eq. (15) increases. The far-field intensity in Fig. 2(b) is also affected: the central minimum is close to zero at small values of ρs, which also provides the widest central dark region.

In Ref. [7] we presented an on-axis encoding technique of the metasurface at plane A, which can be realized in the axially symmetric geometry of Fig. 1. However, our experiments showed that this scheme has axial noise, which proved difficult to avoid. We therefore opted for an alternative, off-axis coding scheme. Our method is effectively similar to that used in Ref. [25], with the difference that we modulate the field amplitude with a more compact metasurface, instead of filtering the diffracted light from a spatial light modulator. In addition, we can obtain extremely high diffraction efficiency with the off-axis encoding. The drawback of our technique is the somewhat trickier alignment of the optical components compared to the on-axis method.

The metasurface is realized using a high-frequency carrier grating as introduced in Ref. [10]. This is a binary surface-relief grating with depth h and fill factor f=c/d as design parameters, as defined in Fig. 3(a). The grating operates at large angles of incidence θ∼30°. Hence the geometry of Fig. 1 is folded at A by twice this angle. The design angle of incidence is the Bragg angle θ=θB given by sinθB=λ02nd.(19)

Assuming that the incident field is linearly polarized in the direction of the grating grooves, we design the carrier grating for maximum efficiency of order l=−1 in TE polarization for λ0=633nm, n=1.49, and d=λ0.

The efficiency of the carrier grating as a function of fill factor f=c/d was calculated using the Fourier modal method (FMM) [26]. The results are shown in Fig. 3(b). The efficiency has two maxima reaching almost unity, at f∼0.16 and f∼0.62. We employ the second since the efficiency drops smoothly to zero when f approaches unity. Since the incident field is not fully collimated, we plot in Fig. 3(c) the efficiency as a function of the angle of incidence θ (Bragg selectivity of the grating) for several fill factors over the range 0°<θ<40°. It turns out that the angular selectivity is not a serious limiting factor for operation of the grating.

The intensity coding is performed as local variations of the fill factor f as in Refs. [10,11]. Fill factors between 0.62 and 0.9, marked with dashed vertical lines in Fig. 3(b), are used for efficiency encoding since they are practically feasible. The target grating depth was chosen to be h=1350nm, giving η−1≈0.95 at the second maximum f=0.62.

The sample was fabricated by a standard electron beam lithography procedure. A 1-inch-diameter and 0.5-mm-thick SiO2 substrate was cleaned and spin-coated with a 1.4-μm-thick positive-tone resist layer (PMMA), which was coated with a conductive 30-nm-thick copper layer via thermal evaporation. After creating the grating patterns with electron beam lithography, the copper layer was removed in HNO3, followed by resist development in MIBK:IPA (1:3 solution) for 75 s and 30 s rinsing with IPA.

Scanning electron microscope images of one of the samples are shown in the Fig. 4. In Fig. 4(a) we show a cross section around a location where the fill factor changes abruptly from 0.8 to 0.62, and in Fig. 4(b) a top view of the region is shown where f is fairly constant. The fill factor is fairly well reproduced, but in view of Fig. 4(a) the groove depth is somewhat reduced especially at large f. This could be improved by increasing the electron beam dose. However, its main effect is to reduce η−1, which should be low at large f according to Fig. 3. In fact the modulation of h helps us to achieve the low values of η−1 required beyond the second vertical line in Fig. 3(b).

The focal lengths of the lenses in Fig. 1(b) before plane A are flexible, but we chose L1=150mm, L2=200mm, and L3=40mm. We measured wD=31μm and wA≈880μm. The GSM was used as in Ref. [5] to estimate the coherence radius at the filter G (Thorlabs NDYR20A) to be σG≈210μm, which means that the field is quasi-homogenous with σ1/wG≈1/21. Measurement with a WFI gives the results in Fig. 5, where we show cross sections of the spectral density and the absolute value of the spatial DOC at plane A, from which wA≈880μm. While |μ(Δx)| is not completely Gaussian but has some side lobes, we can estimate σA≈42μm. The field at plane A is definitely quasi-homogeneous and, since the coherence area is larger than d by a factor of 60, the metasurface works well as a Bragg deflector.

The Gaussian spot at plane A, with wA≈880μm, could be used directly to approximate the Gaussian factor in p(v). Hence, to encode p(v), only the residual RAP factor needed to be realized with the metasurface. Explicitly, this factor is M(v,ϕ)=M0[1+Jm(2πρ0v)cos(mϕ)],(20)for beams with m>0 and M(v,ϕ)=M0[1−J0(2πρ0v)],(21)for m=0. Here, M0 is chosen such that the maximum does not exceed unity. Sampled 100×100 versions of the metasurfaces with spatially varying transmission efficiencies are illustrated in Fig. 6. In particular, we show the m=3, ρ0=3mm solution in Eq. (20) and the m=0, ρ0=2mm solution in Eq. (21), which are demonstrated experimentally. The modulation contrast of p(v) is much weaker in the former case than in the latter. The linear system T, which does the Fourier transform from plane A to plane B, is effectively a positive lens of focal length fH=50mm. Therefore, the spatial frequency v at plane A is connected to the real-space coordinate ρ via v=ρ/λ0fH.

The measured amplitude and phase of the DOC at plane z=z0 are shown in Fig. 7 for the two beams considered. The distributions of μ(Δx,Δy) are evaluated using a mirror-based WFI as in Ref. [27]. We show the measured phases only across a circular area of radius Δρ<75μm for m=3 and Δρ<42μm for m=0, where the absolute value of the DOC is significant and reliable phase measurements could be made.

The DOC of the m=3 solution has an azimuthal periodicity of 60° in the side lobe pattern, whereas the phase in the central region Δρ<35μm is constant. These features, as well as the sharp variations of |μ(Δx,Δy)|, are evident in the measurements even though there are some distortions due to experimental imperfections particularly in the phase measurements within the central region of Fig. 7(e). The abrupt 2π phase transitions in the central region of Fig. 7(f) are due to wrapping into the [−π,π) interval.

The propagation of beams that the DOCs in Fig. 7 generate is illustrated in Figs. 8 and 9, where the experimental and simulated results are shown in top and bottom rows, respectively. The experimental results show some speckle-like granularity that could be contributed by the finite coherence area at plane A, where the field is not strictly speaking globally incoherent as assumed in the model. In particular, the discrepancies caused by the granularity and additional side lobes in the coherence function of the m=3 RAP beam in the Fig. 7 upper panel are visible in the measured intensity profiles in Fig. 8 at z=90mm and z=140mm.

The evolution of the RAP beam with m=3 is shown in Fig. 8 at the source plane and three different propagation distances. We assume parameters R=1/4.5 and q=1/40, which gives zG=567mm and zR=14.2mm. At z=0 the beam shape is Gaussian with a beam waist w0=338μm. As the propagation distance grows the beam acquires its RAP features. The annular m=0 beam in Fig. 9, with w0=210μm, R=1/7.5, q=1/30, zG=218mm, and zR=7.3mm, begins to show its dark center at distances z∼10mm. This property is preserved throughout further propagation.

The RAP beams demonstrated here are not the only possible classes of partially coherent fields that we can generate with the experimental technique discussed in this paper. While Fourier-type kernels allow the demonstration of any Schell-model beam, the use of more general kernels opens the way to produce more general fields. For instance, metasurfaces acting as vortex phase filters with transmittance t(ρ)=t(ρ,ϕ)=exp(ilϕ), where l is an integer, can be placed at plane B to produce orbital angular momentum. Likewise, polarizing optical components can contribute to the linear transformations, allowing the realization of vector beams. Space-variant imaging systems with beam shaping or polarization-control metasurfaces (such as polarization gratings [28,29]) in the intermediate Fourier plane of a 4F imaging system are options, as are transformations involving interferometers. For example, a WFI transforms any Schell-model beam into a corresponding specular beam [30,31]. We have recently demonstrated such specular transformations using a mirror-based WFI [27].

We have presented an implementation of an experimental system for the generation of spatially partially coherent beams with genuine correlation functions, using quasi-planar metasurfaces. The method was applied to particular radially quasi-periodic, azimuthally periodic beams. However, it is applicable to the generation of any physically realizable beams, whether scalar or vectorial.