On-chip interrogator based on Fourier Transform spectroscopy

In this paper, the design and the characterization of a novel interrogator based on integrated Fourier transform (FT) spectroscopy is presented. To the best of our knowledge, this is the first integrated FT spectrometer used for the interrogation of photonic sensors. It consists of a planar spatial heterodyne spectrometer, which is implemented using an array of Mach-Zehnder interferometers (MZIs) with different optical path differences. Each MZI employs a 3$\times$3 multi-mode interferometer, allowing the retrieval of the complex Fourier coefficients. We derive a system of non-linear equations whose solution, which is obtained numerically from Newton's method, gives the modulation of the sensor's resonances as a function of time. By taking one of the sensors as a reference, to which no external excitation is applied and its temperature is kept constant, about 92$\%$ of the thermal induced phase drift of the integrated MZIs has been compensated. The minimum modulation amplitude that is obtained experimentally is 400 fm, which is more than two orders of magnitude smaller than the FT spectrometer resolution.


Introduction
Photonic based sensors find nowadays a wide range of applications. Acoustic and ultrasound sensors [1,2], pressure sensors [3], biochemical and gas sensors [4,5] are examples of sensors based on optical technology. They are low cost, immune to electromagnetic radiation, and operate under a wide range of temperatures. In this paper, we focus our attention on photonic sensors whose transmission or reflection spectra have a peak (or dip) in their lineshape. Examples are sensors based on fiber Bragg gratings (FBGs) [5,6] or on integrated ring resonators [1,2,4]. For these sensors, it is possible to build large and multi-purpose sensor arrays by wavelength multiplexing the spectrum of the sensors [6,7].
The photonic sensors mentioned above are designed in such a way that the signal to be sensed modulates the sensor's resonance wavelength. Interrogation is the technique of demodulating and demultiplexing the response of an array of photonic sensors. Different methods have been proposed in the past. A common approach is to measure the spectrum of the sensor array using a dispersive spectrometer such as an arrayed waveguide grating (AWG) [8,9,10] or an echelle grating [11]. Their sensitivity to the external excitation depends on the spectral resolution of the spectrometer; higher resolution comes at the price of a larger footprint. Another approach is edge filtering, where the output spectra of the photonic sensors is conveyed to an optical filter whose transfer function is linear within certain range. As the spectrum of the sensor shifts due to the sensing signal, the filter converts the resonance wavelength modulation into power modulation which can be obtained by a photodetector. The main drawback is that a high sensitivity may compromise the wavelength operation range [12]. Passaro et al [13] reports the spectral scanning as a possible solution, which features a high sensitivity and a large wavelength operation range. On the other hand, most of these interrogators are based on thermal tuning which limits their interrogation speed to a few kHz. Another approach for interrogation is to use passive interferometers such as Mach-Zehnder interferometers. In combination with a demultiplexing element, such as an AWG, it is possible to interrogate the photonic sensors as demonstrated in [14,15]. Despite the high sensitivity of this interrogator, special care should be taken to match the spectra of the AWG outputs to the sensors spectra. This might be an issue for integrated sensors such as ring resonators [1] since the resonance wavelength, in most of the cases, cannot be predicted during the design due to variations of the fabrication process.
The interrogation method here proposed may be applied to any sensor whose spectrum is finite and is modulated by an external signal. We demonstrate its performance using FBG sensors, but the method is equally suitable to other types of sensors such as ring resonators. To the best of our knowledge, this is the first interrogator based on integrated Fourier Transform (FT) spectroscopy. The technique is promising since it benefits from high flexibility, high sensitivity, and offers a high tolerance to variations of the fabrication process. In the past, FT spectroscopy was applied to demultiplexing FBG sensors [16,17], but at that time, the speed of the method was limited by the mechanical speed of the mirror. Integrated photonics enables the design of new FT spectrometer implementations. The most common one consists of an array of MZIs with different optical path lengths (OPDs) [18,19,20,21]. Thus, the spectrum can be retrieved by calculating the coefficients of the Fourier cosine series from the interferogram. However, since the number of MZIs is finite, the retrieved spectrum is an approximation to the actual one and a large the number of MZIs is required in order to achieve a high spectral resolution.
The design of our integrated FT spectrometer is similar to the one proposed by [22,23], where the complex Fourier coefficients of the system are obtained by using 3×3 multi-mode interferometers (MMIs). In our case, however, instead of retrieving the spectrum, we demonstrate that the complex Fourier coefficients can be written as a sum of the individual contributions of the sensors. We obtain a coupled system of non-linear equations, whose solution gives the modulation of the sensor's resonance wavelength. Since no approximation has been made, the minimum modulation amplitude we experimentally retrieved is 400 fm, more than two orders of magnitude smaller than the spectral resolution of our own FT spectrometer, and limited only by the signal-to-noise ratio of the input signal. Moreover, we demonstrate that the number of interferometers can be as small as the number of sensors, which strongly reduces the device footprint without compromising the interrogator sensitivity. Finally, we propose a novel technique for compensating the slow drift with time of the phases of the MZIs due to temperature fluctuations. [1,24]. This enables the application of this interrogation method for very low speed photonic sensors. Since the speed is only limited by the electronics, our interrogation method is equally suitable for high speed sensors.
2 Design and characterization of the FT Spectrometer Fig. 1a shows a picture of the FT spectrometer. The chip was fabricated by a multi-project wafer run in the Smart Photonics foundry and its dimensions are 4.0 mm × 4.5 mm. The chip has a total of 7 inputs, but inputs #5 and #7 are not used, as indicated in the figure. Following the optical path of main entrance (input #1) the light signal is split into nine beams and guided to nine different Mach-Zehnder interferometers. Other inputs provide access to a limited group of MZIs, allowing the characterization of the sensors using a reduced number of interferometers. For instance, input #6 provides access to MZIs 1-5. The chip is glued to a printed circuit board (PCB), to which the chip pads were wire bonded. Outputs per MZI of this PCB were connected to an other PCB which contain three transimpedance amplifiers (TIAs) for the photo-detectors and a pre-processing module. This module gives a linear combination of the outputs, as indicated in the schematic shown Fig. 1b.
MZIs represent the heart of the on-chip FT spectroscope. All the waveguides have a width of 1.5 µm and were fabricated in the Deep Etch Layer (see Fig. 6 of [25]). The length difference between the arms range from 0.710 mm to 6.39 mm in steps of 0.710 mm. At the end of the MZI, the light signals from the two arms interfere within a 3×3 MMI (360 µm length, 11.4 µm width).
In this section we characterize the MZIs of the FT spectrometer by considering its response to one particular wavelength λ. The transmittance for the given wavelength of l-th output of the m-th MZI is given by: where v ml is the visibility, n ef f,m (λ) is the effective index of waveguides of the m-th MZI, ∆L m the arms length difference of the m-th MZI, and φ l is the MZI phase shift given by (120 • , 0 • , -120 • ) for l = 1,2,3 in case the 3×3 coupler is balanced. In our design, the waveguide effective indexes are all the same except by small deviations caused by variations of the fabrication process. Expanding the term n ef f,m (λ)/λ in Taylor series around λ 0 , we obtain: where δn ef f,m are deviations of the nominal value of the effective index at the m-th MZI and λ 0 a wavelength close to 1550.0 nm. The approximation holds as long as the effect group index (n g ) can be considered constant over the spectrum of interest. Replacing Eq. (2) in Eq. (1) we obtain: where F m = λ 2 0 /(n g ∆L m ) is the free spectral range of the m-th interferometer and In our design, ∆L m is given by ∆L m = m∆L 1 with ∆L 1 = 0.710 mm, leading to F m = F 1 /m, where m is an integer ranging from 1 to 9 and F 1 = 921.7 ± 0.5 pm. Ψ m depends on n ef f (λ 0 ), which might change in case of temperature fluctuations, inducing a phase drift in T ml (λ). Following the schematic of Fig.1b, it is shown that the outputs of the MZIs are connected to integrated photo-detectors (PD). The PD current I ml is given by I ml (λ) = P m R ph T ml (λ), where P m is the optical power delivered at the m-th MZI and R ph is the photodetector responsivity. The outputs of the photodetectors are send to TIAs, whose outputs voltage are given by: where g ml is the transimpedance gain. The 3×3 MMIs were designed to produce interference fringes with similar amplitude and a 120 • shift between each other. Aiming for the interrogation of the photonic sensors, the pre-processing module of the PCB combines the TIA output voltages according to [1]: where V m,x and V m,y are 90 • phase shift voltages, A m,x and A m,y are the voltage amplitudes, x of f,m and y of f,m are voltage offsets, and δφ m is a phase error. If the 3×3 MMI is balanced and the electronic components of the PCB are ideal (ideal operational amplifiers and no variance with respect to the nominal value of the resistors and capacitors), the voltage offsets are zero (x of f,m = y of f,m = 0), δφ m = 0, and A m,x = A m,y = P m R ph gv, where the visibility is v = v m1 = v m2 = v m3 and the TIA gain is g = g m1 = g m2 = g m3 . In this case, the Lissajous curve [V m,x (λ), V m,y (λ)] gives a circle with radius vP m R ph g centred at the origin. The transmission spectrum of each MZI has been measured using a tunable laser (Agilent, 81960A). The laser power is set to 6.0 mW and we performed the laser wavelength sweep ranging from 1550 nm to 1551 nm in steps of 1 pm, while the outputs of the pre-processing module are recorded by the digital acquisition module (DAQ, National Instruments, NI 9220). Fig. 1c shows the measured voltages of the outputs of MZI 1 (∆L = 0.710 mm), as well as a fit of the measured data against to Eq. (6). Since V 1,x and V 1,y have slightly different amplitudes and δφ 1 = 17.9 • , the circle is deformed into a tilted ellipse centred outside of the origin, as shown in Fig. 1d. Non-idealities and variations of electronic components in the PCB are neglected in the following section, whereas the unbalancement of the 3×3 MMI needs to be considered. In Section 3.2 we discuss how to correct for this.

The interrogation method
Here we derive the expressions for determining the resonance wavelengths of the photonic sensors as a function of time. Typically, the spectrum of each sensor has a peaked lineshape, which is modulated by an external signal such as temperature, strain or any other physical or chemical quantity. The photonic sensors are assumed to be wavelength multiplexed. Let there be K sensors with resonance wavelengths λ k (t) at time t, where k = 1, ..., K. The combined spectrum S(λ, λ 1 (t), ..., λ K (t)) received by the interrogator is given by: where s k (λ, λ k (t)) is the spectrum of the k-th sensor. The signals that are to be sensed induce time dependent modulations of the resonance wavelengths. The resonances λ k (t) must be separated so that the curves s k (λ, λ k (t)) do not overlap. In this paper s k (λ, λ k (t)) correspond to the reflection spectra of FBGs sensors. However, the method applies also to integrated photonic sensors as the ones described in [1]. S(λ) is assumed to be a poly-chromatic signal and the values of the TIA output voltages are given by: where the constant G is given by G = (1 − α c )gR ph with α c the coupling losses. The electronic preprocessing module combines the signals from the three outputs of the interferometers according to Eq. (6), resulting in the two 90 • phase shifted voltages V m,x (t) and V m,x (t): As explained in Section 2, the voltage offsets x of f,m and y of f,m are mainly caused by the fact that the 3×3 MMIs are unbalanced. At the end of a calibration process (see Section 3.2), the offsets are removed by averaging and, at this point, they are neglected. By defining a complex voltageV m (t) = V m,x (t) + iV m,y (t) we obtain: The chip is characterized after the MZI phase drift has been stabilized, so Ψ m is constant in time and taken out of the integral in Eq. (11). In Section 3.3, however, a novel method is presented for compensating the environmental phase drift by using one of the sensors as a reference. We assume that S(λ, λ 1 (t), ..., λ K (t)) vanishes outside the interval [λ 0 − F 1 /2, λ 0 + F 1 /2] for all times t, where λ 0 is a wavelength close to 1550.0 nm. Then we have: Eq. (12) are the Fourier coefficients of the function λ → S(λ, λ 1 (t), ..., λ K (t)) when considered as periodic function with period F 1 . This implies that: whereV −m (t) =V * −m (t) since S(λ, λ 1 (t), ..., λ K (t)) is real. The chip contains a finite number of M = 9 interferometers. The retrieved spectrum S M (λ, λ 1 (t), ..., λ K (t)) is given by: Function S(λ, λ 1 (t), ..., λ K (t)) differs from S M (λ, λ 1 (t), ..., λ K (t)) by the fact that the last one features a finite spectral resolution δλ res given by: For M = 9, δλ res = 50 pm. Moreover, S M (λ, λ 1 (t), ..., λ K (t)) is periodic with period F 1 . For a large number of interferometers (M >> K), S M (λ, λ 1 (t), ..., λ K (t)) gives a good approximation to S(λ, λ 1 (t), ..., λ K (t)) and it is possible to obtain the resonance wavelengths by tracking the peaks of S M (λ, λ 1 (t), ..., λ K (t)). However, δλ res represents a limitation to the minimum modulation amplitude to be experimentally obtained.
In order to determine λ k (t) with higher accuracy and using a reduced number of MZIs we derive a non-linear system of equations. We assume in this section that λ(t) is known at t = 0. Let where δ k (t) is the modulation of the resonance wavelength of the k-th sensor that we aim to determine. By substituting Eq . (7) and Eq. (16) into Eq. (11), we obtain: The right-hand side of Eq. (17) represents the Fourier transform of s k (λ − λ k (0) − δ k (t)) evaluated at m/F 1 . Using the shift property of the Fourier transformation, Eq. (17) is rewritten as: We rewrite Eq. (18) as:V for m = 1, ..., M . The coefficients a mk are experimentally determined as explained in Section 3.2. Eq. (20) represents an M×K system of non-linear equations to be solved using Newton's method, where M is the number of interferometers and K is the number of sensors. Hence, the number of interferometers must only be at least as large as the number of sensors (i.e. M >= K), which means that the footprint of the device can be relatively small. In our chip M = 9. The system is explicitly written in the in Eq. (21): It can be show that as long as the phases 2πλ k (t)/F 1 (for k = 1, ..., K) are different and the initial guess for {δ 1 (t), ..., δ K (t)} is close to the actual solution, the Jacobian ∂V m /∂δ k is not singular and the Eqs. Assuming that the FBG sensors spectra have a Lorenzian lineshape, we replace the Fourier transform of s k (λ) into Eq. (19): where s max k is the maximum value of the Lorenzian of the k-th, OP D 1 = n g ∆L 1 is the optical path difference of MZI 1, and L c,k is the cohenrece length given by: where w k is the full width half maxima (FWHM) of the Lorenzian. The coherence length limits the maximum OPD value which allows interferometric fringes to be experimentally resolved. Eq. (22) shows that a mk becomes very small when the MZI free spectral range is comparable or smaller than the FWHM of k-th sensor. As discussed in Section 4, the MZIs with larger OPDs are not used due to the strong attenuation and the reduced signal-to-noise ratio (SNR).

Calibration and experimental determination of the coefficients
The coefficients a mk are experimentally determined via the following calibration procedure. Let t start k be the instant of time when the calibration of k-th sensor starts and t end k be the instant of time when the calibration ends for the same sensor. During the time interval t start k < t < t end k , all sensors are kept at rest, while sensor k is excited. In case sensor k is a temperature sensor, heat is applied (as much as possible) during the calibration. If sensor k is a strain sensor, a large stress is applied (as much as possible). According to Eq. (20), for a balanced 3×3 MMI, the m-th complex voltageV m (t) during the time interval t start k < t < t end k is given by: where δ l (t) = 0 if l = k since no excitation is applied to the other sensors and where c mk = K l =k a ml and θ mk (t) is the complex argument of the term a mk e im2πδ k (t)/F1 , given by: The Lissajous curve {V m (t)}, {V m (t)} for t start k < t < t end k is given by a circular arc: where ( {c mk }, {c mk }) defines the arc centre, |a mk | the radius, and θ mk (t) the instantaneous angle with the real axis. Fig. 2 shows a simulation of the calibration for two sensors. The calibration starts at t = t 0 < 0 and ends at t = 0, when the interrogation procedure starts. During t start 1 < t < t end 1 , sensor 2 is kept at rest, while sensor 1 is excited by moving its resonance wavelength from 1550.50 nm to 1550.16 nm, as shown in Fig. 2a. This induces the oscillations of V 1,x (t) and V 1,y (t) during t start 1 < t < t end 1 as shown in Fig. 2b, which are traced as a circular arc in red shown in Fig. 2c. The procedure is repeated for sensor 2: during t start 2 < t < t end 2 , while sensor 1 is not excited, sensor 2 changes its resonance from 1550.75 nm to 1550.33 nm. This causes the oscillations from t start 2 < t < t end 2 in Fig. 2b which are traced as the circular arc in green shown in Fig. 2c.
As explained in [1], a slight non-ideal behavior of amplitude and phase of 3×3 couplers are not uncommon and result into a deformation of the circle in an ellipse. An ellipse is fitted to the data points (V m,x (t) , V y,m (t) ) during the interval t start k < t < t end k , where V m,x (t) and V m,y (t) are the m-th MZI voltages measured during the calibration. A larger excitation of the k-th sensor results in a larger angular deflection, leading to a more accurate retrieval of geometrical parameters of the ellipse. The fitting gives the ellipse semi-axis r 1,mk and r 2,mk (where r 1,mk > r 2,mk ), the angle α that represents the rotation of the ellipse with respect to the x-axis, and the ellipse centre (x el mk , y el mk ). In order to map the ellipse to an circle, the following transformation is applied: where V m,x and V m,y are the corrected values of the 90 • phase shifted voltages so that the Lissajous curve (V m,x (t), V m,y (t)) for t start k < t < t end k gives a circle arc with radius r 1,mk . The correction of Eq. (27) needs to be performed for all interferometers (m = 1, ..., M ). Although the ellipse semi-axis r 1,mk and r 2,mk , as well as the corrected radius r 1,mk may change according to the sensor (since it depends on its total transmitted or reflected power spectrum) and according to the interferometer (due to the different MZI's coherence lengths), the ellipse eccentricity depends only on the 3×3 MMI, as discussed in Section 2. Thus, for a given interferometer m the ratio r 1,mk /r 2,mk is constant for k = 1, ..., K. The design of the 3×3 MMI is the same for all interferometers, hence the ratio r 1,m /r 2,m is constant for m = 1, ..., M as long as the variations of the fabrication process are negligible.
After calculating the 90 • phase shifted voltages, the modulus of the coefficients a mk can be obtained. Since the radius of the circle arc obtained for the m-th interferometer and the k-th sensor is r 1,mk , the modulus of the coefficients a mk , according to Eq. (26), is given by: Next, the linear transformation of Eq. (27) is applied to the point (x el mk , y el mk ), which gives the centre ( {c m,k }, {c m,k }). The angles θ mk (t) (for m = 1, ..., M and k = 1, ..., K) are given by: where arctan 2 (x, y) is the four quadrant arc tangent. During the final stage of the calibration of sensor k, the angle θ mk (t) remains constant because then no excitation is anymore applied to it. By substituting t = t end k in Eq. (25), we obtain: where the calibration procedure ends at t = 0. According Eq. (16), δ k (0) = 0. Therefore, the argument of a mk is given by: arg(a mk ) = θ mk (t end k ) = θ mk (0).
The values of λ k (t) (for k = 1, ..., K) are in general unknown at the end of the calibration (t = 0), which contradicts the assumption made in Eq. (16). Here, we refine our previous statement by assuming that the values of λ k (t) are known at t = t 0 , before the calibration procedure starts. In most of cases, however, the sensors can be calibrated in such a way that their resonance wavelengths return to their initial value at the end of the calibration (λ k (t 0 ) = λ k (0)). In situations where this is not possible (due to a sensor hysteresis, for instance), the values of λ k (0) can be obtained by following the procedure: (a) determine the value of δ(t start k ) from Eq. After finishing the calibration of all sensors in this way, the offsets are determined by averaging: Finally, the complex voltages are computed as function of time to be used in Eqs. (20) and (21):

Compensation of the phase drift
Since the effective index in Eq.(4) is temperature dependent, local variations of temperature induces the phase Ψ m to drift. In order to account for this effect, we rewrite Eq. (4) according to: where The temperature dependence of the group index n g and to δn nef f have been neglected. Eq. (35) indicates that the phases Ψ m in Eq. (18), (20), and (21) are no longer constant. Eq. (18) can be rewritten as: where The right side of Eq. (36) is identical to Eq. (20) demonstrating that fluctuations of the environmental phase impacts on the solutions of Eq. (20) or Eq. (36). This effect can be corrected by using another sensor as a reference, to which no excitation is applied and its temperature is kept constant. Let δ ref (t) be the solution of Eq. (36) for the reference sensor. The calibration procedure assures that when the interrogation procedure starts (t = 0), the values δ k (0) are zero for all sensors (k = 1, ..., K). Since no excitation is applied to the reference sensor, the function δ ref (t) remains at zero for t > 0. Hence, according to Eq.(37): Thus, the phase drift can be compensated by subtracting the term ∆Ψ(t)F 1 /(2π) in Eq. (37), obtained from Eq.(38):

Experimental setup
The schematics of the experiment is depicted in Fig. 3. Light from a broadband amplified spontaneous emission source (ASE, Optolink, OLS15CGB-20-FA) is sent, through a circulator (OZ Optics, FOC-12N-111-9), to the FBG sensor array (Technicasa, T10). The FBG sensors reflect back to the circulator their combined spectrum, which is amplified by a optical booster amplifier (Thorlabs, S9FC1004P) according to Fig. 3a. The gain is 12 dB and the light is coupled to the chip using lensed fibers (Oz Optics, TSMJ-3A-1550-9). Outputs of the chip are conveyed to a PCB which implements the transimpedance amplifiers for the photodetectors and an pre-processing module in order to implement Eq. (6) electronically (see Fig. 1). The PCB outputs are sampled by the DAQ (National instruments, NI9220), which the maximum sampling speed is 100 kSa/s/channel. The performance of our interrogator is evaluated using four FBG sensors: three as strain sensors one as a reference sensor, used to compensate the environmental phase drift. The calibration is performed in a such way that λ k (t 0 ) = λ k (0). The ends of the fibers containing the FBGs are clamped to the translation stages as shown in Fig. 3b. In order to tune the peak wavelengths λ k (0), stress is applied using the manual positioners, avoiding the angles 2π (λ k (t)) /F 1 to overlap during the experiment. FBG #1 represents the main strain sensor and the translation stage (referred as translation stage 1) to which FBG #1 is attached is controlled by a stepper motor. FBGs #2 and #3 are the secondary strain sensors and they are both attached to translation stage 2 controlled by another stepper motor. FBG #4 is the reference sensor and it is attached only to manual positioners. We programmed the stepper motors to operate in cycles of three steps: (a) the translation stage travels at a constant speed from the position x = 0 to x = ∆ ; (b) The stage rests at x = ∆ ; (c) The stage returns to the original position.
Since FBGs #2 and #3 are secondary strain sensors, we programmed the translation stage to move periodically from the distances x = 0 to x = ∆ (2) = 30µm. In contrast, the translation stage to which FBG #1 is attached, travels to different values of ∆ (1) ranging from 0.5 µm to 200 µm (these values are shown later in Fig. 5). Since the stress to be applied to FBG #1 is much larger compared to FBGs #2 and #3, the translation stage 1 is programmed to move towards −x. Thus, a negative stress applied to FBG #1, avoiding to damage it. Translation stage 1 repeats three times its motion from x = 0 to x = ∆ (1) and from x = ∆ (1) to x = 0. Thus, the travelling distances ∆

Experimental results
As explained in Section 3.4, the performance of our interrogator is evaluated using four FBG sensors: three as strain sensors one as a reference sensor, used to compensate the environmental phase drift. Using manual positioners, a constant stress is applied to all FBGs in such a way that the resonance wavelengths of the sensors are set to λ 1 (0) = 1550.9 nm, λ 2 (0) = 1550.3 nm, λ 3 (0) = 1551.4 nm and λ 4 (0) = 1549.7 nm. The differences of λ k (t) − λ l (t) for l = k can be larger than F 1 (F 1 is the free spectral range of MZI 1) provided that the angles 2πλ k (t)/F 1 = 2πλ l (t)/F 1 for all l, k = 1...K. The light signal is coupled to the chip using input #6 (see Fig. 1a), where the input power is shared among MZIs 1 to 5. Better interrogation results are obtained by sharing the optical power among a reduced number of interferometers since the outputs of the MZIs with larger OPDs are strongly attenuated, according to the discussion in the end of Section 3.1.
In order to retrieve the coefficients a mk , we individually excited the FBG sensors. Following the procedure described in Section 3.2, the complex voltagesV m (t) have been obtained by mapping the ellipse arcs to circle arcs according to Eq. (27), and by removing the voltage offsets according to Eq. (32). Fig. 4a shows the real and imaginary parts ofV 1 (t), to which a low pass filter (cut-off at 45 Hz) has been applied in order to suppress noise. The real and the imaginary parts ofV 1 (t), shown in Fig. 4a, are plotted in Fig. 4b as a Lissajous curve. The figure shows four circular arcs, which correspond to the individual excitation of the sensors, obtained from the outputs of MZI m = 1 during the calibration. The radii and the angles of the arcs at the end of the calibration procedure give the modulus and argument of the coefficients a mk , as described in Section 3.2. Fig. 4b shows, however, that some regions of the Lissajous curve deviate from the expected circular path. This occurs when the resonance wavelengths of two FBGs are about to cross and the spectra of two FBG sensors overlap. This causes that a part of the input optical signal is reflected multiple times in between the FBGs, creating an Fabry-Perot cavity. The interference of the electric field which is reflected multiple times between the FBGs leads to the deviations of the circular arcs. To overcome this issue, we followed the calibration described in Section 3.2 using only the parts of the Lissajous curves that are close to circular. For t > 0 s, the interrogation starts and the three strain sensors are simultaneously excited. As a result, an arbitrary Lissajous curve is obtained. Fig. 4c-f shows the solution of Eq. (20) obtained using the Newton's method. As explained in Section 3.1, the solution obtained at the instant t is used as an initial guess for the Newton's method at the instant t + 1/f s , where f s is the sampling frequency. As a result, the method converges at any t with a maximum of four interactions. For a sampling rate of 10 kSa/s, about one million of systems of equations needs be solved from t =0 s to t = 100 s. Using an Intel i5-3470 processor, the solution is roughly calculated at a rate of a hundred equations per second and the total computational time is about 2h and 45 min.
FBGs #2 and #3 are attached to translation stage 2 which periodically travels from x = 0 to x = ∆ (2) = 30µm. As a result, the functions δ 2 (t) and δ 3 (t) are time periodic, as shown in Figs. 4c and 4d. On the other hand, Fig. 4f shows the solution δ 1 (t), which consists of a succession of dips. The dips are obtained because the stepper motor applies a negative stress to FBG #1, as explained in Section 3.4. Since the translation stage repeats its motion three times to a given distance ∆ (1) , Fig. 4f shows a series of dips grouped by 3 successive ones with approximately the same depth. Fig. 5 shows the modulation amplitude ∆λ (1) for sensor 1 as a function of the strain applied to FBG #1. The strain is assumed to be constant along the fiber and it is defined as: where ε 3j+1 and ∆ (1) 3j+2 are the same, as explained in Section 3.4. On the other hand, the modulation amplitude is defined as: where δ dip 1,3j is the time average of function δ 1 (t) at the three adjacent dips (3j + 1), (3j + 2) and (3j + 3), as indicated in the upper inset of Fig. 5. Similarly, δ max 1,3j is the time average of function δ 1 (t) at its A low pass filter (cut-off at 45 Hz) has been applied to the measured voltages V m,x (t) and V m,y (t). The numbers 1,2,3 and 4 indicate the calibration interval (t start k < t < t end k ) for sensors k =1,...,4. (b) Lissajous plot obtained by plotting the real and imaginary parts ofV 1 (t). During the calibration, the Lissajous curve is a circular arc. During the interrogation, all sensors are simultaneously excited, and an arbitrary Lissajous curve is obtained as shown in orange. (c)-(e) Solutions δ 2 (t), δ 3 (t) and δ 4 (t) of Eq.(36) for t > 0. FBG #4 is the reference sensor. The phase drift was compensated using Eq. (39). (f) Comparison between the solutions δ 1 (t) and δ 1 (t) . The inset shows a zoom of the solution δ 1 (t).
(3j + 1)-th, (3j + 2)-th and (3j + 3)-th maxima, which occur when the translation stage rests around the original position x = 0. The ratio between the amplitude modulation and the strain gives the sensitivity S (1) of FBG #1: By fitting a straight line to the data points (∆λ j ), we retrieved S (1) = 1.217 ± 0.006 pm/µstrain, which agrees with the nominal sensitivity of 1.2 pm/µstrain provided by the manufacturer (Technicasa, T10). The minimum retrieved strain is 365 nanostrain and the corresponding minimum modulation amplitude obtained is ∆λ min = 400±200 fm. This value is more than two orders of magnitude smaller than the resolution of the FT-spectrometer (50 pm). The value of ∆λ min , experimentally retrieved, is not limited by the resolution the FT-spectrometer but only by the SNR of the input signal.
FBG #4 has been taken as a reference sensor and no external excitation is applied to it after the end of the calibration. However, for t < 10 s, Fig. 4d shows small fluctuations of function δ 4 (t) (of the order of a few pm), caused by the cross-talk among sensors. Since the modulation amplitude of FBG #1 is the larger for t < 10 s, its cross-talk with FBG #4 is dominant. The maximum cross talk between FBGs #4 and FBGs #1 is about 1% of the δ 1 (t) value, which is acceptable in most applications.
For t > 60 s, the chip is heated up using a Peltier element. The temperature increases about 0.3 • C in the chip, a value which is comparable to temperature fluctuations in a temperature controlled room. As explained in Section 3.3, the solution δ 4 (t) shown in Fig. 4e is proportional to the drift of the phase Ψ m (t). Fig. 4f shows a comparison between the solutions δ 1 (t) and δ 1 (t), this last one obtained using Eq. (39) (solutions δ 2 (t) and δ 3 (t) are not shown). 92.0% of the phase drift has been compensated. For the sensors presented here, the phase drift influence could have been removed by applying a high pass filter to δ 1 (t) , δ 2 (t) and δ 3 (t) . However, for low speed sensors such as biochemical sensors [4], filtering is not possible since the speed of the sensor is comparable to the phase drift speed.    Although the method can be applied to high speed sensors, its real time implementation is challenging. On one hand, the speed of the FT-spectrometer is limited only by the electronics and the integrated photo-detectors may respond at frequencies higher than 5 GHz. On the other hand, a system of nonlinear equations need to be solved at each instant of time. The computational costs, however, can be reduced by calculating the inverse of the Jacobian ∂V m /∂δ k analytically. Using the transformation z k (t) = 2π(λ k (0)+δ k (t))/F 1 , it can be shown that the Jacobian is given by a product of a diagonal matrix and the Vandermond matrix V (z k ). Since analytic expressions do exist [26] for the inverse of V (z k ), the computational time is mainly governed by the time of calculating product of matrices. Moreover, the reduced number of interactions of Newton's method also contributes in reducing the computational time. Nevertheless, the real time interrogation of high speed sensors may require the usage of an application specific computational solution.

Conclusion
A novel interrogation method based on FT spectroscopy is presented. The technique is promising due to its high flexibility, high sensitivity and reduced interrogator footprint. It can be applied in different situations, in particular, for arrays of integrated sensors where the resonance wavelengths cannot be predicted during the design stage. Three conditions have been identified for the proper interrogation of the sensors: (a) the number of interferometers must only be at least as large as the number of sensors, allowing the interrogator footprint to be relatively small; (b) the MZIs must have different OPDs; (c) the phases 2πλ k (t)/F 1 (for k = 1, ..., K) needs to be different at any time. If the maximum amplitude modulation of the sensors is known, condition (c) is usually not an issue for FBG sensors, since the Bragg wavelengths could be chosen with an accuracy better than 1.0 nm. In case of integrated ring resonators, it is possible in most situations to design rings with a slightly different lengths, assuring a similar free spectral range, but different resonances. Since the phases depend on F 1 , the proper design of the FT spectrometer gives an extra flexibility to avoid the phases 2πλ k (t)/F 1 to overlap.
It has been shown that the minimum modulation amplitude experimentally retrieved is not limited by resolution of the FT-spectrometer, but limited only by the signal-to-noise ratio of the input signal. The minimum modulation amplitude obtained is 400 ± 200 fm and the cross-talk, which also depends on the SNR, is about 1%. Moreover, the phase drift of the interrogator, caused by temperature fluctuations, can be compensated by using one of the sensors as reference sensor to which no external excitation is applied. This is important for low speed sensors where the thermal induced drift of MZI phases is comparable to the speed of the sensors. Our method can also be applied for high speed sensors, but the implementation of real time interrogators require the analytic calculation of the inverse of the Jacobian matrix used in Newton's method. For real time interrogation the the usage of application specific computational solutions may be needed.