1.

Introduction

Anderson localization^1–4 was introduced in condensed matter physics where, e.g., electrons can be localized within a disordered medium such as a semiconductor provided that there is a sufficiently large degree of randomness of the impurities or defects in the semiconductor.

The requirements for Anderson localization of light are difficult to achieve in three-dimensions but are much easier in two-dimensional (2-D) systems.⁵ Experimental results were first reported in a 2-D photonic lattice.⁶ Later, theoretical and experimental results were reported using an optical fiber consisting of a random array of smaller fibers, each having either a higher or lower index of refraction.^7–11 Experimental results agreed well with theoretical simulations.

Previous theoretical simulations were carried out by solving the wave propagation equation using the finite difference beam propagation method (FD-BPM) discussed in Refs. 8,12. However, the computational times for these simulations were reported as extremely long.

The optical fiber community developed a different fast Fourier transform beam propagation method (FFT-BPM).^13–17 In this algorithm, a phase mask is inserted that controls the local phase distortion over a given distance and then the Fresnel diffraction is calculated through a uniform medium over that distance. This process is repeated over multiple cycles. In effect, we are simulating a structured medium through this iterative process.

There are at least three approaches for computing the Fresnel diffraction.¹⁸ We refer to the first as the direct propagation algorithm where each of the $N\times N$ points on the output plane is the sum of the electric fields from each of the $N\times N$ points in the input plane. See Sec. 4.1 in Ref. 18. Even though it is the most accurate approach, it is also the most time consuming.

The second Fresnel diffraction approach incorporates the fast Fourier transform of the product of the input pattern with a diverging lens function in the input plane whose focal length depends on the propagation distance as

This function is then multiplied by another diverging lens function in the output plane whose focal length depends on this same propagation distance. For example, see Eq. (4.13) in Ref. 18. This algorithm is limited to large distances because of Nyquist limitations for encoding the lens functions.

The third and less well-known Fresnel diffraction algorithm uses two fast Fourier transform operations with a converging lens operation (quadratic chirp) placed between them. While notations differ in the references, this one is easily identified because the lens function has the propagation distance in the numerator of the lens function exponent. For example, see Eq. (3-42) in Ref. 18. It is therefore quite useful for short distances because the focal length of the lens increases as the propagation distance decreases. We discuss this algorithm more completely in Sec. 2.

Interestingly, we found that this third algorithm was independently developed within the diffractive optics community, first by Goodman¹⁸ and then expanded by others,^19–21 where it has been labeled as either the “convolution approach” or the “angular spectrum approach.” While some believe that this algorithm is well known, it appears to be constantly re-discovered.²²

In more recent work, we clarified the algorithm using a ray matrix analysis,²³ showing how to relate the propagation distance to the focal length of the lens operation. We have used this algorithm to study the propagation of beams by applying the Fresnel diffraction to a pattern encoded onto a spatial light modulator to avoid moving the detector.²⁴

In this work, we apply this algorithm to the random optical fiber mentioned earlier. Our approach involves two steps: (a) inserting a phase mask that controls the local phase distortion created by a short length of random optical fiber and (b) then computing the Fresnel diffraction through a uniform medium over that distance using the angular spectrum approach. This process is repeated over multiple cycles. In effect, we are simulating a structured medium through this iterative process. We find similar results with much shorter computational times, and then we expand our analysis to cover other applications of the random fiber array.

The paper is organized as follows: after the introduction, Sec. 2 reviews the Fourier transform-based Fresnel diffraction algorithm. Section 3 discusses the application of this algorithm to the random array of fibers. Section 4 shows computer simulation results using this technique, which give similar results to previously published results. In Sec. 5, we examine newer experiments involving the imaging through these fiber arrays and, again, find consistent agreement with previously reported results. Finally, Sect. 6 presents the conclusions of the work.

2.

Review of the Fourier Transform-Based Fresnel Diffraction Algorithm

We begin with a brief review of the “angular spectrum” Fresnel diffraction algorithm and use a one-dimensional notation for convenience. We begin with an input function $g(u)$. The Fresnel diffraction pattern $\widehat{g}(u,d)$ formed at a distance $d$ from the input is given as

Eq. (2) can be easily understood. First, the Fourier transform $G(\xi )$ of the input function $g(u)$ is multiplied by a converging lens function represented as $Z^{*}(\xi ,d)$. The focal length $F$ of this converging lens is related to the propagation distance $d$ by $F=f_{\mathrm{FT}}^2/d$.²² Here, the focal length of the Fourier transform system is given by $f_{\mathrm{FT}}$ and the converging lens function can be written as

Then the inverse Fourier transform of the product is performed. As stated earlier, this algorithm is quite useful for short distances because the focal length of the lens increases as the propagation distance decreases.

Next, we discuss our model and the results of the computer simulations. We provide details that allow easy replication of our results.

3.

Application of the Fast Fourier Transform Beam Propagation Method to the Random Array of Optical Fibers

We applied our approach to the random optical fiber problem as defined in Refs. 78.9.10.–11, in which half of the fibers have one index of refraction while the other half have a different index of refraction. Again following the previous work, we assume that this array has $256\times 256$ fibers, each having a pixel size $\mathrm{\Delta }=1\times 1\mu \mathrm{m}$.

Our model is outlined in Fig. 1. An input Gaussian beam is incident onto the 2-D fiber array where each pixel has either of two phases. The phase difference between the two pixels of the phase mask is given as $\phi =(2\pi \mathrm{\Delta }nd)/{\lambda }_0$, where $\mathrm{\Delta }n$ corresponds to the difference in the indices of refraction between the two fiber arrays, and ${\lambda }_0$ is the vacuum wavelength. We then compute the Fresnel diffraction through some propagation distance $d$ and through a medium having a wavelength $\lambda ={\lambda }_0/n$, which is defined as the vacuum wavelength divided by the average of the two indices of refraction. This process is repeated through numerous iterations. In this case, we are not studying the propagation through the structured random fiber. Instead, we impose the phase variation with the mask and then compute the Fresnel diffraction through a uniform medium with an index of refraction that is the average of the two regions.

Fig. 1

Schematic of iterative process used in computer simultions. Each random mask is separated from the next by a fixed distance over which Fresnel diffraction occurs.

We examine results in which we vary the size of the Gaussian beam and the difference in the indices of refraction (keeping the average value of $n=1.5$). We also examine results using several random masks. At a distance of 0.5 mm, we save the image and compute the RMS radius of the intensity. This is defined by a weighted integral of the transmitted energy divided by the integral of the transmitted intensity as

As mentioned earlier, the key to the process is the calculation of the focal length $F=f_{\mathrm{FT}}^2/d$ of the lens required in Eq. (3), and we use our diffractive optics experience as a guide. With an array of $N$ pixels, each having a size $\mathrm{\Delta }$ and a wavelength $\lambda ={\lambda }_0/n$, the effective focal length for the FFT operation is given as $f_{\mathrm{FT}}=N{\mathrm{\Delta }}^2/\lambda$.²³ Using this, the equation defining the focal length in the Fresnel diffraction algorithm is given as

Using the values above ($N=256$, $\mathrm{\Delta }=1\mu \mathrm{m}$, $\lambda =0.5\mu \mathrm{m}$, an average index of refraction of $n=1.5$, and a step size of $d=1\mu \mathrm{m}$), this focal length is $F=2621\mathrm{m}$ and is easily encoded onto the $N\times N$ array.

Next, we show some simulation results that compare our approach with previously reported results.

4.

Computer Simulation Results for the Application of the Fast Fourier Transform Beam Propagation Method to the Random Array of Optical Fibers

As stated earlier, our model assumes a square fiber array having $256\times 256$ smaller fibers, each having a size of $\mathrm{\Delta }=1\mu \mathrm{m}$. Figure 2 shows the RMS intensity versus propagation distance where the step size was $d=1\mu \mathrm{m}$ and $\mathrm{\Delta }n=0.1$ and $n=1.5$ for 10 different random patterns. Here we use an input Gaussian waist of $3.5\mu \mathrm{m}$ with an RMS value for the radius of the input intensity beam of $2.5\mu \mathrm{m}$. Our results are consistent with the computational and experimental results reported earlier.^7,8 The beam broadens as it propagates but reaches a constant value that depends on the random array. The variation in the widths of the beams can be explained since the size of the input beam is roughly the same as the pixel size and the broadening effect depends strongly on the array of fibers that are locally excited by the input beam.

Fig. 2

RMS intensity radius as a function of propagation for 10 different random fiber arrays.

We had some concerns about the effect of the step size on the results, as well as the effects of different sized input beams. Results are shown in Fig. 3. We used input beams having Gaussian electric field waists of 3.5, 14.1, and $35.3\mu \mathrm{m}$ giving RMS intensity values of 2.5, 10, and $25\mu \mathrm{m}$, respectively. For each case, we compare step distances of $d=1.0,0.5,0.25\mu \mathrm{m}$. Each curve shows the overlap between the results using different step sizes, and the results of changing the step sizes are indistinguishable. We also see that the large size beams do not spread as much (by a factor of about 1.4) as the smaller size input beams (by a factor of about 10). These results also are consistent with previously reported results and will have ramifications with experiments regarding imaging through random fiber arrays as discussed in Sec. 5.

Fig. 3

RMS intensity radius as a function of propagation distance for different sized input beams. For each case, we compare step distances of 1, 0.5, and $0.25\mu \mathrm{m}$. The results when we vary the step sizes overlap and are indistingishable.

Next, we examined the effects where the difference between the indices of refraction for the two kinds of fiber was changed. Results are shown in Fig. 4 for a single random mask—again for step size of $d=1\mu \mathrm{m}$ and where the difference in the indices of refraction changed from $\mathrm{\Delta }n=0.005$ (top) to $\mathrm{\Delta }n=0.4$ (bottom) again with an average value of $n=1.5$. As expected, when the difference in indices was small, the beam is not confined very well. However, as the index differences increase, the confinement increases. These results are consistent with previous results^7,8 and will also have implications on imaging results, shown in Sec. 5.

Fig. 4

RMS intensity radius as a function of propagation distance for a given random array where the index of refraction difference varies from (top to bottom) 0.005, 0.01, 0.02, 0.025, 0.05, 0.1, and 0.4.

In an attempt to understand the spreading phenomena, we decided to examine the effects when we selectively illuminated the fast fibers (having the lower index of refraction) and when we selectively illuminated the slow fibers (having the higher index of refraction), again where the difference in the indices of refraction is $\mathrm{\Delta }n=0.1$.

Figure 5(a) shows the results for an input beam having an intensity RMS value of $10\mu \mathrm{m}$. The middle curve in Fig. 5(d) shows the RMS radius as a function of propagation distance and is the same as in Fig. 3. The beam width increases and saturates. Next, we multiplied the input Gaussian by a binary pattern that restricted it to the higher index fibers as shown in Fig. 5(b). Now the beam width saturates as a lower RMS radius as shown in the lower curve in Fig. 5(d). Alternatively, we multiplied the input Gaussian by a binary pattern that restricted it to the lower index fibers as shown in Fig. 5(c). Now the beam width saturated as a higher radius as shown by the upper curve in Fig. 5(d).

Fig. 5

(a) Incident Gaussian beam having an intensity RMS radius value of $10\mu \mathrm{m}$. (b) Input Gaussian beam multiplied by 1,0 random mask. (c) Input Gaussian beam multiplied by 0,1 random mask, and (d) RMS intensity radius as a function of propagation distance for the three cases above. Middle curve corresponds to Fig. 3. Lower and higher curves correspond to the cases where the incident beam is restricted to the higher or lower index fibers, respectively.

We can understand this effect physically by realizing that light is trapped in the higher index core of an optical fiber. However, there is an evanescent electric field that extends into the lower index region. If this evanescent electric field can reach another high index region, tunneling will cause the electric field to pass into the neighboring region. So the beam size will begin to expand through this process. However, the degree of the spreading depends on the randomness of the array in the region where the incident beam is located, as well as the difference between the two indices of refraction. Therese parameters affect the extent of the evanescent field and will determine the degree of the spreading. Surprisingly, the Fresnel diffraction process accomplishes the same effect by allowing the electric field to spread to neighboring fibers.

5.

Imaging

In more recent studies,^10–11 the imaging of patterns through these random fiber arrays has been studied. As would be expected, the broadening of the beams as shown in Fig. 3 would suggest that larger patterns with coarser features might be more faithfully transmitted compared with smaller beams having finer features. In addition, we would expect that the broadening can be reduced by increasing the difference in the indices of refraction as shown in Fig. 4. Our simulations are consistent with those reported earlier. Figures 6Fig. 7–8 show our results using the input number “2” having heights of roughly 120, 60, and 30 pixels.

Fig. 6

Imaging of an input pattern over 10 cm of random fiber. (a) Input pattern having a height of 120 pixels, (b) image after propagation over 10 cm of random fiber with an index difference of $\mathrm{\Delta }n=0.1$, and (c) image after propagation over 10 cm of random fiber with an index difference of $\mathrm{\Delta }n=0.4$.

Fig. 7

Imaging of an input pattern over 10 cm of random fiber. (a) Input pattern having a height of 60 pixels, (b) image after propagation over 10 cm of random fiber with an index difference of $\mathrm{\Delta }n=0.1$, and (c) image after propagation over 10 cm of random fiber with an index difference of $\mathrm{\Delta }n=0.4$.

Fig. 8

Imaging of an input pattern over 10 cm of random fiber. (a) Input pattern having a height of 30 pixels, (b) image after propagation over 10 cm of random fiber with an index difference of $\mathrm{\Delta }n=0.1$, and (c) image after propagation over 10 cm of random fiber with an index difference of $\mathrm{\Delta }n=0.4$.

Figure 6(a) shows a magnified view of an input pattern having a vertical height of 120 pixels, in rough agreement with the size from Fig. 3 in Ref. 11. Figure 6(b) shows the computed results for the image after propagation through 10 cm in a random fiber where the difference in the indices of refraction is $\mathrm{\Delta }n=0.1$. The image is very recognizable. Figure 6(c) shows considerable improvement when the index difference increases to $\mathrm{\Delta }n=0.4$. Note that these propagation distances are twice as long as in Ref. 11.

Figure 7(a) shows a more magnified view of an input pattern having a smaller vertical height of 60 pixels, in rough agreement with the size from Fig. 4 in Ref. 11. Figure 7(b) shows the computed results for the image after propagation through 10 cm in a random fiber where the difference in the indices of refraction is $\mathrm{\Delta }n=0.1$. The image is much more recognizable than the results shown in Ref. 11. Figure 7(c) again shows considerable improvement when the index difference increases to $\mathrm{\Delta }n=0.4$.

Figure 8(a) shows a much more magnified view of an input pattern having a much smaller vertical height of only 30 pixels, in rough agreement with the size from Fig. 5 in Ref. 11. Figure 8(b) shows the computed results for the image after propagation through 10 cm in a random fiber where the difference in the indices of refraction is $\mathrm{\Delta }n=0.1$. The image is again much more recognizable than the results shown in Ref. 11. However, the detail size of this object is expected to be more distorted since each feature of the input pattern would spread as seen from the results in Fig. 3. Figure 8(c) again shows some improvement when the index difference increases to $\mathrm{\Delta }n=0.4$.

Overall, we were surprised by these results. Again, note that we propagated the images through 10-cm distances compared with Ref. 11 where the propagation distances were 5 cm. However, we expect that the results would depend on such factors as the kind of image being transmitted, including the font type. Nevertheless, we draw two conclusions. First, the broadening of the beams as shown in Fig. 3 will ultimately limit the image sizes that are propagated. Second, as mentioned in Ref. 11, larger index differences are achievable with glass–air disordered fibers but would be difficult to manufacture.

6.

Conclusions

In this work, we have used an FFT-BPM algorithm to examine the propagation of light beams in an optical fiber consisting of a random array of smaller fibers having two different indices of refraction. This algorithm consists of having the beam pass through a binary phase mask representing the phase shifts of the two fibers and then calculating the Fresnel diffraction pattern through some propagation distance. This iterative process is repeated over long intervals. This algorithm was developed in two separate optical communities, with, as far as we can see, no overlapping references. However, the approach allows extremely rapid computation of difficult propagation processes. Our application of this approach to the random optical fiber problem produces results in excellent agreement with experimental studies and computational results using a much more time intensive (and more accurate) approach.

One can get an informal comparison of the computational times required between the two approaches. Reference 8 reports computation times of ${10}^5$ CPU hours using a transverse area of ${10}^5$ to ${10}^6$ ${\lambda }^2$, ${10}^6$ to ${10}^7$ longitudinal steps, and 100 random fiber arrays using a large HPC cluster with over 1000 nodes and 24 gigabytes of memory per node for the FD-BPM method. By contrast, using an iMac with a 2.66-GHz Intel Core 2 Duo processor, we examine an array of $256\times 256$ fibers, each having one pixel, and ${10}^5$ steps in about 15 min! Nevertheless, we find results that are consistent with the much more accurate, but also much more computationally intensive, approach. In addition, we were able to examine the effects of selectively exciting either the high or low index fibers in the random array.

In conclusion, this approach would allow the capability to examine a large number of variables before applying the much more accurate, but more time consuming, approach.