Modal decomposition of fiber modes based on direct far-field measurements at two different distances with a multi-variable optimization algorithm

Byungho Kim; Jeongkyun Na; Juhwan Kim; Hansol Kim; Yoonchan Jeong; Yoonchan Jeong

doi:10.1364/OE.430161

1. Introduction

Following the remarkable development of high-power fiber laser technology in recent years [1,2], the needs for large-mode-area (LMA) fibers have been increasing [3]. Since the core radius of an LMA fiber is large, its V-number is normally greater than the cut-off value for single-mode excitation (i.e., $\textrm{V} = 2.405$) [4]. Thus, higher-order modes (HOMs) sometimes become predominant in an LMA fiber for various reasons, including the transverse mode instability, some unwanted mode couplings, etc. [5,6]. Excitation of HOMs often becomes a critical issue in LMA-fiber-based devices and systems [7,8], so that it is of great importance to quantify the fractions of HOMs with high accuracy and reliability.

To date, a number of modal decomposition (MD) methods have been introduced and investigated in order to figure out how HOMs are excited in LMA fibers [9–18]. In fact, most MD methods utilize the near-field measurement in that the near-field pattern measured just in front of the fiber is nearly identical to that inside the fiber, because diffraction effects are negligible under such conditions. In addition, to obtain the optimal MD, various numerical algorithms, such as the genetic algorithm [10], the stochastic parallel gradient descent (SPGD) algorithm [11], and deep learning [15–18], are often utilized. Whilst the existing MD methods based on the near-field measurement are effective in some specific circumstances, there is room for improvement for a couple of reasons: (1) Their measurement resolution is strongly affected by the detector resolution, which is normally determined by the charge-coupled device (CCD)’s unit pixel size. In general, the pixel size of a standard high-resolution CCD camera is as large as ∼ 2 $\mu \textrm{m}$, which is sometimes not fine enough for resolving the near-field patterns of LMA-fiber modes if they are considerably multimoded to result in rather complicated field patterns out of a core of 10 to 20 $\textrm{}\mu \textrm{m}$ in diameter. To overcome this limitation, Huang et al. and An et al. have proposed to use an additional 4f-imaging system [12,15] although the use of extra imaging systems may increase complexity and inaccuracy in terms of the measurement and alignment. (2) MD performed only by the near-field measurement can hit a conjugate-pair issue [11], so that the solution obtained by the near-field measurement should always be double-checked by additional far-field measurements.

In order to overcome these limitations, one can instead utilize direct far-field measurements of fiber modes without relying on an additional 4f-imaging system for MD as long as one can properly deal with the mode conversion from fiber modes to free-space modes in the first place. In this paper, we present a novel MD method, exploiting free-space Hermite Gaussian (HG) modes converted from individual LMA-fiber modes as the basis modes for reconstructing the far-field pattern of a multimode beam from an LMA fiber [19,20]. We also implement an SPGD-based multi-variable optimization algorithm equipped with the $\textrm{D}4\mathrm{\sigma }$ method to obtain the optimal MD regardless of the existence of an offset of the beam pattern at the CCD image plane. The proposed method directly measures the far-field beam pattern at two different distances, so that the double boundary conditions ensure the unique solution to the MD. We demonstrate the feasibility and effectiveness of the proposed method by means of both numerical simulations and experimental measurements.

The detailed discussion is thus given as the following: In Section 2, we explain the theoretical background of the proposed MD method and verify its feasibility and effectiveness regarding the far-field beam patterns numerically generated from a typical fiber under LMA condition. In Section 3 and Section 4, we carry out numerical simulations and experimental measurements regarding a multimode beam at 654 nm generated through an optical fiber of $\textrm{V} = 5.5$, respectively. Based on the measurements, we carry out the MD by the proposed method and characterize how accurately and reliably it can determine the modal fractions and figure out the further evolution of the beam pattern at an even farther distance. In Section 5, we extend our discussion to the potential of the proposed method for scaling to an LMA fiber of a larger core dimension that guides more modes. In Section 6, we finally draw the conclusion of our investigation on the background that the proposed far-field-based MD method is capable of providing the unique solution to the MD of LMA-fiber modes, readily overcoming the resolution issue caused by the limited pixel size of a CCD camera even with a simple experimental set-up.

2. Theory and the numerical method

2.1 Composition of a multimode beam as the superposition of linearly polarized modes

The eigen-modes of a weakly-guiding optical fiber can be represented by linearly polarized (LP) modes [21]. We note that whilst the eigen-modes and the corresponding field profiles are determined individually, the total field of an actual beam guided by an LMA fiber under a certain condition is represented by a composite field ${E_{\textrm{Fiber}}}({x,y,z} )$ through the superposition of the entire guided modes:

(1a)$${E_{\textrm{Fiber}}}({x,y,z} )= \mathop \sum \nolimits_{j = 1}^N {\rho _j}{\psi _j}({x,y,z} )= \mathop \sum \nolimits_{j = 1}^N {\rho _j}{\psi _{0,j}}({x,y} ){e^{ - i{\beta _j}z}}$$

with

(1b)$${\rho _j} = |{{\rho_j}} |{e^{i{\phi _j}}},$$

where ${\psi _j}$ and ${\psi _{0,j}}$ denote the LP-mode field for the $j$-th guided mode and the corresponding LP-mode-field profile at $z = 0$, respectively, the latter of which is represented by a real-valued, even or odd function in terms of x or y, normalized as $\mathop {\smallint\!\!\!\smallint }\nolimits_{ - \infty }^\infty {|{{\psi_{0,j}}} |^2}dxdy = 1$ [21]; ${\beta _j}$ denotes the wave number for the $j$-th guided mode; ${|{{\rho_j}} |^2}$ denotes the relative power of the $j$-th guided mode, normalized via $\mathop {\smallint\!\!\!\smallint }\nolimits_{ - \infty }^\infty {|{{E_{\textrm{Fiber}}}} |^2}dxdy = 1$; ${\phi _j}$ denotes the relative phase of the $j$-th guided mode to that of the lowest-order guided mode, i.e., the fundamental mode (FM), the latter of which can be set to 0 without loss of generality; N denotes the number of total guided modes in the fiber. We also note that the longitudinal position at $z = 0$ coincides with the fiber’s output end for simplicity, so that ${E_{\textrm{Fiber}}}({x,y,z} )$ is justified in the region of $z \le 0$. In fact, the number of the total guided modes depends on the V-number of the fiber [19]. We here consider that the V-number is given by 5.5 as a typical value for a step-index LMA fiber [8–10]. In this case, total 8 modes, including $\textrm{L}{\textrm{P}_{01}}$, $\textrm{L}{\textrm{P}_{11,\;\textrm{odd}}}$, $\textrm{L}{\textrm{P}_{11,\; \textrm{even}}}$, $\textrm{L}{\textrm{P}_{02}}$, $\textrm{L}{\textrm{P}_{21,\textrm{ odd}}}$, $\textrm{L}{\textrm{P}_{21,\;\textrm{even}}}$, $\textrm{L}{\textrm{P}_{31,\; \textrm{odd}}}$, and $\textrm{L}{\textrm{P}_{31,\;\textrm{even}}}$, are guided by the fiber if the guided beam is linearly polarized. For the MD, what we need to do is to determine the relative power ${|{{\rho_j}} |^2}$ and relative phase ${\phi _j}$ for each mode which has a non-zero contribution to the composite beam.

2.2 Conversion into Hermite Gaussian modes

Whilst we have figured out how the fiber-optic guided modes are determined in the preceding section, the evolution of the output beam transmitted from the fiber end into free space cannot be described directly by the guided modes of the fiber, because it propagates in unbounded, free space undergoing diffraction immediately after exiting the fiber end. This is really true for a beam to be measured under far-field conditions. Thus, in order to describe the diffraction effects in free space properly, the out beam should be represented by the “eigen-modes of free space”, i.e., Gaussian modes, which can be in the form of either HG modes or Laguerre Gaussian modes [19,20]. Since we represent the fiber modes with LP modes, it should be sensible to decompose the free-space beam in terms of the former rather than the latter. The HG modes in free space are given by the following form [21]:

(2a)$${F_{l,m}}({x,y,z} )= {A_{l,m}}\frac{{{W_0}}}{{W(z )}}{G_l}\left( {\frac{{\sqrt 2 x}}{{W(z )}}} \right){G_m}\left( {\frac{{\sqrt 2 y}}{{W(z )}}} \right)exp\left[ { - ikz - ik\frac{{{x^2} + {y^2}}}{{2R(z )}} + i({l + m + 1} )\zeta (z )} \right]$$

with

(2b)$${G_l}(u )= {H_l}(u )\textrm{exp}\left( { - \frac{{{u^2}}}{2}} \right),\; {m},l = 0,1,2,\; \textrm{} \ldots ,$$

where ${A_{l,m}}$ denotes a normalizing factor to yield $\mathop {\smallint\!\!\!\smallint }\nolimits_{ - \infty }^\infty {|{{F_{l,m}}} |^2}dxdy = 1$; ${H_l}$ denotes the Hermite polynomial of order l; ${W_0}$ denotes the beam waist, $ k$ denotes the wave number of the beam in free space; $W(z )$, $R(z )$, and $\zeta (z )$ denote the beam size, the radius of curvature of the wavefront, and the phase shift, respectively. We note that the z axis given here coincides with the center axis of the given fiber, assuming that the fiber end is flat-cleaved. Thus, we hereafter call these unprimed $xyz$ coordinates as the “HGM” coordinates, the origin of which is located at the center of the fiber end.

The mode conversion across the fiber-air interface is determined by the field-continuity condition at $z = 0$. For simplicity, the reflection by the fiber facet can be ignored in that it has no influence on the field-continuity condition in the transverse directions and can also be taken into account via simply adjusting the normalizing factors of the fields if necessary. This suggests that at $z = 0$ the field of each LP mode should be matched to a certain individual field constructed in free space. Thus, if the free-space field justified in the region of $z \ge 0$ is represented by ${\xi _j}({x,y,z} )$ through the superposition of the eigen-mode fields in free space, i.e., a complete set of mutually orthogonal HG modes, this should satisfy the field-continuity condition at $z = 0$ with the corresponding LP-mode field, thereby yielding [22,23]

(3a)$${\xi _j}({x,y,z} )= \mathop \sum \nolimits_{l,m} {a_{j,l,m}}\; \textrm{}{F_{j,l,m}}({x,y,z} )$$

with

(3b)$${a_{j,l,m}} = \mathop {\smallint\!\!\!\smallint }\nolimits_{ - \infty }^\infty {\psi _{0,j}}({x,y} ){F_{j,l,m}}({x,y,0} )dxdy,$$

where ${F_{j,l,m}}$ denotes the normalized HG-mode field of the beam waist ${W_{0,j}}$. It is noteworthy that we have introduced an additional index j to ${F_{l,m}}$ to clarify that the corresponding HG-mode field is dedicated to a specific LP-mode field ${\psi _{0,j}}$. In fact, Eq. (3b) comes from the field-continuity condition at $z = 0,$ i.e., ${\psi _{0,j}}({x,y} )= {\xi _j}({x,y,0} )$. Thus, the coefficient ${a_{j,l,m}}$, which is related with the magnitude and initial phase of the corresponding HG-mode field, is readily be determined by the field-overlap integral as given by Eq. (3b), because ${F_{j,l,m}}$’s are all normalized and mutually orthogonal for the given j [22,23]. By the uniqueness theorem [24], the free-space field ${\xi _j}({x,y,z} )$ that satisfies both the wave equation in free space and the given boundary condition, must be the general solution justified in the region of $z \ge 0$, which is, in fact, the free-space field transmitted from the fiber end on account of the LP-mode field ${\psi _j}({x,y,z} )$ in the fiber. In addition, it is noteworthy that since both ${\psi _{0,j}}({x,y} )$ and ${F_{j,l,m}}({x,y,0} )$ are defined as real quantities [see Eqs. (1) and (2) as well as the description given below Eq. (1)], the field-overlap integral between them, i.e., the corresponding coefficient ${a_{j,l,m}}$, invariably results in a real quantity signified with a plus or minus sign.

After all, the total field of the output beam ${E_{\textrm{FS}}}({x,y,z} )$ transmitted from the fiber end into free space, i.e., into the region of $z \ge 0$, is given by the total sum of ${\xi _j}({x,y,z} )$’s along with the corresponding modal coefficients ${\rho _j}$’s for the whole individual LP modes within the fiber:

(4)$${E_{\textrm{FS}}}({x,y,z} )= \mathop \sum \nolimits_j {\rho _j}{\xi _j}({x,y,z} )= \mathop \sum \nolimits_{j,l,m} |{{\rho_j}} |{e^{i{\phi _j}}}{a_{j,l,m}}\; \textrm{}{F_{j,l,m}}({x,y,z} ),$$

where we note that ${a_{j,l,m}}$’s and ${F_{j,l,m}}({x,y,z} )$’s are known coefficients and functions as defined by Eqs. (3b) and (2a), respectively; however, $|{{\rho_j}} |$’s and ${\phi _j}$’s remain unknown for the moment, which will be determined via the SPGD-based MD to be discussed in the sections that follow. In fact, ${E_{\textrm{FS}}}({x,y,z} )$, which already satisfies the wave equation in free space in the region of $z \ge 0$, can uniquely be determined as long as it also satisfies the field-continuity condition with ${E_{\textrm{Fiber}}}({x,y,z} )$ at $z = 0$, i.e., ${E_{\textrm{FS}}}({x,y,0} )= {E_{\textrm{Fiber}}}({x,y,0} )$ [19], for which $|{{\rho_j}} |$’s and ${\phi _j}$’s should correctly be determined.

In principle, as discussed above, the superposition of HG modes of an arbitrarily large number can represent any field pattern if an appropriate boundary condition is given, because they form a complete basis set [23]; however, in practice one may need to restrict the number of modes to a certain finite number in order that they effectively represent a specific LP-mode-field pattern from a perspective of accuracy and efficiency. It is important to use an appropriate number of HG modes in that decomposition with too many HG modes would require too much computing power, whereas decomposition with too few HG modes would yield too low accuracy. We thus implement a computational strategy such that we maximize and make the total sum of the HG-mode powers, i.e., $\sum {|{{a_{j,l,m}}} |^2}$, for each LP mode greater than at least 99% relative to the individual LP-mode power via iteratively optimizing ${W_{0,j}}$ as well as minimizing the number of HG modes to take.

We present an example calculation result for the LMA fiber specified in the preceding section [$\textrm{V}\, = \; \textrm{}\,5.5$ at 654 nm: 8.2-$\mu $m diameter core with 0.14 numerical aperture (NA)] in Table 1. We have found that HG modes of ${F_{00}}$ to ${F_{60}}$ were sufficient to represent each LP mode such that we could readily obtain over 99% in terms of the ratio of the total sum of the HG-mode powers relative to the corresponding LP-mode power with the given combinations of HG modes.

Table 1. MD of 8 LP modes with HG modes

View Table | View all tables in this article

2.3 Off-axis reconstruction

In practice, when we measure a beam intensity profile at a specific distance by using a CCD camera, the center axis of the fiber or the center axis of the output beam transmitted from the fiber end is not necessarily matched to the center axis of the image plane of the CCD camera. Moreover, if the fraction of HOMs is significantly larger than that of the FM, the x and y centroids for the composite beam, i.e., $\bar{x}$ and $\bar{y}$, can be different from those for the lowest HG mode, i.e., the fundamental Gaussian (FG) mode. This is due to the interference effects amongst the multiple modes, although the centroids of all individual higher-order HG modes should coincide with that of the FG mode. To resolve this issue, we thus devise an iterative technique based on the $\textrm{D}4\mathrm{\sigma }$ method as the following.

For ease of description, based on Eq. (4) we represent the total field of the output beam at a specific distance of ${z_d}$ from the fiber end as

(5)$${E_d}({x,y} )= {E_{\textrm{FS}}}({x,y,{z_d}} )= \mathop \sum \nolimits_{j,l,m} |{{\rho_j}} |{e^{i{\phi _j}}}{a_{j,l,m}}\; \textrm{}{F_{j,l,m}}({x,y,{z_d}} ).$$

We note that since the CCD camera cannot measure the relative phases directly, we instead use the intensity distribution as defined by

(6)$${I_{\textrm{HGM},d}}({x,y} )= {I_{\textrm{HGM}}}({x,y,{z_d}} )= {|{{E_d}({x,y} )} |^2}.$$

In the following, the primed $x^{\prime}y^{\prime}z^{\prime}$ coordinates denote the coordinates for the CCD camera, i.e., the “CCD” coordinates, assuming that $x^{\prime}$, $y^{\prime}$, and $z^{\prime}$ may not necessarily coincide with x, y, and z of the HGM coordinates. The iterative technique for determining the centroids of the free-space beam is based on the following three steps:

(1) Obtain the initial centroid values based on Eq. (7), using the beam intensity profile ${I_{\textrm{CCD}}}$ measured by the CCD camera. Assuming that these initial centroids represent the virtual center of the output beam, reconstruct the beam intensity profile ${I_{\textrm{HGM}}}$ in the HGM coordinates. $(7)$$\overline {{x_{\textrm{CCD},d}}} = \frac{{{\smallint\!\!\!\smallint }{I_{\textrm{CCD},d}}({x{^{\prime}},y{^{\prime}}} )x{^{\prime}}dx{^{\prime}}dy{^{\prime}}}}{{{\smallint\!\!\!\smallint }{I_{\textrm{CCD},d}}({x{^{\prime}},y{^{\prime}}} )dx{^{\prime}}dy{^{\prime}}}},\; \textrm{}\overline {{y_{\textrm{CCD},d}}} = \frac{{{\smallint\!\!\!\smallint }{I_{\textrm{CCD},d}}({x{^{\prime}},y{^{\prime}}} )y{^{\prime}}dx{^{\prime}}dy{^{\prime}}}}{{{\smallint\!\!\!\smallint }{I_{\textrm{CCD},d}}({x{^{\prime}},y{^{\prime}}} )dx{^{\prime}}dy{^{\prime}}}}\; \textrm{}({d = 1,\; \textrm{}2} ).$$$
(2) Obtain new centroids based on the reconstructed beam intensity profile ${I_{\textrm{HGM}}}$ in the HGM coordinates, which are set to be $\overline {{x_{\textrm{HGM},d}}} $ and $\overline {\; {y_{\textrm{HGM},d}}} $ as the following: $(8)$$\overline {{x_{\textrm{HGM},d}}} = \frac{{{\smallint\!\!\!\smallint }{I_{\textrm{HGM},d}}({x,y} )xdxdy}}{{{\smallint\!\!\!\smallint }{I_{\textrm{HGM},d}}({x,y} )dxdy}},\; \textrm{}\overline {{y_{\textrm{HGM},d}}} = \frac{{{\smallint\!\!\!\smallint }{I_{\textrm{HGM},d}}({x,y} )ydxdy}}{{{\smallint\!\!\!\smallint }{I_{\textrm{HGM},d}}({x,y} )dxdy}}\; \textrm{}({d = 1,\; \textrm{}2} ).$$$
(3) Perform the coordinate transformation from ${F_{j,l,m}}({x,y,{z_d}} )$ to ${F_{j,l,m}}({x_d^{^{\prime}},y_d^{{\prime}},z_d^{{\prime}}} )$ or vice versa, using the relations $x_d^{{\prime}} = \; \textrm{}x - ({\overline {{x_{\textrm{CCD},d}}} - \overline {{x_{\textrm{HGM},d}}} } )$ and $\textrm{}y_d^{{\prime}} = \; \textrm{}y - ({\overline {{y_{\textrm{CCD},d}}} - \overline {{y_{\textrm{HGM},d}}} } )$. Assume $z_d^{{\prime}} = {z_d}$ for simplicity. Repeat this procedure while updating the reconstructed beam intensity profiles at the same time, the latter procedure of which is to be discussed in Section 2.4, until obtaining the final outcome that satisfies the optimization tolerance. (If both the reconstructed beam in the HGM coordinates and the measured beam in the CCD coordinates are fully matched, the differences of their centroid points must be identical to the differences of their origins regardless of the existence of an offset of the beam pattern measured at the CCD image plane, whilst they may not necessarily be matched in the beginning of the iteration steps.)

In Figs. 1(a) and (b), we show the beam intensity profiles to be reconstructed in the HGM coordinates and to be measured by the CCD camera, respectively, for example, in which the black and red dots denote the calculated centroids of the corresponding beam intensity profiles and the centers of the HG modes in the corresponding coordinates (i.e., the HGM coordinates and the CCD coordinates), respectively. One can see that the differences between the two x and y centroids in Figs. 1(a) and (b) become identical to the differences between the HGM center and the CCD center in Fig. 1(b) once the beam intensity profile reconstructed by the SPGD algorithm has been made exactly the same as the beam intensity profile measured by the CCD camera.

Fig. 1. Schematic for finding the transverse offset between the HGM and CCD coordinates: (a) Beam intensity profile to be reconstructed by the SPGD algorithm and (b) beam intensity profile to be measured by a CCD camera. The black and red dots denote the calculated centroids of the corresponding beam intensity profiles and the centers of the HG modes in the corresponding coordinates (i.e., the HGM coordinates and the CCD coordinates), respectively.

Download Full Size | PDF

2.4 Field optimization procedure

The field optimization is carried out based on the SPGD algorithm twice at $z = {z_1}$ and $z = {z_2}$, which are the two different distances from the fiber end to the CCD camera. Thus, the SPGD algorithm for the MD of 8 LP modes (1 FM and 7 HOMs) is eventually run with 16 unknown variables, including 7 relative powers (i.e., ${|{{\rho_j}} |^2}$), 7 relative phases (i.e., ${\phi _j}$), and 2 nearly known constraints of ${z_1}$ and ${z_2}$, the last two of which are optional. It is noteworthy that 2 variables (one each for the relative powers and the relative phases) become redundant due to the fact that the total sum of the relative powers is normalized to unity and the relative phase of the FM is set to 0, as assumed in Section 2.1. For testing the reconstructed intensity profiles during the SPGD algorithm, we use the following error or cost function, which is commonly used in the modal analysis [25]:

(9)$${\mathrm{\Delta }_d} = 1 - \frac{{{\smallint\!\!\!\smallint }{I_{\textrm{CCD},d}}({x{^{\prime}},y{^{\prime}}} )\; \textrm{}{I_{\textrm{HGM},d}}({x_d{^{\prime}},y_d{^{\prime}}} )dx{^{\prime}}dy{^{\prime}}}}{{\sqrt {{\smallint\!\!\!\smallint }I_{\textrm{CCD},d}^2({x{^{\prime}},y{^{\prime}}} )dx{^{\prime}}dy{^{\prime}\; }{\smallint\!\!\!\smallint }I_{\textrm{HGM},d}^2({x_d^{{\prime}},y_d^{{\prime}}} )dx{^{\prime}}dy{^{\prime}}} }}.$$

It noteworthy that since the composite beam is reconstructed based on dual procedures of testing at two different distances, the conjugate-pair solution issue can effectively be avoided [11]. The detailed numerical procedures are summarized in the following along with the flow chart illustrated in Fig. 2. The solution is iteratively determined such that it eventually results in the best matches with the intensity profiles measured by the CCD camera at both $z = {z_1}$ and $z = {z_2}$ through implementing SPGD-based random perturbations to the unknown variables. In addition, it is noteworthy that while the size of the random perturbations is basically made proportional to the magnitudes of its individual components, it is also dynamically adjusted during the iteration loop, with taking account of the rate of convergence if necessary. In fact, the additional dynamic adjustment significantly helps increase the speed of convergence and avoid local extremum solutions [11,26]. For example, in the case when the rate of convergence is too slowed down while the cost function is yet far off the target tolerance, one can refresh and enlarge the size of the random perturbations.

1) Calculate the centroids $\{{({\overline {{x_{\textrm{CCD},d}}} ,\; \overline {{y_{\textrm{CCD},d}}} } )\textrm{|}d = 1,2} \}$ of the measured beam intensity profiles via Eq. (7) at $z = \; \textrm{}{z_1}$ and $z = \; \textrm{}{z_2}$.
2) Set the initial values for the unknown variables, i.e. the variable set for optimization, which include $|{{\rho_j}} |$’s, ${\phi _j}$’s, $ \textrm{}{z_1}$, and ${z_2}$ and are iteratively determined through the SPGD algorithm. Note that fine adjustment of the distance parameters, i.e., ${z_1}$ and ${z_2}$, is additionally carried out for taking account of the uncertainties that may arise during the experimental measurement, which may be skipped if unnecessary.
3) Generate random perturbations and apply them to the variable set, the size of which is made proportional to the magnitudes of its individual components or is dynamically adjusted with taking account of the rate of convergence if necessary.
4) Generate the HG intensity profiles ${I_{\textrm{HGM},1}}$ and ${I_{\textrm{HGM},2}}$ based on the refreshed variable set, referring to Table 1 and calculate the centroids of the generated beam intensity profiles $\{{({\overline {{x_{\textrm{HGM},d}}} ,\; \overline {{y_{\textrm{HGM},d}}} } )\textrm{|}d = 1,2} \}$ via Eq. (8).
5) Reconstruct a pair of new beam intensity profiles out of the beam intensity profiles obtained from Step 4) by translating them via $\{{x_d^{\prime} = \; x - ({\overline {{x_{\textrm{CCD},d}}} - \overline {{x_{\textrm{HGM},d}}} } ),\; \textrm{}y_d^{\prime}} =$ $y - ({\overline {{y_{\textrm{CCD},d}}} - \overline {{y_{\textrm{HGM},d}}} } ) {|{d = 1,2} } \} $.
6) Calculate the total error between the reconstructed and measured beam intensity profiles via the cost function of Eq. (8). If the error is “sufficiently” lower than that of the previous step, which means the overall rate of its convergence is also within the acceptable range, update the variable set such that it includes the perturbations applied. Otherwise, simply proceed without updating the variable set.
7) Terminate the SPGD algorithm, keeping the updated variable set if the total error is lower than the tolerance value. Otherwise, return to Step 3).

Fig. 2. Flow chart of the SPGD algorithm for MD.

Download Full Size | PDF

3. Modal decomposition of a theoretically constructed composite beam

To verify the feasibility of the proposed method, we test it to a known composite beam arbitrarily generated by theory. We consider the same conditions for the wavelength of the incident radiation and the LMA fiber as given in the preceding sections. The relative powers and phases are arbitrarily determined, and ${z_1}$ and ${z_2}$ are set to $7$ and $16\; \textrm{mm}$, as a typical example. In addition, we assume that the CCD camera is “deliberately” positioned with an x-offset of 0.5 mm and a y-offset of 1 mm to the fiber axis, i.e., to the origin of the HGM coordinates in the transverse directions.

In Fig. 3, we illustrate the generated and reconstructed beam intensity profiles for a typical case of the initial conditions. The errors defined by Eq. (9) were estimated to $9.66 \times {10^{ - 8}}$ and $9.64 \times {10^{ - 8}}$ for ${z_1}$ and ${z_2}$, respectively, which were primarily incurred by the numerical tolerance set to $1 \times {10^{ - 7}}$. The numerical results are both qualitatively and quantitatively in good match.

Fig. 3. SPGD simulation results. (a) and (b): generated and reconstructed beam intensity profiles at ${z_1} = 7\; \textrm{mm}$, respectively; (c) and (d): generated and reconstructed beam intensity profiles ${z_2} = 16\; \textrm{mm}$, respectively.

Download Full Size | PDF

In Fig. 4(a), we illustrate the number of iterations required for obtaining the reconstructed beam intensity profiles for 15 different cases of the initial conditions. We note that on average it took ${\sim} \; \textrm{}3.42 \times {10^4}$ iterations to complete the MD via the SPGD algorithm. In Fig. 4(b), we also illustrate how the magnitude of the cost function given by Eq. (9) evolved with iterations for a typical case of the initial conditions. The spikes in Fig. 4(b) were due to the enlarged perturbations made through the dynamic adjustment when the overall rate of convergence was too slowed down. In more detail, the dynamic adjustment was made when the convergence rate of the cost function failed to yield at least a 10-dB reduction over 1,000 iterations, the condition of which will hereafter hold unless stated otherwise. Even though we took account of 8 LP modes and additional parameters for representing the off-axis and longitudinal information, the average iteration number required for completion of the MD is comparable with that of the existing MD method applied to an LMA fiber of 6 modes [11]. All the detailed information, including the mean relative powers and phases, and root-mean-squared (RMS) errors over the 15 different cases of the initial conditions, is summarized in Table 2.

Table 2. SPGD based MD results

View Table | View all tables in this article

Fig. 4. Progression of the SPGD simulation: (a) Iteration numbers required for completion of the MD, depending on how the initial conditions are set and (b) evolution of the magnitude of the cost function given by Eq. (9) with iterations for a typical case of the initial conditions.

Download Full Size | PDF

It is noteworthy that the RMS errors for the relative powers are represented in percentage, being normalized by the relative powers of the individual modes. We stress that the RMS percentage errors of relative powers are well below $1{\%}$ for the entire individual modes. In particular, the RMS percentage error of the FM’s relative power is less than $0.05{\%}$, which demonstrates the high accuracy of the proposed MD method.

In Fig. 5, we illustrate the generated and reconstructed beam intensity profiles at $z = {z_3}({ = 25\; \textrm{mm}} )$ based on the results shown in Fig. 3, in which we have translated the generated beam to the center of the image plane just to compare them side by side. One can see that both beam intensity profiles are still in good match even at a farther distance. The overall error defined by Eq. (9) is estimated to $9.62 \times {10^{ - 8}}$, remaining below the numerical tolerance set to $1 \times {10^{ - 7}}$, which clearly indicates that there is no conjugate-pair issue in the solution obtained.

Fig. 5. SPGD simulation results for far-field intensity profiles at ${z_3}({ = 25\; mm} ):$ (a) Generated and (b) reconstructed beam intensity profiles.

Download Full Size | PDF

In Fig. 6, we also illustrate the near-field reconstruction result at $z = 0$ based on the results shown in Fig. 3, in which we have again translated the generated beam to the center of the image plane just to compare them side by side. One can see that both beam intensity profiles look nearly identical. It is noteworthy that over the all 15 simulations, there has been no conjugate-pair issue with in the solution obtained. This result also confirms that the mode conversion from LP modes into HG modes summarized in Table 1 have correctly been carried out.

Fig. 6. SPGD simulation results for near-field intensity profiles at $z = 0$: (a) Generated and (b) reconstructed beam intensity profiles.

Download Full Size | PDF

In addition, we emphasize that the transverse offset between the generated and reconstructed beam intensity profiles was well figured out via Eqs. (7) and (8), so that there was no issue at all in performing the MD even when the generated beam intensity profiles had considerable transverse and longitudinal offsets. This is a great merit of the proposed method, because the beam intensity profiles measured by a CCD camera can invariably have a variety of transverse and longitudinal offsets depending on the experimental conditions.

4. Modal decomposition of experimentally constructed composite beams

4.1 Experimental arrangement

For demonstrating the full feasibility of the proposed method, we perform an experiment, generating typical composite-mode beam intensity profiles through an LMA fiber and reconstructing the corresponding beam intensity profiles by the proposed method. In Fig. 7, we illustrate the schematic of the experiment arrangement. In the experiment, an SMF-28 was used as a typical LMA fiber for the incident laser at 654 $\textrm{nm}$ with a full-width-at-half-maximum (FWHM) bandwidth of ∼ 0.9 $\textrm{nm}$ and an average power of 0.1 $\textrm{mW}$. It is noteworthy that although an SMF-28 is supposed to operate in single mode at ∼ 1550 $\textrm{nm}$, it operates as an LMA fiber for the incident light at 654 $\textrm{nm}$, because it has a step-index core with a diameter of 8.2 $\mu \textrm{m}$ and an NA of 0.14, thereby yielding the V-number of 5.5, which definitely guides 8 LP modes for the incident light at 654 $\textrm{nm}$. In addition, to ensure the excitation of LP modes, we collimate the incident beam through a polarizing beam splitter, so that only linearly polarized light was incident on the LMA fiber. The fiber had the length of 10 $\textrm{cm}$ and was kept straight in order for it to preserve the incident polarization state without having any bend-induced birefringence [27]. The beam intensity profiles of the transmitted laser light were measured by a CCD camera initially at ${z_1} = 7\; \textrm{mm}$ and ${z_2} = 16\; \textrm{mm}$, which are determined by the distance from the fiber end. Another beam intensity profile was measured at ${z_3} = 25\; \textrm{mm}$ for further comparison between the measured and reconstructed beams intensity profiles. These experimental conditions are exactly the same as those investigated in the preceding sections. In addition, the CCD camera had $2752 \times 2192$ pixels with a unit pixel size of $4.54\; \textrm{}\mu \textrm{m}$ ${\times} $ $4.54\; \textrm{}\mu \textrm{m}$. We removed all the glassware mounted in front of the CCD camera sensor array, so that there was nothing in front of the senor array of the CCD camera.

Fig. 7. Experiment setup to measure the beam intensity profiles at two different distances. To ensure the excitation of LP modes, the beam was incident on the LMA fiber through a polarizing beam splitter. The spectrum of the incident laser beam is shown on the left-hand side (resolution: 0.05 nm).

Download Full Size | PDF

We note that the $\textrm{D}4\mathrm{\sigma }$ method and the numerical MD itself can be sensitive to camera noise. Thus, we paid extra attention to minimizing camera noise and environmental disturbances: We carried out the experiment in the standard laboratory environment at room temperature, keeping all the optics and CCD camera built upon an optical table, in order that undesirable environmental disturbance issues were appropriately controlled. We kept the CCD camera’s exposure time at 200 $\mu \textrm{s}$, so that any short-term fluctuations in the beam intensity profile were averaged out within the given time scale. We also took 10 consecutive CCD images and averaged them out to determine the beam intensity profile at a fixed position in the z direction although we did not notice any considerable fluctuations among them. In addition, it is noteworthy that the CCD camera we used could operate at 27 fps at full resolution, so that any long-term fluctuations in the beam intensity profile could be resolved out within that rate.

4.2 Modal decomposition of the measured beams

In Fig. 8, we illustrate the measured and reconstructed beam intensity profiles for two different cases of composite modes excited through the LMA fiber under slightly different conditions. We set the error tolerance of the SPGD optimization to less than $3.5 \times {10^{ - 3}}$ in both SPGD optimizations. The beam intensity profiles shown in Figs. 8(a) to (d) denote Case 1, and those shown in Figs. 8(e) to (f) denote Case 2. In Fig. 9, we also illustrate how the cost function evolved with iterations for Case 1 and Case 2, respectively. One can see that the SPGD-based optimizations progressed in a similar manner as shown in Fig. 4.

Fig. 8. Experimental measurements and the corresponding MD results by the proposed SPGD algorithm. (a) to (d): for Case 1; (e) to (f): for Case 2. (a), (c), (e), and (g): beam intensity profiles measured by the CCD camera; (b), (d), (f), and (h): reconstructed beam intensity profiles reconstructed by the proposed SPGD algorithm.

Download Full Size | PDF

Fig. 9. Evolutions of the magnitude of the cost function given by Eq. (9) with respect to iteration number for (a) Case 1 and (b) Case 2, respectively.

Download Full Size | PDF

We note that unlike the case discussed in Section 3, in which we had the exact MD information on the theoretically generated composite beams, for the experimental results it is impossible to estimate the modal mismatches individually, because the individual mode fractions are entirely unknown. Nevertheless, one can see that the reconstructed beam intensity profiles are qualitatively and quantitatively in good match with those measured by the CCD camera. In addition, one can also see that the proposed method worked fine even when there were considerable amount of the transverse offsets in the beam intensity profiles measured by the CCD camera, which were completely resolved out by means of the additional algorithm via Eqs. (7) and (8). It is noteworthy that the distance mismatches were determined to be well below $1\; \textrm{mm}$, which might be due to the experimental errors. The quantitative mismatches between the measured and reconstructed beam intensity profiles based on Eq. (9) were estimated to $3.49 \times {10^{ - 3}}$ and $3.00 \times {10^{ - 3}}$ at ${z_1}$ and $3.48 \times {10^{ - 3}}$ and $3.00 \times {10^{ - 3}}$ at ${z_2}$ for Case 1 and Case 2, respectively. All the MD results are summarized in Table 3.

Table 3. SPGD-based MD results

View Table | View all tables in this article

4.3 Justification of the uniqueness of the modal decomposition

As discussed in Section 3, if the MD results obtained by the SPGD algorithm were uniquely determined and all correct, it should be able to reconstruct the beam intensity profile of the transmitted beam at any distance with preserved accuracy. Thus, to justify the uniqueness and validity of the MD results obtained in the preceding section, we further performed an extra measurement for both Case 1 and Case 2 at a farther distance of ${z_3} = 25\; \textrm{mm}$ and compared them with the corresponding beam intensity profiles reconstructed based on the MD results obtained from the preceding measurements as shown in Fig. 10.

Fig. 10. Experimental measurements and the corresponding numerical reconstructions based on the preceding MD results. (a) and (c): beam intensity profiles at ${z_3} = 25\; \textrm{mm}$ measured by the CCD camera; (b) and (d): beam intensity profiles at ${z_3} ={\sim} \; 25\; \textrm{mm}$ reconstructed by the preceding MD results.

Download Full Size | PDF

One can see that the reconstructed beam intensity profiles are still qualitatively in good match with the measured beam intensity profiles. We note that the quantitative mismatches between them based on Eq. (9) were estimated to $3.49 \times {10^{ - 3}}$ and $2.44 \times {10^{ - 3}}$ at ${z_3}$, which still remain below the original error tolerance set to $3.5 \times {10^{ - 3}}$.

5. Discussion

The simulations and experiments presented above are for the case when 8 LP modes are guided in the LMA fiber. In general, 6 to 8 LP modes are guided in a typical LMA fiber; however, in principle, the proposed method is not completely limited to such a typical case. In this section, we further discuss the potential of the proposed method for scaling to an LMA fiber that guides more than 10 LP modes, for example.

We thus consider an additional case that the fiber’s V-number is scaled up to 6.7, which can be realized by increasing the core diameter of the LMA fiber investigated in the preceding sections to 10 $\mu \textrm{m}$ while keeping all the other parameters the same. Then, this fiber becomes to guide 12 LP modes in total. In this case, we have found that HG modes of ${F_{00}}$ to ${F_{71}}$ were sufficient to represent the given 12 LP modes such that we could obtain over 99% of the modal power fraction relative to each LP-mode power as summarized in Table 4.

Table 4. MD of 12 LP modes with HG modes

View Table | View all tables in this article

With the MD results for the case of 12 LP modes, we have repeated the same numerical procedures as we took for the case of 8 LP modes. In Fig. 11, we illustrate the generated and reconstructed beam intensity profiles. In this trial, the numerical errors defined by Eq. (9) were estimated to $9.89 \times {10^{ - 8}}$ and $9.87 \times {10^{ - 8}}$ at ${z_1}$ and ${z_2}$, respectively, which were primarily incurred by the numerical tolerance set to $1 \times {10^{ - 7}}$.

Fig. 11. SPGD simulation results for an extended case of 12 LP modes. (a) and (b): generated and reconstructed beam intensity profiles at ${z_1} = 7\; \textrm{mm}$, respectively; (c) and (d): generated and reconstructed beam intensity profiles ${z_2} = 16\; \textrm{mm}$, respectively.

Download Full Size | PDF

In Fig. 12(a), we illustrate the number of iterations required for obtaining reconstructed beam intensity profiles for 15 different cases of the initial conditions in the similar manner with the case of 8 LP modes. We note that on average it took ${\sim} \; \textrm{}5.59 \times {10^4}$ iterations to complete the MD via the SPGD algorithm, which is of a comparable order of magnitude with that for the case of 8 LP modes [see Fig. 4(a) for comparison]. In Fig. 12(b), we also illustrate how the magnitude of the cost function given by Eq. (9) evolved with iterations for a typical case of the initial conditions.

Fig. 12. Progression of the SPGD simulation for an extended case of 12 LP modes: (a) Iteration numbers required for completion of the MD, depending on how the initial conditions are set and (b) evolution of the magnitude of the cost function given by Eq. (9) with iterations for a typical case of the initial conditions.

Download Full Size | PDF

All the detailed information on the simulation results shown above, including the mean relative powers and phases, and RMS errors over the 15 different cases of the initial conditions, is summarized in Table 5. In fact, one can still see that the reconstructed beam intensity profiles are both qualitatively and quantitatively in good match with the generated beam intensity profiles. We stress that the RMS percentage errors of relative powers are still below $1{\%}$ for the entire individual modes. It is noteworthy that the mismatch of the relative power of the FM remain as low as ${\sim} \; \textrm{}0.06{\%}$, which is slightly higher than, but still comparable with the case of 8 LP modes. In addition, one can also see that the convergence characteristics of the SPGD optimization were well maintained even when the number of modes was extended to 12 LP modes.

Table 5. SPGD-based MD results with 12 LP modes

View Table | View all tables in this article

Whilst we have verified that the SPGD algorithm utilized to determine unknown multi-variables of modes works nicely with high accuracy and reliability for most typical cases of LMA fibers that guide up to 12 LP modes [10–17], the algorithm’s scalability in terms of the number of modes that it can handle, really depends on the efficacy of the algorithm itself as well as the computing power. Since the accuracy and speed of convergence of the SPGD algorithm inherently tend to diminish with the size of the variable set for optimization, one may consider exploiting other multi-variable optimization methods, including deep-learning techniques [16–18], in conjunction with the proposed direct far-field measurement procedure without relying on an additional 4f-imaging system. In particular, when the MD of even larger numbers of modes, e.g., tens of modes, is required, one can consider exploiting them as a viable alternative.

In addition, whilst we have also verified that the error tolerance could maintain as low as ${10^{ - 7}}$ with the proposed method for the MD of theoretically constructed composite beams of up to 12 LP modes, this was for an ideal situation without including any noise or defects in imaging. In fact, in a real situation when we experimentally carried out the MD of the composite beams, the error tolerance was increased to $3.5 \times {10^{ - 3}}$. This could be explained in two folds: One is that in the real situation the measured intensity beam profiles might contain various inherent CCD noise components, although temporal fluctuations were rectified as much as possible by keeping the CCD camera’s exposure time as long as 200 $\mu \textrm{s}$ as well as taking the averaged image over 10 consecutive CCD measurements. The other is that there might be some uncertainties in determining intrinsic fiber parameters and polarization states of the composite beams. In order to mitigate the inherent CCD noise or random noise further, one can consider implementing low-noise readout techniques for CCD devices [28] or exploiting extended numerical techniques as discussed in [16–18]. Moreover, to obtain the more rigorous MD, one can also consider using full-vectorial fiber modes rather than using LP modes. We thus think that the proposed direct far-field measurement procedure can further be investigated in conjunction with a variety of other numerical MD methodologies.

6. Conclusion

In order to overcome the resolution and conjugate-pair issues of the existing MD methods, we have proposed to exploit the direct far-field measurements at two different distances without relying on an additional 4f-imaging system along with the MD based on an SPGD-based multi-variable optimization algorithm equipped with the $\textrm{D}4\mathrm{\sigma }$ technique. We convert the fiber modes (LP modes) into free-space HG modes, the latter of which play as the basis modes for decomposing the far-field beam intensity profiles measured at a distance. In particular, the proposed method is capable of obtaining the optimal MD solution regardless of the amount of the offsets of the beam intensity profiles measured at the CCD image plane. We have verified both numerically and experimentally that the proposed method can determine the power fractions of all the excided modes in an LMA fiber having V-number of 5.5 with high accuracy, and satisfactorily figure out the beam intensity profiles at a farther distance. In the numerical verification, the error tolerance of the beam intensity mismatch between the theoretically generated and reconstructed beams could go down to less than $1 \times {10^{ - 7}}$. In the experimental demonstration, the error tolerance of the beam intensity mismatch between the measured and reconstructed beams could go down to less than $3.5 \times {10^{ - 3}}$. We note that this level of error tolerance could maintain for beam intensity profiles measured at a father distance, which confirms that the MD results obtained by the proposed method are sufficiently unique and valid. In addition, we have also looked into the potential of the proposed method for scaling to a multimode fiber that guides more than 10 LP modes, for example. We have verified that when the fiber’s V-number was scaled up to 6.7, thereby guiding 12 LP modes in total, the accuracy of the MD results and the convergence of the SPGD optimization could maintain as comparable as the case of 8 LP modes. We expect that the proposed and demonstrated MD method, including the direct far-field measurement procedure, will effectively be utilized in determining and analyzing the characteristics of the composite modes formed in various LMA fibers.

Funding

Agency for Defense Development (UD180040ID).

Disclosures

The authors declare no conflicts of interest.

Data availability

Data underlying the results presented in this paper are not publicly available.

References

1. Y. Jeong, J. K. Sahu, D. N. Payne, and J. Nilsson, “Ytterbium-doped large-core fiber laser with 1.36 kW continuous-wave output power,” Opt. Express 12(25), 6088–6092 (2004). [CrossRef]

2. A. Liem, J. Limpert, H. Zellmer, A. Tünnermann, V. Reichel, K. Mörl, S. Jetschke, S. Unger, H.-R. Müller, J. Kirchhof, T. Sandrock, and A. Harschak, “1.3 kW Yb-doped fiber laser with excellent beam quality,” in Conference on Lasers and Electro-Optics/International Quantum Electronics Conference and Photonic Applications Systems Technologies, Technical Digest (CD) (Optical Society of America, 2004), paper CPDD2.

3. N. G. R. Broderick, H. L. Offerhaus, D. J. Richardson, R. A. Sammut, J. Caplen, and L. Dong, “Large mode area fibers for high power applications,” Opt. Fiber Technol. 5(2), 185–196 (1999). [CrossRef]

4. J. P. Koplow, D. A. V. Kliner, and L. Goldberg, “Single-mode operation of a coiled multimode fiber amplifier,” Opt. Lett. 25(7), 442–444 (2000). [CrossRef]

5. A. V. Smith and J. J. Smith, “Mode instability in high power fiber amplifiers,” Opt. Express 19(11), 10180–10192 (2011). [CrossRef]

6. H. Lü, P. Zhou, X. Wang, and Z. Jiang, “General analysis of the mode interaction in multimode active fiber,” Opt. Express 20(6), 6456–6471 (2012). [CrossRef]

7. H. Yoda, P. Polynkin, and M. Mansuripur, “Beam quality factor of higher order modes in a step-index fiber,” J. Lightwave Technol. 24(3), 1350 (2006). [CrossRef]

8. R. Ismaeel, T. Lee, B. Oduro, Y. Jung, and G. Brambilla, “All-fiber fused directional coupler for highly efficient spatial mode conversion,” Opt. Express 22(10), 11610–11619 (2014). [CrossRef]

9. 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]

10. L. Li, J. Leng, P. Zhou, and J. Chen, “Multimode fiber modal decomposition based on hybrid genetic global optimization algorithm,” Opt. Express 25(17), 19680–19690 (2017). [CrossRef]

11. K. Choi and C. Jun, “High-precision modal decomposition of laser beams based on globally optimized SPGD algorithm,” IEEE Photonics J. 11(5), 1–10 (2019). [CrossRef]

12. L. Huang, S. Guo, J. Leng, H. Lü, P. Zhou, and X. Cheng, “Real-time mode decomposition for few-mode fiber based on numerical method,” Opt. Express 23, 4620–4629 (2015). [CrossRef]

13. D. B. S. Soh, J. Nilsson, S. Baek, C. Codemard, Y. Jeong, and V. Philippov, “Modal power decomposition of beam intensity profiles into linearly polarized modes of multimode optical fibers,” J. Opt. Soc. Am. A 21(7), 1241–1250 (2004). [CrossRef]

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. Select. Topics Quantum Electron. 15(1), 61–70 (2009). [CrossRef]

15. Y. An, L. Huang, J. Li, J. Leng, L. Yang, and P. Zhou, “Deep learning-based real-time mode decomposition for multimode fibers,” IEEE J. Select. Topics Quantum Electron. 26(4), 1–6 (2020). [CrossRef]

16. X. Fan, F. Ren, Y. Xie, Y. Zhang, J. Niu, J. Zhang, and J. Wang, “Mitigating ambiguity by deep-learning-based modal decomposition method,” Opt. Commun. 471, 125845 (2020). [CrossRef]

17. Y. An, J. Li, L. Huang, L. Li, J. Leng, L. Yang, and P. Zhou, “Numerical mode decomposition for multimode fiber: From multi-variable optimization to deep learning,” Opt. Fiber Technol. 52, 101960 (2019). [CrossRef]

18. E. S. Manuylovich, V. V. Dvoyrin, and S. K. Turitsyn, “Fast mode decomposition in few-mode fibers,” Nat. Commun. 11(1), 5507 (2020). [CrossRef]

19. R. Brüning, Y. Zhang, M. McLaren, M. Duparré, and A. Forbes, “Overlap relation between free-space Laguerre Gaussian modes and step-index fiber modes,” J. Opt. Soc. Am. A 32(9), 1678–1682 (2015). [CrossRef]

20. B. M. Heffernan, R. D. Niederriter, M. E. Siemens, and J. T. Gopinath, “Tunable higher-order orbital angular momentum using polarization-maintaining fiber,” Opt. Lett. 42(14), 2683–2686 (2017). [CrossRef]

21. B. E. A. Saleh and M. C. Teich, Fundamentals of Photonics, 2nd ed. (Wiley, 2007).

22. P. Horak and F. Poletti, “Multimode nonlinear fibre optics: theory and applications,” in Recent progress in optical fiber research, M. Yasin, ed. (Books on Demand, 2012), Chap. 1, pp. 3–24.

23. X. Xue, H. Wei, and A. G. Kirk, “Intensity-based modal decomposition of optical beams in terms of Hermite–Gaussian functions,” J. Opt. Soc. Am. A 17, 108–1091 (2000). [CrossRef]

24. D. K. Cheng, Field and Wave Electromagnetics, 2nd ed. (Addison-Wesley, 1989).

25. R. Bruning, P. Gelszinnis, C. Schulze, D. Flamm, and M. Duparré, “Comparative analysis of numerical methods for the mode analysis of laser beams,” Appl. Opt. 52(32), 7769–7777 (2013). [CrossRef]

26. Q. Hu, L. Zhen, Y. Mao, S. Zhu, X. Zhou, and G. Zhou, “Adaptive stochastic parallel gradient descent approach for efficient fiber coupling,” Opt. Express 28(9), 13141–13154 (2020). [CrossRef]

27. R. T. Schermer, “Mode scalability in bent optical fibers,” Opt. Express 15(24), 15674–15701 (2007). [CrossRef]

28. G. F. Moroni, F. Chierchie, L. Stefanazzi, E. E. Paolini, M. S. Haro, G. Cancelo, J. Tiffenberg, and J. Estrada, “Interleaved readout of charge-coupled devices (CCDs) for correlated noise reduction,” IEEE Trans. Instrum. Meas. 69(10), 7580–7587 (2020). [CrossRef]