Analysis of the use of tapered graded-index polymer optical fibers for refractive-index sensors

The behavior of tapered graded-index polymer optical fibers is analyzed computationally for different refractive indices of the surrounding medium. This serves to clarify the main parameters affecting their possible performance as refractive-index sensors and extends an existing study of similar structures in glass fibers. The ray-tracing method is employed, its specific implementation is explained, and its results are compared with experimental ones, both from our laboratory and from the literature. The results show that the current commercial graded-index polymer optical fibers can be used to measure a large range of refractive indices with several advantages over glass fibers. ©2008 Optical Society of America OCIS codes: (060.2370) Fiber optics sensors; (250.5460) Polymer waveguides; (080.2720) Mathematical methods (general). References and links 1. J. Villatoro, D. Monzón-Hernández, and D. Luna-Moreno, “In-line optical fiber sensors based on cladded multimode tapered fibers,” Appl. Opt. 43, 5933-5938 (2004). 2. G. Shangping and S. Albin, “Transmission property and evanescent wave absorption of cladded multimode fiber tapers,” Opt. Express 11, 215-223 (2003). 3. S. Xue, M.A. van Eijkelenborg, G.W. Barton, and P. Hambley, “Theoretical, Numerical, and Experimental Analysis of Optical Fiber Tapering,” J. Lightwave Technol. 25, 1169-1176 (2007). 4. O. Ziemann, H. Poisel, A. Bachmann, and J. Vinogradov, “Special problems measuring POF,” in POF Modelling: Theory, Measurement and Application, C.-A. Bunge and H. Poisel, ed. (Books on Demand GmbH, Norderstedt, Germany, 2008). 5. J. Zubia and J. Arrue, “Plastic Optical Fibers: An Introduction to Their Technological Processes and Applications,” Opt. Fiber Technol. 7, 101-140 (2001). 6. J. Munisami and D. Kalymnios, “High NA POF performance versus the requirements of the recent standard ISO/IEC JTC 1 FDIS 24702,” in Proceedings of 15 International Conference on Plastic Optical Fibers and Applications POF’2006, (Korea, 2006), pp. 102-109. 7. M. Lomer, J. Arrue, C. Jáuregui, P. Aiestaran, J. Zubia, and J.M. López Higuera, “Lateral polishing of bends in plastic optical fibres applied to a multipoint liquid-level measurement sensor,” Sens. Actuators A 137, 68-73 (2007). 8. J. Arrue, G. Aldabaldetreku, G. Durana, J. Zubia, and F. Jiménez, “Computational research on the behaviour of bent plastic optical fibres in communications links and sensing applications,” in Recent Research Developments in Optics, S.G. Pandalai, ed. (Research Signpost, Kerala, India, 2005), Chap. 5. 9. A. W. Snyder and J. D. Love, Optical waveguide theory (Chapman & Hall, New York, 1983). 10. Chromis FiberOptics Co., “Chromis Fiberoptics,” (Head office in 6 Powder Horn Dr., Warren, NJ 07059, USA). http://www.chromisfiber.com. 11. FiberFin Inc., “FiberFin,” (Head office in 201 Beaver Street, Yorkville, Illinois, USA.). http://www.fiberfin.com. 12. Asahi Glass Co. Ltd., “Asahi Glass,” (Lucina Division; Head office in 1-12-1, Yurakucho, Chiyoda-ku, Tokyo 100-8405, Japan). http://www.lucina.jp/eg_lucina/productsengf2.htm. 13. A. Sharma, D. Vizia, and A. K. Ghatak, “Tracing rays through graded-index media: a new method,” Appl. Opt. 21, 984-987 (1982). #99997 $15.00 USD Received 11 Aug 2008; revised 21 Sep 2008; accepted 22 Sep 2008; published 2 Oct 2008 (C) 2008 OSA 13 October 2008 / Vol. 16, No. 21 / OPTICS EXPRESS 16616 14. F. Jiménez, J. Arrue, G. Aldabaldetreku, and J. Zubia, “Numerical Simulation of Light Propagation in Plastic Optical Fibres of Arbitrary 3D Geometry,” WSEAS Trans. Math. 3, 824-829 (2004). 15. T.A. Birks and Y.W. Li, “The shape of fiber tapers,” J. Lightwave Technol. 10, 432-438 (1992). 16. J. Arrue, F. Jimémez, M. Lomer, G. Aldabaldetreku, G. Durana, and J. Zubia “Characterization of tapered, polished or uncladded SI and GI POF geometries for use in tapers and multipoint sensors,” in Proceedings of 15 International Conference on Plastic Optical Fibers and Applications POF’2006, (Korea, 2006), pp. 187-192. 17. C. McAtamney, A. Cronin, R. Sherlock, G.M. O’Connor, and T. J. Glynn, “Reproducible method for fabricating fused biconical tapered couplers using a CO2 laser based process,” in Proceedings of 3 rd International WLT Conference on Lasers in Manufacturing, (Munich, 2005), pp. 673-678. 18. J. Zubia, J. Arrue, G. Fuster, and D. Kalymnios, “Light power behavior when bending plastic optical fibers,” IEE P-Optoelectron. 145, 313-318 (1998). 19. J. Mateo, I. Garces, and A. Losada, “A novel technique to fabricate low loss POF tapers,” in Proceedings of 9 International Conference on Plastic Optical Fibers and Applications POF’2000, (Boston, 2000), pp. 7276. 20. M. Kezmah and D. Donlagic, “Multimode all-fiber quasi-distributed refractometer sensor array and crosstalk mitigation,” Appl. Opt. 46, 4081-4091 (2007).


Introduction
Graded-index (GI) optical fibers with appropriate geometries can be used as refractive-index sensors, even without removing the cladding.For sensing purposes, light attenuation should vary when the fiber is immersed in fluids of different refractive indices, which only requires using fibers without protective jacket (no more than core and cladding).Since there is no step between the core and the cladding refractive indices, light rays are not reflected at the corecladding interface.However, there is a step between the cladding and the external medium.Light rays follow curved paths in the core, which, in certain cases, reach the interface with the cladding, the trajectory becoming a straight line inside it.At the interface with the external medium, rays can get lost by refraction, if the angle of incidence on this surface is lower than the critical one.If there are refracted rays, the attenuation in the fiber will depend on the refractive index of the medium surrounding it.To ensure a noticeable loss of power, fiber optic tapers can be employed [1][2][3].
The main novelty in this paper is that we analyze this sensing principle for the case of GI polymer optical fibers (GI POFs) instead of glass ones.POFs have many practical advantages over glass optical fibers: on the one hand, the material is much more flexible, which allows much larger diameters.This facilitates handling and installation, with low-cost connectors and components and less strict geometrical tolerances.On the other hand, they have high numerical apertures (which facilitates the use of inexpensive light sources) and also low attenuation in the visible region and, in some cases, even in the second window (1300 nm).Specifically, with the introduction of perfluorinated poly(methyl methacrylate) (PMMA), the maximum transmission distance due to attenuation has been extended from 100 m to over 1 km [4].However unimportant this may be in the sensing head itself, an already installed POF-based local area network may include remote sensing among many other services, so POFs of low attenuation are still desirable.
For these and other reasons, POFs have been attracting a considerable amount of attention.For example, they are often used as transducers in sensing applications requiring no more than the POF itself, an inexpensive LED and a photodiode, the sensing principle often relying on the modulation of light intensity, so sophisticated signal interrogation techniques are not necessary.Besides, the set-up is often easier [5][6].One advantage of sensors based on GI POFs, in comparison with those achieved with other solutions with step-index (SI) POFs, is that the former do not require removing part of the cladding, which may be difficult in practice [7,8], while GI POFs retain the advantages already commented for POFs in general over glass fibers.
In this paper we also describe a new implementation of the ray-tracing method for the analysis of tapered GI fibers.A previous method [1] for the specific case of glass GI tapers consisted in using a numerical integration technique proposed in [2] for the case of step-index tapers.This may be admissible for GI fibers of low numerical aperture, but it could not be employed for the high-numerical-aperture GI fibers analyzed by us, as will be explained later.
This method serves us to analyze tapers with the presently three existing types of commercial GI POFs.Experimental and theoretical results available with sensors using tapered graded-index glass fibers are employed to validate the numerical calculations.
In theory, the maximum outer refractive index that can be measured in all cases coincides with that of the cladding (n 2 ), since all rays escaping from the core into the cladding would also escape from the cladding if the outer medium's refractive index is greater than or equal to n 2 .Since n 2 depends on the material used to manufacture the fiber, which cannot be arbitrarily changed, each fiber will have its own upper bound of measurable refractive indices.In Table 1 we summarize the main characteristics of the commercial fibers studied, where n is the refractive index of the outer medium, ρ 0 is the core radius, ρ 0,clad is the radius of the cladding, n 1 and n 2 are the core and the cladding refractive indices respectively, and g is the refractive index profile exponent [9].Note that the OM-Giga and the GigaPOF fibers have the same characteristics, despite being from different manufacturers [10,11].The characteristics of the Lucina fiber (the other commercial GI-POF) can be found in [12].

Implementation of the ray-tracing method
The use of a computer model based on the ray-tracing method is one of the most flexible ways to calculate the propagation of light along a highly multimode optical fiber.It is based on a geometric-optics approach that uses light rays in combination with Fresnel transmission coefficients for the calculation of radiation losses along the fiber, which, in the case of gradedindex fibers, may occur either when the ray reaches the outer interface of the cladding or at turning points, by means of a tunneling mechanism [9].A good simplified approach consists in neglecting tunneling effects, considering only refraction losses and using a transmission coefficient that is 1 whenever refraction occurs and 0 otherwise [8,9].At reflection points, both refraction and tunneling loss due to the evanescent field can occur, but the latter is negligible unless the incidence angle α i is only slightly greater than the critical angle α c .On the other hand, when α i < α c , refraction occurs, the Fresnel transmission coefficient T being approximately 1 if α i /α c < 0.9 [8].When α i is closer to α c , a precise value of T is given by a generalized Fresnel transmission coefficient that takes into account the local curvature at the interface, T being less than 1.Therefore, by using T = 1 whenever α i < α c and T = 0 otherwise, we are obtaining a value of attenuation slightly greater than the real one.Thus, our results provide an upper limit for the attenuation of the tapers considered, serving this value as a conservative reference, in the sense that in practice slightly smaller attenuations can be expected.
The ray paths in the case of GI fibers are usually calculated numerically, since, very often, it is not possible to find an analytical expression for the ray path considered.
The mathematical expressions for the definition of the geometrical and optical parameters for the numerical resolution the eikonal equation in graded-index tapered multimode fibers, not previously available in the literature, will also be provided in this section.

Discretization of the eikonal equation
In this section we will explain how to apply the ray-tracing method inside a graded-index tapered fiber with cladding.The path of a light ray propagating in a graded-index medium is a curve governed by the eikonal equation [9]: where r is the position vector of a generic point in the curve, s is its natural parameter (distance from any given point r 0 along the ray path), n(r) is the refractive index as a function of the point in space, and ∇n(r) is its gradient.Note also that dr/ds is the unit tangent t unit to the ray path at r. Figure 1 illustrates this.The eikonal equation is valid for any geometry as long as n(r) is sufficiently smooth.At points where n is discontinuous, reflection or refraction occurs as governed by Snell's law.When an analytic or exact solution for the ray path does not exist, as usually happens, the eikonal equation must be solved numerically in order to obtain the trajectory of the light ray approximately.
Furthermore, in many practical implementations of the ray-tracing method it is critical that the approximate solution of the eikonal equation be obtained in a very efficient way from the numerical point of view, i.e. with a computational effort as small as possible for any given accuracy.This is due to the fact that the ray-tracing computational simulation of light propagation along a POF typically involves tracing the paths of several thousands of individual rays [8], the trajectory of each one of them being discretized into many thousands of individual steps.Some very efficient numerical methods to solve the eikonal equation are provided in [13].By means of the following change of variable, the differential equation is first rewritten in a form that is more suitable for numerical integration: Defining a tangent t to the ray path as t = dr/dt = n t unit , we finally write: This set of differential equations determines the ray path from the initial conditions r 0 (position vector of the ray at the beginning of its trajectory) and t 0 , which is n(r 0 ) times the unit tangent t unit to the ray at r 0 .Note that, unlike the eikonal equation, Eq. ( 3) is very similar in structure to a set of ordinary differential equations u" = f(t,u) with initial conditions (t 0 ,u 0 ), but without t appearing on the right-hand side of the equations, and with three components in each equation (x, y and z), which can be handled using matrix notation.Therefore, starting with (t 0 ,r 0 ), a succession of points (t 1 ,r 1 ), (t 2 ,r 2 )… can be generated by means of any of the standard algorithms for solving sets of ordinary differential equations with initial conditions.For our calculations, we have used Runge-Kutta's third-order algorithm [13], with which we have obtained an excellent trade-off between accuracy and computational cost: where Here, h is the step size chosen to run the algorithm.The smaller h is, the higher the accuracy is, but also the computational cost.Since h cannot be arbitrarily low, it must be selected carefully in order for the global (accumulated) error at the end of ray paths to be acceptable.
The selection was made by comparing the numerical solutions obtained with different step sizes for a given ray path.By using a very exact ray path as reference, i.e. one obtained numerically with a very small step size h (close to 10 −6 m, with which rounding errors are still negligible using double-precision arithmetic), we found that a value of h on the order of 10 −5 m was sufficient for the global error to be negligible in the typical lengths of a few centimeters needed in our studies.For the calculation of power loss, only radiation loss through refraction will be taken into account, since tunneling losses are much smaller [9] and material absorption is very small in a few centimeters (in addition, it does not affect the sensitivity -slope of the curves-of the tapers analyzed).A number of rays between 5 000 and 15 000 has been used in our computational simulations, depending on the taper geometry considered.This number was chosen by comparing a few results with those obtained with many more rays.In this way, we saw that the number of significant digits that we can reach in the quotient between output and input powers is two, which is accurate enough for our purposes (e.g., in one of the geometries studied we obtained 0.316 with 5 000 rays, 0.333 with 15 000 rays and 0.336 with 30 000 rays).

Particular case of a tapered GI fiber with core and cladding
For the sake of programming flexibility, it is convenient to define the geometry of the fiber in a separate software module, in such a way that the ray-tracing module be valid for any geometry.The procedure we have used is very similar to that described in a previous work presented by some of us about the simulation of very flexible three-dimensional geometries with step-index fibers [14].The main difference now is that, since we are dealing with gradedindex fibers, the geometry module must also return the refractive index and its gradient for any point r of the ray path.The position along the path is specified by means of a parameter p that determines the fiber cross section where the point of the ray path lies.For the type of geometries studied here, p was chosen to be the coordinate along the fiber symmetry axis (i.e.p = z in Fig. 2).The geometry module also returns some other geometrical parameters, such as the radii of the cross sections of the core and the cladding as functions of p.The geometry that will be used as sensing probe is shown in Fig. 2, which is the same as that described in [15] for the specific case of glass fibers, but we will show that the appropriate dimensions of the sensor differ from one type of fiber to another, and also that they depend on the range of refractive indices to be measured.The cladding of the fiber has not been depicted in the figure, since we will focus our attention on the core before considering the cladding as well.The radius of the cladding decreases or increases maintaining the same proportion with that of the core.Firstly, we will obtain mathematical expressions for the refractive-index profile n(r) and its gradient ∇n(r) (grad n) as functions of the x, y and z coordinates (see also Fig. 1).The surfaces of constant refractive index are represented in Fig. 2 by dotted lines, which are their intersections with the meridional plane.There are three parts that should be considered for the calculation of grad n: the narrowing section of length L, the central part of constant diameter (waist) of length L 0 , and the expanding section following it, of length L.
Regarding the narrowing section, if we take the origin of the z axis at the point where the fiber begins to be tapered, the core radius (interface with the cladding) follows a decreasing exponential given by 0 0 ( ) where ρ 0 and ρ min are the core radii of the undeformed fiber and of the waist, respectively, and L is the length of the section [1].This exponential shape is obtained when the heat source moves over a constant length L 0 (if the width of the heating zone were varied, a variety of taper shapes could be obtained).The density of the core can be assumed to remain constant after deformation, and the refractive index n to depend only on the dopant concentration.Under these assumptions, the surfaces of constant refractive index n are, in the 2-D view, exponentials proportional to that of the interface with the cladding.Therefore, if d is the variable distance from a point on the dotted line to the fiber symmetry axis, and r is its initial value at the entrance to the taper, we have: The explicit expression of n(d,z) is readily obtainable, since n(r) is a known function (a quasiparabolic profile with exponent g) given by where n 1 is the maximum refractive index, which takes place at the fiber symmetry axis, and n 2 is the minimum one, i.e. that of the cladding [9].The result is: Since d = (x 2 +y 2 ) 1/2 , the three-dimensional gradient of n is obtained as: , , , , which, after substituting d for (x 2 +y 2 ) 1/2 , yields grad n as an explicit function of x, y and z, as needed to solve the eikonal equation numerically.
For the waist we have , and the procedure is the same as for the narrowing section.
Regarding the expanding section following the waist, the result is nearly the same as that obtained for the narrowing section.There are only two differences: the third component of the gradient has the opposite sign, and z should be substituted by z − (L+L 0 ).
In the cladding, the gradient of the refractive index is null.Therefore, if a ray reaches the core-cladding interface its trajectory becomes a straight segment (Fig. 3).For such a ray, it is necessary to calculate the normal n to the surface at the point where it reaches the interface with the external medium, so as to be able to determine the incidence angle and to compare it with the critical angle of refraction [14].The unit normal n to the interface can be obtained analytically by taking advantage of the existing axial symmetry.For example, for the narrowing section the result can be obtained in two steps, by first considering the unit normal in two dimensions n 2D , i.e. in the meridional plane, which is given by (see Fig. 3): ( ) where the denominator normalizes the vector, and where the minus sign in the first component of n 2D is necessary to take the absolute value of the derivative dρ(z) / dz.The three-dimensional vector n can be obtained by rotating n 2D around the fiber axis, since its axial and radial components are equal to those of n 2D , resulting: This expression is also valid for the expanding section following the waist, since the minus sign in the first component of n 2D also appears in this section, this time in order to take the correct orientation backwards.

Tapers with sinusoidal shape
Depending on the manufacturing process (as will be seen later), the shape of a taper may be not exponential, but approximately sinusoidal, as shown in Fig. 4: L tot Fig. 4. Core of a taper with a sinusoidal shape.
The tapers produced in our laboratory in the way described below have a shape which is very approximately of this kind, as can be seen in Fig. 5: Fig. 5. Tapers manufactured by the authors.Proportions have been altered slightly, specially on the right (with 0.6-mm taper at the top, 0.4-mm one at the bottom) to increase cross-dimensions over axial ones, so that the taper shape can be appreciated better.
The implementation of the ray-tracing method for this particular case is analogous to the one described for exponential tapers.Some details are provided in [16].The pictures in Fig. 5 (right) show the tapers used for the experimental measurements of subsection 3.1.

Validation of the numerical algorithms and experimental results
The numerical algorithms described in the previous section have been implemented in a computer program.In order to validate its results, we have compared them with experimental measurements carried out in our laboratory with sinusoidal POF tapers manufactured by us (subsection 3.1) and also with experimental results reported in the literature for glass tapers of exponential shape (subsection 3.2).

Validation of the numerical algorithms by comparison with our own experiments
We have carried out experimental measurements in our laboratory with the two tapers on the right of Fig. 5.In this subsection we compare them with the corresponding attenuations calculated numerically for these geometries.
The waist diameters are 0.6 mm and 0.4 mm respectively, with an estimated dimensional tolerance of ±0.05 mm.If we call narrowing ratio (NR) the quotient ρ min /ρ max (which will be an important parameter in the rest of the paper), we have NR∈[0.55,0.65]for the 0.6-mm taper and NR∈[0.35,0.55]for the 0.4-mm one.The measurements of their attenuations have been carried out using a low-numerical-aperture light source (NA ≈ 0.16).Table 2 summarizes the experimental and the computational results obtained for both tapers: The computational and the experimental results agree satisfactorily well.The discrepancies may be attributed to a number of different reasons, like the practical impossibility of knowing the exact illumination conditions at the entrance of the taper, tolerances in experimental measurements, etc.

Comparison of our numerical results with experimental ones reported in the literature
Let us now consider the experimental results provided in [1] for a sensor based on a tapered GI glass optical fiber.The geometry of that sensor is the exponential one described in the previous section, with the manufacturing process imposing a relationship between L and L 0 , namely L = L 0 ln(1/NR) [1,15,17].In the present paper all the calculations with exponential tapers are carried out taking this relationship into account, e.g.whenever the narrowing ratio is changed for a given L 0 , the length L is also conveniently adjusted.
The exponential shape occurs whenever the heat source moves over a constant length L 0 , whatever its speed, inasmuch it is much larger than the elongation speed.The length L becomes larger when the narrowing ratio NR = ρ min,clad / ρ 0,clad decreases [1], ρ min,clad being the radius of the cladding at the waist (the narrowing ratio applies both to the core and to the cladding).
Figure 6 shows the influence of the outer refractive index n on the fraction of power reaching the output end of the sensor, which has been either calculated by using two different methods based on ray-tracing, or measured in the laboratory [1].In all cases, the optical fiber is the standard graded-index glass one indicated in Table 1
In Fig. 6 we have included the experimental results reported in [1] corresponding to NR = 40/62.5= 0.64, which is the highest value used in the measurements, together with our own computational results for the same value of NR, obtained as described in the previous section.Although our aim is to analyze GI-POF tapers instead of glass ones, in this way we have had an experimental reference with which to be able to compare some of our results.
Although there were other experimental results available for lower values of NR, such as NR = 0.48 or lower, we have not used them in order to ensure a sufficiently high number of modes propagating along the fiber for a reliable application of the ray-tracing method [5,9].Since the number of modes M inside the core can be approximated by [9], M ≈ 160 at the waist of the taper in the case of the experimental results in [1] with NR = 0.64 and λ = 850 nm.This number is on the verge of validity of the ray-tracing method; however, the real number of modes propagating inside the fiber at the waist of the taper increases considerably if n < n 2 = 1.459, since total internal reflections at the external interface of the cladding serve to guide modes along the fiber.In the case of a POF we will deal with thousands or millions of modes, so the applicability of the method will be out of doubt.In Fig. 6 we have also included the results calculated in [1] for a graded-index glass taper with the narrowing ratios closest to the experimental one of 0.64 (namely, for NR = 0.7 and NR = 0.6).These results were calculated by using a numerical integration technique proposed in [2] for the case of step-index tapered fibers, in which ray paths are straight (as in Fig. 1 of [2]).This employment of straight rays may be admissible for GI fibers of low numerical aperture, in which light rays, even if they are curved, propagate almost parallely to the fiber symmetry axis, but it is not admissible for high-numerical-aperture fibers, like some of the ones analyzed by us.Our technique approximates the exact, curved ray paths in the taper, which requires both more elaborate numerical methods and higher computational costs, but it will serve us to deal with POFs.In the case of Fig. 6 (glass taper) both methods yield reasonably similar results (compare the gray shaded region limited by NR = 0.7 and NR = 0.6 with the squared blue line).When n decreases and approaches 1.43, our results agree better with the experimental ones.One of the possible reasons for the agreement not to be complete is that it is impossible in practice to reproduce the exact illumination conditions at the entrance of the device, since these conditions change from one LED to another and depend on the distance from the LED to the fiber, etc.Our results have been obtained with a light source of 10 000 rays initially parallel to the fiber axis and uniformly illuminating its whole input cross section.The length of the waist is L 0 = 5 mm.
In order to facilitate the design of sensors with this type of geometry employing any type of fiber, let us first consider the trajectory of light rays graphically.Figure 7 shows, for the same launched ray, the values of some incident angles at the external surface of the taper, illustrating the influence of its length L. In Fig. 7(b) L is larger than in Fig. 7(a) (5 mm against 2 mm), but the waist radius is the same in both cases (ρ min,clad = 20 μm), as well as the optical fiber employed (the glass fiber of Table 1).We can notice that the ray path from the entrance of the taper to the point where it enters the cladding has more oscillations for L = 5 mm (Fig. 7(b)) than for L = 2 mm (in both cases it is approximately sinusoidal with a decreasing amplitude).Moreover, in a hypothetical case in which the narrowing section of the taper were much shorter than 2 mm, the ray would tend to reach the cladding without oscillations, because, in such a case, the surfaces of constant refractive index would be almost perpendicular to the ray, and the radial component of the gradient of the refractive index would decrease considerably.As a consequence, the light power attenuation would be much higher.In contrast, for moderate lengths the total refraction loss becomes greater when L increases.This indicates that the amplitude of the oscillations of the ray path does not decrease as rapidly as the diameter of the fiber.As for the expanding section following the waist of the taper, the attenuation tends to be small, because the diameter of the fiber increases with z.Anyway, there can be refracting rays in this section as well.For example, this could be the case of the ray plotted in Fig. 3 for a sufficiently high outer refractive index n.Moreover, a significant amount of attenuation can occur in the waist.As was shown in Fig. 3, rays can enter the cladding in the waist, and our computations have shown them not to be negligible.As soon as a ray has crossed the core-cladding interface, it is susceptible of escaping if the refractive index of the outer medium is high enough.
It is clear that very short tapered sections (values of L close to 0) tend to yield high losses, because in such small distances many rays cannot bend sharply enough to avoid the corecladding interface standing almost perpendicularly to their trajectories.As a consequence, losses decrease as L increases for small values of L (and a given narrowing ratio).However, this tendency does not continue forever, and losses in the tapered section start to increase again when L becomes large enough.Figure 8 helps to explain this effect.It shows two tapered sections with the same narrowing ratio (NR = 0.3) but with different lengths L, namely 2 mm and 10 mm.In both cases, the core of the tapered fiber has been depicted together with the trajectories of some light rays.These rays are identical at the entrances of both geometries (the 2-mm taper and the 10-mm one).Therefore, we can compare the number of rays that reach the core-cladding interface in both cases.Let us consider, to facilitate the explanation, that any ray entering the cladding will escape from the fiber (i.e. that the refractive index of the outer medium is high enough).The ray paths have been obtained numerically (considering that the fiber is the Lucina of Table 1).In this way, we can notice that 5 of the rays considered escape when L = 2 mm, whereas this number is greater (7) when L has been increased to 10 mm.Although the depicted rays have been chosen arbitrarily, the fact that losses can increase when L increases is not a coincidence.In fact, we can notice that the rays that escape are those whose initial distance to the fiber symmetry axis is greater than a certain threshold value for each taper, which establishes the limit between escaping or not.This distance is smaller in the 10-mm taper considered than in the 2-mm one, and it would decrease slowly if we chose larger values of L (although the rate of change would tend to decrease).The dashed lines in Fig. 8 also help to understand what is happening.
On the other hand, we have seen that the length L 0 is not very important if it is large enough.This can be easily explained: because of the translational symmetry in the straight section, the trajectory of each ray becomes periodical in shape [9], so the ray will be lost by refraction within one period, or never.

Sensor design with a commercial GI POF
We will use the same type of geometry as that employed for glass graded-index fibers, but we will analyze the dimensions that would be appropriate for the sensor in the case of utilizing any of the commercial GI POFs available nowadays (namely the GigaPOF from Chromis [10], the OM-Giga from Optimedia [11], and the Lucina from Asahi Glass [12]).Let us begin with the GigaPOF / OM-Giga fibers (it does not matter which, since they both have the same optical and geometrical characteristics: see Table 1).We have also carried out simulations for the Lucina fiber.We will choose L 0 = 5 mm for the length of the waist, which is a value that has already been reported in the manufacturing process [1]. Figure 9  We have found that this higher loss of the Lucina fiber is due to its much thicker cladding (the core-to-cladding diameter ratio is 0.24 in the Lucina and 0.9 in the OM-Giga / GigaPOF).On one hand, we have calculated the power loss in a hypothetical fiber with the same characteristics as the Lucina but with a smaller cladding diameter.Specifically, we have taken a diameter of 66 μm instead of the real value of 245 μm, in order to have the same 0.9 core-tocladding diameter ratio.For such a fiber, we have obtained P out /P in = 0.31, which almost coincides with the result corresponding to the OM-Giga / GigaPOF.On the other hand, the fact that a thicker cladding yields higher losses has a simple geometrical explanation, which is illustrated in Fig. 11.This figure shows the trajectory of a light ray from the core-cladding interface (point P 1 on the first contour line) until it reaches the outer surface of a thin cladding (point P 2 on the second contour line) or, alternatively, of a thick one (point P 3 on the third contour line).Since a ray within the cladding propagates following a straight line, its angle of incidence onto the interface between the cladding and the outer medium depends only on the normal to it (n 2 or n 3 ).Figure 11 shows an angle of incidence α 3 that is typically smaller than α 2 (for the geometry represented, α 2 = 71.3ºand α 3 = 67.4º).Therefore, a ray that would reflect at P 2 if the cladding were narrow, could refract and get lost at P 3 if the cladding were thicker.This turns out to be the dominant effect explaining the higher losses for thicker claddings. )and a thick one (P 3 ), with the corresponding angles of incidence α 2 and α 3 .
The choice of NR influences the range of measurable refractive indices with each type of fiber.For example, if NR = 0.4 we can measure values from n = 1.34 (α c = 90º) to n = 1.10 (α c = 55º), whereas, in the case of the OM-Giga / GigaPOF fibers, the range is smaller and extends from n = 1.492 (α c = 90º) to n = 1.352 (α c = 65º).If the narrowing ratio were 0.7 instead of 0.4, we could only measure refractive indices down to 1.469 (α c = 80º) in the case of the OM-Giga / GigaPOF fibers, and down to 1.259 (α c = 70º) in the case of the Lucina fiber.The maximum measurable refractive index is always n = n 2 .Therefore, the range of measurable refractive indices can be adjusted at will up to n 2 , although the narrowing ratio should not be too close to 1 for losses to occur.Note that, if the range of attenuations does not change much, any increase in the range of the measurable refractive indices must be achieved at the expense of the sensitivity of the sensor to variations of n, and vice-versa.Another advantage of this type of sensor over other solutions (e.g.bends with SI fibers of thick cladding [18]) is that it is not necessary to remove the cladding, since light rays enter the cladding without refraction.For tapered structures, we can now compare the performance of GI POFs with that of GI glass fibers.We can observe that, for a given narrowing ratio, the range of refractive indices that can be measured is larger in the case of a Lucina POF than in the case of the glass fiber reported in [1].For example, in Fig. 10 we can notice that, for NR = 0.7 and a Lucina fiber, the range of refractive indices extends from 1.34 (α c = 90º) to 1.259 (α c = 70º), so the difference is Δn (Lucina,NR=0.7)= 0.081, whereas, in the case of a glass fiber with NR = 0.64 only, the difference is Δn (glass,NR=0.64)= 0.029 (see the computational curve in Fig. 6).As for the OM-Giga / GigaPOF fibers, the range is not so large because the cladding is much thinner, as explained before.In any case, the OM-Giga / GigaPOF could be preferred over the Lucina fiber when the advantage of having a large fiber diameter is an important issue, since their diameter is larger than that of the Lucina and, to a greater extent, larger than that of the conventional glass GI fibers.Notice that, even in the hypothetical case of having two fibers with the same geometry, profile exponent and numerical aperture, the maximum measurable refractive index would be different if the cladding is of a different material (n 2 is not a parameter that can be changed for a given material).
Manufacture of PMMA POF tapers by drawing and heating requires taking care of the fact that this material is thermoplastic with a low softening temperature, but the technique has already been put into practice satisfactorily.For example, a method reported in [19] and illustrated graphically in Fig. 12 has served to manufacture a taper at any position in the fiber, both in step-index POFs and in GI POFs.The manufacture of PMMA tapers has also been reported recently in [3].Finally, it is worth mentioning that glass tapers have also been made by etching the fiber instead of heating and pulling it [20].In this way, refraction losses are greater, since ray paths inside the core are unaltered by the etching, which just brings the outer interface closer to the fiber symmetry axis.In consequence, the range of refractive indices that can be measured is larger for the same value of NR.For example, for NR = 40/62.5= 0.64 the difference between the maximum and the minimum measurable values of n in an etched glass taper is Δn (glass,etched) = 0.035, whereas this difference is smaller with a stretched glass taper, namely, Δn (glass,stretched) = 0.029.All in all, in both cases, the range is much smaller than that obtained by us for the Lucina fiber with the same narrowing ratio, which is Δn (Lucina,NR=0.64)= 0.13.

Conclusions
We have shown how to implement a method based on ray tracing to analyze the propagation of light and its power attenuation in tapered multimode graded-index polymer optical fibers.With this implementation, tested with the aid of experimental results existing for glass fibers, we have shown that similar geometries can serve to construct refractive-index sensors using commercial GI POFs (known to have distinct advantages over glass fibers).In the case of the Lucina POF, the narrowing ratio does not need to be very small for a large range of refractive indices to be measurable, but the maximum refractive index that can be measured is smaller than in the case of the OM-Giga / GigaPOF fibers.For the same narrowing ratio, the range of refractive indices that can be measured with a tapered Lucina POF is greater than with a GI glass fiber.

Fig. 1 .
Fig. 1.Variables in the eikonal equation of the ray path.

Fig. 2 .
Fig. 2. Cross section of the sensing probe.The surfaces of constant refractive index are represented by dotted lines.The undeformed core radius is ρ 0 .

Fig. 3 .
Fig. 3. Rays follow curved paths in the core and straight ones in the cladding.Calculating the normal n at interfaces is needed in all three sections of the taper.

Fig. 7 .
Fig. 7. Ray paths inside a taper, with some angles of incidence on the external surface.Angles appear distorted because of the different scales employed in the horizontal and vertical axes for a better visualization.The internal(external) exponential curves represent the contour of the core(cladding) in the meridional plane.The oscillating curve illustrates the trajectory of a meridional ray for the same launching conditions in two cases: (a) Length L = 2 mm; (b) L = 5 mm.

Fig. 8 .
Fig. 8. Example to illustrate that, for certain values of the taper length L, power loss can increase as L increases for a given narrowing ratio.The fiber is a Lucina one with ρ 0 = 60 μm.The rays entering the taper and the narrowing ratio (NR = 0.3) are identical in (a) and (b), but in (a), with L = 2 mm, 5 rays escape, while in (b), with L = 10 mm, 7 rays escape.

is a 3 - 9 (Fig. 10 .
Fig.10.Influence of NR on the sensitivity (slope) of the sensor and on the range of measurable refractive indices (critical angles), for the OM-Giga / POF and the Lucina fibers.

Fig. 12 .
Fig. 12. Manufacturing of a taper in a plastic optical fiber by heating and pulling.

Table 1 .
Different possible fiber candidates for a refractive index sensor.

Table 2 .
Comparison between author's experimental and computational results.