Analytical description of optical vortex beams generated by discretized vortex producing lenses

Discretized vortex-producing lenses programmed on low performance spatial light modulators have been used for the generation of optical vortices. However, the description of these vortices has been supported only by numerical simulations. In this work, a general analytical treatment (any topological charge - any discretization levels) for the propagation of a Gaussian beam through a discretized vortex-producing lens is presented. The resulting field could be expressed as a sum of Kummer beams with different amplitudes and topological charges focalized at different planes, whose characteristics of formation can be modified by tuning the parameters of the setup. Likewise, vortex lines are analyzed to understand the mechanism of formation of observed topological charges, which appear in specific planes. Conservation of the topological charge is demonstrated. Theoretical predictions are supported by experiments.


I. INTRODUCTION
An optical vortex is a singularity in the optical field. They have been widely studied because there exist beams presenting optical vortices which are also endowed with orbital angular momentum [1]. This particularity has allowed the use of this type of light structures in a wide variety of fields, such as image processing [2,3], microscopy [4], optical tweezers [5,6], optical communications [7], optical metrology [8], integrated optics [9], nanophotonics [10], among others.
Several methods have been proposed for optical vortex generation, among them, the most used are those based on fork-type holograms [11] and spiral phase plates (SPP) [12,13]. SPPs are widely used, mainly due to its simplicity and high energy efficiency. More recently, SPPs have been combined with radial dependence structures such as Fractal zone plates, Dammann zone plates, and lens-phase masks, allowing the creation of multifoci systems and the relaxation of phase modulation SPP conditions for the creation of high quality optical vortices [14][15][16] and, although these diffractive elements have been studied, the analysis is limited to computational propagation algorithms, which are not enough to offer a deep explanation of the associated phenomena.
Guo et al. [17] developed an interesting analytical procedure for discretized SPPs (DSPP), where the angular periodicity of the DSPP was expressed as a Fourier series, interpreting therefore, the discrete mask as a linear combination of ideal continuous SPPs of different topological charge, whose associated coefficients depend on the beam principal topological charge and the number of phase levels N . In this sense, the field at the observation plane can be understood as a superposition of * Corresponding author: alencina@faa.unicen.edu.ar optical vortices of different topological charges converging to the same plane.
A particular case of multifoci systems is a discretized vortex-producing lens (DVPL), which corresponds to the discretization of a phase mask formed by the superposition of a SPP and a Fresnel lens, therefore, for its analytical treatment, it is necessary to combine angular and radial effects. In this case, the diffraction properties will depend not only on the lens-phase parameters, but also on the N discretization levels. In this paper, we present a general analytical treatment (any topological charge-any discretization level) for the propagation of a Gaussian beam through a DVPL. In the procedure, the DVPL phase mask is decomposed as a combination of continuous SPPs. However, unlike the DSPP case, the coefficients of the series include an additional quadratic-phase term that explains why optical vortices of certain topological charges are observed at different planes along the optical axis. A complete analysis of the field propagation is developed to define the exact location, topological charge and energy contribution of each component. The results allow describing the intrinsic mechanism associated to the generation of the multifoci system. A study of the mechanism associated to the formation of the optical vortices is performed by studying the vortex lines and the conservation of the topological charge. Likewise, the formalism allows identifying relevant parameters of the setup to conveniently manipulate the diffraction conditions and, therefore, the energy structure along the optical axis. All theoretical predictions are supported by experimental results. These results not only provide insight on the basic physical mechanisms involved in the generation of optical vortices by DVPLs, but can also be useful in engineering applications such as: optical angular momentum channels and compact vortex inline metrological applications.

II. DVPL FOURIER EXPANSION
A vortex-producing lens phase mask is expressed as where (ρ,φ) are polar coordinates, is the topological charge (from now on it will be called principal topological charge), k is the wave number, and f F R is the diffractive-lens phase focal length. If the continuous phase profile is discretized in N phase steps and each step has a constant phase increment of φ = 2π N , the discretized complex transmittance of CON (ρ,φ) will correspond to a DVPL phase mask, and is mathematically expressed as Floor[x] being the function which takes the nearest integer smaller than or equal to x. Since the transmittance of Eq. (2) is a periodic function of the azimuthal angle φ with a period of 2π , it can be expanded into a Fourier series, where the coefficients t m (ρ) depend on the radial variable ρ and are obtained from Solving the integral of Eq. (4) (see Appendix A) gives the following coefficients: with sinc(x) = sin(x)/x. Then, the phase mask of Eq. (2) can be interpreted as a superposition of SPPs with topological charges The weights t m have two phase factors: one of them is constant, whereas the other one has a radial quadratic dependency analogous to a diffractive-lens phase focal length m f F R . The presence of this quadratic phase means that each SPP focuses at a different plane, depending on the number of levels N used and the diffractive-lens phase focal length. The number of levels N also plays a key role in the SPPs superposition, making the principal topological charge more dominant as N increases. On the other hand, notice that the integer number j refers to the different terms of the linear expansion without mentioning any specific topological charge. Thus, we call each term the j th order of the expansion. Figure 1 shows the weights for some orders of the expansion, for four different discretization levels. From the figure, some aspects can be emphasized: there is no weight symmetry with respect to the zeroth order, only for two-level discretization, orders −1 and 0, will have the same weight in the superposition and the same topological charge with different sign and, for four or more discretization levels, the higher orders are of negligible energy.
It is interesting to compare Eq. (5) with a very similar one that appears in Ref. [17] for DSPPs (it corresponds to Eq. (5) but with ρ = 0). Because in Ref. [17] the terms of the expansion are not accompanied by a quadratic phase factor, all the vortices of the linear combination are focused at the same z plane, whereas, for DVPLs, they are focused at different planes, as it is shown in the next section.

III. FIELD PROPAGATION AFTER THE DVPL
In this section, the analytical expression for the field after crossing the DVPL is developed. The input plane is composed of the illumination beam A(ρ) and the DVPL phase mask (ρ,φ), and can be written as This field is propagated a distance f in free-space towards a thin lens of focal distance f . Then, it is propagated again an extra distance z 0 . The complete analytical procedure is detailed in Appendix B. Here we report the result of the output field: where (r,θ ) are the output-plane coordinates and u m (r,θ ; z 0 ) has the form Equation (8) where m = , ± N , ± 2N , . . . . For the case of the principal topological charge, m = , z m corresponds to the focus of the optical system formed by the physical lens and the Fresnel lens, while the other topological charges are focused at distances equal to multiples of N f 2 f F R , with respect to the principal topological charge. This can be easily seen replacing the possible values of m in Eq. (10) by the jth orders, Thus, for a given optical setup, it is possible to modify the distance between the focused vortex beams by changing the number N of levels or by changing the focal distance of the lenses. Figure 2 presents schematically these results. There, the position and charge of each order is displayed as a function of the discretization level N . Note that the vertical axis is centered at the focus of the optical system f − f 2 /f F R and the scale is in units of f 2 /f F R . From the figure, it is apparent that at the zeroth order, the principal vortex is obtained, irrespective the value of the discretization level N . As N increases higher orders are more distant. This behavior, combined with the fact that higher orders carry less energy (see Fig. 1), allows us to understand why high quality vortex can be obtained with DVPL of lower N . Moreover, as the discretization increases the topological charge of the higher order increases, then increasing the dark disk contribution at the optical axis of the system. This situation, combined with the defocusing of each vortex away of the plane of its order, also contributes to improve the quality of the generated vortex. The latter will be evident in the next section.
An interesting point to highlight is that, because at each z m a different topological charge is focused, it seems that the topological charge is not conserved under propagation. However, as is discussed later, this is not the case.

IV. GAUSSIAN BEAM INPUT
For the special case of a Gaussian beam, with beam waist ω 0 , and amplitude by solving the Hankel transform [18, Sec. 8.6, p. 29, Eq. (9)], the optical field at a distance z 0 is where I n (x) is the modified Bessel function, m = + jN with j = 0, ± 1, ± 2,..., and Each term of Eq. (12) resembles a Kummer beam [8], with topological charge m, focused at z m . Note that the Gaussian factor before the Bessel functions is responsible for focusing each term at its corresponding z m . From this view, it is clear that the beam waist ω 0 plays also an important role in defining the quality of the vortex at the zeroth order because it dominates the defocusing of the ±1th orders through b m . It is interesting also to note that, although each term is focused at a corresponding z m , the other terms are present and act as a background which could degrade the quality of the vortex by splitting the charge m into m unitary charges [19]. This point is addressed next. In Fig. 3, we present simulations of the resulting intensity and phase for the case of a DVPL illuminated by a Gaussian beam, for different discretization levels, orders and principal topological charges. The parameters employed in the simulation are λ = 532 nm, ω 0 = 5 mm, f = 20 cm, and f F R = 1.6 m. The first thing to notice from Fig. 3 is that at a scale of the order of the doughnut-shaped intensity, the topological charge is as expected (see Fig. 2). Yet, a closer look towards the optical axis shows that, except for the zero order (and the ±1th order for N = 2), in all the other cases there is a splitting of the topological charge m into a bundle of vortices with lower topological charge. However, the principal charge prevails at the optical axis. This result suggests that the principal charge is robust against the background field corresponding to the remaining terms, telling us that the topological charge of the principal vortex is conserved in the vicinity of the optical axis (see the next section for a more detailed analysis).
The case of N = 2 is out of the previous analysis because the topological-charge at the optical axis reverses its sign for negative orders. These reversed charges are also robust against the background field corresponding to the remaining terms (not shown). This sign change could lead us to think that topological charge is not conserved under propagation. However, as it is well known, the dynamical inversion of the topological charge of an optical vortex occurs for noncanonical vortices in presence of astigmatic transformations [20,21]. It was verified (not shown) that in between the planes of the ±1th and zeroth order there is a place where a Freund's critical foliation appears. It is expected, as in the astigmatic case, that another inversion plane exists where the topological charge sign reverts again but its finding is out of the focus of the present paper. In these critical foliations vortices go off to (or come in from) far (but finite) distances from the optical axis of the system (in a scale of the beam waist). Nevertheless, owing to its solenoidal nature, they must close on themselves. This fact is explained in detail in the next section where it is shown that the overall (not local) topological charge is conserved (as expected) under propagation.

A. Vortex lines
As it is well known, optical vortices are lines in space [22,23]. In analogy with fluid mechanics [24], the direction of these lines is described by the vorticity = 1 2 ∇ × j, j being the current associated with the field, i.e., the Poyting vector [25,26]. Since vorticity is defined by a curl, it is solenoidal, i.e., it has null divergence, which implies that vortex lines are closed curves. If a finite volume is considered, vortex lines close on themselves or the number of lines entering the volume is the same as the number of lines that leave it. Figure 4 shows the numerically obtained oriented trajectories [27] of vortices for a beam generated with a DVPL with N = 3 and = 3. The same parameters as in Fig. 3 are used. Three views of the same situation are depicted. In each view, a plane is highlighted. From the figure, the presence of a principal vortex trajectory at the optical axis (winding number 3) and some vortex trajectories that concentrate around the ±1th order planes (winding number 1) are evident. It is clear that all vortex trajectories that enter the displayed volume also come out of it, in accordance with the definition of vorticity. For the j = −1 order, nine vortex trajectories concentrate towards the principal vortex line. Owing the directions and winding number of these trajectories, this gives rise to the m = −6 vortex observed at this plane, as predicted by Eq. (11) (see Fig. 2) and shown in Fig. 3. For j = 0, only the principal vortex line is present, corresponding to an = 3 vortex. In the case of j = 1, again nine vortex trajectories concentrate close to the principal vortex line, but this time their directions coincide, thus contributing to the topological charge m = 12 observed at this plane. The previous analysis encompasses the concept of topological charge conservation from a geometrical point of view: If the only singularity introduced to the field is the principal topological charge and this generates a principal vortex trajectory with its corresponding winding number then, because of the solenoidal property, the other trajectories must belong to closed curves so that the total topological charge evaluated over the entire field equals the principal one. This implies that lines other than the principal line can close far away from the optical axis but at finite distances. Otherwise, some charge "would leave" the beam and the total charge would change, which contradicts the topological charge conservation. This statement is supported in the next subsection.

B. Field at large r
The topological charge of a field is defined as the times of 2π its phase φ varies in a closed path, i.e., by calculating the Burgers vector as dφ = 2π . As it is observed from Figs. 3 and 4 at z m the corresponding order dominates and the field strongly resembles a vortex field of charge m. This fact could suggest that the topological charge changes as the field propagates. However, as discussed in the previous section, the vortex trajectories other than the principal ones, are closed curves whose segments group together around the focus of an order. If this is what happens, the topological charge calculated with a path enclosing all loops must return the value of the topological charge imprinted in the DVPL. To ensure all loops are included the path must be taken at a large radius r. By considering the first and second term of the asymptotic expansion of the Bessel functions [28, Sec 10.40.1, p. 255], for large r, Eq. (12) can be written as Note that, despite an overall phase term depending on z 0 , the field at large r does not change on propagation, i.e., does not depends on z 0 . This is a first evidence that any property determined from this expression is conserved as the field propagates in space. After some algebraic manipulation, the field in the previous expression can be written as From this expression, the phase φ can be easily obtained. However, owing to the periodicity of the field in the azimuthal coordinate, care must be taken to unwrap the phase before any calculation. Then, having unwrapped the phase, the Burgers 063813-5 FIG. 5. Unwrapped phase variation for large r as a function of θ, for different values of discretization level and principal topological charge. At the top-left corner of each subfigure, the discretization level is indicated. The corresponding principal topological charge of each curve is indicated in the right side of each subfigure. Because all curves were shifted to satisfy φ(0) = 0, the unwrapped phase value for θ = 2π is a direct reading of the Burgers vector for that case.
vector can be calculated. If we select a circumference centered at the optical axis as the path to evaluate the integral, the Burgers vector reduces to 2π 0 1 r where φ uw is the unwrapped phase. From this result, it is evident that by only evaluating the phase at θ = 2π and θ = 0 the Burgers vector can be calculated. Then, the Burgers vector value can be inferred by simple inspection of the functional behavior of the phase as a function of θ .
Unwrapping is a numerical operation, so the Burgers vector is calculated for particular cases of interest. Figure 5 shows the unwrapped phase for different values of discretization levels and principal charges. The discretization level is shown at the top-left corner of each subfigure. Also, along the right-side axis the principal topological charge corresponding to the contiguous phase curve is indicated. All curves were shifted to fulfill φ(0) = 0. Then, by reading the phase value at θ = 2π , the resulting Burgers vector is obtained for each case. Thus, by dividing these values by 2π , the total topological charge of the field for a given principal topological charge and discretization level is obtained. From the figures, it is apparent that phase variation is smooth for = 1,3 whereas there are jumps of π for the other cases. The number of jumps is given by the principal topological charge times the discretized levels. Regardless of whether the curves have jumps or not, it can be observed that given a principal topological charge the phase value for θ = 2π is the same for all discretized levels. This result tells us that the Burgers vector does not depend on the discretization level (besides that, it did not depend on z 0 , the propagation coordinate). By simple inspection, it is evident that the value of the unwrapped phase for θ = 2π , i.e., the Burgers vector, divided by 2π , coincides with the principal topological charge. Therefore, the calculated total topological charge matches the principal topological charge imprinted to the DVPL. These results are important evidence for the conservation of the topological charge under propagation, irrespectively of any principal charge and discretization level considered.

VI. EXPERIMENTAL VERIFICATION
To experimentally validate the results obtained in the previous sections, a DVPL was built. The mounted experimental setup is presented in Fig. 6. A laser of wavelength λ = 532 nm is filtered and collimated by a spatial filter SF and a lens L.1. This beam passes through a Mach-Zehnder interferometer formed by beam splitters BS.1, BS.2, and mirrors M.1 and M.2. One arm of the interferometer is used as a reference wave for phase retrieval, while, on the other, the DVPL was built up through a phase-mostly spatial light modulator, composed by a polarizer P.1, a quarter wave-plate QW P , a Holoeye LC2002 twisted-nematic liquid crystal display (TN-LCD), and an analyzer P.2. A lens L.2 is placed at a distance equal to the focal length f 2 = 20 cm from the TN-LCD. Finally, by using a 40x microscope objective, the intensity and the interference patterns of the complex field are registered with a CMOS camera (DCC1545M Thorlabs). Both, microscope objective and camera, could be moved longitudinally, to explore the z 0 dependence of the beam. A shutter is employed to block the reference beam when intensity images are recorded. The phase was recovered by using a five-step phase-shifting technique by encoding the phase delays directly in the LC-SLM [29].
To build up the DVPLs, the Holoeye LC2002 TN-LCD was characterized in a phase-mostly configuration using the procedure developed by Amaya et al. [30]. A maximum phasemodulation close to 1.5π , with 5% of coupled amplitude was obtained. This performance allows us to implement DVPLs  Fig. 3. For example, at the zeroth order, central column (j = 0 view), the vortices obtained were registered at the focus of the optical system f − (f 2 /f F R ) and it can be observed that their topological charges coincide with the ones programed in the DVPL independently of the number of the discretization levels N and the value of the topological charge .

VII. CONCLUSIONS
The beams resulting of the system proposed by Rueda et al [16], i.e., by using DVPLs, can be considered as a superposition of Kummer beams of different topological charges, each of which is focused at a different plane. The additional quadratic term in Eq. (5) generates a longitudinal separation of the focal plane of each Kummer. Simple expressions relating the relevant parameters of the system, the position of the foci, and the appearing topological charges are derived. This allows us to understand the presence of a high-quality vortex at z reported by Rueda et al.: the most important term of Eq. (12) is focused at that plane, while all other components are out of focus (their contribution to the intensity and the total phase of the beam are negligible at that plane). The presence of the additional foci and the topological charge of the phase singularity present at each of them was verified experimentally. This complete understanding of the system helps us see that the role of the physical lens in the separation of the different components of the beam is irrelevant (although it ensures they be equidistant), and the quality of the system could be improved by omitting it. An interesting point is that, although this kind of system generates optical vortices of different topological charges in different planes along the propagation axis, we were able to obtain an expression that allows the verification of the topological charge conservation at any particular point, any charge, and any discretization levels. This work also highlights the importance of the beams resulting from this system, given that they possess a complex topological structure. Finally, the understanding of vortex beam generation could be of importance in different applications, such as optical angular momentum channels and compact vortex inline metrological applications.
where the coefficients t m (ρ) depend on the radial variable ρ and are given by t m (ρ) = 1 2π 2π 0 (ρ,φ) exp (−imφ)dφ By making the change of variableφ = 1 φ ( φ − kρ 2 2f F R ), replacing φ by 2π N , and taking out the terms that do not depend on the integration variable, Eq. (A2) can be rewritten as Calling the term in square brackets c m , for 0 p < N, p ∈ N, Floor(φ) takes constant values in the interval p φ < p + 1, thus Defining κ = 2πm N , each of the terms in the sum of Eq. (A4) can be expressed in the form: which is the Fourier transform of the function rect(φ − 2p+1 p ), describing a rectangle of unit height. The solution is equal to exp − i 2p + 1 p κ sinc(κ/2).
Replacing this result in Eq. (A4), taking into account that sinc(z) = sin(z)/z, and using the definition of κ, c m can be written as If m− N is an integer n, all the terms in the sum of Eq. (A7) are equal to exp(2πpn) = 1, and so the sum is equal to N. In any other case, the sum is null. Finally, the coefficients of the expansion can be written as in Eq. (5).