Identification and correction of Sagnac frequency variations: an implementation for the GINGERINO data analysis

Ring laser gyroscopes are top sensitivity inertial sensors used in the measurement of angular rotation. It is well known that the response of such remarkable instruments can in principle access the very low frequency band, but the occurrence of nonlinear effects in the laser dynamics imposes severe limitations in terms of sensitivity and stability. We report here general relationships aimed at evaluating corrections able to effectively account for nonlinear laser dynamics. The so-derived corrections are applied to analyse thirty days of continuous operation of the large area ring laser gyroscope GINGERINO leading to duly reconstruct the Sagnac frequency ωs\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\omega _s$$\end{document}. The analysis shows that the evaluated corrections affect the measurement of the Earth rotation rate Ω⊕\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\varOmega _{\oplus }$$\end{document} at the level of 1 part in 1.5×103\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$1.5\times 10^{3}$$\end{document}. The null shift term ωns\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\omega _{ns}$$\end{document} plays a non negligible role. It turns out proportional to the optical losses μ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mu $$\end{document} of the ring cavity, which are changing in time at the level of 10%\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$10\%$$\end{document} within the considered period of thirty days. The Allan deviation of estimated Ω⊕\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\varOmega _{\oplus }$$\end{document} shows a remarkable long term stability, leading to a sensitivity better than 10-10\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$10^{-10}$$\end{document} rad/s with more than 10 s of integration time, and approaching (8.5±0.5)×10-12\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(8.5\pm 0.5)\times 10^{-12}$$\end{document} rad/s with 4.5×105\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$4.5\times 10^{5}$$\end{document} s of integration time.


Introduction
Ring laser gyroscopes (RLGs) are inertial sensors based on the Sagnac effect [1][2][3]. They are largely used for inertial navigation, and applications in geodesy, geophysics and even for General Relativity, where tests are foreseen [4]. Since 2011 we are studying the feasibility of the test of Lense-Thirring dragging of the rotating Earth at the level of 1% with an array of large frame RLGs [5][6][7]. For that purpose it a e-mail: angela.divirgilio@pi.infn.it (corresponding author) is necessary to push the relative accuracy of the Earth rotation rate Ω ⊕ measurement in the range from 1 part in 10 9 up to 1 part in 10 12 . RLG consists of a laser with a cavity comprising of three or four mirrors, rigidly attached to a frame; large frame RLGs are utilised to measure the Earth rotation rate, being attached to the Earth crust. Because of the Sagnac effect, the two counter-propagating cavity beams have slightly different frequency, and the beat note of the two beams is proportional to the angular rotation rate of the ring cavity. Large frame RLGs are the most sensitive instruments for inertial angular rotation measurements. The Sagnac frequency of a RLG is in fact proportional to the component of the angular velocity − → Ω felt by the instrument along the normal to the cavity plane: where A is the area of the ring cavity, L is its perimeter, λ the wavelength of the light, and θ is the angle between the area vector of the ring and − → Ω . For RLGs lying horizontally (area vector vertical) θ is the co-latitude angle, while for RLGs aligned at the maximum Sagnac signal θ = 0. Equation 1 defines the scale factor S, which is a function of the geometry and of λ, quantities than can be measured with a very high accuracy.
Further to sensitivity, other key points of such instruments rely on their broad bandwidth, which can span from kHz down to DC, and their very large dynamical range. In fact the same device can record microseismic events and high magnitude nearby earthquakes [8], being the signal based on the measurement of the beat note between the two counterpropagating beams. It has been proven that large size RLGs, equipped with state of the art mirrors, can reach the relative precision of 3 parts in 10 9 with one day of integration time, in the measurement of Ω ⊕ [1]. If shot noise limited, sensitivity level scales with the square of the size of the ring cavity [1]. However, other limitations can affect the measurement.
The laser dynamics is non linear and plays a role in determining the RLG signal. We have recently developed a model to reconstruct the Sagnac frequency ω s starting from the measured beat note ω m and using the mono-beam signals. 1 The Sagnac angular frequency can be expressed as the linear sum of several terms. As we have discussed in a recent paper [9], reconstruction of the Sagnac angular frequency is mostly affected by the correction ω s0 associated with the so called back-scattering. In this paper we provide the necessary information to evaluate the other correction terms. Thirty days of continuous operation of GINGERINO are analysed taking into account the set of first order expansion terms, and results discussed.
Results demonstrate that the so called null shift term, the first order correction ω ns1 , plays a non negligible role in affecting both accuracy and sensitivity of the apparatus. The finding is in contrast with conventional treatments of RLG data, where null shift effects are typically neglected or considered simply proportional to the difference of the two mono-beam signals. On the contrary we found that including such correction leads to improve the measurement accuracy. Furthermore, since null shift effects can be related to the optical losses of the system, which inherently depend on operation time, a time dependent effect is found.

The analysis scheme to take into account laser dynamics
Our approach is based on the model of laser dynamics we recently developed in details in [9] based on the Aranowitz model of RLG [10,11], where laser dynamics is described in terms of several dimensionless parameters (Lamb parameters). The polarization of the laser plasma is described up to the third order expansion in powers of the field [10,11]. In the case of large frame RLGs, to avoid mode competition, a 50:50 mixture of two isotopes, 20 N e and 22 N e, is utilised, and the laser operation is set close to threshold to guarantee singlemode operation. This particular choice allows for some simplifications: the Lamb parameters of cross-saturation (θ 12 and θ 21 ) can be neglected, and we can assume the self-saturation terms equal to each other (i.e., β 1 = β 2 = β). The general goal of the analysis is to evaluate the Sagnac frequency correcting systematics due to the laser nonlinear dynamics by using the available data: the measured beat note and the mono-beam signals [9,12]. Figure 1 shows a scheme of the experimental apparatus, where photodiodes collecting the relevant signals are indicated.
The model provides two analytical expressions for the Sagnac angular frequency ω s = 2π f s : a non-linear relationship (Eq. 8 of Ref. [9]). This analytical expression was expanded up to the second order in term of two parameters δ ns and K , and can be calculated from the Lamb parameters. Accordingly, ω s is expressed as the linear sum of a total of six terms: As demonstrated in [9] ω s0 , accounting for back-scatter effects, is the dominant term in the expansion. The term ω s0 does not depend on Lamb parameters and can be expressed starting from the knowledge of experimental parameters according to: where ω m is the measured beat note (expressed as angular frequency), I 1,2 , I S1,S2 , and are the DC level of the mono-beam signals, their amplitude demodulated at the beat note frequency, and the dephasing between the two signals, respectively.
In order to better define ω s0 and take into account experimental issues, such as dark currents in the photodiode signals, the term ω sξ , see Eq. 4, has been introduced, where ξ is a proportionality constant which could be evaluated by statistical means.
The other terms in Eq. 2 give smaller contributions. They can be evaluated as detailed in the Appendix, starting from the knowledge of the Lamb parameters [10], a careful estimation of all operating parameters of the system, and the evaluation of the total optical losses μ. The latter can be measured with a 1% accuracy via the ring down method [11], i.e., by suddenly switching off the laser discharge and measuring the time behaviour of the intracavity signals through photodiodes P H1 and P H2 in Fig. 1. The method cannot be applied continuously during the RLG operation, since it requires to switch off the laser. Moreover, optical losses can change during operation, therefore their evaluation should be repeated frequently.
Because of that, and also to avoid uncertainties and fluctuations related to estimation of operating parameters such as, laser beam size, pressure and polarization conditions of the discharge (see Appendix for their use in the relevant expressions), in the present paper we have adopted a statistical approach for evaluating μ. As shown in the next subsection, this parameter enters as a linear proportionality factor Typical scheme of RLG with a square ring cavity. Photodiodes used to measure the beat note, necessary to evaluate the ω m , the two mono-beams (P H1 and P H2) and the gain monitor (DM) are shown. PMT indicates an additional photodetector not used for the purposes of the present paper. Note that the signal produced by the photodiode P H1 is used also in a feedback loop to control the power of the laser into ω ns1 . Therefore, we can express ω ns1 = μ×ω ns1 , wherẽ ω ns1 is evaluated starting from known data and μ determined through a linear regression procedure. The terms in Eq. 2 containing the subscript K require some further discussion. As reported in the appendix, the parameter K is proportional to oscillating terms at the frequency ω s . The evaluation of the K terms is not straightforward, and several approaches are in principle feasible. In the current analysis, we assume that ω s is nearly constant in time, and as a consequence for frequencies much below ω s the contribution of these terms should average to zero. As reported in the appendix, assuming sin ω s t = cos ω s t = M K the equations is simplified, emerging a proportionality term M K , that is evaluated via the statistical approach.
In summary, the analysis procedure is as follows: in the first level, different termsω ns1,2 are evaluated and stored on disk; in the second level, updated ξ , μ and M K coefficients are found with linear regression models, minimising the sum of Eq. 2. It is important to remark that μ can be stable with time, or varying with time. In this second case it is necessary to add to the model the information necessary to follow the time behaviour. Section 2.3 presents the method used for GINGERINO. We remind that the available signals from the RLG are sampled at high frequency (typically 5 kHz), and the whole calculation is carried out at sufficiently high frequency, in our case 2.5 kHz, while data are stored after decimation, usually corresponding to a sampling frequency of 2, 1, or 0.1 Hz. The diagram in Fig. 2 shows the first level involved in the method, while Fig. 3 presents a diagram of the second level of analysis. In the following, 30 days of data produced by GIN-GERINO will be analysed. In the analysis, corrections are included only up to the first order: ω ns1 = ω ns , ω K 1 = ω K . Second order terms are in fact smaller and do not contribute to improve the results compared to the present state of accuracy. The Appendix contains however all information to evaluate the different terms in the general case.

Some details on the evaluation of first order expansion terms
As already mentioned, we start from the DC components, I 1 and I 2 , of the mono-beam signals, their AC components at the beat frequency ω m /2π , I S1 and I S2 , and the relative phase between the two AC components. Dissipative processes (like diffusion or absorption from the mirrors) can produce non reciprocal losses between the two counter-propagating beams, so that two distinct loss parameters μ 1 = μ 2 must be used. Without loss of generality, it is possible to take μ 1 = μ and μ 2 = μ + δμ, where μ represents the reciprocal losses term and δμ the non reciprocal ones. The plasma dispersion function depends on several parameters: gas pressure, line width of the excited isotopes, resonance frequency of the two counter propagating beams, area of the beam profile at the discharge, and temperature of the plasma. Assuming that β 1 = β 2 = β (a rough evaluation reported in [13] gives β 1 − β 2 10 −14 ), combining the information of the mono-beam DC signals I 1 and I 2 with the plasma dispersion function, the gains of the two beams Gain1 and Gain2 (see Appendix) are evaluated. Since β 1 = β 2 the non reciprocal loss term δμ can be hence determined. Using other relationships reported in the literature [14,15] it is possible to write ω ns1,2 , ω K 1,2 and ω ns K as a function of the known quantities ω m , I 1 , I 2 , I S1 , I S2 , , and μ. As shown in the Appendix, losses appear as a multiplicative factor in ω ns1,2 .
Further to be dependent on losses, terms correcting the null shift depend also on , as well as on geometry and size of the cavity, in a rather complicated way. Therefore, minimization of ω ns1,2 is not straightforward. To give an example using the parameters of GINGERINO, the minimum is found for δμ = 0; assuming typical values in physical units I 1 = I 2 = 1 V and I S1 = I S2 = 0.01 V (see Appendix for the calculation), the null shift correction is of the order of 0.4 μHz, with a modulation induced by variations of of the order of 5 × 10 −6 μHz. In general δμ = 0, as can be inferred from the difference in mono-beam signals. While keeping all other parameters fixed as in the example above, a difference in the demodulated photodiode signals I S1 and I S2 of 1 mV leads to remarkable values of the null shift correction. In such conditions, the null shift contributes to around 500 μHz, with modulations induced by of the order of 0.013 μHz.
The above examples suggest that, in general, variations with are rather small and can be neglected. Therefore, in the ideal case of equal losses the effect of the null shift is a systematic term affecting the accuracy of the measurement at the level of 1.5 × 10 −9 , assuming variations of the losses with time below 10%, while when δμ is not negligible the effect is orders of magnitude higher. In the considered case the requirement on the stability of μ with time is of the order of 1 part in 10 4 in order to achieve a precision of 1 part in 10 9 . As far as accuracy is concerned, the systematic effect is of the order of 1 part in 10 6 , one thousand times the goal of GINGER (1 part in 10 9 at least).
While being less relevant than those related to null shift effects, also the terms with the K suffix deserve discussion. Equation 2 shows that such terms are linearly summed to ω s . As a consequence, their contribution to low frequency signals, i.e., those stimulating the major interest in typical applications of RLGs, can be reduced by low pass filtering. Despite of the possibility to remove related contributions through low pass filtering, the evaluation of ω K requires some further assumptions. ω K oscillates at the Sagnac frequency, containing sin(ω s t) and cos(ω s t) terms. However, due to nonlinear dynamics, the long time average of such oscillating terms can deviate from zero. Therefore, we can express ω K as the product of a small parameter M K with a given function of known parameters. As already mentioned, the value of M K and its very slow variations can be found via linear regression.
A Mathematica notebook is reported in the Appendix to show details of the symbolic calculations leading to ω ns and ω K functions, to demonstrate how losses are a scale factor for the ω ns terms, and to provide numerical examples of the involved values calculated for the GINGERINO setup. We note that, according to our model, ω K can be written as the ratio of two polynomials where the parameter μ enters in one term, preventing direct evaluation of such a parameter through a simple fit procedure. However, we have checked that, with typical μ values of 100 − 1 ppm, variations in calculated results are not significant.
In the following, we focus onto calculations for a general RLG dynamics, and their implementation for the specific case of GINGERINO.

The selected playground and the first step of the analysis
Thirty days of continuous operation between June 16 and July 15 2018 have been selected for this test. After June 21, heavy operations took place in the underground laboratories where GINGERINO is placed for the building of another experiment. Disturbances induced in GINGERINO operation, consisting of mostly mechanical effects, are well visible in the second half of the data stream. We included also these data to show that, even in non ideal conditions, very high sensitivity and stability can be obtained. The data are acquired with our DAQ system at 5 kHz sampling rate [12]. The whole set of data has been analysed on a hourly base, and the terms ω m , ω ns /μ, and ω K 1 /M K (assuming μ = 10 −4 ) were evaluated and stored for offline analysis. In the calculation, average losses μ = 10 −4 are used to evaluate ω K , since it has been checked that the dependence on μ is negligible at the present sensitivity limit. The analysis, based on the Hilbert transform of interferograms and mono-beam signals, leads to the beat note frequency ω m /2π , the amplitudes of the mono beams at the beat frequency I S1 , I S2 , and the relative phase . In order to avoid spurious oscillations due to the boundaries between contiguous hours in the hourlybased analysis, for each hour 6 minutes of data are added at the beginning and at the end; these extra samples are removed after the analysis (overlap-save method). The analysed data are stored after decimation at 20 and 0.1 Hz; the decimation is implemented via the standard Matlab function decimate, which applies 8th order Chebyshev Type I lowpass filter with cutoff frequency 0.8 × (Fs/2), where Fs is 20 or 0.1 Hz.

The second step of the analysis
A second routine collects offline different days for longer analysis. The first operation is to identify the portions of data which are not at normal operation, typically less than 5% of the data are removed. To this aim, we make use of the fringe contrast of the interferogram acquired at 20 Hz rate.
The choice of such a sampling rate minimizes the discarded sample. This operation associates the flag 1 to the good samples and the flag 0 to the bad ones; from this routine we get the mask to select data decimated at lower frequency.
In the present analysis, decimated at half an hour rate, above 90% of data are selected. Both the intensity of the mono-beam signal (I 2 , collected by the photodiode indicated as P H2 in Fig. 1

Reconstruction of the Sagnac frequency by a linear regression model
The vectors ω ns , ω ns × F F DM , ω ns × F F P H2 , ω ξ , and ω K are collected, and the linear regression method is utilised to determine the unknown parameters μ, M K and ξ , and evaluate ω s accordingly. Let us remind that ω s0 accounts for back-scatter noise, and an additional term ω ξ has been implemented to account for inaccuracies on the measured quantities I S1,2 , I 1,2 and [9]. In summary, in the procedure μ and M K are physical quantities which could be evaluated independently, while ξ compensates noises in the measured quantities. Figure 4 shows data before (top panel) and after (bottom panel) the application of this procedure. The dispersion across the mean value is clearly suppressed using ω s , meaning that the new regressive analysis, which includes the terms ω ns and ω K , duly accounts for it. The Matlab function fitlm has been used to perform linear regression. Results indicate that the p-value, which tests the null hypothesis, is approximately zero for all vectors excluding the one related to ω K , whose p-value is 0.038. On the other hand, the relative error in the estimation of the M K term involved in the ω K expres- Fig. 4 Top: time variations expressed as Ω and θ, see Eq. 1, of the analysed data [utilising ω s0 with the mean value subtracted; the mean value is 2π × (280.208 ± 0.001 Hz), compatible with a RLG with area vector vertical within a few mrad error]. Bottom: as above, but using ω s evaluated with the linear regression model. Note: around day 20 some data were lost because of a failure in the data acquisition system sion is around 50%, whereas it amounts to 1 − 12% for all other vectors. 2 Figure 5 shows the contribution of the laser systematics. A dominant role is played by ω ns , which produces a slowly variable level ranging, in units of frequency, between 180-200 mHz, with standard deviation 2.3 mHz. ω ξ gives a small correction with mean value −5.7 mHz and standard deviation 2.3 mHz. These two terms affect residuals, whereas ω K is smaller by more than a factor 10; it has a mean value compatible with zero and standard deviation 0.05 mHz, see Fig. 6. Figure 7 shows the evaluated losses. During the 30 day period considered, a variation of 12% is observed in μ, with a visible trend towards higher losses for increasing time. Figure 8 compares the distributions of the beat frequency ω m /2π , of the ω s0 /2π reconstructed in the first stage of our approach and that of the finally determined ω s /2π . It should be noticed that the distribution of ω m is highly non-Gaussian while that of both reconstructed frequencies becomes close Fig. 5 Top: ω ns evaluated through the best-fit procedure and the one calculated assuming average losses. Bottom: the contribution using the DM and P H2 signals to follow the time behaviour of the losses μ Fig. 6 The first step of the analysis evaluates ω s0 , which accounts for back-scatter noise; ω ξ has been developed to further improve backscatter cancellation, taking into account inaccuracies in the signals produced by photodiodes P H1, P H2, their demodulated intensity I S1,2 and dephasing to Gaussian, indicating that the nonlinear terms of dynamics are correctly accounted for. Moreover, ω s is shifted towards lower frequencies as a consequence of the null shift term ω ns . Figure 9 reports the modified Allan deviation, expressed in angular rotation rate. We observe that, after the application of the procedure (red line), the variance is decreasing with the integration time up to more than 2 days. This suggests that with our procedure we are effectively correcting the long term laser dynamic effects. It is well known that RLGs are sensitive to Chandler Wobble effect, which is typically below Fig. 8 Distribution of the beat frequency ω m /2π , of the first stage of the reconstruction ω s0 /2π and of the final ω s /2π Fig. 9 Modified Allan deviation of the measured angular velocity Ω from the beat note (ω m ), ω s0 and evaluated ω s (blue, green and red lines, respectively), relative to the mean value, expressed in angular velocity [rad/s] 100 nrad. We have checked that the level reached in the long term stability is above the limit imposed by the Chandler Wobble effect. The so far obtained long term instability is about 10 times larger than the scale factor changes induced by temperature variations considering 5 × 10 −6 / o C the thermal expansion coefficient of granite, the material used for the RLG frame, and 0.02 o C as typical temperature variation in the 30 day period. The obtained modified Allan deviation represents a fair improvement compared to the analysis in [12] (the improvement is around a factor 7) carried out by using available signals on a pure phenomenological basis, for instance the long time behaviour was cancelled using the DM signal.

Discussion and conclusions
The Sagnac frequency ω s of the GINGERINO RLG has been evaluated taking into account the laser dynamics. The model leads to a linear sum of several terms, which can be determined with the available signals of beat frequency and mono-beams. The so-called back-scatter noise, which has been so far considered the most severe limitation of RLG, is accounted for by terms ω s0 and ω ξ . It turns out that the null shift term ω ns is the second dominant contribution. It depends linearly on cavity losses μ; this term definitely affects accuracy of the apparatus, since in the best case, if the losses are constant, it represents a non-negligible DC contribution. Furthermore, it is likely that losses are not constant with time, hence ω ns affects stability. If a long term stability of 1 part in 10 9 is required for the Ω ⊕ measurement, the requirement for the stability of μ with time, or equivalently for the relative accuracy in the evaluation of μ, is ω ns ω s × 10 9 . According to the model, it is straightforward to see that the null shift level is minimised when δμ 0.
Changes of μ at the level of 10% are evident in GIN-GERINO. In order to follow those changes we have utilised the DC mono-beam signal of the uncontrolled beam, that is the signal which is not used in feedback loops controlling the laser operation (beam 2 in our case), and the discharge monitor DM, also called gain monitor. ω s is evaluated with a linear regression model.
Remarkably, during the measurements considered in this paper, in particular after the twentieth day of the 30 day record, heavy activity was present for the construction of a new experiment inside the underground laboratory hosting the RLG. Despite that, and the fact that GINGERINO is not equipped with control of the geometrical scale factor (accordingly its wavelength is not fixed), in 5.4 days of integration time it has reached the relative stability of (1.7±0.07)×10 −7 [corresponding to a resolution in the angular velocity measurement of (8.5 ± 0.5) × 10 −12 rad/s].
The analysis poses the problem of the accuracy of angular velocity estimations; however, it is not straightforward to understand the accuracy of the final measurement, an issue to be addressed by future Monte Carlo studies to be carried out following the approach we have already employed in [13].
Main goal of the present analysis was to demonstrate that null shift effects can be evaluated in analytical terms and that they are not negligible at the present sensitivity level of our instrumentation. In order for the linear regression method to give reliable results, the model must be complete. Our data are affected by seismic contributions. They are obviously not considered in the model, but, due to their typical frequency, their role is suppressed thanks to low pass filtering used in our procedure. Furthermore, environmental conditions, in particular temperature, can affect RLG operation. We have verified that their role is negligible for the present data set. In order to check that the procedure does not artificially affect shape and properties of the input signals, we have added a spurious signal synthetically produced to the data, with amplitude of the order of fractions of nrad/s, and verified that it was not modified due to the procedure.
The analysis has also shown that losses are the main limitation of GINGERINO. The question is now which part of the apparatus is the main responsible for losses and for their behaviour as a function of time. The hypothesis is that the main contribution comes from the gain tube of the laser. In the standard high sensitivity RLG scheme, it acts both as gain tube and spatial filter. The laser discharge is in fact produced within a pyrex capillary, with inner diameter of 4 mm, which is a bit small compared to the waist of the beams.
In the near future, a new capillary tube with slightly larger inner diameter will be installed in GINGERINO to limit losses. Furthermore, additional data, such as current of the discharge tube and laser wavelength, will be continuously recorded in order to better characterise the instrument operation in the long term.
As far as the analysis is concerned, the effect of terms θ 1,2 will be eventually included in the model to account for possible deviations of the laser mixture behaviour from the expected one. Furthermore, the data analysis will be extended in order to evaluate the null shift second order terms ω ns2 based on the recipes given in the present paper and to determine the associated correction in data from GINGERINO. Finally, we plan to devote efforts aimed to better understand the influence of each parameter, in particular of the change of wavelength, in the RLG operation, which for the present set of data matters, since GINGERINO is free running and the wavelength can change with even small temperature variations.

Data Availability Statement
This manuscript has no associated data or the data will not be deposited. [Authors' comment: The data of this manuscript will not be deposited. The data can be made available from the corresponding author on reasonable request.] Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecomm ons.org/licenses/by/4.0/. Funded by SCOAP 3 .   In the following a few relations are explicitly reported using the parameters of GINGERINO In this analysis μ indicates the losses of mode 2, while δμ is such that the losses of beam 1 is μ + δμ δμ is perfectly deducible: δμ = DMU