{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Import supporting package"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import xarray as xr\n",
    "import numpy as np\n",
    "import copy\n",
    "\n",
    "from uncertainties import ufloat\n",
    "from uncertainties import unumpy as unp\n",
    "from uncertainties import umath\n",
    "import random\n",
    "import matplotlib.pyplot as plt\n",
    "plt.rcParams['font.size'] = 12\n",
    "\n",
    "from DataContainer.ReadData import read_hdf5_file\n",
    "from Analyser.ImagingAnalyser import ImageAnalyser\n",
    "from Analyser.FitAnalyser import FitAnalyser\n",
    "from Analyser.FitAnalyser import NewFitModel, DensityProfileBEC2dModel\n",
    "from ToolFunction.ToolFunction import *\n",
    "\n",
    "from scipy.optimize import curve_fit\n",
    "\n",
    "from ToolFunction.HomeMadeXarrayFunction import errorbar, dataarray_plot_errorbar\n",
    "xr.plot.dataarray_plot.errorbar = errorbar\n",
    "xr.plot.accessor.DataArrayPlotAccessor.errorbar = dataarray_plot_errorbar\n",
    "\n",
    "imageAnalyser = ImageAnalyser()\n",
    "\n",
    "import Analyser.FitAnalyser as fa\n",
    "# %matplotlib notebook"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Start a client for parallel computing"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "dict_keys(['camera_3'])\n",
      "The detected scaning axes and values are: \n",
      "\n",
      "{}\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.collections.QuadMesh at 0x2729cb61be0>"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "[encoded data, 184140 characters]",
      "text/plain": [
       "<Figure size 640x480 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "groupList = [\n",
    "    \"images/MOT_3D_Camera/in_situ_absorption\",\n",
    "    \"images/ODT_1_Axis_Camera/in_situ_absorption\",\n",
    "    \"images/ODT_2_Axis_Camera/in_situ_absorption\",\n",
    "    \"images/Vertical_Axis_Camera/in_situ_absorption\",\n",
    "]\n",
    "\n",
    "dskey = {\n",
    "    \"images/MOT_3D_Camera/in_situ_absorption\": \"camera_0\",\n",
    "    \"images/ODT_1_Axis_Camera/in_situ_absorption\": \"camera_1\",\n",
    "    \"images/ODT_2_Axis_Camera/in_situ_absorption\": \"camera_2\",\n",
    "    \"images/Vertical_Axis_Camera/in_situ_absorption\": \"camera_3\"\n",
    "}\n",
    "\n",
    "img_dir = '//DyLabNAS/Data/'\n",
    "SequenceName = \"Evaporative_Cooling\" + \"/\"\n",
    "folderPath = img_dir + SequenceName + '2023/08/30'# get_date()\n",
    "\n",
    "shotNum = f\"0077\"\n",
    "filePath = folderPath + \"/\" + shotNum + \"/*.h5\"\n",
    "\n",
    "dataSetDict = {\n",
    "    dskey[groupList[i]]: read_hdf5_file(filePath, groupList[i])\n",
    "    for i in [3]\n",
    "}\n",
    "print(dataSetDict.keys())\n",
    "\n",
    "dataSet = dataSetDict[\"camera_3\"]\n",
    "\n",
    "\n",
    "print_scanAxis(dataSet)\n",
    "\n",
    "scanAxis = get_scanAxis(dataSet)\n",
    "\n",
    "dataSet = auto_rechunk(dataSet)\n",
    "\n",
    "dataSet = imageAnalyser.get_absorption_images(dataSet)\n",
    "\n",
    "imageAnalyser.center = (550, 750)\n",
    "imageAnalyser.span = (900, 900)\n",
    "imageAnalyser.fraction = (0.1, 0.1)\n",
    "\n",
    "dataSet_cropOD = imageAnalyser.crop_image(dataSet.OD)\n",
    "dataSet_cropOD = imageAnalyser.substract_offset(dataSet_cropOD).load()\n",
    "\n",
    "dataSet_cropOD.plot.pcolormesh(cmap='jet', col=scanAxis[0], row=scanAxis[1], vmin=0, vmax=3)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.collections.QuadMesh at 0x20902035310>"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "[encoded data, 187112 characters]",
      "text/plain": [
       "<Figure size 640x480 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "dataSet_cropOD.plot.pcolormesh(cmap='jet', col=scanAxis[0], row=scanAxis[1], vmin=0, vmax=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "from scipy.ndimage import rotate"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "data_rot = copy.deepcopy(dataSet_cropOD)\n",
    "data_rot.data = rotate(data_rot.data, 36, reshape=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "x = dataSet_cropOD.x.data\n",
    "y = dataSet_cropOD.y.data\n",
    "\n",
    "X,Y = np.meshgrid(x,y)\n",
    "X_ = X.flatten()\n",
    "Y_ = Y.flatten()\n",
    "# print(y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "ename": "type",
     "evalue": "only integers, slices (`:`), ellipsis (`...`), numpy.newaxis (`None`) and integer or boolean arrays are valid indices",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mIndexError\u001b[0m                                Traceback (most recent call last)",
      "Cell \u001b[1;32mIn[25], line 1\u001b[0m\n\u001b[1;32m----> 1\u001b[0m params\u001b[39m.\u001b[39;49mdata[\u001b[39m'\u001b[39;49m\u001b[39mx0_bec\u001b[39;49m\u001b[39m'\u001b[39;49m]\n",
      "\u001b[1;31mIndexError\u001b[0m: only integers, slices (`:`), ellipsis (`...`), numpy.newaxis (`None`) and integer or boolean arrays are valid indices"
     ]
    }
   ],
   "source": [
    "params.data['x0_bec']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "shape: (900, 900)\n"
     ]
    },
    {
     "data": {
      "image/png": "[encoded data, 42752 characters]",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x smaller y, 1d fit along x\n",
      "\n",
      "1d fit initialization\n",
      "center = [419. 428.]\n",
      "BEC widths: [ 80 339]\n",
      "\n",
      "1d init fit values\n",
      "Name          Value      Min      Max   Stderr     Vary     Expr Brute_Step\n",
      "amp_bec       1.145        0    2.977     None     True     None     None\n",
      "amp_th        1.145        0    2.977     None     True     None     None\n",
      "deltax          240        0      900     None     True     None     None\n",
      "sigma_bec     65.57        0      160     None     True     None     None\n",
      "sigma_th      165.8        0      inf     None    False 0.632*sigma_bec + 0.518*deltax     None\n",
      "x0_bec          419      409      429     None     True     None     None\n",
      "x0_th           419      409      429     None     True     None     None\n",
      "1d fitted values\n",
      "x0_bec: 421.007, (init = 419.000), bounds = [409.00 : 429.00] \n",
      "x0_th: 409.000, (init = 419.000), bounds = [409.00 : 429.00] \n",
      "amp_bec: 1.923, (init = 1.145), bounds = [0.00 : 2.98] \n",
      "amp_th: 0.351, (init = 1.145), bounds = [0.00 : 2.98] \n",
      "sigma_bec: 63.902, (init = 65.574), bounds = [0.00 : 160.00] \n",
      "sigma_th: 40.386, (init = 165.763), bounds = [0.00 : inf] \n",
      "\n"
     ]
    },
    {
     "data": {
      "image/png": "[encoded data, 44276 characters]",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Init Params\n",
      "Name                    Value      Min      Max   Stderr     Vary     Expr Brute_Step\n",
      "amp_bec                 2.348        0    3.812     None     True     None     None\n",
      "amp_th                 0.4282        0    3.812     None     True     None     None\n",
      "atom_number_bec      5.239e+04     -inf      inf     None    False amp_bec / 5 * 2 * 3.14159265359 * sigmax_bec * sigmay_bec     None\n",
      "atom_number_th           3211     -inf      inf     None    False amp_th * 2 * 3.14159265359 * 1.20206 / 1.643 * sigma_th * sigma_th     None\n",
      "condensate_fraction    0.9423     -inf      inf     None    False atom_number_bec / (atom_number_bec + atom_number_th)     None\n",
      "rot_angle                 -36      -66       -6     None    False     None     None\n",
      "sigma_th                40.39        0      900     None     True     None     None\n",
      "sigmax_bec               63.9        0      160     None     True     None     None\n",
      "sigmay_bec              277.9        0      678     None     True     None     None\n",
      "x0_bec                    419      409      429     None     True     None     None\n",
      "x0_th                     419      409      429     None     True     None     None\n",
      "y0_bec                    428      418      438     None     True     None     None\n",
      "y0_th                     428      418      438     None     True     None     None\n",
      "\n",
      "bval first fit\n",
      "amp_bec: 2.632, (init = 2.348), bounds = [0.00 : 3.81] \n",
      "amp_th: 0.000, (init = 0.428), bounds = [0.00 : 3.81] \n",
      "x0_bec: 409.000, (init = 419.000), bounds = [409.00 : 429.00] \n",
      "y0_bec: 438.000, (init = 428.000), bounds = [418.00 : 438.00] \n",
      "x0_th: 409.002, (init = 419.000), bounds = [409.00 : 429.00] \n",
      "y0_th: 438.000, (init = 428.000), bounds = [418.00 : 438.00] \n",
      "sigmax_bec: 69.914, (init = 63.902), bounds = [0.00 : 160.00] \n",
      "sigmay_bec: 283.514, (init = 277.869), bounds = [0.00 : 678.00] \n",
      "sigma_th: 65.907, (init = 40.386), bounds = [0.00 : 900.00] \n",
      "rot_angle: -36.000, (init = -36.000), bounds = [-66.00 : -6.00] \n",
      "\n"
     ]
    }
   ],
   "source": [
    "data =  dataSet_cropOD\n",
    "fitModel = DensityProfileBEC2dModel(is_debug=True)\n",
    "fitAnalyser_1 = FitAnalyser(fitModel, fitDim=2)\n",
    "params = fitAnalyser_1.guess(data, guess_kwargs=dict(rot_angle=36, vary_rot=False, pre_check=False, post_check=False))\n",
    "result = fitAnalyser_1.fit(data, paramsArray=params)\n",
    "# result = fitModel.fit(data.data.flatten('F'), x= X_, y=Y_)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'plt' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "Cell \u001b[1;32mIn[1], line 2\u001b[0m\n\u001b[0;32m      1\u001b[0m \u001b[39m# fig, axs = plt.subplots(shape[0] , 5, figsize=(14, 4 * shape[0] ),dpi = 100)\u001b[39;00m\n\u001b[1;32m----> 2\u001b[0m fig, axs \u001b[39m=\u001b[39m plt\u001b[39m.\u001b[39msubplots(\u001b[39m1\u001b[39m , \u001b[39m3\u001b[39m, figsize\u001b[39m=\u001b[39m(\u001b[39m13\u001b[39m,\u001b[39m4\u001b[39m  ),dpi \u001b[39m=\u001b[39m \u001b[39m100\u001b[39m)\n\u001b[0;32m      4\u001b[0m vmax\u001b[39m=\u001b[39m\u001b[39m2\u001b[39m\n\u001b[0;32m      5\u001b[0m cropOD \u001b[39m=\u001b[39m data\u001b[39m.\u001b[39mto_numpy()\n",
      "\u001b[1;31mNameError\u001b[0m: name 'plt' is not defined"
     ]
    }
   ],
   "source": [
    "# fig, axs = plt.subplots(shape[0] , 5, figsize=(14, 4 * shape[0] ),dpi = 100)\n",
    "fig, axs = plt.subplots(1 , 3, figsize=(13,4  ),dpi = 100)\n",
    "\n",
    "vmax=2\n",
    "cropOD = data.to_numpy()\n",
    "\n",
    "\n",
    "# axs[i,0].set_title(f'image {i}, cond. frac = {Fitmodel.cond_frac(result[i], X, Y) *1e2:.2f} %', loc='left', fontsize=18)\n",
    "man_fit = fa.density_profile_BEC_2d(X,Y, amp_bec=2.348, amp_th=0.4282, sigma_th=40.39, sigmax_bec=64, sigmay_bec=278, \n",
    "                                x0_bec=419,x0_th=419,y0_bec=428, y0_th=428, rot_angle=-36)\n",
    "\n",
    "\n",
    "ax = axs[0]\n",
    "ax.pcolormesh(X, Y, cropOD.T, vmin=0, vmax=3, cmap='jet', shading='auto')\n",
    "#plt.colorbar(art, ax=ax, label='z')\n",
    "\n",
    "\n",
    "# Plot gaussian 2d Fit + legend including Width parameters\n",
    "ax = axs[1]\n",
    "\n",
    "ax.pcolormesh(X, Y, man_fit, cmap='jet', shading='auto')\n",
    "#plt.colorbar(art, ax=ax, label='z')\n",
    "\n",
    "ax = axs[2]\n",
    "\n",
    "ax.pcolormesh(X, Y, man_fit-cropOD.T, vmin=0, vmax=0.2, cmap='jet', shading='auto')\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "0\n",
      "0.5\n",
      "0.2\n",
      "0.5\n",
      "0.2\n",
      "1.3375883632180532\n",
      "\n",
      "30\n",
      "0.5\n",
      "0.2\n",
      "0.5330127018922194\n",
      "-0.07679491924311221\n",
      "1.3375883632180532\n",
      "\n",
      "60\n",
      "0.5\n",
      "0.2\n",
      "0.4232050807568878\n",
      "-0.33301270189221926\n",
      "1.3375883632180534\n",
      "\n",
      "90\n",
      "0.5\n",
      "0.2\n",
      "0.20000000000000004\n",
      "-0.5\n",
      "1.3375883632180532\n",
      "\n",
      "120\n",
      "0.5\n",
      "0.2\n",
      "-0.07679491924311213\n",
      "-0.5330127018922193\n",
      "1.3375883632180536\n",
      "\n",
      "150\n",
      "0.5\n",
      "0.2\n",
      "-0.3330127018922194\n",
      "-0.42320508075688773\n",
      "1.3375883632180532\n"
     ]
    }
   ],
   "source": [
    "import Analyser.FitAnalyser as fa\n",
    "import numpy as np\n",
    "for angle in np.arange(0,180,30):\n",
    "    print()\n",
    "    print(angle)\n",
    "    print(fa.density_profile_BEC_2d(0.5,0.2, rot_angle=angle))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "x = dataSet_cropOD.x.data\n",
    "y = dataSet_cropOD.y.data\n",
    "\n",
    "X,Y = np.meshgrid(x,y)\n",
    "X_ = X.flatten()\n",
    "Y_ = Y.flatten()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "keys = ['amp_bec', 'amp_th','x0_bec','y0_bec','x0_th','y0_th','sigmax_bec','sigmay_bec','sigma_th','rot_angle']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'amp_bec': 2.63233857911176, 'amp_th': 2.2068623923891248e-06, 'x0_bec': 409.00000000030843, 'y0_bec': 437.9999999999945, 'x0_th': 409.00151243134604, 'y0_th': 437.99999780448076, 'sigmax_bec': 69.91365732143392, 'sigmay_bec': 283.51382854536473, 'sigma_th': 65.90653264163306, 'rot_angle': -36.0}\n"
     ]
    }
   ],
   "source": [
    "res = fitAnalyser_1.get_fit_value(result)\n",
    "bval = {}\n",
    "for key in keys:\n",
    "    bval[key] = float(res[key].data)\n",
    "print(bval)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "[encoded data, 298460 characters]",
      "text/plain": [
       "<Figure size 2000x400 with 5 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# fig, axs = plt.subplots(shape[0] , 5, figsize=(14, 4 * shape[0] ),dpi = 100)\n",
    "shape = [1]\n",
    "fig, axs = plt.subplots(shape[0] , 5, figsize=(20,4 * shape[0] ),dpi = 100)\n",
    "# fig.suptitle(f'{file}_{saveprefix}: ν = ({data_x.trap_f[0]}, {data_x.trap_f[1]}, {data_x.trap_f[2]}) Hz, T = {data_x.temp[0]*1e9:.0f}nK, tof = {data_x.tof * 1e3:.0f}ms, a = {data_x.a_fac :.0f}a$_0$, N = {data_x.N:.0e} ', y=-0.005*shape[0]+0.965, fontsize=24)\n",
    "# fig.suptitle('test', va='bottom')\n",
    "cropOD = data.to_numpy()\n",
    "for i in range(0,shape[0]):\n",
    "\n",
    "    # axs[i,0].set_title(f'image {i}, cond. frac = {Fitmodel.cond_frac(result[i], X, Y) *1e2:.2f} %', loc='left', fontsize=18)\n",
    "    man_fit = fa.density_profile_BEC_2d(X,Y, amp_bec=2.348, amp_th=0.4282, sigma_th=40.39, sigmax_bec=64, sigmay_bec=278, \n",
    "                                    x0_bec=419,x0_th=419,y0_bec=428, y0_th=428, rot_angle=-36)\n",
    "    \n",
    "    # bval = result.best_values\n",
    "    fit = fa.density_profile_BEC_2d(X,Y, **bval)\n",
    "    vmax = np.max(cropOD)\n",
    "    # vmax = 1\n",
    "\n",
    "    if bval['amp_bec'] > bval['amp_th']:\n",
    "        cen_str = 'bec'\n",
    "    else:\n",
    "        cen_str = 'th'\n",
    "    cen_x = round(bval[f'x0_{cen_str}'])\n",
    "    cen_y = round(bval[f'y0_{cen_str}'])\n",
    "\n",
    "    ax = axs[0]\n",
    "    ax.pcolormesh(X, Y, cropOD.T, vmin=0, vmax=vmax, cmap='jet', shading='auto')\n",
    "    #plt.colorbar(art, ax=ax, label='z')\n",
    "\n",
    "\n",
    "    # Plot gaussian 2d Fit + legend including Width parameters\n",
    "    ax = axs[1]\n",
    "\n",
    "    ax.pcolormesh(X, Y, fit, cmap='jet', shading='auto')\n",
    "    #plt.colorbar(art, ax=ax, label='z')\n",
    "\n",
    "    ax = axs[2]\n",
    "\n",
    "    ax.pcolormesh(X, Y, fit-cropOD.T, vmin=0, vmax=0.2, cmap='jet', shading='auto')\n",
    "\n",
    "    for i in range(0,2):\n",
    "        axs[i].hlines(bval['y0_bec'], 0,900)\n",
    "        axs[i].hlines(bval['y0_th'], 0,900)\n",
    "\n",
    "        axs[i].vlines(bval['x0_bec'], 0,900)\n",
    "        axs[i].vlines(bval['x0_th'], 0,900)\n",
    "\n",
    "        axs[i].vlines(419, 0,900, color='red')\n",
    "        axs[i].hlines(428, 0,900, color='red')\n",
    "\n",
    "    ax = axs[3]\n",
    "\n",
    "    ax.plot(x, cropOD.T[cen_y, :])\n",
    "    ax.plot(x, fit[cen_y, :], label='bimodal')\n",
    "    ax.plot(x, fa.thermal(x, bval['x0_th'], bval['amp_th'], bval['sigma_th']), label='thermal')\n",
    "    ax.legend(fontsize=8)\n",
    "\n",
    "    ax = axs[4]\n",
    "\n",
    "    ax.plot(y, cropOD.T[:, cen_x])\n",
    "    ax.plot(y, fit[:, cen_x], label='bimodal')\n",
    "    ax.plot(y, fa.thermal(y, bval['y0_th'], bval['amp_th'], bval['sigma_th']), label='thermal')\n",
    "    ax.legend(fontsize=8)\n",
    "\n",
    "\n",
    "t_fonts=20\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2.527851655733842e-109"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "fit[300, 400]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[  0,   1,   2, ..., 897, 898, 899],\n",
       "       [  0,   1,   2, ..., 897, 898, 899],\n",
       "       [  0,   1,   2, ..., 897, 898, 899],\n",
       "       ...,\n",
       "       [  0,   1,   2, ..., 897, 898, 899],\n",
       "       [  0,   1,   2, ..., 897, 898, 899],\n",
       "       [  0,   1,   2, ..., 897, 898, 899]], dtype=int64)"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib notebook\n",
    "shotNum = \"0024\"\n",
    "filePath = folderPath + \"/\" + shotNum + \"/*.h5\"\n",
    "\n",
    "dataSetDict = {\n",
    "    dskey[groupList[i]]: read_hdf5_file(filePath, groupList[i], excludeAxis = ['sweep_start_freq', 'sweep_stop_freq'])\n",
    "    for i in [0]\n",
    "}\n",
    "\n",
    "dataSet = dataSetDict[\"camera_0\"]\n",
    "\n",
    "print_scanAxis(dataSet)\n",
    "\n",
    "scanAxis = get_scanAxis(dataSet)\n",
    "\n",
    "dataSet = auto_rechunk(dataSet)\n",
    "\n",
    "dataSet = imageAnalyser.get_absorption_images(dataSet)\n",
    "\n",
    "imageAnalyser.center = (135, 990)\n",
    "imageAnalyser.span = (250, 250)\n",
    "imageAnalyser.fraction = (0.1, 0.1)\n",
    "\n",
    "dataSet_cropOD = imageAnalyser.crop_image(dataSet.OD)\n",
    "dataSet_cropOD = imageAnalyser.substract_offset(dataSet_cropOD).load()\n",
    "\n",
    "Ncount = imageAnalyser.get_Ncount(dataSet_cropOD)\n",
    "Ncount_mean = calculate_mean(Ncount)\n",
    "Ncount_std = calculate_std(Ncount)\n",
    "\n",
    "fig = plt.figure()\n",
    "ax = fig.gca()\n",
    "Ncount_mean.plot.errorbar(ax=ax, yerr = Ncount_std, fmt='ob')\n",
    "\n",
    "plt.ylabel('NCount')\n",
    "plt.tight_layout()\n",
    "#plt.ylim([0, 3500])\n",
    "plt.grid(visible=1)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "l = list(np.arange(0.15, 0.25, 0.005))\n",
    "# l = np.logspace(np.log10(250e-6), np.log10(500e-3), num=15)\n",
    "\n",
    "l = [round(item, 7) for item in l]\n",
    "random.shuffle(l)\n",
    "\n",
    "print(l)\n",
    "print(len(l))\n",
    "np.mean(l)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "[10.25, 10.255, 10.26, 10.265, 10.27, 10.275, 10.28, 10.285, 10.29, 10.295, 10.3, 10.305, 10.31, 10.315, 10.32, 10.325, 10.33, 10.335, 10.34, 10.345, 10.35, 10.355]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "pixel = 5.86e-6\n",
    "M     = 0.6827\n",
    "F = (1/(0.3725*8.4743e-14)) * (pixel / M)**2\n",
    "NCount = 85000\n",
    "AtomNumber = NCount * F / 1e8\n",
    "print(AtomNumber)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "muB = 9.274e-24\n",
    "hbar = 6.626e-34 / (2 * np.pi)\n",
    "gJ = 1.24\n",
    "Delta = 2 * np.pi * 100 * 1e3\n",
    "\n",
    "Bz = (Delta*hbar) / (muB*gJ)\n",
    "print(Bz * 1e4)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## ODT 1 Calibration"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "v_high = 2.7\n",
    "\"\"\"High Power\"\"\"\n",
    "P_arm1_high = 5.776 * v_high - 0.683\n",
    "\n",
    "v_mid = 0.2076\n",
    "\"\"\"Intermediate Power\"\"\"\n",
    "P_arm1_mid = 5.815 * v_mid - 0.03651\n",
    "\n",
    "v_low = 0.0587\n",
    "\"\"\"Low Power\"\"\"\n",
    "P_arm1_low = 5271 * v_low - 27.5\n",
    "\n",
    "print(round(P_arm1_high, 3))\n",
    "print(round(P_arm1_mid, 3))\n",
    "print(round(P_arm1_low, 3))"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## ODT 2 Power Calibration"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "v = 0.7607\n",
    "P_arm2 = 2.302 * v - 0.06452\n",
    "print(round(P_arm2, 3))"
   ]
  }
 ],
 "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.9.13"
  },
  "vscode": {
   "interpreter": {
    "hash": "c05913ad4f24fdc6b2418069394dc5835b1981849b107c9ba6df693aafd66650"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
