kalman_filters_tests.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Simple notebook that loads test data and applies the ensemble Kalman Filters.\n",
    "\n",
    "We compare the mean of the analysis ensembles and the covariance matrices.\n",
    "Actually only 1000 entries of the covariance matrix are computed due to the large size of the prior ensemble. Ther perturbations don't have to be equal"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[Errno 2] No such file or directory: 'kalmanfilters'\n",
      "/home/ldap-server/draco/notebooks/kalmanfilters\n"
     ]
    }
   ],
   "source": [
    "%cd kalmanfilters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "from ensrf import *\n",
    "from ensrf_direct import *\n",
    "from estkf import *\n",
    "from etkf import *\n",
    "from etkf_livings import *\n",
    "from ensrf_serial import *\n",
    "import numpy as np\n",
    "import scipy\n",
    "from time import time"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The autoreload extension is already loaded. To reload it, use:\n",
      "  %reload_ext autoreload\n"
     ]
    }
   ],
   "source": [
    "%load_ext autoreload\n",
    "%autoreload 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [],
   "source": [
    "#You can look into the imported files using\n",
    "\n",
    "%less estkf.py"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 63,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Y shape: (293,)\n",
      "R shape: (293,)\n",
      "Xf shape: (55296, 100)\n",
      "HXf shape: (293, 100)\n"
     ]
    }
   ],
   "source": [
    "###LOAD TEST DATA\n",
    "\n",
    "Y=np.load('../testdata/Y.npy',allow_pickle=True)\n",
    "R=np.load('../testdata/R.npy',allow_pickle=True)\n",
    "Xf=np.load('../testdata/Xf.npz',allow_pickle=True)['arr_0']\n",
    "HXf=np.load('../testdata/HXf.npy',allow_pickle=True)\n",
    "\n",
    "print('Y shape:',np.shape(Y))\n",
    "print('R shape:',np.shape(R))\n",
    "print('Xf shape:',np.shape(Xf))\n",
    "print('HXf shape:',np.shape(HXf))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "-------------------------------\n",
      "ESTKF  executed in  0.056397438049316406 seconds\n",
      "-------------------------------\n",
      "-------------------------------\n",
      "EnSRF  executed in  0.1246793270111084 seconds\n",
      "-------------------------------\n",
      "-------------------------------\n",
      "EnSRF_serial  executed in  10.726128101348877 seconds\n",
      "-------------------------------\n",
      "-------------------------------\n",
      "ENSRF_direct  executed in  0.41527724266052246 seconds\n",
      "-------------------------------\n",
      "-------------------------------\n",
      "ETKF  executed in  0.05635523796081543 seconds\n",
      "-------------------------------\n",
      "-------------------------------\n",
      "ETKF_livings  executed in  0.06243419647216797 seconds\n",
      "-------------------------------\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "\n",
    "variables=[Xf, HXf, Y, R]\n",
    "\n",
    "mean={}\n",
    "cov={}\n",
    "\n",
    "funcs=[ESTKF, EnSRF, EnSRF_serial,ENSRF_direct, ETKF, ETKF_livings]\n",
    "\n",
    "for i,f in enumerate(funcs):\n",
    "    name=str(f.__name__)\n",
    "    begin = time()\n",
    "    full=f(*variables)\n",
    "    end=time ()\n",
    "    print('-------------------------------')\n",
    "    print(name,' executed in ',end-begin, 'seconds')\n",
    "    print('-------------------------------')\n",
    "    mean[name]=np.mean(full,axis=1)\n",
    "    cov[name]=np.cov(full[:1000,:],ddof=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'ESTKF': array([225.77951468, 225.75510173, 225.75508779, ..., 252.08298075,\n",
       "        252.08490147, 252.08658989]),\n",
       " 'EnSRF': array([225.77951468+0.j, 225.75510173+0.j, 225.75508779+0.j, ...,\n",
       "        252.08298075+0.j, 252.08490147+0.j, 252.08658989+0.j]),\n",
       " 'EnSRF_serial': array([225.77951468, 225.75510173, 225.75508779, ..., 252.08298075,\n",
       "        252.08490147, 252.08658989]),\n",
       " 'ENSRF_direct': array([225.77951468+1.02921696e-32j, 225.75510173+9.65516829e-33j,\n",
       "        225.75508779+1.08391337e-32j, ..., 252.08298075+5.84616035e-33j,\n",
       "        252.08490147+1.06922816e-32j, 252.08658989+1.66809607e-32j]),\n",
       " 'ETKF': array([225.77951468, 225.75510173, 225.75508779, ..., 252.08298075,\n",
       "        252.08490147, 252.08658989])}"
      ]
     },
     "execution_count": 55,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mean"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 62,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "True\n",
      "True\n",
      "True\n",
      "True\n",
      "True\n"
     ]
    }
   ],
   "source": [
    "#check equality of means\n",
    "m=list(mean.values()) \n",
    "for i,me in enumerate(m):\n",
    "    if i<len(m)-1:\n",
    "        print(np.allclose(m[i],m[i+1]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 61,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "True\n",
      "True\n",
      "True\n",
      "True\n",
      "True\n"
     ]
    }
   ],
   "source": [
    "#check equality of covariance matrix entries\n",
    "c=list(cov.values()) \n",
    "for i,me in enumerate(c):\n",
    "    if i<len(c)-1:\n",
    "        print(np.allclose(c[i],c[i+1]))"
   ]
  }
 ],
 "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.10.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}