notebooks/stage-one-analysis.ipynb

{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "85abd34c",
   "metadata": {},
   "source": [
    "# TGP Stage 1 analysis\n",
    "\n",
    "This notebook performs an interactive analysis of one of the **Stage 1** datasets.\n",
    "It starts by loading the raw data and proceeds with pre-processing and analysis steps, where each step is visualized.\n",
    "The goal of this analysis is to extract an interesting region in the magnetic field ($B$) and plunger voltage ($V$) space where the device is both gapped and has zero-bias peaks (ZBPs).\n",
    "The extracted regions with gapped ZPBs are where Stage 2 measurements will be performed."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9fc86e8a",
   "metadata": {},
   "outputs": [],
   "source": [
    "from pathlib import Path\n",
    "\n",
    "import tgp\n",
    "import xarray as xr\n",
    "\n",
    "tgp.plot.set_mpl_rc()\n",
    "\n",
    "folder = Path(\"..\") / \"data\"\n",
    "\n",
    "print(\"Using tgp code version:\", tgp.__version__)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "163127f2",
   "metadata": {},
   "source": [
    "# Analysis on experimental data\n",
    "\n",
    "We load the raw data and prepare it for the `tgp` code's analysis.\n",
    "The `prepare` function performs renames and joins the left and right cutter arrays into a single index (if required)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "35db4d21",
   "metadata": {},
   "outputs": [],
   "source": [
    "ds = xr.load_dataset(folder / \"experimental\" / \"deviceA_stage1.nc\", engine=\"h5netcdf\")\n",
    "ds = tgp.prepare.prepare(ds)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "fad9e4de",
   "metadata": {},
   "source": [
    "The dataset we have just loaded contains many data variables.\n",
    "\n",
    "The following are the *raw* measured complex arrays\n",
    "- `I_1w_L`, `I_1w_R`: the local $1\\omega$ component of the current (corresponds to the local differential conductance, which the 1st derivative of current with respect to voltage) on the left and right\n",
    "- `I_2w_L` and `I_2w_R`: the local $2\\omega$ componentof the  current (2nd derivative of current with respect to voltage) on the left and right\n",
    "- `I_2w_RL`, `I_2w_LR`: the non-local $2\\omega$ component of the current on the left and right\n",
    "- `I_3w_L`, `I_3w_R`: the local $3\\omega$ component of the current (3rd derivative of current with respect to voltage) on the left and right\n",
    "\n",
    "These are converted to the following arrays (*used in the analysis*) by only keeping the real part\n",
    "- `L_2w_nl`, `R_2w_nl`: the non-local $2\\omega$ component of the current on the left and right\n",
    "- `L_3w`, `R_3w`: the local $3\\omega$ component of the current on the left and right\n",
    "\n",
    "We use $3\\omega$ local signals to identify ZBP peaks and $2\\omega$ non-local signals to distinguish gapped and gapless states. "
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "51909125",
   "metadata": {
    "lines_to_next_cell": 0
   },
   "source": [
    "## Visualize the data\n",
    "\n",
    "We can view the $1\\omega$ signals where we plot the real parts of the conductance with magnetic field on the x-axis and plunger voltage on the y-axis.\n",
    "Every separate plot is for a different cutter value pair.\n",
    "First we show the conductance on the left then on the right."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "91336770",
   "metadata": {},
   "outputs": [],
   "source": [
    "lim_1w = 1.0\n",
    "\n",
    "tgp.one.print_lim(ds.L_1w, lim_1w)\n",
    "tgp.one.print_lim(ds.R_1w, lim_1w)\n",
    "\n",
    "# With no labels and titles:\n",
    "# tgp.plot.one.plot_1w(ds, lim_1w, plot_kwargs=dict(figsize=(12, 8)), minimal=True)\n",
    "\n",
    "tgp.plot.one.plot_1w(ds, lim_1w)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "1ed04e3e",
   "metadata": {},
   "source": [
    "Similar to the plots above, we plot the $2\\omega$ nonlocal signals."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3afd236f",
   "metadata": {},
   "outputs": [],
   "source": [
    "lim_2w = 1e3\n",
    "\n",
    "tgp.one.print_lim(ds.L_2w_nl, lim_1w)\n",
    "tgp.one.print_lim(ds.R_2w_nl, lim_1w)\n",
    "\n",
    "tgp.plot.one.plot_2w(ds, lim_2w)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "40f5a3fa",
   "metadata": {},
   "outputs": [],
   "source": [
    "tgp.plot.one.plot_2w_at(ds, cutter_pair_index=0, vmax=lim_2w)  # change the cutter_pair_index to see different plots"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "ea0bff68",
   "metadata": {
    "lines_to_next_cell": 0
   },
   "source": [
    "## Extract the gap from $2\\omega$ non-local signal\n",
    "\n",
    "After we have inspected the data, we are ready to perform the first analysis step.\n",
    "The $2\\omega$ component of the conductance informs us about the energy gap of the spectrum.\n",
    "We threshold the non-local $2\\omega$ component above its noisy background to get a boolean array of the gap.\n",
    "This works because the non-local conductance is suppressed below the gap."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "80bfb4e3",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Manually set the threshold region\n",
    "# tgp.one.set_2w_th(ds, B_max=0.7, V_max=-1.410, verbose=True)\n",
    "\n",
    "# Automatically set\n",
    "tgp.one.set_2w_th(ds)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "677a938f",
   "metadata": {},
   "source": [
    "The result is an array with the dimensions `(cutter_pair_index, B, V)`.\n",
    "We can visualize this by plotting the gap per cutter pair value."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "00cc62b8",
   "metadata": {},
   "outputs": [],
   "source": [
    "tgp.plot.one.plot_2w_th(ds)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "5d63fea4",
   "metadata": {},
   "source": [
    "The `cutter_pair_index` dimension is averaged out, which results in the plot below. This shows which fraction of the cutter pairs is gapped for each point in $(B, V)$ space."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3739ef3d",
   "metadata": {},
   "outputs": [],
   "source": [
    "tgp.plot.one.plot_2w_th_avg(ds, vmax=None)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "7d3066c6",
   "metadata": {},
   "source": [
    "To convert this array of fraction to a boolean array, were `True` is gapped and `False` is gapless, we threshold the array again.\n",
    "We require that the system is gapped for at least 50\\% of the cutter value pairs."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3213ade9",
   "metadata": {},
   "outputs": [],
   "source": [
    "tgp.one.set_gapped(ds, th_2w_p=0.5)\n",
    "tgp.plot.one.plot_gapped(ds)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "5109f7f2",
   "metadata": {},
   "source": [
    "## Extract the ZBPs from $3\\omega$ signal\n",
    "\n",
    "We extract whether ZBPs occur based on the $3\\omega$ component.\n",
    "This works because the 3rd derivative of the conductance informs us about the curvature of the conductance with respect to bias voltage.\n",
    "We first plot the data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5ece3fe5",
   "metadata": {},
   "outputs": [],
   "source": [
    "lim_3w = 1e9\n",
    "\n",
    "tgp.one.print_lim(ds.L_3w, lim_1w)\n",
    "tgp.one.print_lim(ds.R_3w, lim_1w)\n",
    "\n",
    "tgp.plot.one.plot_3w(ds, vmin=-lim_3w, vmax=lim_3w)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3e2f5501",
   "metadata": {},
   "outputs": [],
   "source": [
    "tgp.plot.one.plot_3w_at(ds, cutter_pair_index=0, lim=lim_3w)  # change the cutter_pair_index to see different plots"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "90b7afb6",
   "metadata": {},
   "source": [
    "To go to a boolean array of ZPBs, we threshold the $3\\omega$ data below a negative number, this corresponds with the blue areas in the plots above.\n",
    "The result of the thresholding operation, and its resulting ZBPs, is plotted below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4062c3ec",
   "metadata": {},
   "outputs": [],
   "source": [
    "tgp.one.set_3w_th(ds, th_3w=1e7)\n",
    "tgp.plot.one.plot_3w_th(ds)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "5d75744f",
   "metadata": {},
   "source": [
    "The `cutter_pair_index` dimension is averaged out, which results in the plot below. This shows which fraction of the cutter pairs for each point in $(B, V)$ space we get ZBPs."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7a9f8ca9",
   "metadata": {},
   "outputs": [],
   "source": [
    "tgp.plot.one.plot_zbp(ds)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "72050e18",
   "metadata": {},
   "source": [
    "We require that at least 50\\% of the cutter pair values show ZBPs."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "33be7cad",
   "metadata": {},
   "outputs": [],
   "source": [
    "tgp.one.set_3w_tat(ds, th_3w_tat=0.5)  # try, 0.3, 0.5, 0.7\n",
    "tgp.plot.one.plot_3w_tat(ds)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "bb3c7b1d",
   "metadata": {},
   "source": [
    "We cluster the ZPB regions and plot its result."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cab9ede5",
   "metadata": {},
   "outputs": [],
   "source": [
    "tgp.one.set_clusters(ds)\n",
    "tgp.plot.one.plot_clusters(ds)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "42987cdd",
   "metadata": {},
   "source": [
    "## Combining gap and ZBP\n",
    "\n",
    "Finally, we combine the gap and ZBP arrays to get the following."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "88910df4",
   "metadata": {},
   "outputs": [],
   "source": [
    "# This variable indicates we include `clusters_n` of the largest clusters in the plots\n",
    "clusters_n = 11\n",
    "tgp.plot.one.plot_stage_1(ds, clusters_n)\n",
    "tgp.plot.one.plot_stage_1_gapped_clusters(ds, clusters_n)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "4dae2c69",
   "metadata": {},
   "source": [
    "We can plot the interesting regions from most promising to least promising (based on its size.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5367ffb7",
   "metadata": {},
   "outputs": [],
   "source": [
    "zoomin_ranges = tgp.one.get_zoomin_ranges(ds, clusters_n, zoomin_V_height=0.01)\n",
    "tgp.plot.one.plot_zoomed_clusters(ds, zoomin_ranges, clusters_n)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "87af59b3",
   "metadata": {},
   "source": [
    "Alternatively, we can directly see the numerical values of the regions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1e8f393b",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Get bounding boxes and area of the clusters\n",
    "tgp.common.cluster_infos(\n",
    "    ds[\"clusters\"],\n",
    "    pct_box=5,  # Margin around box of at least 5%\n",
    "    min_margin_box=(0.1, 0.002), # and at least 0.1 T and 0.002 V\n",
    ")"
   ]
  }
 ],
 "metadata": {
  "language_info": {
   "name": "python",
   "pygments_lexer": "ipython3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}