23_DL_data_expanded.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Data pre-processing for u-net model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We create a new dataset for the u-net model with the additional variable IR_016."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "## You are using the Python ARM Radar Toolkit (Py-ART), an open source\n",
      "## library for working with weather radar data. Py-ART is partly\n",
      "## supported by the U.S. Department of Energy as part of the Atmospheric\n",
      "## Radiation Measurement (ARM) Climate Research Facility, an Office of\n",
      "## Science user facility.\n",
      "##\n",
      "## If you use this software to prepare a publication, please cite:\n",
      "##\n",
      "##     JJ Helmus and SM Collis, JORS 2016, doi: 10.5334/jors.119\n",
      "\n"
     ]
    }
   ],
   "source": [
    "import glob\n",
    "import os\n",
    "import pyart\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "from copy import deepcopy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/opt/anaconda3/envs/exts-ml2/lib/python3.7/site-packages/ipykernel/ipkernel.py:287: DeprecationWarning: `should_run_async` will not call `transform_cell` automatically in the future. Please pass the result to `transformed_cell` argument and any exception that happen during thetransform in `preprocessing_exc_tuple` in IPython 7.17 and above.\n",
      "  and should_run_async(code)\n"
     ]
    }
   ],
   "source": [
    "# suppress anoying iypthon warnings. Not ideal since we suppress also potentially relevant warnings\n",
    "import warnings\n",
    "warnings.filterwarnings('ignore')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Auxiliary functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Function to read original dataset\n",
    "# data is stored as (nz, ny, nx), we return (nx, ny)\n",
    "def read_nc(fname):\n",
    "    sat_grid = pyart.io.read_grid(fname)\n",
    "    for field_name in sat_grid.fields.keys():\n",
    "        data = np.transpose(np.squeeze(sat_grid.fields[field_name]['data']))\n",
    "    return data        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Function for minmax scaling\n",
    "def minmax_scaling(data, vmin, vmax):\n",
    "    data2 = deepcopy(data)\n",
    "    data2[data2>vmax] = vmax\n",
    "    data2[data2<vmin] = vmin\n",
    "    return (data2-vmin)/(vmax-vmin)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Some global variables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "fbasepath = '/data/pyrad_products/MSG_ML/'\n",
    "dl_basepath = '/data/ml_course/05_Capstone_project/dl_data_expanded/'\n",
    "features = ['HRV_norm', 'IR_108', 'WV_062-IR_108', 'IR_016']\n",
    "nfeatures = len(features)\n",
    "target = 'POH90'\n",
    "\n",
    "vmins = [0., 200., -78., 0.]\n",
    "vmaxs = [100., 311., 9., 67.]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We use minmax normalization to put all variables within the 0-1 range. The min, max values for each variable have been obtained from the EDA. The features matrix has shape nx, ny, n channels (HRV_norm, IR_108, WV_062-IR_108 and IR_016). The target matrix has shape nx, ny, n classes. The classes are no hail (0) and hail (1). We transform the POH90 to 0 (no hail or not computed) and 1 probabilty of hail above 90%). The shape of those matrices is the one required by the u-net as implemented in the unet package. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "20200731173000\r"
     ]
    }
   ],
   "source": [
    "years = ['2018', '2019', '2020']\n",
    "months = ['04', '05', '06', '07', '08', '09']\n",
    "for year in years:\n",
    "    for month in months:\n",
    "        # Get list of files and data size\n",
    "        flist = glob.glob(fbasepath+'*/NETCDF/'+features[0]+'/'+year+month+'*.nc')\n",
    "        if len(flist) == 0:\n",
    "            continue\n",
    "        flist.sort()\n",
    "        img_size = read_nc(flist[0]).shape\n",
    "        data_size = img_size[0]*img_size[1]\n",
    "        \n",
    "        for fname in flist:\n",
    "            # Get time step\n",
    "            bfile = os.path.basename(fname)\n",
    "            dt_str = bfile[0:14]\n",
    "            print(dt_str, end=\"\\r\", flush=True)\n",
    "            \n",
    "            # Read all files corresponding to a time step\n",
    "            # Put them in features and target matrices\n",
    "            X = np.empty((img_size[0], img_size[1], nfeatures), dtype=np.float32)\n",
    "            for i, (vmin, vmax, feature) in enumerate(zip(vmins, vmaxs, features)):\n",
    "                flist_aux = glob.glob(fbasepath+'*/NETCDF/'+feature+'/'+dt_str+'*.nc')\n",
    "                data = read_nc(flist_aux[0])\n",
    "                data = minmax_scaling(data, vmin, vmax)  \n",
    "                X[:, :, i] = data\n",
    "               \n",
    "            flist_aux = glob.glob(fbasepath+'*/NETCDF/'+target+'/'+dt_str+'*.nc')\n",
    "            y = read_nc(flist_aux[0])\n",
    "            \n",
    "            # Only hail/no hail\n",
    "            y[y == 1] = 0\n",
    "            y[y == 2] = 1\n",
    "            \n",
    "            # onehot encoding\n",
    "            y_onehot = np.eye(2)[y]\n",
    "             \n",
    "            # Save data into a .npz file\n",
    "            np.savez(dl_basepath+dt_str+'_data.npz', features=X, targets=y_onehot)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.7.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}