## Abstract

We propose a measurement protocol and parameter estimation algorithm to recover the powers and relative phases of each of the vector modes present at the output of an optical fiber that supports the HE_{11}, TE_{01}, HE_{21}, and TM_{01} modes. The measurements consist of polarization filtered near-field intensity images that are easily implemented with standard off-shelf components. We demonstrate the accuracy of the method on both simulated and measured data from a recently demonstrated fiber that supports stable orbital angular momentum states.

© 2013 Optical Society of America

## 1. Introduction

Recent years have seen a surge in interest in few-mode optical fibers for many applications, including mode-division optical multiplexing [1] and generation and transmission of optical vortices [2–4]. In contrast with standard graded-index multimode fiber, applications of such fibers depend on properties of individual modes, as opposed to averages over many modes. Therefore, standard mode power distribution measurements, such as those based on encircled flux [5], do not apply, and new methods are required.

In this paper, we demonstrate such a method that is applicable to six-mode fibers that support the pair of fundamental HE_{11} modes, along with the higher order TM_{01}, TE_{01}, and pair of HE_{21} modes. We propose a measurement protocol and statistical parameter estimation algorithm that separates both the powers and relative phases of all six modes. The measurements consist of polarization-filtered near-field intensity (NFI) images, which are easily implemented with standard off-shelf components. The complex mode amplitude recovery is, in general, only unique up to a two-fold ambiguity in mode labeling, but under most situations of interest a small amount of prior information allows a unique recovery.

We demonstrate our technique with simulated and measured data from a recently demonstrated six-mode “vortex” fiber [6] that supports the stable propagation of orbital angular momentum (OAM) states. The index profile, radial mode fields, and guided mode spectrum of this fiber is presented in Fig. 1. The feature that allows the stable propagation of OAM in superpositions of the two HE_{21} modes, as well as that of the azimuthally- and radially-polarized TE_{01} and TM_{01} modes, is the enhancement of the effective index splitting among these modes due to the vector term in the wave equation, as first exploited for microbend grating applications [7, 8]. Information about both the magnitude and phase relationships between the modes is valuable for studying superpositions of the orbital angular momentum states [4]. One of the contributions of the current paper is to generalize the analysis method of [4] to recover phase information and to relax the dominant mode assumption made in that paper.

A variety of other methods for modal content analysis of few-mode fibers have been developed in recent years. Using only the geometry of NFI patterns, the recovery of a subset of the complex mode amplitude parameters of few-mode fibers is described in [9]. Modal powers and phases in the scalar approximation are recovered from far-field intensity patterns in [10], although the modes within a degenerate group are not distinguished. A similar method, using NFI patterns, was used to retrieve partial information about the coherent superposition of vector modes in [11]. A full coherent vector reconstruction is given in [12] by use of a specially designed holographic optical element and corresponding reconstruction algorithm. In [13, 14] the mode powers of the scalar modes are recovered by a self-interference method that does not require a priori knowledge of the mode structure. Low-coherence interferometry with an external reference beam, as done in [15], enables a similar, but more versatile, measurement. The modal content of a six-mode fiber, including the vortex content, was analyzed via a different interferometric measurement employing a reference beam in [16]. As part of a MIMO demonstration, the coherent mode content of a six-mode fiber was extracted with a specially designed free-space multiplexer and coherent detection in [1].

The primary features of the method we present that distinguish it from existing methods are its simplicity of implementation combined with the wealth of information it provides. Indeed, in most situations of interest, complete information about the coherent superposition of modes present at the fiber output are recovered from NFI images alone, as filtered through standard quarter wave plates and polarizers. In particular, no special-purpose optical elements nor interferometry with reference beams are required. We discuss our recovery algorithm in Section 2, and show that it retrieves all possible information given the postulated set of measurements. Experimental details are presented in Section 3, and examples on simulated and measured data in Section 4.

## 2. Algorithms

#### 2.1. Mode field representations

Under monochromatic illumination, the transverse electric field **E*** _{t}*(

*r*,

*θ*,

*t*) at the fiber end-face may be represented as a complex superposition of the six guided mode fields

**e**

*(*

_{i}*r*,

*θ*)

*i*. In the full vectorial treatment of the Maxwell equations, the six fields

**e**

*represent $\mathbf{H}{\mathbf{E}}_{11}^{\left(e\right)}$, $\mathbf{H}{\mathbf{E}}_{11}^{\left(o\right)}$,*

_{i}**TE**

_{01},

**TM**

_{01}, $\mathbf{H}{\mathbf{E}}_{21}^{\left(e\right)}$, and $\mathbf{H}{\mathbf{E}}_{21}^{\left(o\right)}$. The representation (2.1) also applies to unchirped pulses as long as the coherent relationship between the modes is maintained, which will be true if the temporal pulse width is larger than the group delay differences of the fiber modes. For a weakly guided fiber, it is a very good approximation to express these vector mode fields as linear combinations of the LP

_{01}and LP

_{11}degenerate sets of solutions to the scalar wave equation [17]. Denoting the radial wave functions of the LP states by

*F*(

_{ℓm}*r*), we have a basis set that we refer to as the

*Vector*set

*Vortex*basis set [18–20] is given by

*ħ*/photon each of orbital and spin angular momentum, which are aligned, for a total angular momentum of 2

*ħ*/photon. The ${\mathbf{V}}_{01}^{\left(T\pm \right)}$ elements are not true fiber modes, as the TE and TM modes are separated in effective index. However, they are valid elements of a basis that may be used to expand the transverse electric fields at a given point in the fiber, which is the only way in which we will use them. Also, the ${\mathbf{V}}_{11}^{\left(HE\pm \right)}$ elements are not vortex states, as they do not vanish on the axis; we label them with the symbol

**V**for consistency of notation. We note that, unlike the Vector basis, the polarization of each element of the Vortex set does not depend on position; all are uniformly circularly polarized. The same is true of the

*LP*basis set, given by

#### 2.2. Vortex mode power estimation

For reasons that will soon be apparent, the Vortex basis is most convenient for analysis. Referring to the expansion (2.1), we define

The intensity of the transverse electric field at the fiber end-face (2.1), after passing through a positive helicity polarizer *P*_{+}, is given by

*r*,

_{i}*θ*}. By linear regression on the two circularly polarized NFI’s, we may estimate the ten coefficients

_{i}_{35}and Ξ

_{46}, of which both the real and imaginary parts may be independently recovered. We observe that the fundamental mode powers |

*γ*

_{1}|

^{2}and |

*γ*

_{2}|

^{2}are directly observed, while those of the LP

_{11}states must be calculated from |Ξ

_{3+5}|, |Ξ

_{4+6}|, |Ξ

_{35}|, and |Ξ

_{46}|. From these four quantities we may compute

*γ*

_{3}|, |

*γ*

_{5}|, and similarly for the |

*γ*

_{4}|, |

*γ*

_{6}| pair. In many cases of interest, this ambiguity may be removed with prior knowledge of the mode powers. For instance, if most power is localized in one of the transverse modes, we expect that |

*γ*

_{5}| ≫ |

*γ*

_{3}| and |

*γ*

_{4}| ≫ |

*γ*

_{6}|. If most power is localized in one of the HE

_{21}modes, we expect that |

*γ*

_{3}| ≫ |

*γ*

_{5}| and |

*γ*

_{4}| ≫ |

*γ*

_{6}|. If the power is localized in ${\mathbf{V}}_{21}^{\left(HE+\right)}$, then again |

*γ*

_{3}| ≫ |

*γ*

_{5}|. When the launch is chosen to maximize ${\mathbf{V}}_{21}^{\left(HE-\right)}$, and the

*P*

_{+}projection is taken (or vice versa), then we encounter the ambiguity of how to assign the values of |

*γ*

_{4}| and

*γ*

_{6}. In this case, we remark that all six mode

*powers*are estimated (in the vortex basis); but the

*complex amplitudes*are not, and the powers of two of the six states are ambiguous. In most cases of practical interest, this ambiguity is resolved with the full complex mode amplitude estimation, as discussed in Section 2.3.

#### 2.3. Complex amplitude estimation

We now proceed to estimate the relative phases *ϕ _{i}* in (2.6). Since there will be an ambiguous overall phase, we define

*ϕ*

_{1}= 0. We denote the relative phases

_{35}and and Ξ

_{46}we may read off

*δ*

_{35}and

*δ*

_{46}. To go further we need more information.

The next step is to project the field pattern onto two orthogonal LP subspaces via orthogonal linear polarizers, which yields the NFI images

In order to recover the phases in (2.6) we form

*δ*

_{43}, which we find as the solution to

*γ*

_{3},

*γ*

_{4},

*γ*

_{5}, and

*γ*

_{6}, but not their absolute phases (i.e. relative to

*γ*

_{1}), nor that of

*γ*

_{2}. We obtain these from Ξ

_{13+5}, Ξ

_{13−5}, Ξ

_{26+4}and Ξ

_{26−4}. From the definitions,

*ϕ*

_{3}, and, immediately,

*ϕ*

_{4},

*ϕ*

_{5}, and

*ϕ*

_{6}. Similarly,

*δ*

_{24}, and therefore

*ϕ*

_{2}, as

*ϕ*

_{4}is already known.

Inspection of (2.9) and (2.14) reveals that all of the regression coefficients are invariant under the mapping ** γ** →

**, with**

*γ̃**γ*

_{3}| > |

*γ*

_{5}| or |

*γ*

_{4}| > |

*γ*

_{6}|, and which pair of modes satisfies the inequality. Fortunately, in most cases of practical interest this ancillary information is available.

Suppose we are given that |*γ*_{3}| > |*γ*_{5}| (the other case works analogously). Then the assignment of
${\mathrm{\Gamma}}_{46}^{\left(\text{min}\right)}$ and
${\mathrm{\Gamma}}_{46}^{\left(\text{max}\right)}$ to |*γ*_{4}| and |*γ*_{6}|, as done at the end of Section 2.2, is ambiguous. We may proceed through the regression analysis of the current section to estimate the phases *ϕ _{i}* under both hypotheses |

*γ*

_{4}| ≷ |

*γ*

_{6}|. In light of (2.19), because the values of |

*γ*

_{3}| and |

*γ*

_{5}| are, by assumption, correct, only the correct assignment of powers to |

*γ*

_{4}| and |

*γ*

_{6}| will result in values of all regression coefficients that are consistent with those observed. Therefore, a consistency check may be performed by comparing the estimated value of, say, Ξ̂

_{36+54}= |

*γ̂*

_{3}+

*γ̂*

_{6}|

^{2}+ |

*γ̂*

_{5}+

*γ̂*

_{4}|

^{2}with the observed value. The hypothesized assignment of powers to |

*γ*

_{4}| and |

*γ*

_{6}| that yields the best match should be chosen. By including more regression coefficients in the consistency metric a more robust choice can be made. We remark that only the coefficients ${\mathrm{\Xi}}_{*}^{\left(x\right)}$ from the

**x**̂ linearly polarized image were employed to uniquely reconstruct the complex mode amplitudes. The

**ŷ**polarized image is, in principle, unnecessary, but its inclusion in the non-linear regression discussed in Section 3 will result in better estimates.

Finally, we comment on the conditions under which the ambiguity (2.19) can be avoided. Typically, we are interested in putting most power into either the HE_{21} or the TE/TM_{01} pairs, and it is exactly these situations in which the ambiguity does not arise. If it is one of the transverse modes that is preferentially excited, then either |*γ*_{5}| > |*γ*_{3}| or |*γ*_{6}| > |*γ*_{4}|. The propagation constants of the transverse modes are split both from one another and from the HE_{21} pair, so there is no ambiguity as to which mode is present. The HE_{21} pair, by contrast, is degenerate, so mode coupling may, in principle, spoil the state that is prepared at the input end of the fiber. However, when most power is contained in the degenerate pair, then one of the vortex states
${\mathbf{V}}_{21}^{\left(HE+\right)}$ or
${\mathbf{V}}_{21}^{\left(HE-\right)}$ must contain at least nearly half of the total power, and we may determine which one is dominant by observing the circularly polarized NFI images. If the dominant state is
${\mathbf{V}}_{21}^{\left(HE+\right)}$, say, then we have |*γ*_{3}| > |*γ*_{5}|.

## 3. Implementation

Two experimental setups that have been used to obtain the measurements necessary for the recovery algorithms are presented in Fig. 3. In both, the vortex states were excited by the methods described in [4]. The output of a narrowband tunable continuous-wave (CW) laser is coupled into a standard single-mode fiber (SMF), with its polarization state determined by a polarization controller before coupling into the vortex fiber. A microbend grating converts power in the fundamental mode into a vortex state that is determined by the input polarization state. In the experiments we discuss, the conversion efficiency is better than 25 dB at the resonant wavelength of the grating. The output of the vortex fiber is imaged onto a camera. The appropriate polarization projections may be obtained, in series, with appropriately oriented polarizers and quarter-wave plates. In long fibers, temperature fluctuations in the vortex fiber can change the phase relationships between the modes over the time it takes to obtain the images in series. To avoid this problem, the scheme shown in Fig. 3(b) was implemented to capture both linear and circular polarization projections simultaneously. The results presented below use data from this latter setup.

Once the requisite images are obtained, the analyses of Section 2 can, in principle, be implemented directly via standard linear regression algorithms. However, in practice there are additional steps that must be taken to obtain reliable results. The first issue involves the radial mode fields *F*_{01}(*r*) and *F*_{11}(*r*). These can, in principle, be obtained from a mode solver applied to a measured radial index profile of the fiber. In practice, we found that these calculated mode fields did not agree with the NFI measurements, rendering the regression results unstable. We attribute the differences to a combination of inaccuracies in the fiber profile measurement, uncertainties in translating the measured index profile to wavelength used in the experiments, aberrations and vignetting in the optical train in front of the camera, scattering at the fiber end-face; and camera calibration error. We did not attempt to unravel the relative importance of these possible causes, as the problem is easily circumvented.

There are two possible resolutions. The approach taken in [4] was to confine attention to a narrow annulus of pixels *r* ≃ *r*_{0}, for which *F*_{01}(*r*_{0}) ≃ *F*_{11}(*r*_{0}). The information at other radii is not used. When it is assumed that most power in the fiber resides in the HE_{21} mode pair, then it is possible to uniquely recover all mode powers (although not their relative phases). In order to relax the assumption of a dominant HE_{21} mode, as we do in this paper, we must also consider the power in pixels near the fiber center *r* ≃ 0, and have knowledge of the ratio *F*_{01}(*r*_{0})/*F*_{01}(0). This ratio is easily measured by imaging the vortex fiber output with the microbend grating removed. In this paper we take a slightly different approach, in that we use the entire image but with measured, as opposed to calculated, radial mode fields. We therefore require NFI images of one of the LP_{11} states without the fundamental mode present. Since the microbend grating can transfer power from the fundamental to the LP_{11} mode with better than 25 dB efficiency, such images are easy to obtain. With the two images in hand, the mode field center and radial mode field are obtained with a non-linear fit.

Once the preprocessing steps are completed, the complex mode amplitude recovery algorithm consists of two additional steps. The first is the linear regression analysis of Section 2. This analysis requires only three NFI images: those that project onto the two circularly polarized subspaces, and one of the two linearly polarized images. The results are estimates of all mode powers and relative phases, up to the two-fold ambiguity that was previously discussed. The second step is to use the results of the linear regression as a starting point for a non-linear regression to refine the complex mode amplitude estimates. This non-linear regression uses all four available images: both, instead of just one, linear projections, along with the two circularly polarized projections. The non-linear programming problem is

*p*∈ {+, −,

**x**̂,

**ŷ**} are the polarization states on which we project,

**x**is the pixel index,

*I*are the measured polarization filtered images,

_{p}*P*are the polarization projections, and

_{p}**E**

*is the model transverse electric field given in Eq. (2.1), with basis set given by (2.5). In addition to providing a more efficient estimator, and using the additional data contained in the extra linearly polarized projection, an additional advantage of the non-linear technique is that it enforces physically meaningful values of the parameters, which need not be the case in the linear regression approach. We have observed that the results of the non-linear regression can, rarely, violate the assumptions used to disambiguate the linear parameter estimates with respect to the invariance (2.19). The disambiguation is therefore repeated after the non-linear regression.*

_{t}Finally, we comment on sources of error in the complex mode amplitude recovery. In Section 4.1 we show results for simulated images that include sources of noise similar to those in the experiments. With a single camera shot for each image, we find that the mode powers and phases are estimated to within a few percent at worst. There are also several sources of bias that must be minimized. The beam splitters employed may exhibit polarization dependent loss and birefringence, which will tend to underestimate mode purity. The preparation of the fiber end-face can distort the beams; while we measure the radial wavefunction, we assume the angular dependence is that of a flat interface. Camera miscalibration will bias the results, and dark counts can impose a floor on sensitivity.

## 4. Results

#### 4.1. Simulations

In this section we present the results of the analysis of Sections 2.2 and 2.3 applied to two simulated data sets, one in which most power is concentrated in the the positive helicity vortex state
${\mathbf{V}}_{21}^{\left(HE+\right)}$, the other in which the **TM**_{01} state dominates. In both cases the “nuisance” mode powers are more than 20 dB lower than the dominant mode. In each case shot noise, dark noise, and quantization noise was simulated, with a twelve-bit camera with mean dark count of 250 assumed. The four polarization filtered images were assumed to have equal exposure time, set so that the brightest pixel across all images would have a count of 3200 in the absence of dark current. These assumptions correspond to the experimental conditions in Section 4.2. The images, including all noise sources, for the vortex mode are presented in Fig. 4, and those for the **TM**_{01} are shown in Fig. 5.

The linear regression analysis was performed in each case for 100 noise realizations; the resulting mode power and phase estimates are presented in Figs. 6 and 7. We observe that the reconstructions are quite stable in both cases, with some variability present in a few of the parameters.

#### 4.2. Experimental

The measurement setup of Fig. 3 was applied to a short length (3 m) of the fiber shown in Fig. 1. In the top row of Fig. 8 we show the three measured polarization-filtered images necessary to carry out the analysis, applied to a fiber in which the ${\mathbf{V}}_{21}^{\left(HE-\right)}$ was preferentially excited. The fitted values obtained after the complex mode amplitude recovery are shown in the bottom row of Fig. 8. The close agreement between the two rows is evidence for the accuracy of the reconstructed amplitudes, since, as described in Section 2, the model is identifiable under these launch conditions.

Figure 9(a) shows the powers estimated in each of the three pairs of modes HE_{11}, HE_{21}, and TE/TM_{01}, at a number of wavelengths surrounding the microbend grating resonance. A quantitative measure of the accuracy of the reconstruction can be obtained by looking at the power in the fundamental mode. An independent measurement of this quantity was obtained by coupling the vortex fiber into a single mode fiber, thereby stripping the power in the higher order vortex modes. The comparison beween this measurement and that estimated from the NFI images is shown in Fig. 9(b), with excellent agreement. During this measurement the two OAM states that comprise HE_{21} were excited with roughly equal power and stable relative phase, as shown in Fig. 10. We note that the variance of the phase estimate increases dramatically at the extreme ends of the wavelength range, and some apparent errors in mode assignment are visible in the power estimates in this region as well. Such estimation error is expected in these regions, as the assumptions that allow us to resolve the two-fold ambiguity (2.19) do not hold as the OAM-state powers tend to zero.

## 5. Conclusion

We proposed a novel technique for recovering the mode powers and relative phases at the output of any six-mode optical fiber that supports the HE_{11}, TE_{01}, HE_{21}, and TM_{01} modes. The method relies only on polarization filtered NFI measurements that are easily implemented with standard off-the-shelf components. The reconstruction algorithm is based on a statistical regression analysis and is unique up to a two-fold ambiguity that is resolvable in most cases of practical interest. We demonstrated the method on simulated and measured data from a recently developed fiber that supports stable OAM propagation. The experiments we have described employed a short length of fiber, but the method is equally applicable to longer fiber lengths, and such studies will be presented elsewhere. The analysis of fibers supporting larger numbers of modes is complicated by the high dimensionality of the corresponding vector spaces, but we believe that generalizations of our method can give useful information about these fibers as well. In addition to probing the properties of few-mode fibers themselves, this method is also valuable for characterizing launching devices.

## Acknowledgments

This work is sponsored by the Defense Advanced Research Projects Agency Information in a Photon (InPho) program under Air Force Contract FA8721-05-C-0002 and grant HR0011-11-1-0004. Opinions, interpretations, conclusions, and recommendations are those of the authors and not necessarily endorsed by the United States Government.

## References and links

**1. **R. Ryf, S. Randel, A. H. Gnauck, C. Bolle, A. Sierra, S. Mumtaz, M. Esmaeelpour, E. C. Burrows, R.-J. Essiambre, P. J. Winzer, D. W. Peckham, A. H. McCurdy, and R. Lingle, “Mode-Division Multiplexing Over 96 km of Few-Mode Fiber Using Coherent 6 × 6 MIMO Processing,” J. Lightwave Technol. **30**(4), 521–531 (2012). [CrossRef]

**2. **C. N. Alexeyev, A. V. Volyar, and M. A. Yavorsky, “Vortex-preserving weakly guiding anisotropic twisted fibres,” J. Opt. A-Pure Appl. Op. **6**(5), S162–S165 (2004). [CrossRef]

**3. **N. K. Viswanathan and V. V. G. Inavalli, “Generation of optical vector beams using a two-mode fiber,” Opt. Lett. **34**(8), 1189–1191 (2009). [CrossRef] [PubMed]

**4. **N. Bozinovic, S. Golowich, P. Kristensen, and S. Ramachandran, “Control of orbital angular momentum of light with optical fibers,” Opt. Lett. **37**(13), 2451–2453 (2012). [CrossRef] [PubMed]

**5. **“Launched Power Distribution Measurement Procedure for Graded-Index Multimode Fiber Transmitters (FOTP-203),” ANSI/TIA/EIA-455-203-2001, (2001).

**6. **S. Ramachandran, P. Kristensen, and M. F. Yan, “Generation and propagation of radially polarized beams in optical fibers,” Opt. Lett. **34**(16), 2525–2527 (2009). [CrossRef] [PubMed]

**7. **S. Golowich and S. Ramachandran, “Impact of fiber design on polarization dependence in microbend gratings,” Opt. Express **13**(18), 6870–6877 (2005). [CrossRef] [PubMed]

**8. **S. Ramachandran, S. Golowich, M. F. Yan, E. Monberg, F. V. Dimarcello, J. Fleming, S. Ghalmi, and P. Wisk, “Lifting polarization degeneracy of modes by fiber design: a platform for polarization-insensitive microbend fiber gratings,” Opt. Lett. **30**(21), 2864–2866 (2005). [CrossRef] [PubMed]

**9. **A. Volyar, “Algebra and geometry of optical vortices: measuring of a mode spectrum in optical fibres,” in First International Conference on Advanced Optoelectronics and Lasers, 2003. Proceedings of CAOL 2003 , **1**, 30–35 (2003).

**10. **M. Skorobogatiy, C. Anastassiou, S. G. Johnson, O. Weisberg, T. D. Engeness, S. A. Jacobs, R. U. Ahmad, and Y. Fink, “Quantitative characterization of higher-order mode converters in weakly multimoded fibers,” Opt. Express **11**(22), 2838–2847 (2003). [CrossRef] [PubMed]

**11. **B. A. Alvi, M. Asif, and A. Israr, “Measurement of complex mode amplitudes in multimode fiber,” in SPIE OPTO , 82840K–82840K (2012). [CrossRef]

**12. **T. Kaiser, D. Flamm, S. Schrter, and M. Duparr, “Complete modal decomposition for optical fibers using CGH-based correlation filters,” Opt. Express **17**(11), 9347–9356 (2009). [CrossRef] [PubMed]

**13. **J. W. Nicholson, A. D. Yablon, S. Ramachandran, and S. Ghalmi, “Spatially and spectrally resolved imaging of modal content in large-mode-area fibers,” Opt. Express **16**(10), 7233–7243 (2008). [CrossRef] [PubMed]

**14. **J. W. Nicholson, A. D. Yablon, J. M. Fini, and M. D. Mermelstein, “Measuring the modal content of large-mode-area fibers,” IEEE J. Sel. Top. Quantum Electron. **15**(1), 61–70 (2009). [CrossRef]

**15. **D. N. Schimpf, R. A. Barankov, and S. Ramachandran, “Cross-correlated (C^{2}) imaging of fiber and waveguide modes,” Opt. Express **19**(14), 13,008–13,019 (2011). [CrossRef]

**16. **V. Inavalli and N. Viswanathan, “Switchable vector vortex beam generation using an optical fiber,” Opt. Commun. **283**(6), 861–864 (2010). [CrossRef]

**17. **A. Snyder and J. Love, *Optical waveguide theory*, vol. 190 (Springer, 1983).

**18. **J. J. Kapany and N. S. Burke, *Optical Waveguides* (Academic Press, 1972).

**19. **A. V. Volyar and T. A. Fadeeva, “Vortex nature of optical fiber modes: I. structure of the natural modes,” Tech. Phys. Lett. **22**, 330–332 (1996).

**20. **R. J. Black and L. Gagnon, *Optical Waveguide Modes: Polarization, Coupling and Symmetry*, 1st ed. (McGraw-Hill Professional, 2010).