analyseScript/Joschka/Crop_data.ipynb

408 lines
1.0 MiB
Plaintext
Raw Normal View History

2023-07-27 17:16:08 +02:00
{
"cells": [
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"collapsed": true,
"ExecuteTime": {
"end_time": "2023-07-26T13:29:26.579448700Z",
"start_time": "2023-07-26T13:29:25.300714900Z"
}
},
"outputs": [],
"source": [
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"from scipy.ndimage import gaussian_filter"
]
},
{
"cell_type": "code",
"execution_count": 18,
"outputs": [],
"source": [
"data = np.zeros((2,11, 1200, 1920))\n",
"data[0] = np.load('Data_Britta/OD_ft_flatfield.npy')\n",
"data[1] = np.load('Data_Britta/OD_ft_manual.npy')\n",
"\n",
"shape = np.shape(data)"
],
"metadata": {
"collapsed": false,
"ExecuteTime": {
"end_time": "2023-07-26T14:04:42.011001700Z",
"start_time": "2023-07-26T14:04:41.128955800Z"
}
}
},
{
"cell_type": "code",
"execution_count": 24,
"outputs": [],
"source": [
"def calc_thresh(data):\n",
" shape = np.shape(data)\n",
" thresh = np.zeros(shape)\n",
" sigma = 0.4\n",
"\n",
" if len(shape) == 4:\n",
" blurred = gaussian_filter(data, sigma=sigma)\n",
" for i in range(0,shape[0]):\n",
" for j in range(0, shape[1]):\n",
" thresh[i,j] = np.where(blurred[i,j] < np.max(blurred[i,j])*0.5,0,1)\n",
"\n",
" elif len(shape) == 3:\n",
" blurred = gaussian_filter(data, sigma=sigma)\n",
" for i in range(0,shape[0]):\n",
" thresh[i] = np.where(blurred[i] < np.max(blurred[i])*0.5,0,1)\n",
"\n",
" else:\n",
" print(\"Shape of data is wrong, output is empty\")\n",
"\n",
" return thresh\n",
"\n",
"def calc_cen(thresh1):\n",
" \"\"\"\n",
" returns array: [Y_center,X_center]\n",
" \"\"\"\n",
" cen = np.zeros(2)\n",
" (Y,X) = np.shape(thresh1)\n",
"\n",
"\n",
" thresh1 = thresh1 /np.sum(thresh1)\n",
"\n",
" # marginal distributions\n",
" dx = np.sum(thresh1, 0)\n",
" dy = np.sum(thresh1, 1)\n",
"\n",
" # expected values\n",
" cen[0] = np.sum(dx * np.arange(X))\n",
" cen[1] = np.sum(dy * np.arange(Y))\n",
" return cen\n",
"\n",
"def calc_cen_bulk(thresh):\n",
" \"\"\"\n",
" returns array in shape of input, containing array with [Y_center,X_center] for each image\n",
" \"\"\"\n",
" shape = np.shape(thresh)\n",
" cen = np.zeros((shape[0], shape[1], 2))\n",
" for i in range(0, shape[0]):\n",
" for j in range(0, shape[1]):\n",
" cen[i,j] = calc_cen(thresh[i,j])\n",
" return cen\n",
"\n",
"def guess_BEC_width(thresh, center):\n",
" \"\"\"\n",
" returns width of thresholded area along both axis through the center with shape of thresh and [X_width, Y_width] for each image\n",
" \"\"\"\n",
" shape = np.shape(thresh)\n",
" BEC_width_guess = np.zeros((shape[0], shape[1], 2))\n",
"\n",
" for i in range(0, shape[0]):\n",
" for j in range(0, shape[1]):\n",
" BEC_width_guess[i, j, 0] = np.sum(thresh[i, j, round(center[i,j,1]), :])\n",
" BEC_width_guess[i, j, 1] = np.sum(thresh[i, j, :, round(center[i,j,0])])\n",
"\n",
" return BEC_width_guess\n",
"\n",
"\n",
"\n",
"def gaussian(x, x0, sigma, A):\n",
" return A * np.exp(-0.5 * (x-x0)**2 / sigma**2)\n",
"\n",
"# def polylog(power, numerator, order = 15):\n",
"#\n",
"# dataShape = numerator.shape\n",
"# numerator = np.tile(numerator, (order, 1))\n",
"# numerator = np.power(numerator.T, np.arange(1, order+1)).T\n",
"#\n",
"# denominator = np.arange(1, order+1)\n",
"# denominator = np.tile(denominator, (dataShape[0], 1))\n",
"# denominator = denominator.T\n",
"#\n",
"# data = numerator/ np.power(denominator, power)\n",
"#\n",
"# return np.sum(data, axis=0)\n",
"\n",
"def thermal(x, x0, amp, sigma):\n",
" res = np.exp(-0.5 * (x-x0)**2 / sigma**2)\n",
" return amp/1.643 * polylog(2, res)\n",
"\n",
"def Thomas_Fermi_1d(x, x0, amp, sigma):\n",
" res = (1-(( x - x0 ) / sigma) **2) **(3/2)\n",
" return amp * np.where(res > 0, res, 0)\n",
"\n",
"def density_1d(x, x0_bec, x0_th, amp_bec, amp_th, sigma_bec, sigma_th):\n",
" return thermal(x, x0_th, amp_th, sigma_th) + Thomas_Fermi_1d(x, x0_bec, amp_bec, sigma_bec)\n",
"\n",
"\n",
"def polylog(pow, x):\n",
" order = 15\n",
" sum = 0\n",
" for k in range(1,order):\n",
" sum += x ** k /k **pow\n",
" return sum\n",
"\n",
"\n",
"def ThomasFermi_2d(x, y=0.0, centerx=0.0, centery=0.0, amplitude=1.0, sigmax=1.0, sigmay=1.0):\n",
" res = (1- ((x-centerx)/(sigmax))**2 - ((y-centery)/(sigmay))**2)**(3 / 2)\n",
" return amplitude * np.where(res > 0, res, 0)\n",
" # return amplitude * 5 / 2 / np.pi / max(tiny, sigmax * sigmay) * np.where(res > 0, res, 0)\n",
"\n",
"\n",
" # return amplitude / 2 / np.pi / 1.20206 / max(tiny, sigmax * sigmay) * polylog(2, np.exp( -((x-centerx)**2/(2 * (sigmax)**2))-((y-centery)**2/( 2 * (sigmay)**2)) ))\n",
"\n",
"\n",
"def density_profile_BEC_2d(x, y=0.0, amp_bec=1.0, amp_th=1.0, x0_bec=0.0, y0_bec=0.0, x0_th=0.0, y0_th=0.0,\n",
" sigmax_bec=1.0, sigmay_bec=1.0, sigmax_th=1.0, sigmay_th=1.0):\n",
" return ThomasFermi_2d(x=x, y=y, centerx=x0_bec, centery=y0_bec,\n",
" amplitude=amp_bec, sigmax=sigmax_bec, sigmay=sigmay_bec\n",
" ) + polylog2_2d(x=x, y=y, centerx=x0_th, centery=y0_th,\n",
" amplitude=amp_th, sigmax=sigmax_th, sigmay=sigmay_th)"
],
"metadata": {
"collapsed": false,
"ExecuteTime": {
"end_time": "2023-07-27T09:17:14.540399600Z",
"start_time": "2023-07-27T09:17:14.391125700Z"
}
}
},
{
"cell_type": "code",
"execution_count": 17,
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"2304000.0\n",
"(2, 11, 1200, 1920)\n"
]
}
],
"source": [
"sigma = 0.4\n",
"blurred = gaussian_filter(data, sigma=sigma)\n",
"\n",
"thresh = np.zeros(shape)\n",
"for i in range(0,shape[0]):\n",
" for j in range(0, shape[1]):\n",
" thresh[i,j] = np.where(blurred[i,j] < np.max(blurred[i,j])*0.5,0,1)\n",
" print(np.sum(thresh[i,j]))\n",
"center = calc_cen_bulk(thresh)\n",
"print(np.shape(thresh))\n",
"BEC_width_guess = guess_BEC_width(thresh, center)"
],
"metadata": {
"collapsed": false,
"ExecuteTime": {
"end_time": "2023-07-26T14:04:12.103846400Z",
"start_time": "2023-07-26T14:04:09.168960300Z"
}
}
},
{
"cell_type": "code",
"execution_count": 32,
"outputs": [
{
"data": {
"text/plain": "<Figure size 2000x500 with 22 Axes>",
"image/png": "iVBORw0KGgoAAAANSUhEUgAABkYAAAGyCAYAAACx73h8AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAAEAAElEQVR4nOz9eXxVRbb3j78z52Q4CSSQAAkECUMQ6CBBQMALEiYFBEFFGQRFBQWnxgZncIQWFRUUFRVFFBQFBQU0SFpAQEAiUxiCBBJIgknIBCEQ2L8/PrWzYz/dz+3rvfB97q/3er1OTs4ealjrs1atqlpV5WNZloVLLrnkkksuueSSSy655JJLLrnkkksuueSSSy655NK/Afn+f10Al1xyySWXXHLJJZdccskll1xyySWXXHLJJZdccsmlS0XuxIhLLrnkkksuueSSSy655JJLLrnkkksuueSSSy659G9D7sSISy655JJLLrnkkksuueSSSy655JJLLrnkkksuufRvQ+7EiEsuueSSSy655JJLLrnkkksuueSSSy655JJLLrn0b0PuxIhLLrnkkksuueSSSy655JJLLrnkkksuueSSSy659G9D7sSISy655JJLLrnkkksuueSSSy655JJLLrnkkksuufRvQ+7EiEsuueSSSy655JJLLrnkkksuueSSSy655JJLLrn0b0PuxIhLLrnkkksuueSSSy655JJLLrnkkksuueSSSy659G9D7sSISy655JJLLrnkkksuueSSSy655JJLLrnkkksuufRvQ+7EiEsuueSSSy655JJLLrnkkksuueSSSy655JJLLrn0b0P/5YmRH374gYEDB9KwYUN8fHxYvnz57+6PGTMGHx+f33369ev3u2eKi4sZMWIEXq+XyMhI7rjjDioqKv5bFXHp/31ysePSHyUXOy79UXKx49IfJRc7Lv1RcrHj0h8lFzsu/VFysePSHyUXOy79UXKx49IfJRc7Lv2/RP/liZFTp07xpz/9iblz5/7TZ/r160deXl7N55NPPvnd/REjRrBnzx6+++47Vq5cyQ8//MBdd931Xy+9S/+ryMWOS3+UXOy49EfJxY5Lf5Rc7Lj0R8nFjkt/lFzsuPRHycWOS3+UXOy49EfJxY5Lf5Rc7Lj0/xRZ/w0CrGXLlv3u2m233WZdf/31//SdvXv3WoC1devWmmurVq2yfHx8rGPHjv13iuPS/yJysePSHyUXOy79UXKx49IfJRc7Lv1RcrHj0h8lFzsu/VFysePSHyUXOy79UXKx49IfJRc7Lv1/Tf4XY7IlPT2d+vXrU6dOHa655hqeffZZoqKiANi0aRORkZGkpKTUPJ+amoqvry9btmxhyJAh/0d6VVVVVFVV1fy+cOECxcXFREVF4ePjczGq4NIloNOnT1NWVlbz+9y5c6xbt47o6GiioqJc7Lj0T8nFjkt/lFzsuPRH6R9hJz09nXr16uH1eunduzfPPfecix2X/g9ysePSH6VLgR1w8fP/j+Rix6U/Si52XPqj5GLHpT9KLnZc+p8my7IoLy+nYcOG+Pr+J5tl/XdmVfgHM3uffPKJ9eWXX1o7d+60li1bZiUlJVkdO3a0qqurLcuyrOeee85q0aLF/5FWvXr1rDfeeOMf5vPUU09ZgPv5N/q88847Lnbcj4sd9+Nix/38r/k0b97cxY77cbHjfv7XYcfFz7/nx8WO+3Gx435c7Lif/y0fFzvu549+cnJy/ikebPofXzEyfPjwmv/btm1Lu3btaNasGenp6fTq1esPpfnII4/w0EMP1fwuLS2lcePG3J8zkVenPcd1Ly+liLpszuvK0gY3MCzvC95uMI5fuYyRfMRaetGK/fQ78Q1f1L+Bv9GDnfyJPbTmBr4gglJm5D1K6Yr6nLrNh3ZBu/BQSQrbiKCUKIppy04q8fAjV3GYplzPl5QQQQhnuLd4Po3DDnBL4GLe+vR+Ft40lFHbPsen8SmerP8MHzGS+3iNb+lDCttoyy6Kqcsb3Ms9zCWRg+ynFa8xiV6s5blDzxFdJ4fudTdwI58CcDl72EInVjKQH05046f6nRjNB4Rxije4h+d5lGU5N/Fq/D3M5V56ks639CGvLJZvvP35jJu4gB++nOd23uMZnqAXaxlTuISnoqcyr2w8ALO8D1NJMPMYz6G1bfi+V2euydjMqOS36cMajpIAwGMZL/FJ8hDuuTCXN3zvZR8t2Uk7zuNHBKUEcZYr2cJ2UnjrL/fz578+x2vFkxhRdxGd2MKEiGXw+CKefvgH1tKLv+X1YFCD5bQo28as+Le58cYbueaaay4adpJzlpKxozc3/sdCfuUytt/TFWLMQ2eA1kAhUAVsBXoDBcD1gD8ENCzj3CQvLy6ayOfcwOZvr4G3gAFQ57bjnFzUkODrTxIYfJayG2MIeL+Mcw97afbubg5ltIFdQBdgNXjvKqAsNoY6hcdp57uLv73TF9oBk6Hf+mWs3j0ECqBhr0NEUELmNx2gCJX3XRi0ZAmrivsTGHyOUzui4FdgE3R+43s2R1xD09JMUtjGZ/eMgodgUOISttCJgm8TIA3oYertD537fM8x4siZ1YKkydvJnNcBrqvW/XR/OAp0BHaa/Ash4LYyzlUF4Y0qoWxRDDRW2Tvu+IGtj1+Nz32niKt/jNwTjbBmhdLsr7sBOPR4GxgJd7d6lbcG3M8DK19gds7DcMAf1gD54icfwqsr7+L+998Gf2BiBLz8EQ3vaE0UxezKawc7guB8GYyMv+jY4ZUcUm9PJ23tICiBOkOOc/LJhnqoP/AU0A6Cp53kTHwxbGkmfjUCroUbSxfy2cZRsAuix+cQQSmH7mkD3aHtLVv59fRlnNoUBaOBFGAExNyUTcEjCdACksZuJ/P9DsLjzRZc7wP3ANXCbvzkA+TMagERkmnDsYc4vrsZHAduBL4C7oTgbSc580odADo/9T2lRHAefw688yfoDqRD8K0n6e9dxdqzvSgriIJp/hACo15/m4X33CU9aQGUAy0RrntUwdEg8aE7BEwoo7raj671f2TDK73Fh1+BF4qhR128nxVQtj4GdsKgB5dwGb8ye8AjMNTUKRO4HYiqgkVB8CHQD9WvArinGv7qD1cClwGngPrAZiAUCIL4mw6QM6cFPBYB7RbBDQMgDHjPPHtDGdx3CbDzbg54vbAO6fBlyMZ0BBYh+wJwDOgJ7EA8Xgd8DtwDwUNPcuZUCKwMghbg7VJA2YQYeNyiQ9Mf2T6vK3QDvkQ61M6kGWbyuQV4E9gGvIHkdDswG2hrrqdAaOsiYkPyOPRJG+V/GXSe+j2bp18DDYHJJ2BefYgyaTxo3m0Lwdef5Mx3whafQsMlhzj+TjN4FrgOGAGsAu6uhl3+NLt2N805yI9nr6JsZ4zS6mbk293w47ThVwRwVNg883UdQq8v4tTnURCEMDvjJGe61YH/AEoRLjchO+dneHoe4bAbwksRsBfZ5OWGZ43NvVLDxw8j4OZF0HcA0UNzKJwUr+faAT+XwdR4PvvsM5KTky8Odq7PgZZe6cRqVN/eSBa7gGVAcTGwB+7oDu8i+9ka+Kup7yzzf5jhw83AISDYXNtkMk1C9qIhsvNfAw+bez8BV6H0TxkeDkAYGIxkFQasNLxrDBwohpvrQhMgBHgfaA5cYfJ8GHgB6AtMA54D3gYeRzp6BhiOcNgEmGDB1z5K/xjSocYIH4sQRhrh6JR/rf9bG341NnyMMmkeM2V5yOTTxdTnRyAA6cqdCCtxwFLghEmno8l7IrKHp0x605Ct+zAC+iyCygGyXaVID54Bostg8UXGTpccGOqFyWXwiFe+Snfgb0DJNngiRXWJMOWuQPIH6Iza4+5GDhUId0tMPSoQjlYDRxA2SoFeQDrStSKk84uQntn/j0CyaGryjgA+RbiLMO83ASYBz5u0q5FctiIfrbXJz8+k1R/oa0G6DzyA2rvXTBlOmXevQ/bxZ1PH3sB3COsDgA3m+cctmOCjdiQU4XoVcDfSiZ3AX0ydbEzVBYbh6EM35GcBHDDlT0f4wuT3IjAK8GBsDbJXxcBO02Y9PEBlPI4w1QI4WAaH/me
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": "<Figure size 2000x500 with 22 Axes>",
"image/png": "iVBORw0KGgoAAAANSUhEUgAABkYAAAGyCAYAAACx73h8AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAABwYElEQVR4nO3df3RUhZ3//1ciJfgrSfkZcopCqVa3ClJ/pHzWuipZkfb4o6W2RVbxR6Xuqnsq3VPhnK2EXXvI0f3QfmypbHe3uvupWOvno1DZih8EEW0RFeTg0UrBrxEUEreyJPwoEeR+/wgZk8m9c++8596Z++P54IyQmTuTO5mnM3fuOzNT5TiOIwAAAAAAAAAAgAyorvQKAAAAAAAAAAAAlAuDEQAAAAAAAAAAkBkMRgAAAAAAAAAAQGYwGAEAAAAAAAAAAJnBYAQAAAAAAAAAAGQGgxEAAAAAAAAAAJAZDEYAAAAAAAAAAEBmMBgBAAAAAAAAAACZwWAEAAAAAAAAAABkBoMRAAAAAAAAAACQGUUPRtatW6crrrhCjY2Nqqqq0rJly/qdfsMNN6iqqqrf4fLLL++3zJ49ezRz5kzV1taqvr5eN998s/bv31/SFUH80Q6saAdWtAMr2oEV7cCKdmBFO7CiHVjRDqxoB3FS9GDkwIEDmjhxohYvXuy5zOWXX67du3fnDo888ki/02fOnKnXX39dq1at0ooVK7Ru3TrNnj27+LVHotAOrGgHVrQDK9qBFe3AinZgRTuwoh1Y0Q6saAex4pRAkvPEE0/0O27WrFnOVVdd5XmeN954w5HkvPzyy7njnnrqKaeqqsp57733SlkdJAjtwIp2YEU7sKIdWNEOrGgHVrQDK9qBFe3AinZQaYOiGLasXbtWI0eO1Cc/+UldeumluueeezRs2DBJ0vr161VfX6/zzjsvt3xzc7Oqq6u1YcMGfeUrXxlwed3d3eru7s59ffToUe3Zs0fDhg1TVVVVFFcBZXDw4EF1dXXlvj58+LCeffZZDR8+XMOGDaMdeKIdWNEOrNzaWbt2rUaMGKHa2lr95V/+pX7wgx/QDgagHViVox2JftKIdmBFO7CiHVjRDsLmOI727dunxsZGVVf7vFlWKVMVuUz2HnnkEWf58uXOli1bnCeeeMI588wznfPPP985cuSI4ziO84Mf/MA5/fTTB1zWiBEjnJ/+9Keu32f+/PmOJA4ZOvzLv/wL7XCgHQ60wyExh9NOO412ONAOh8S1Qz/ZPNAOB9rhQDscknKgHQ7Ww86dOz176BX6K0a++c1v5v599tlna8KECRo/frzWrl2rKVOmmC5z3rx5mjNnTu7rzs5OnXLKKZLulFRT4hp/bK5a1aq5oV0eCmmV9FVJp+cd3y3ph7rmmmt06aWXJqYdlBPtwIp2YOXVjtTbz2OPPaZzzjmHdpCHdmBVnnYk+kkf2oEV7cCKdmBFO4hCTzsnn3yy75KRvJVWX5/+9Kc1fPhwbd++XVOmTFFDQ4Pef//9fsscOXJEe/bsUUNDg+tl1NTUqKbGLdIaSUNCW9chff6LcviEvH7eVVVViWoH5UY7sKIdWHm3I0njxo2jHXigHVhF345EP+lEO7CiHVjRDqxoB9EI8hZpPm+0Vbp3331XH3zwgUaPHi1Jmjx5svbu3auNGzfmllmzZo2OHj2qpqamkr5Xi1oqev6oLy9rytkO0oV2YEU7sHrvvfdoBya0AyvagRXtwIp2YEU7sKIdRKnKcRynmDPs379f27dvlyRNmjRJixYt0iWXXKKhQ4dq6NChWrBggaZPn66Ghga99dZb+t73vqd9+/bptddey03mpk2bpo6ODi1ZskSHDx/WjTfeqPPOO09Lly4NtA5dXV2qq6uTNFdRTPVajv1B2Lol7Tn273+WNFXSWEnHHzs8J2m8pP+t5cuXa8GCBYlrB1GhHVjRDqyCtHOmen7DaYkmTpyogwcP0g5EO7CrfDsS/SQT7cCKdmBFO7CiHUTtkKRWdXZ2qra2tuCSRb9i5JVXXtGkSZM0adIkSdKcOXM0adIk3X333TruuOO0ZcsWXXnllTr99NN1880369xzz9Xzzz/f7+VKDz/8sM444wxNmTJFX/rSl3ThhRfqZz/7WbGrEhmGIlHZpZ47vX8+9vXTx/79rHpS7JD0fyVJd9xxRyLbQVRoB1a0A6sg7TwiqaeFc845h3ZwDO3AinZgRTuwoh1Y0Q6saAfxUfQrRuIgyqleOV4tEuX3SPr6FzPVs2AinGa0AyvaQSmi64d20o52YMXjFqxoB1a0AyvagRXtwCrCV4ykXTleLeL2PcL6vpVa/zheJgAAAAAAAAAA+TI3GInrDvi+6xXXdYxSFq8zAAAAAAAAAKD8MjcYAQAAAAAAAAAA2ZW5wUgSXpmQhHUEAAAAAAAAACCJMjcYSSuGKQAAAAAAAAAA+GMwkhJZ/4wSAAAAAAAAAACCYDCSQgxGAAAAAAAAAABwl5rBSBjDgLgMFOKyHgAAAAAAAAAApA2DkZAvAwAAAAAAAAAAxFdqBiMAAAAAAAAAAAB+GIzEUNSvXOGVMQAAAAAAAACArGIwkjEtx/4AAAAAAAAAAJBFmRmMlGsYEPehQznWL+4/AwAAAAAAAABAdjEYqdD3yV8uTcOENF0XAAAAAAAAAEC6ZGYwEjdpHowAAAAAAAAAABBXDEYAAAAAAAAAAEBmMBgBAAAAAAAAAACZwWCkRLwFFgAAAAAAAAAAyVH0YGTdunW64oor1NjYqKqqKi1btix32uHDh3XXXXfp7LPP1oknnqjGxkZdf/312rVrV7/LGDt2rKqqqvodWltbTVeg0oOJpH3/yq5vm6Slkv5JUouk3/c57SNJqyT9myTps5/9bOTtIEnaRDuwaRPtwKZN/u38VNL/lCR9+9vfph0c0ybagU2baAc2baId2LSJdmDTJtqBTZtoB3FR9GDkwIEDmjhxohYvXjzgtIMHD2rTpk36/ve/r02bNunxxx/X1q1bdeWVVw5Y9h/+4R+0e/fu3OGOO+4wXYFKDyYqLVmDkcOSRkn6ssdpuyX9D0nSL37xi8jbQZLQDqxoB1ZB2rlI0g2SpG3bttEOjqEdWNEOrGgHVrQDK9qBFe0gPgYVe4Zp06Zp2rRprqfV1dVp1apV/Y77yU9+ogsuuEA7duzQKaeckjv+5JNPVkNDQ7HfPjFajv1J+vcI12nHDm6GSLpe0iFJy3X++ednth24oR1Y0Q6sgrQj9fQj3Xfffbr00ktpB6Id2NEOrGgHVrQDK9qBFe0gPiL/jJHOzk5VVVWpvr6+3/Gtra0aNmyYJk2apPvuu09HjhzxvIzu7m51dXX1O+SL24AgivXJv8y4XeewlasdpA/twIp2YNXV1UU7MKEdWIXRjkQ/WUQ7sKIdWNEOrGgHUSr6FSPFOHTokO666y7NmDFDtbW1ueP/9m//Vp///Oc1dOhQ/e53v9O8efO0e/duLVq0yPVyFi5cqAULFhT8XmkfEkju1zHqV41U6lUp5WwH6UI7sKIdlGL+/Pm0AxPagVUY7Uj0k0W0AyvagRXtwIp2EKUqx3Ec85mrqvTEE0/o6quvHnDa4cOHNX36dL377rtau3Ztv4Dz/fznP9e3v/1t7d+/XzU1NQNO7+7uVnd3d+7rrq4ujRkzRtJc9bzMKruSOxhpkfQNSWfmHX9IUqumTZum9vZ22oGLFtEObFpEO7BpkXs7knRA0n2aMGGCnn/+edpBnhbRDmxaVI52JPpJnxbRDmxaRDuwaRHtwKZFtIPw9ezj6ezsLNiNFNFbaR0+fFhf//rX9c4772jVqlW+K9HU1KQjR46ora3N9fSamhrV1tb2O2SR1ytGilk+jO8ZrY8kSTt37qQdFIl2YEU7sPpI0jJJ0vLly2kHRaAdWIXbjkQ/2UE7sKIdWNEOrGgH5RH6YKR3KLJt2zY988wzGjZsmO95Nm/erOrqao0cOTLs1UmVYocUQZYPeplRfO+B+t/x0Q6Cox1Y0Q6sPpL0mKT/liQNHTrU9xy0gx60AyvagRXtwIp2YEU7sKIdlE/RnzGyf/9+bd++Pff122+
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"fig, ax = plt.subplots(shape[0],shape[1], figsize=(20,5))\n",
"\n",
"for i in range(0,shape[0]):\n",
" for j in range(0,shape[1]):\n",
" ax[i,j].pcolormesh(cropOD[i,j], cmap='jet')\n",
" ax[i,j].plot(center[i,j,0], center[i,j,1], marker='x', color='b', markersize=12)\n",
"plt.show()\n",
"\n",
"fig, ax = plt.subplots(shape[0],shape[1], figsize=(20,5))\n",
"\n",
"for i in range(0,shape[0]):\n",
" for j in range(0,shape[1]):\n",
" ax[i,j].pcolormesh(thresh[i,j], cmap='jet')\n",
"plt.show()\n"
],
"metadata": {
"collapsed": false,
"ExecuteTime": {
"end_time": "2023-07-27T09:21:25.397250500Z",
"start_time": "2023-07-27T09:21:20.771713200Z"
}
}
},
{
"cell_type": "code",
"execution_count": 30,
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"[[[75.43843844 74.33233233]\n",
" [74.49700599 75.66467066]\n",
" [74.20763723 75.66109785]\n",
" [74.70685579 75.38534279]\n",
" [74.58390805 75.63678161]\n",
" [74.44421488 75.28099174]\n",
" [74.56997972 75.31845842]\n",
" [74.41902834 75.50404858]\n",
" [74.36344538 75.48739496]\n",
" [74.44730679 75.54332553]\n",
" [74.6125 75.37321429]]\n",
"\n",
" [[70.24350649 81.52272727]\n",
" [75.33898305 77.08474576]\n",
" [74.57142857 75.00571429]\n",
" [75.12068966 75.64655172]\n",
" [74.47509579 75.85823755]\n",
" [74.67948718 75.84615385]\n",
" [74.9796748 75.65447154]\n",
" [75.30769231 75.05494505]\n",
" [74.69072165 75.39175258]\n",
" [75.11167513 75.88324873]\n",
" [74.90547264 75.09452736]]]\n"
]
}
],
"source": [
"thresh = calc_thresh(data)\n",
"center = calc_cen_bulk(thresh)\n",
"\n",
"shape = np.shape(data)\n",
"cropOD = np.zeros((shape[0], shape[1], 150, 150))\n",
"\n",
"for i in range(0,shape[0]):\n",
" for j in range(0, shape[1]):\n",
" cropOD[i,j] = data[i,j, round(center[i,j,1]-75):round(center[i,j,1]+75), round(center[i,j,0]-75):round(center[i,j,0]+75)]\n",
"\n",
"thresh = calc_thresh(cropOD)\n",
"center = calc_cen_bulk(thresh)\n",
"print(center)"
],
"metadata": {
"collapsed": false,
"ExecuteTime": {
"end_time": "2023-07-27T09:20:10.349363100Z",
"start_time": "2023-07-27T09:20:07.775299200Z"
}
}
},
{
"cell_type": "code",
"execution_count": 33,
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"[[[17. 18.]\n",
" [19. 18.]\n",
" [24. 20.]\n",
" [26. 23.]\n",
" [26. 19.]\n",
" [28. 23.]\n",
" [26. 22.]\n",
" [26. 25.]\n",
" [25. 22.]\n",
" [25. 20.]\n",
" [30. 25.]]\n",
"\n",
" [[ 4. 3.]\n",
" [ 7. 8.]\n",
" [ 7. 14.]\n",
" [10. 4.]\n",
" [15. 13.]\n",
" [22. 17.]\n",
" [15. 15.]\n",
" [ 7. 4.]\n",
" [13. 12.]\n",
" [19. 15.]\n",
" [14. 12.]]]\n"
]
}
],
"source": [
"BEC_width_guess = guess_BEC_width(thresh, center)\n",
"\n",
"print(BEC_width_guess)"
],
"metadata": {
"collapsed": false,
"ExecuteTime": {
"end_time": "2023-07-27T09:26:35.539786700Z",
"start_time": "2023-07-27T09:26:35.510570900Z"
}
}
},
{
"cell_type": "code",
"execution_count": null,
"outputs": [],
"source": [],
"metadata": {
"collapsed": false
}
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 2
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython2",
"version": "2.7.6"
}
},
"nbformat": 4,
"nbformat_minor": 0
}