<jupyter_start><jupyter_text>4. Calculate variance of Stokes and anti-Stokes measurements The goal of this notebook is to estimate the variance of the noise of the Stokes measurement. The measured Stokes and anti-Stokes signals contain noise that is distributed approximately normal. We need to estimate the variance of the noise to:- Perform a weighted calibration - Construct confidence intervals<jupyter_code>import os

from dtscalibration import read_silixa_files
from matplotlib import pyplot as plt

%matplotlib inline
filepath = os.path.join('..', '..', 'tests', 'data', 'double_ended2')

ds = read_silixa_files(
    directory=filepath,
    timezone_netcdf='UTC',
    file_ext='*.xml')<jupyter_output>6 files were found, each representing a single timestep
6 recorded vars were found: LAF, ST, AST, REV-ST, REV-AST, TMP
Recorded at 1693 points along the cable
The measurement is double ended
Reading the data from disk<jupyter_text>And we define the sections as we learned from the previous notebook. Sections are required to calculate the variance in the Stokes.<jupyter_code>sections = {
    'probe1Temperature': [slice(7.5, 17.), slice(70., 80.)],  # cold bath
    'probe2Temperature': [slice(24., 34.), slice(85., 95.)],  # warm bath
    }
ds.sections = sections<jupyter_output><empty_output><jupyter_text>Lets first read the documentation about the `ds.variance_stokes` method.<jupyter_code>print(ds.variance_stokes.__doc__) 
I_var, residuals = ds.variance_stokes(st_label='ST')
print("The variance of the Stokes signal along the reference sections "
      "is approximately {} on a {} sec acquisition time".format(I_var, ds.userAcquisitionTimeFW.data[0]))
from dtscalibration import plot

fig_handle = plot.plot_residuals_reference_sections(
        residuals,
        sections,
        title='Distribution of the noise in the Stokes signal',
        plot_avg_std=I_var ** 0.5,
        plot_names=True,
        robust=True,
        units='',
        method='single')<jupyter_output>/Users/bfdestombe/Projects/dts-calibration/python-dts-calibration-dev/src/dtscalibration/plot.py:301: FutureWarning: Conversion of the second argument of issubdtype from `float` to `np.floating` is deprecated. In future, it will be treated as `np.float64 == np.dtype(float).type`.
  if (np.issubdtype(resid[time_dim].dtype, np.float) or
/Users/bfdestombe/Projects/dts-calibration/python-dts-calibration-dev/src/dtscalibration/plot.py:302: FutureWarning: Conversion of the second argument of issubdtype from `int` to `np.signedinteger` is deprecated. In future, it will be treated as `np.int64 == np.dtype(int).type`.
  np.issubdtype(resid[time_dim].dtype, np.int)):
/Users/bfdestombe/Projects/dts-calibration/python-dts-calibration-dev/.tox/docs/lib/python3.7/site-packages/numpy/lib/nanfunctions.py:1667: RuntimeWarning: Degrees of freedom <= 0 for slice.
  keepdims=keepdims)
/Users/bfdestombe/Projects/dts-calibration/python-dts-calibration-dev/.tox/docs/lib/python3.7/site-packages/xarray/core/na[...]<jupyter_text>The residuals should be normally distributed and independent from previous time steps and other points along the cable. If you observe patterns in the residuals plot (above), it might be caused by:- The temperature in the calibration bath is not uniform- Attenuation caused by coils/sharp bends in cable- Attenuation caused by a splice<jupyter_code>import scipy
import numpy as np

sigma = residuals.std()
mean = residuals.mean()
x = np.linspace(mean - 3*sigma, mean + 3*sigma, 100)
approximated_normal_fit = scipy.stats.norm.pdf(x, mean, sigma)
residuals.plot.hist(bins=50, figsize=(12, 8), density=True)
plt.plot(x, approximated_normal_fit);<jupyter_output><empty_output>
StarcoderdataJupyterStructuredCleanDedup