{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 15. Calibration using matching sections"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In notebook 14 we showed how you can take splices or connectors within your calibration into account. To then calibrate the cable we used reference sections on both sides of the splice. If these are not available, or in other cases where you have a lack of reference sections, matching sections can be used to improve the calibration.\n",
    "\n",
    "For matching sections you need two sections of fiber than you know will be the exact same temperature. This can be, for example, in duplex cables or twisted pairs of cable."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Demonstration\n",
    "To demonstrate matching sections, we'll load the same dataset that was used in previous notebooks, and modify the data to simulate a lossy splice, just as in notebook 14."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-04-06T08:11:42.175912Z",
     "iopub.status.busy": "2022-04-06T08:11:42.175058Z",
     "iopub.status.idle": "2022-04-06T08:11:43.633971Z",
     "shell.execute_reply": "2022-04-06T08:11:43.633285Z"
    }
   },
   "outputs": [],
   "source": [
    "import os\n",
    "\n",
    "from dtscalibration import read_silixa_files\n",
    "\n",
    "# The following line introduces the .dts accessor for xarray datasets\n",
    "import dtscalibration  # noqa: E401  # noqa: E401\n",
    "from dtscalibration.variance_stokes import variance_stokes_constant\n",
    "\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-04-06T08:11:43.636636Z",
     "iopub.status.busy": "2022-04-06T08:11:43.636444Z",
     "iopub.status.idle": "2022-04-06T08:11:43.916858Z",
     "shell.execute_reply": "2022-04-06T08:11:43.916339Z"
    }
   },
   "outputs": [],
   "source": [
    "filepath = os.path.join(\"..\", \"..\", \"tests\", \"data\", \"double_ended2\")\n",
    "\n",
    "ds_ = read_silixa_files(directory=filepath, timezone_netcdf=\"UTC\", file_ext=\"*.xml\")\n",
    "\n",
    "ds = ds_.sel(x=slice(0, 110))  # only calibrate parts of the fiber\n",
    "\n",
    "\n",
    "sections = {\n",
    "    \"probe1Temperature\": [slice(7.5, 17.0)],  # cold bath\n",
    "    \"probe2Temperature\": [slice(24.0, 34.0)],  # warm bath\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Again, we introduce a step loss in the signal strength at x = 50 m. For the forward channel, this means all data beyond 50 meters is reduced with a 'random' factor. For the backward channel, this means all data up to 50 meters is reduced with a 'random' factor."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-04-06T08:11:43.919445Z",
     "iopub.status.busy": "2022-04-06T08:11:43.919271Z",
     "iopub.status.idle": "2022-04-06T08:11:43.940161Z",
     "shell.execute_reply": "2022-04-06T08:11:43.939718Z"
    }
   },
   "outputs": [],
   "source": [
    "ds[\"st\"] = ds.st.where(ds.x < 50, ds.st * 0.8)\n",
    "ds[\"ast\"] = ds.ast.where(ds.x < 50, ds.ast * 0.82)\n",
    "\n",
    "ds[\"rst\"] = ds.rst.where(ds.x > 50, ds.rst * 0.85)\n",
    "ds[\"rast\"] = ds.rast.where(ds.x > 50, ds.rast * 0.81)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will first run a calibration without adding the transient attenuation location or matching sections. A big jump in the calibrated temperature is visible at x = 50. \n",
    "\n",
    "As all calibration sections are before 50 meters, the first 50 m will be calibrated correctly."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-04-06T08:11:43.942589Z",
     "iopub.status.busy": "2022-04-06T08:11:43.942387Z",
     "iopub.status.idle": "2022-04-06T08:11:44.853767Z",
     "shell.execute_reply": "2022-04-06T08:11:44.853280Z"
    }
   },
   "outputs": [],
   "source": [
    "st_var, resid = variance_stokes_constant(\n",
    "    ds.dts.st, sections, ds.dts.acquisitiontime_fw, reshape_residuals=True\n",
    ")\n",
    "ast_var, _ = variance_stokes_constant(\n",
    "    ds.dts.ast, sections, ds.dts.acquisitiontime_fw, reshape_residuals=False\n",
    ")\n",
    "rst_var, _ = variance_stokes_constant(\n",
    "    ds.dts.rst, sections, ds.dts.acquisitiontime_bw, reshape_residuals=False\n",
    ")\n",
    "rast_var, _ = variance_stokes_constant(\n",
    "    ds.dts.rast, sections, ds.dts.acquisitiontime_bw, reshape_residuals=False\n",
    ")\n",
    "\n",
    "out = ds.dts.calibrate_double_ended(\n",
    "    sections=sections,\n",
    "    st_var=st_var,\n",
    "    ast_var=ast_var,\n",
    "    rst_var=rst_var,\n",
    "    rast_var=rast_var,\n",
    ")\n",
    "\n",
    "out.isel(time=0).tmpw.plot(label=\"calibrated\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we run a calibration, adding the keyword argument '**trans_att**', and provide a list of floats containing the locations of the splices. In this case we only add a single one at x = 50 m.\n",
    "\n",
    "We will also define the matching sections of cable. The matching sections have to be provided as a list of tuples. A tuple per matching section. Each tuple has three items, the first two items are the slices of the sections that are matching. The third item is a bool and is True if the two sections have a reverse direction (as in the \"J-configuration\").\n",
    "\n",
    "In this example we match the two cold baths to each other.\n",
    "\n",
    "After running the calibration you will see that by adding the transient attenuation and matching sections the calibration returns the correct temperature, without the big jump.\n",
    "\n",
    "*In single-ended calibration the keyword is called '**trans_att**'.*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-04-06T08:11:44.856090Z",
     "iopub.status.busy": "2022-04-06T08:11:44.855889Z",
     "iopub.status.idle": "2022-04-06T08:11:46.067901Z",
     "shell.execute_reply": "2022-04-06T08:11:46.067290Z"
    }
   },
   "outputs": [],
   "source": [
    "matching_sections = [(slice(7.5, 17.6), slice(69, 79.1), False)]\n",
    "\n",
    "out2 = ds.dts.calibrate_double_ended(\n",
    "    sections=sections,\n",
    "    st_var=st_var,\n",
    "    ast_var=ast_var,\n",
    "    rst_var=rst_var,\n",
    "    rast_var=rast_var,\n",
    "    trans_att=[50.0],\n",
    "    matching_sections=matching_sections,\n",
    ")\n",
    "\n",
    "out2.isel(time=0).tmpw.plot(label=\"calibrated\")\n",
    "\n",
    "out.isel(time=0).tmpw.plot(label=\"normal calibration\")\n",
    "out2.isel(time=0).tmpw.plot(label=\"matching sections\")\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}