Testing a novel μ-Raman method to retrieve refractive index and density field from femtosecond laser-written optical waveguides

In this work, we analysed the refractive index field and other local material properties of femtosecond laser-written waveguides properly combining a novel and direct Raman strategy, waveguides coupling and instrumented nano-indentation. At first, we measured a 2D Raman map within the cross section of femtosecond laser-written waveguides. Then, the Raman Phonon Deformation Constants of lithium niobate were employed to retrieve the strain and density variation from the A1(TO) phonon shifts in the analysed region. We test the results obtained with different combinations of phonons by computing the numerical guided modes and comparing them with those experimentally measured. As a relevant finding, we found that the combination of the A1(TO)1 and A1(TO)4 phonon shifting is the most proper one to compute strain, density and refractive index variation, almost in this kind of waveguide. Finally, a linear path across waveguides cross section was explored with instrumented nano-indentation and the expected variation of local density was detected through a softening of the elastic module observed in the region directly modified by the ultra-fast laser.


Introduction
Femtosecond laser writing is an important tool to fabricate 2D and 3D optical circuits; this is a relative new technique that can be applied in almost every transparent material. It consists of producing a tightly focusing and a precise displacement of the femtosecond laser inside of a transparent sample to generate a permanent refractive index increase that allows optical waveguide formation [1][2][3]. Employing an adequate level of pulse energy (around 1 lJ) and 200 fs of pulse duration, a phenomenon that includes a plasma generation and a permanent breakdown of the material is triggered in the focal volume region. As a consequence, an expansion that affects the leading regions takes place and the refractive index suffers a permanent increase within them; hence the light can be guided in these regions. In the whole physical process of waveguide formation, many complex ultra-fast sub-processes occur in the focal region and the resulting material structure inside the femtosecond laser interaction zone and within the leading regions it is not well known yet. In the last years, much effort has been done to retrieve the strain and refractive index field using different simulations and measurements. As a result, residual strain was determined as a relevant cause of wave guiding formation [4] and l-Raman and l-luminescence maps were widely employed to analyse femtosecond laserwritten structures [5][6][7][8][9][10]. In particular, iterative approaches using an elastic and numerical model were employed to compute the residual strain and the refractive index field [11]. Also, the empirical relation between phonons shifts and the strain field, allowed the computation of the Raman Phonon Deformation Potential Constants (PDPC) in a previous work [12].
A particular disadvantage of using numerical elastic models and making successive iterations to retrieve the refractive index field, is to be time consuming. To overcome this fact, we present and evaluate in this work, a direct and simple method to compute the strain field, and so the refractive index field using only the PDPC of the crystalline material and a Raman map of the studied region. To evaluate the resultant refractive index profile, we computed the guided intensities and compared these with the experimental guided modes. Complementary, the expected density variation in the analysed region was detected by a softening of the elastic modulus. This was measured using a nano-indentation technique in the modified region.

Experimental setup
Double track optical waveguides were fabricated in an xcut lithium niobate crystal doped with Nd (0.3%) and Mg (5%) by femtosecond laser inscription. The Ti:sapphire ultra-fast laser with a pulse duration of 120 fs and a repetition rate of 1 kHz was centred at a wavelength of 800 nm. The energy of the laser pulse was 1 lJ and its linear polarization was aligned with the crystal z-axis. The pulsed laser was focused 350 lm below the crystal surface with a 10x (NA = 0.25) microscope objective and it was used a writing speed of 30 lm/s. A scheme of the fs laser writing process and an image of the performed laser tracks taken with an optical microscope can be seen in Fig. 1.
After waveguides fabrication, the end faces were carefully polished with diamond powder with controlled size up to 0.25 lm. After that, a 1550 nm solid state continuous laser was coupled to the waveguides by using the wellknown ''end fire'' coupling method. This laser was launched in the front edge of the waveguides using an optical fibre attached to a micrometre translation stage. In the ending edge of the waveguides, the near field was amplified using a 10x (NA = 0.25) microscope objective and collected, after a Glam Thomson Prism (GTP) polarizer, with a standard Newport beam profiler analyser. We present a scheme of this method in Fig. 2.
Following, a 2D l-Raman map was measured in the generated waveguides cross section by using the y(zz)y back-scattering configuration. For this scanning, a 200 mW Argon laser, centred at a wavelength of 532 nm, and a Confocal Olympus BX-41 with an objective of 50x (NA = 0.25) were used. We focused the laser 7 lm below the crystal surface. The steps of the map were 0.7 lm in zand x-direction. The motorized translators had a spatial resolution of 0.1 lm. The well-known Raman spectrum of lithium niobate for the mentioned configuration ( Fig. 3a) was collected in each position within the studied region. The A 1 (TO) 1 , A 1 (TO) 2 , A 1 (TO) 3 and A 1 (TO) 4 phonons were processed with LABSPEC5 software: by fitting these peaks with Lorentzian functions, the central wavenumber shifts were retrieved for each point of the scanned area. A scheme of this procedure is shown in Fig. 3a, b. No signal was detected in the wavenumber shifting of A 1 (TO) 3 phonon, this should be related with the low signal to noise ratio corresponding to this phonon (Fig. 3a).
Finally, a linear map of indentations was performed in the path ab across the waveguides (as scheme is presented in Fig. 3c) using a TriboIndenter Hysitron TI 950 Nanoindentation equipment. This system allows a controlled indentation process to obtain a precise load versus normal displacement curve (p vs. h, Fig. 3c). From this information, the elastic modulus was computed by the TriboIndenter software in each point of the path ab using the well-known Oliver-Pharr method [13]. The resolution of the nano-indentation equipment in normal displacement is \0.2 nm and in load is \30 nN. The translation resolution of the motorized platform in the plane axes is 0.5 lm. The indentations were performed using a spherical indenter of 2 lm of radius with a max load of 9000 lN and with zero holding time. The indentations were separated by a distance of 1.3 lm.

Method and results
In the studied method, we obtained the refractive index field from the wavenumber shift of the Raman A 1 (TO) phonons. To achieve this end, two important steps were followed. At first, we computed the local strain components-by using the Raman Linear Potential Theory-from 20 m Fig. 1 a Scheme of waveguides fabrication by using fs laser writing and b images of the cross section of waveguides surface taken by an optical microscope. The arrows at the top indicate the propagation direction of writing laser. The circular dotted curves indicate the main waveguiding regions the Raman phonon shifts measured in each point of the studied region. Secondly, the local variation of the refractive index field-by using piezo-optic relation-was retrieved from previously obtained strain field.
In lithium niobate crystals, by assuming the Raman Linear Potential Theory [11], the wavenumber shifting (Dx i ) and the normal strain components (e x , e y , e z ) are related between each other by the following linear equation, valid for each A 1 (TO)i phonon: where e i and f i are the lithium niobate piezo-spectroscopic constants of the four A 1 (TO 1 )i phonons and they are: e 1 = -270 ± 70, e 2 = -230 ± 60, e 4 = -550 ± 130, and f 1 = -13 ± 5, f 2 = -215 ± 50, f 4 = -160 ± 40 [12] (e 3 and f 3 are still unknown).
It can be seen that Eq. 1 assumes that the phonon wavenumber shifting is equal to a weighted linear combination of the deformation components e x and e z in each point of the map under plane strain assumption (e y = 0). As the Dx 1,2,3 (z, x) were experimentally measured, a possible deformation state was directly computed solving each of the three possible combinations of these equations. For example, by using the relation presented in Eq. 1 for A 1 (TO) 1 and A 1 (TO) 4 phonons (i.e. A 1 (TO) 1&4 ), we can obtain the strain components as follows:  After this, the Dn x and Dn z fields were computed from the deformation components fields e x (z, x) and e z (z, x) by means of the piezo-optic relation presented in Eqs. 4 and 5, where the plane strain condition is assumed (e y = 0).
From the obtained n x (z, x) and n z (z, x), the supported guided modes for the different optical polarizations (E x and E z ) can be computed by solving the Helmholtz equation with the finite element method to compare them with experimentally guided modes [15,16].
As mentioned above, we measured a map of the y(zz)y Raman spectrum in the waveguides cross section. After fitting the four A 1 (TO) peaks, the central wavenumber shift of each of them was obtained and this is presented in Fig. 4, with the exception of phonon A 1 (TO) 3 , for which no acceptable signal was observed. The range of wavenumber shift for the A 1 (TO) 1 , A 1 (TO) 2 and A 1 (TO) 4 phonons was between -1 and 1 cm -1 .
In Fig. 4, it can be seen that the shifting of A 1 (TO) 1 and A 1 (TO) 4 phonons has similar shapes outside the laser track regions (enclosed by a dashed line), but inside these, the mentioned phonons have qualitatively different behaviours: while the A 1 (TO) 1 phonon has a positive shift, the A 1 (TO) 4 phonon shift is negative.
On the other hand, the A 1 (TO) 2 phonon, outside the modified laser track region, has an opposite behaviour in comparison with that of A 1 (TO) 1 and A 1 (TO) 4 phonons: in the regions where the A 1 (TO) 2 has a positive wavenumber shift, the A 1 (TO) 1 and A 1 (TO) 4 have a negative shift, and vice versa.
In Fig. 4d, the intensity of the A 1 (TO) 4 phonon is plotted in grey scale (for the other phonons, the same effect was observed). This map depicts with dark grey colour the regions where the material structure is considerably modified by the incident femtosecond laser interaction; this strong reduction in the phonon intensity is associated with a partial amorphization of the crystalline structure [17].
By taking into account the measured Dx 1 (z, x), Dx 2 (z, x) and Dx 4 (z, x) surfaces and given e i and f i values, there are three different linear equations (from Eq. 1 with i = 1, 2, 4) with two unknowns surfaces: e z (z, x) and e x (z, x). Hence, we used the three possible combinations (i.e. A 1 (TO) 1&2 , A 1 (TO) 2&4 and A 1 (TO) 1&4 ) of these linear equations to retrieve the strain fields e z (z, x) and e x (z, x), as shown in Eqs. 2 and 3. For the studied kind of waveguide, x-polarized guided modes are reported in the literature as more intense than zpolarized guided modes [15]. Then, the refractive index field Dn x was firstly computed using Eq. 5 for each of the mentioned combinations A 1 (TO) i&j . These Dn x fields are plotted in Fig. 5a, b, c. From these, we can see that the Dn x obtained with the different combinations of phonons are considerably different within the laser track region. On one hand, the resulting Dn x for A 1 (TO) 1&2 and A 1 (TO) 2&4 has a positive increase in this region, which is not physically possible; in type II waveguides, a decrease in the refractive index of about Dn x * -1 9 10 -2 takes place in the direct interaction region [4]. However, the refractive index field obtained with mentioned combinations was also used to compute the guided intensity because the relevant region to support the guided modes is in the surrounding of the laser track region, not inside it. To this end, a decrease in refractive index was set within the laser track regions as a boundary condition. On the other hand, the Dn x retrieved from combination A 1 (TO) 1&4 has a reasonable change of about Dn x * -7 9 10 -3 within the laser track region.
Then, by comparing the refractive index fields obtained with the different combinations of phonons in the whole analysed region, we classified three different regions within the studied area. This can be clearly seen in Fig. 5d, where the refractive index variation computed across the path ab is shown. Each line corresponds to a different combination of phonons. The different regions are labelled and shown with a different background colour: Region R1 where a strong change in the material structure arisen due to the laser interaction (usually called ''laser track'' region) and thus, refractive index computed with the different combinations of phonons is considerably different. This region can be seen as a dark region through an optical microscope (Fig. 1b).
Region R2 where there is a lower modification of the material structure. Although the refractive index fields obtained from the different combinations of phonons are different, there is less difference between them than that observed in Region R1. This region cannot be defined with an optical microscope.
Region R3 where the refractive index fields obtained with the different pair of phonons matched and so, it can be concluded that the crystalline structure of lithium niobate is maintained.
In the mentioned regions (R1, R2, R3), we associated a difference between the different retrieved Dn x with the grade of crystallinity of the material and we could see, as it is mentioned in other works, that inside the laser tracks and around these there are different grades of crystallinity [17].
Then, to compare the performance of the different combinations of phonons to predict the refractive index Dn x field in the studied waveguides, we computed the guided intensity for the refractive index fields shown in Fig. 5. As mentioned before, the computed refractive index variations obtained from A 1 (TO) 1&2 and A 1 (TO) 2&4 combinations in the laser track region (R1) are not physically possible for this kind of waveguide structure. Hence, for these two combinations, we set a step of Dn x = -1 9 10 -2 in the laser track regions to perform the guided modes simulations. The perimeter of these laser track regions is shown in Fig. 5a, b with a dashed line. At any rate, as mentioned before, the most relevant refractive index region for guided modes simulation is within the guided region. On the other hand, in the Dn x field obtained from combination A 1 (TO) 1&4 , no additionally refractive index decrease was set in the laser track region to perform the guided modes computation, as previously explained. In Fig. 6a, b, c, we plotted the x-polarized guided intensity (i.e. the square sum of all the guided modes) computed for each combination of phonons and for a propagation wavelength of 1550 nm. In Fig. 6d, the experimentally acquired guided intensities polarized in xaxis for the same wavelength is plotted.
After comparing the computed intensities with the experimental guided modes, it was found that combinations A 1 (TO) 1&2 and A 1 (TO) 2&4 (shown in Fig. 6a, b) did not agree with that experimentally measured, mainly because they super-estimate the refractive index variation in the mentioned Region R2; for example, by taking into account the computed guided intensities in the central waveguide (enclosed by dashed curves in Fig. 6), these combinations clearly did not agree with those experimental measured and presented in Fig. 6d: in those two combinations, the distribution of intensity has two maximums while the experimentally guided intensity shows only one.
On the other hand, comparing the guided intensity within the central waveguide in Fig. 6c with that in Fig. 6d, we can see that the guided intensity expected from combination A 1 (TO) 1&4 is the most representative of the experimental data (Fig. 6d). Then, this combination of phonons is the most proper to compute refractive index filed, even in the modified Region R1, where this combination predicted a refractive index decrease.
Referring to z-polarized guided modes, it is reported in the literature that they are supported by this kind of femtosecond laser-written waveguides by launching 650, 800 and 1085 nm wavelengths [15,18]. In the current analysis, we could not measure z-polarized guided modes at 1550 nm wavelength, and accordingly no guided modes were obtained in simulations using retrieved Dn z field. However, in Fig. 7 it can be seen that the expected Dn z (z, x) field reveals an increase mainly in regions enclosed by curves I and II, where z-polarized light was measured in previous works, though generally with less intensity than x-polarized modes. In the waveguides analysed in this work, these regions have a predicted Dn z maximum of 3 9 10 -4 and this value was not enough to obtain simulated guided modes from the retrieved refractive index field. This discrepancy between currently and previously published works can be attributed to differences in the writing parameters employed.
Complementary to the preceding analysis, a linear scan of nano-indentations was carried out across the path ab in the waveguides cross section in order to explore the mechanical properties variation in these guiding structures. An image taken by a microscope of the indented path ab is shown in Fig. 8a. The elastic modulus was computed in each point of the map following the process previously mentioned. In order to understand the measured elastic modulus variation in the waveguides cross section, we obtained the density variation expected from the strain field retrieved above (with combination A 1 (TO 1 ) 1&4 ). This variation was computed using the well-known Eqs. 6 and 7. Then, the obtained density map was fitted to the measured elastic modulus variation using the Gibson-Ashby law presented in Eq. 8 [19]. To this end, only the exponent b in this equation was used as a fitting parameter. The obtained value of b after this fitting was equal to 9, which is not an empirically expected value for solid materials. According to previous works, this value should be between 0.9 and 4 [20]. That non-expected value obtained from our fitting might be produced due to the fact that the Gibson-Ashby law does not take into account changes in the material structure which indeed occur in the analysed region. where DV/V is the relative volume change in a certain point of the path where q/q o is the variation of relative density The relative elastic modulus (E/E o ) measured along path ab is plotted in Fig. 8b with blue points. Also, the relative density variation obtained from l-Raman map (using combination A 1 (TO) 1&4 ) and fitted using the Gibson-Ashby law to the elastic modulus variation is included in this figure (red solid curve). Although the implemented mechanical method has considerable noise, we found a clear correlation between the expected density variation and the elastic modulus measurement along the ab path: a considerable decrease in the elastic modulus (of about a 10%) agreed with a decrease in the local density in the laser track regions, where the material density is lower because there exists a more disordered structure.

Conclusion and future work
In this work, we studied femtosecond laser-written waveguides in lithium niobate by three different complementary techniques: end-fire coupling method, micro-Raman map and nano-indentation map. We computed the refractive index field directly from the Raman map by means of the Deformation Potential Theory, and we found that the most proper combination of A 1 (TO) phonons to compute refractive index field in type II waveguides is formed by A 1 (TO) 1 and A 1 (TO) 4 phonons. This combination allowed the computation of a physically expected refractive index variation Dn x and Dn z , even in the laser modified region. We also measured an elastic modulus decrease in the laser modified region, and this behaviour verifies the local density variation computed from the l-Raman map.
In conclusion, we presented a direct method to compute the strain and refractive index variation in femtosecond waveguides type II in lithium niobate, and we tested this method comparing the experimental x-polarized guided modes with those expected by the retrieved refractive index variations. Finally, we tested with instrumented nano-indentation the decrease in density occurring in the region directly modified by the ultrafast laser.
In future works, we will apply the method presented in this work to study other crystalline materials like LiCAF and LiSAF processed by femtosecond laser.  Fig. 8 a Image taken by a microscope of structural modified regions and b elastic modulus ratio E/Eo (blue points) and density ratio fitting (q/q 0 ) 9 (red line)