dataset.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "ename": "ModuleNotFoundError",
     "evalue": "No module named '_gdal'",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mImportError\u001b[0m                               Traceback (most recent call last)",
      "File \u001b[0;32m~/anaconda3/envs/SARLens/lib/python3.10/site-packages/osgeo/__init__.py:30\u001b[0m, in \u001b[0;36mswig_import_helper\u001b[0;34m()\u001b[0m\n\u001b[1;32m     29\u001b[0m \u001b[39mtry\u001b[39;00m:\n\u001b[0;32m---> 30\u001b[0m     \u001b[39mreturn\u001b[39;00m importlib\u001b[39m.\u001b[39;49mimport_module(mname)\n\u001b[1;32m     31\u001b[0m \u001b[39mexcept\u001b[39;00m \u001b[39mImportError\u001b[39;00m:\n",
      "File \u001b[0;32m~/anaconda3/envs/SARLens/lib/python3.10/importlib/__init__.py:126\u001b[0m, in \u001b[0;36mimport_module\u001b[0;34m(name, package)\u001b[0m\n\u001b[1;32m    125\u001b[0m         level \u001b[39m+\u001b[39m\u001b[39m=\u001b[39m \u001b[39m1\u001b[39m\n\u001b[0;32m--> 126\u001b[0m \u001b[39mreturn\u001b[39;00m _bootstrap\u001b[39m.\u001b[39;49m_gcd_import(name[level:], package, level)\n",
      "File \u001b[0;32m<frozen importlib._bootstrap>:1050\u001b[0m, in \u001b[0;36m_gcd_import\u001b[0;34m(name, package, level)\u001b[0m\n",
      "File \u001b[0;32m<frozen importlib._bootstrap>:1027\u001b[0m, in \u001b[0;36m_find_and_load\u001b[0;34m(name, import_)\u001b[0m\n",
      "File \u001b[0;32m<frozen importlib._bootstrap>:1006\u001b[0m, in \u001b[0;36m_find_and_load_unlocked\u001b[0;34m(name, import_)\u001b[0m\n",
      "File \u001b[0;32m<frozen importlib._bootstrap>:674\u001b[0m, in \u001b[0;36m_load_unlocked\u001b[0;34m(spec)\u001b[0m\n",
      "File \u001b[0;32m<frozen importlib._bootstrap>:571\u001b[0m, in \u001b[0;36mmodule_from_spec\u001b[0;34m(spec)\u001b[0m\n",
      "File \u001b[0;32m<frozen importlib._bootstrap_external>:1176\u001b[0m, in \u001b[0;36mcreate_module\u001b[0;34m(self, spec)\u001b[0m\n",
      "File \u001b[0;32m<frozen importlib._bootstrap>:241\u001b[0m, in \u001b[0;36m_call_with_frames_removed\u001b[0;34m(f, *args, **kwds)\u001b[0m\n",
      "\u001b[0;31mImportError\u001b[0m: libpoppler.so.126: cannot open shared object file: No such file or directory",
      "\nDuring handling of the above exception, another exception occurred:\n",
      "\u001b[0;31mModuleNotFoundError\u001b[0m                       Traceback (most recent call last)",
      "\u001b[1;32m/home/roberto/PythonProjects/SSFocus/dataset.ipynb Cella 1\u001b[0m line \u001b[0;36m1\n\u001b[0;32m----> <a href='vscode-notebook-cell://ssh-remote%2Besago/home/roberto/PythonProjects/SSFocus/dataset.ipynb#W0sdnNjb2RlLXJlbW90ZQ%3D%3D?line=0'>1</a>\u001b[0m \u001b[39mfrom\u001b[39;00m \u001b[39mosgeo\u001b[39;00m \u001b[39mimport\u001b[39;00m gdal\n\u001b[1;32m      <a href='vscode-notebook-cell://ssh-remote%2Besago/home/roberto/PythonProjects/SSFocus/dataset.ipynb#W0sdnNjb2RlLXJlbW90ZQ%3D%3D?line=1'>2</a>\u001b[0m \u001b[39mimport\u001b[39;00m \u001b[39mos\u001b[39;00m\n\u001b[1;32m      <a href='vscode-notebook-cell://ssh-remote%2Besago/home/roberto/PythonProjects/SSFocus/dataset.ipynb#W0sdnNjb2RlLXJlbW90ZQ%3D%3D?line=2'>3</a>\u001b[0m \u001b[39mimport\u001b[39;00m \u001b[39mnumpy\u001b[39;00m \u001b[39mas\u001b[39;00m \u001b[39mnp\u001b[39;00m\n",
      "File \u001b[0;32m~/anaconda3/envs/SARLens/lib/python3.10/site-packages/osgeo/__init__.py:46\u001b[0m\n\u001b[1;32m     42\u001b[0m                 \u001b[39mraise\u001b[39;00m \u001b[39mImportError\u001b[39;00m(traceback_string \u001b[39m+\u001b[39m \u001b[39m'\u001b[39m\u001b[39m\\n\u001b[39;00m\u001b[39m'\u001b[39m \u001b[39m+\u001b[39m msg)\n\u001b[1;32m     43\u001b[0m         \u001b[39mreturn\u001b[39;00m importlib\u001b[39m.\u001b[39mimport_module(\u001b[39m'\u001b[39m\u001b[39m_gdal\u001b[39m\u001b[39m'\u001b[39m)\n\u001b[0;32m---> 46\u001b[0m _gdal \u001b[39m=\u001b[39m swig_import_helper()\n\u001b[1;32m     47\u001b[0m \u001b[39mdel\u001b[39;00m swig_import_helper\n\u001b[1;32m     49\u001b[0m __version__ \u001b[39m=\u001b[39m _gdal\u001b[39m.\u001b[39m__version__ \u001b[39m=\u001b[39m _gdal\u001b[39m.\u001b[39mVersionInfo(\u001b[39m\"\u001b[39m\u001b[39mRELEASE_NAME\u001b[39m\u001b[39m\"\u001b[39m)\n",
      "File \u001b[0;32m~/anaconda3/envs/SARLens/lib/python3.10/site-packages/osgeo/__init__.py:43\u001b[0m, in \u001b[0;36mswig_import_helper\u001b[0;34m()\u001b[0m\n\u001b[1;32m     41\u001b[0m         traceback_string \u001b[39m=\u001b[39m \u001b[39m'\u001b[39m\u001b[39m'\u001b[39m\u001b[39m.\u001b[39mjoin(traceback\u001b[39m.\u001b[39mformat_exception(\u001b[39m*\u001b[39msys\u001b[39m.\u001b[39mexc_info()))\n\u001b[1;32m     42\u001b[0m         \u001b[39mraise\u001b[39;00m \u001b[39mImportError\u001b[39;00m(traceback_string \u001b[39m+\u001b[39m \u001b[39m'\u001b[39m\u001b[39m\\n\u001b[39;00m\u001b[39m'\u001b[39m \u001b[39m+\u001b[39m msg)\n\u001b[0;32m---> 43\u001b[0m \u001b[39mreturn\u001b[39;00m importlib\u001b[39m.\u001b[39;49mimport_module(\u001b[39m'\u001b[39;49m\u001b[39m_gdal\u001b[39;49m\u001b[39m'\u001b[39;49m)\n",
      "File \u001b[0;32m~/anaconda3/envs/SARLens/lib/python3.10/importlib/__init__.py:126\u001b[0m, in \u001b[0;36mimport_module\u001b[0;34m(name, package)\u001b[0m\n\u001b[1;32m    124\u001b[0m             \u001b[39mbreak\u001b[39;00m\n\u001b[1;32m    125\u001b[0m         level \u001b[39m+\u001b[39m\u001b[39m=\u001b[39m \u001b[39m1\u001b[39m\n\u001b[0;32m--> 126\u001b[0m \u001b[39mreturn\u001b[39;00m _bootstrap\u001b[39m.\u001b[39;49m_gcd_import(name[level:], package, level)\n",
      "\u001b[0;31mModuleNotFoundError\u001b[0m: No module named '_gdal'"
     ]
    }
   ],
   "source": [
    "from osgeo import gdal\n",
    "import os\n",
    "import numpy as np\n",
    "import cv2\n",
    "import re\n",
    "import matplotlib.pyplot as plt\n",
    "import sentinel1decoder"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import osgeo"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def start_points(size, split_size, overlap=0): # overlap = % of overlap between tiles, both rows and columns \n",
    "    points = [0]\n",
    "    stride = int(split_size * (1-overlap))\n",
    "    counter = 1\n",
    "    while True:\n",
    "        pt = stride * counter\n",
    "        if pt + split_size >= size:\n",
    "            points.append(size - split_size)\n",
    "            break\n",
    "        else:\n",
    "            points.append(pt)\n",
    "        counter += 1\n",
    "    return points\n",
    "\n",
    "\n",
    "def create_patches_labels(L1, folder, i):     \n",
    "    patch_dir = os.path.join(folder, 'L1_tile')\n",
    "    if not os.path.isdir(patch_dir):\n",
    "        os.mkdir(patch_dir)\n",
    "        \n",
    "    img_h, img_w = L1.shape\n",
    "    split_width = 3034 # define tile size \n",
    "    split_height = 916\n",
    "\n",
    "    X_points = start_points(img_w, split_width, 0)\n",
    "    Y_points = start_points(img_h, split_height, 0)\n",
    "    name = 'splitted'\n",
    "    \n",
    "    tot = []\n",
    "    for s in Y_points:\n",
    "        for j in X_points:\n",
    "            split = L1[s:s+split_height, j:j+split_width]\n",
    "            patch_name = f'patch_{name}_{i}_{s:02}_{j:02}.npy'\n",
    "            split = split.astype(np.float32)\n",
    "            np.save(os.path.join(patch_dir, patch_name), split)\n",
    "            tot.append(split)        \n",
    "    return tot\n",
    "\n",
    "\n",
    "\n",
    "def create_patches(re, im, folder, i):      \n",
    "    patch_dir = os.path.join(folder, 'L0_tile')\n",
    "    if not os.path.isdir(patch_dir):\n",
    "        os.mkdir(patch_dir)\n",
    "        \n",
    "    img_h, img_w = re.shape #or im\n",
    "    split_width = 3034\n",
    "    split_height = 916\n",
    "\n",
    "    X_points = start_points(img_w, split_width, 0)\n",
    "    Y_points = start_points(img_h, split_height, 0)\n",
    "    name = 'splitted'\n",
    "    \n",
    "    tot = []\n",
    "    for s in Y_points:\n",
    "        for j in X_points:\n",
    "            split_re = re[s:s+split_height, j:j+split_width]\n",
    "\n",
    "            split_im = im[s:s+split_height, j:j+split_width]\n",
    "            \n",
    "            split = cv2.merge([split_re, split_im])\n",
    "            split = split.astype(np.float32)\n",
    "            \n",
    "            patch_name = f'patch_{name}_{i}_{s:02}_{j:02}.npy'\n",
    "            np.save(os.path.join(patch_dir, patch_name), split)\n",
    "            tot.append(split)\n",
    "        \n",
    "    return tot"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "path = \"L0\"\n",
    "L0_files =[]\n",
    "for root, dirs, files in os.walk(\"L0/\"):\n",
    "    for file in files:\n",
    "        strin1 = 'annot'\n",
    "        string2 = 'index'\n",
    "        if file.endswith(\".dat\") and strin1 not in file and string2 not in file:\n",
    "            L0_files.append(root + '/' + file)\n",
    "print(L0_files)\n",
    "\n",
    "path = \"L1\"\n",
    "L1_files = []\n",
    "for root, dirs, files in os.walk(path):\n",
    "    for file in files:\n",
    "        check = str(os.path.join(root, file))\n",
    "        if(file.endswith(\"001.tiff\")):\n",
    "            L1_files.append(os.path.join(root,file))\n",
    "print(L1_files)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "path = \"L1\"\n",
    "i = 0\n",
    "L1_size = []\n",
    "for root, dirs, files in os.walk(path):\n",
    "    for file in files:\n",
    "        check = str(os.path.join(root, file))\n",
    "        if(file.endswith(\"001.tiff\")):\n",
    "\n",
    "            print(\"Current file: \", os.path.join(root, file))\n",
    "            #Load L1 to obtain same dimension for L0 and L1\n",
    "            raster = gdal.Open(os.path.join(root,file))\n",
    "\n",
    "            rasterArray = raster.ReadAsArray()\n",
    "            print('SLC data type: ', rasterArray.dtype)\n",
    "            rasterArray = rasterArray.reshape(len(rasterArray),len(rasterArray[0]))\n",
    "            print('SLC shape: ', rasterArray.shape)\n",
    "\n",
    "            image_SLC = abs(rasterArray)\n",
    "            print('SLC abs data type: ',image_SLC.dtype)\n",
    "\n",
    "            del rasterArray\n",
    "            #Delate L1 raws and columns to have azimuth and range proportional to 916 and 3034\n",
    "            Range = 916\n",
    "            Azimuth = 3034\n",
    "            num_row = len(image_SLC)//Range\n",
    "            num_col = len(image_SLC[0])//Azimuth\n",
    "            print(\"New number of row: \", num_row*Range)\n",
    "            print(\"New number of columns: \", num_col*Azimuth)\n",
    "            SLC_r = len(image_SLC)\n",
    "            SLC_c = len(image_SLC[0])\n",
    "            r = SLC_r - num_row*Range #to be subtracted\n",
    "            c = SLC_c - num_col*Azimuth\n",
    "            L1 = image_SLC[0:SLC_r-r, 0:SLC_c-c]\n",
    "            L1_size.append(L1.shape)\n",
    "            print(\"New L1 dimension: \", L1.shape)\n",
    "\n",
    "            del image_SLC\n",
    "            #plt.figure(figsize = (20,20))\n",
    "            #plt.imshow(L1[0:916,0:3034], cmap = 'gray', vmin = 0, vmax = 255)\n",
    "            #plt.colorbar()\n",
    "\n",
    "            folder = ''\n",
    "            L1_patches = create_patches_labels(L1, folder, i)\n",
    "            print('Tot number of patches for the current image: ', len(L1_patches))\n",
    "            i += 1\n",
    "            print(\"Dim of a single patch: \", L1_patches[len(L1_patches)-1].shape)\n",
    "            del L1_patches, L1\n",
    "\n",
    "print(L1_size)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "path = 'L0'\n",
    "i = 0 # per la size di L1 corrispondente\n",
    "j = 0\n",
    "path = \"L0\"\n",
    "for file in L0_files[0:len(L0_files):2]:\n",
    "\n",
    "    print(\"Current file: \", file)\n",
    "    l0file = sentinel1decoder.Level0File(file)\n",
    "    # print(l0file.packet_metadata)\n",
    "    # print(l0file.ephemeris)\n",
    "\n",
    "    selected_burst = 9\n",
    "    selection = l0file.get_burst_metadata(selected_burst)\n",
    "    print(selection)\n",
    "    # Decode the IQ data\n",
    "    radar_data = l0file.get_burst_data(selected_burst)\n",
    "\n",
    "    plt.figure(figsize=(12, 12))\n",
    "    plt.title(\"Sentinel-1 Raw I/Q Sensor Output\")\n",
    "    plt.imshow(abs(radar_data[:200,:200]), vmin=0, vmax=15, origin='lower')\n",
    "    plt.xlabel(\"Fast Time (down range)\")\n",
    "    plt.ylabel(\"Slow Time (cross range)\")\n",
    "    plt.show()\n",
    "\n",
    "    raw_L1, col_L1 = L1_size[i]\n",
    "    raw_L0, col_L0 = radar_data.shape\n",
    "    print(\"RAW data size: \", radar_data.shape)\n",
    "    print(\"SLC data size: \", L1_size[i])\n",
    "    raw = raw_L0 - raw_L1\n",
    "    col = col_L0 - col_L1\n",
    "    print(\"Raw and Columns to be subtraced: \", raw, col)\n",
    "    L0 = radar_data[0:raw_L0-raw,0:col_L0-col]\n",
    "    print(\"Check final L0 shape: \", L0.shape)\n",
    "    \n",
    "    del radar_data, raw, col\n",
    "    #plt.figure(figsize = (20,20))\n",
    "    #plt.imshow(L0_re[0:917,0:3034], cmap = 'gray', vmin = L0_re.min(), vmax = L0_re.max())\n",
    "    #plt.colorbar()\n",
    "    #plt.close()\n",
    "    \n",
    "    #Cut patches\n",
    "    folder = ''\n",
    "    L0_patches = create_patches(L0.real, L0.imag, folder, j)\n",
    "    print('Tot number of patches for the current image: ', len(L0_patches))\n",
    "    print(\"Dim of a single patch: \", L0_patches[len(L0_patches)-1].shape)\n",
    "    i = i + 1\n",
    "    j = j + 1\n",
    "    print(i)\n",
    "    del L0, L0_patches\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "sarglasses",
   "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"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}