Retrieving cell area with xESMF (calling ESMPy)

xESMF_area.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "- Install the development branch by `pip install --upgrade git+https://github.com/JiaweiZhuang/xESMF.git@area`\n",
    "- Reference code: https://gist.github.com/rokuingh/058e1bc6e08c4d5e33f37657bcfb13d6"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt\n",
    "import xesmf as xe"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<xarray.Dataset>\n",
       "Dimensions:  (x: 144, x_b: 145, y: 90, y_b: 91)\n",
       "Coordinates:\n",
       "    lon      (y, x) float64 -178.8 -176.2 -173.8 -171.2 -168.8 -166.2 -163.8 ...\n",
       "    lat      (y, x) float64 -89.0 -89.0 -89.0 -89.0 -89.0 -89.0 -89.0 -89.0 ...\n",
       "    lon_b    (y_b, x_b) float64 -180.0 -177.5 -175.0 -172.5 -170.0 -167.5 ...\n",
       "    lat_b    (y_b, x_b) int64 -90 -90 -90 -90 -90 -90 -90 -90 -90 -90 -90 ...\n",
       "Dimensions without coordinates: x, x_b, y, y_b\n",
       "Data variables:\n",
       "    *empty*"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ds = xe.util.grid_global(2.5, 2)\n",
    "ds"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Help on function cell_area in module xesmf.util:\n",
      "\n",
      "cell_area(ds)\n",
      "    Get cell area of a grid.\n",
      "    Assume unit sphere (radius is 1, total area is 4*pi)\n",
      "    \n",
      "    Parameters\n",
      "    ----------\n",
      "    ds : xarray DataSet or dictionary\n",
      "        Contains variables ``lon``, ``lat``, ``lon_b``, ``lat_b``\n",
      "        Note that boundary is required for computing cell area\n",
      "    \n",
      "    Returns\n",
      "    -------\n",
      "    area : 2D numpy array for cell area\n",
      "\n"
     ]
    }
   ],
   "source": [
    "help(xe.util.cell_area)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "area = xe.util.cell_area(ds)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(90, 144)"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "area.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "12.566370614359172"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "area.sum() # 4*pi, unit sphere"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.collections.QuadMesh at 0x2ac4dfdd0198>"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD8CAYAAABn919SAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAFBNJREFUeJzt3X+QXfV53/H3Z1dCQgKKBIYqggbcMq5dpthUQyBuMy7YjXE9QDrQwfW0mpap/nEbHCcTRDxTN512RkwzsdOZFEcTHGs6lB+WoWKYiVNFxe3kjxILG9tgTIUxwTIKMg4Y8VurffrHPYo3Ylf37I97dvfs+zWzc+4595zzffiy99lH3+8956SqkCQtf2OLHYAkaWGY0CWpJ0zoktQTJnRJ6gkTuiT1hAldknqiVUJPcnOSx5I8nuSTzbaNSfYmOdAsN4w2VEnSyQxN6EkuBv41cBlwCfDRJBcB24F9VXURsK9ZlyQtkjYV+ruB/1tVr1XVBPC/gV8CrgV2NfvsAq4bTYiSpDZWtdjnMeA/JTkLeB34CLAfOLeqDgFU1aEk50x3cJJtwDaAccb/3jrOWJDAJWmlOMKLL1TVO4btNzShV9UTSW4D9gKvAN8EJtoGUlU7gZ0AZ2Rj/VyuanuoJAn449r9Z232azUpWlV3VNWlVfULwF8AB4Dnk2wCaJaH5xqsJGn+2gy5kOScqjqc5G8A/wS4ArgQ2ArsaJZ7hp5nbIzx006fR7iStAK93G63Vgkd+HIzhn4U+ERVvZhkB3BvkpuAZ4Eb5hKnJGlhtEroVfUPptn2Y8ABcUlaItpW6AtjfJycflqnTUrSstdyyMVL/yWpJ7qt0MfGYP26TpuUpJXCCl2SesKELkk90fGkaJg8bU2nTUrSSmGFLkk90WmFXmPh2PpTumxSklYMK3RJ6onOK/SJdd0O20vSSmGFLkk9YUKXpJ7o+ErRcOzU8U6blKSVwgpdknqi40lROLouXTYpSStGqwo9ya8keTzJY0nuSrI2yYVJHk5yIMk9SfyCuSQtoqEJPclm4JeBLVV1MTAO3AjcBny2qi4CXgRuGmWgkqSTazvksgo4NclRYB1wCLgS+GfN+7uAfw/cfrKT1BhMnOqQiySNwtAKvap+CPwWg+eGHgJ+AjwCvFRVE81uB4HNowpSkjTc0Ao9yQbgWuBC4CXgS8DV0+xaMxy/DdgGsPq0DUystUKXpFFoMyn6QeD7VfWjqjoK3Af8PHBmkuN/EM4Dnpvu4KraWVVbqmrLqlPXL0jQkqS3azOG/ixweZJ1wOvAVcB+4CHgeuBuYCuwZ9iJagyOrZ17sJKkmbUZQ38Y2A18Hfh2c8xO4BbgU0meAs4C7hhhnJKkIVp9y6WqPgN85oTNTwOXLXhEkqQ56fZeLnHIRZJGxXu5SFJPdH4vl2M+I1qSRsIKXZJ6woQuST2xCJOi015QKkmaJyt0SeqJ7idFrdAlaSSs0CWpJ0zoktQTHU+KFrVmstMmJWmlsEKXpJ7otkIfA6zQJWkkrNAlqSc6H0MfX3Os0yYlaaWwQpeknhia0JO8K8mjU35eTvLJJBuT7E1yoFlu6CJgSdL0hg65VNWTwHsBkowDPwTuB7YD+6pqR5LtzfotJztXxmD1KRPzDlqS9HazHXK5CvheVf0ZcC2wq9m+C7huIQOTJM3ObCdFbwTual6fW1WHAKrqUJJzpjsgyTZgG8Dqd5zB2jVH5xqrJOkkWlfoSU4BrgG+NJsGqmpnVW2pqi2r/tr62cYnSWppNkMuVwNfr6rnm/Xnk2wCaJaHFzo4SVJ7sxly+Rg/HW4BeADYCuxolnuGnWAsk6w/5a1ZBShJaqdVhZ5kHfAh4L4pm3cAH0pyoHlvx8KHJ0lqq1WFXlWvAWedsO3HDL710tpYinWrrdAlaRS8UlSSeqLTe7mMpzjNCl2SRsIKXZJ6woQuST3R8ZDLJGesfqPLJiVpxbBCl6Se6LRCH0uxftWbXTYpSSuGFbok9YQJXZJ6ovtJ0VVOikrSKFihS1JPdDspSnH6uBW6JI2CFbok9YQJXZJ6ovNJ0dMccpGkkWj7gIszk+xO8t0kTyS5IsnGJHuTHGiWG0YdrCRpZm0r9N8BvlJV1zcPi14H/Aawr6p2JNkObAduOdlJxpjk9DErdEkahaEVepIzgF8A7gCoqreq6iXgWmBXs9su4LpRBSlJGq5Nhf5O4EfAHyS5BHgEuBk4t6oOAVTVoSTnTHdwkm3ANoCzf+YUTh9/fUEClyT9VW3G0FcBlwK3V9X7gFcZDK+0UlU7q2pLVW05Y2Onc7CStKK0SegHgYNV9XCzvptBgn8+ySaAZnl4NCFKktoYWjJX1Z8n+UGSd1XVk8BVwHean63Ajma5Z9i5xphk/Zi3z5WkUWg7BvJvgTubb7g8DfxLBtX9vUluAp4FbhhNiJKkNlol9Kp6FNgyzVtXzaaxsRTr89ZsDpEkteSl/5LUEyZ0SeqJzm+fu27MIRdJGgUrdEnqic4r9LWZ6LJJSVoxrNAlqSc6rtBhbSa7bFKSVgwrdEnqCRO6JPVEp0MuSbE21WWTkrRiWKFLUk90PCka1sa/IZI0CmZXSeoJE7ok9US3k6LAav+GSNJItEroSZ4BjgDHgImq2pJkI3APcAHwDPBPq+rF0YQpSRpmNhX6P6yqF6asbwf2VdWOJNub9VtOdoIQVjspKkkjMZ/sei2wq3m9C7hu/uFIkuaqbYVewP9MUsDvVdVO4NyqOgRQVYeSnDPdgUm2AdsAzt88zphj6JI0Em0T+vur6rkmae9N8t22DTTJfyfApZes8TJRSRqRVuVyVT3XLA8D9wOXAc8n2QTQLA+PKkhJ0nBDK/Qk64GxqjrSvP5HwH8AHgC2Ajua5Z6h5yKszvj8IpYkTavNkMu5wP1Jju//36vqK0m+Btyb5CbgWeCG0YUpSRpmaEKvqqeBS6bZ/mPgqtk0FnBSVJJGxOwqST1hQpeknjChS1JPmNAlqSdM6JLUEyZ0SeoJE7ok9YQJXZJ6woQuST1hQpeknjChS1JPmNAlqSdM6JLUEyZ0SeqJto+gWxAFTDLZZZOStGK0rtCTjCf5RpIHm/ULkzyc5ECSe5KcMrowJUnDzKZCvxl4AjijWb8N+GxV3Z3k88BNwO0nO0FRHK1jcwpUknRyrSr0JOcB/xj4/WY9wJXA7maXXcB1owhQktRO2yGXzwG/Dn85AH4W8FJVTTTrB4HN0x2YZFuS/Un2v/Bjq3NJGpWhQy5JPgocrqpHknzg+OZpdq3pjq+qncBOgEsvWVNOikrSaLQZQ38/cE2SjwBrGYyhfw44M8mqpko/D3hudGFKkoYZmtCr6lbgVoCmQv+1qvp4ki8B1wN3A1uBPUPPRXG0rNAlaRTmc2HRLcCnkjzFYEz9joUJSZI0F7O6sKiqvgp8tXn9NHDZrI4HjjqGLkkj4aX/ktQTJnRJ6olO7+UySfGGk6KSNBJW6JLUE93ebbHCGzXdNUmSpPmyQpeknjChS1JPdDwpCm+Uf0MkaRTMrpLUEx1X6OGN6rRJSVoxrNAlqSdM6JLUE50Pubw26bOkJWkUrNAlqSe6rdArvFpW6JI0CkMr9CRrk/xpkm8meTzJbzbbL0zycJIDSe5JYqaWpEXUpkJ/E7iyql5Jshr4kyR/CHwK+GxV3Z3k88BNwO0nO9EkY7w6uWbeQUuS3m5ohV4DrzSrq5ufAq4EdjfbdwHXjSRCSVIrrSZFk4wneRQ4DOwFvge8VFUTzS4Hgc0zHLstyf4k+1/+i4npdpEkLYBWk6JVdQx4b5IzgfuBd0+32wzH7gR2Alxw8Wl15NipcwxVknQys/raYlW9xOAh0ZcDZyY5/gfhPOC5hQ1NkjQbQyv0JO8AjlbVS0lOBT4I3AY8BFwP3A1sBfYMO9ckYxyZXDu/iCVJ02oz5LIJ2JVknEFFf29VPZjkO8DdSf4j8A3gjhHGKUkaYmhCr6pvAe+bZvvTwGWjCEqSNHudXil6rMZ45ZhDLpI0Ct7LRZJ6ovO7LR6xQpekkbBCl6Se6HwM/eUJK3RJGgUrdEnqCRO6JPVE9w+4mPD2uZI0ClboktQT3U+KHnVSVJJGwQpdknrChC5JPdHxkEt45ajPkpakUbBCl6Se6Pxri69ZoUvSSAyt0JOcn+ShJE8keTzJzc32jUn2JjnQLDeMPlxJ0kzaDLlMAL9aVe9m8CzRTyR5D7Ad2FdVFwH7mnVJ0iJp88SiQ8Ch5vWRJE8Am4FrgQ80u+1i8PDoW052rska49W3HHKRpFGY1aRokgsYPI7uYeDcJtkfT/rnzHDMtiT7k+yf+Mmr84tWkjSj1pOiSU4Dvgx8sqpeTtLquKraCewEWPu3Ntcbb66eS5ySpCFaVehJVjNI5ndW1X3N5ueTbGre3wQcHk2IkqQ2hlboGZTidwBPVNVvT3nrAWArsKNZ7hl2rpqEo291+k1JSVox2mTX9wP/HPh2kkebbb/BIJHfm+Qm4FnghtGEKElqo823XP4EmGnA/KqFDUeSNFfdjn9UOPbmeKdNStJK4b1cJKknuq3QJ4E3/RsiSaNgdpWknjChS1JPdD4pGodcJGkkzK6S1BOdVuiZhPE32t0DRpI0O1boktQTHY+hW6FL0qhYoUtST5jQJaknup8UfbPLFiVp5bBCl6SeWIRJ0U5blKQVY2iFnuQLSQ4neWzKto1J9iY50Cw3jDZMSdIwbYZcvgh8+IRt24F9VXURsK9ZlyQtojZPLPo/SS44YfO1wAea17uArwK3DDvX4ErRWcUnSWpprpOi51bVIYBmec7ChSRJmouRT4om2QZsAzhl/QZWvVGjblKSVqS5VujPJ9kE0CwPz7RjVe2sqi1VtWXVqevn2JwkaZi5JvQHgK3N663AnoUJR5I0V0OHXJLcxWAC9OwkB4HPADuAe5PcBDwL3NCmsUzCqtcdcpGkUWjzLZePzfDWVQsciyRpHjq/l8vq16zQJWkUvJeLJPVEt/dymSzGXz/WaZOStFJYoUtST5jQJaknOp4ULVa9NtFlk5K0YlihS1JPdF6hj7/6VpdNStKKYYUuST1hQpeknuj2e+jHirFX3uy0SUlaKazQJaknOr5SdBJefa3TJiVppbBCl6Se6HgM/Rh15JVOm5SklcIKXZJ6Yl4JPcmHkzyZ5Kkk2xcqKEnS7M15yCXJOPC7wIeAg8DXkjxQVd+Z6ZianOTYkSNzbVKSdBLzqdAvA56qqqer6i3gbuDahQlLkjRb85kU3Qz8YMr6QeDnTtwpyTZgW7P65h/X7sfm0eZiORt4YbGDmKPlGvtyjRuWb+zLNW5YvrG3jftn25xsPgk902x72wNDq2onsBMgyf6q2jKPNhfFco0blm/syzVuWL6xL9e4YfnGvtBxz2fI5SBw/pT184Dn5heOJGmu5pPQvwZclOTCJKcANwIPLExYkqTZmvOQS1VNJPk3wB8B48AXqurxIYftnGt7i2y5xg3LN/blGjcs39iXa9ywfGNf0LhT9bZhb0nSMuSVopLUEyZ0SeqJThL6crpFQJLzkzyU5Ikkjye5udm+McneJAea5YbFjnU6ScaTfCPJg836hUkebuK+p5nAXnKSnJlkd5LvNn1/xXLo8yS/0vyePJbkriRrl2qfJ/lCksNJHpuybdo+zsB/aT6z30py6RKL+z83vyvfSnJ/kjOnvHdrE/eTSX5xcaL+y1jeFvuU934tSSU5u1mfd5+PPKFPuUXA1cB7gI8lec+o252HCeBXq+rdwOXAJ5p4twP7quoiYF+zvhTdDDwxZf024LNN3C8CNy1KVMP9DvCVqvrbwCUM/huWdJ8n2Qz8MrClqi5m8OWAG1m6ff5F4MMnbJupj68GLmp+tgG3dxTjdL7I2+PeC1xcVX8X+H/ArQDNZ/VG4O80x/zXJgctli/y9thJcj6D26Y8O2Xz/Pu8qkb6A1wB/NGU9VuBW0fd7gLGv6fp+CeBTc22TcCTix3bNLGex+BDeSXwIIOLv14AVk33/2Kp/ABnAN+nmaSfsn1J9zk/vVp6I4NvjD0I/OJS7nPgAuCxYX0M/B7wsen2Wwpxn/DeLwF3Nq//Sn5h8C28K5ZSnzfbdjMoXJ4Bzl6oPu9iyGW6WwRs7qDdeUtyAfA+4GHg3Ko6BNAsz1m8yGb0OeDXgclm/SzgpaqaaNaXat+/E/gR8AfNcNHvJ1nPEu/zqvoh8FsMqqxDwE+AR1gefX7cTH28nD63/wr4w+b1ko87yTXAD6vqmye8Ne/Yu0jorW4RsNQkOQ34MvDJqnp5seMZJslHgcNV9cjUzdPsuhT7fhVwKXB7Vb0PeJUlNrwynWa8+VrgQuBngPUM/tl8oqXY58Msi9+dJJ9mMEx65/FN0+y2ZOJOsg74NPDvpnt7mm2zir2LhL7sbhGQZDWDZH5nVd3XbH4+yabm/U3A4cWKbwbvB65J8gyDO19eyaBiPzPJ8QvIlmrfHwQOVtXDzfpuBgl+qff5B4HvV9WPquoocB/w8yyPPj9upj5e8p/bJFuBjwIfr2aMgqUf999kUAB8s/msngd8PclfZwFi7yKhL6tbBCQJcAfwRFX99pS3HgC2Nq+3MhhbXzKq6taqOq+qLmDQx/+rqj4OPARc3+y25OIGqKo/B36Q5F3NpquA77DE+5zBUMvlSdY1vzfH417yfT7FTH38APAvmm9eXA785PjQzFKQ5MPALcA1VTX1yfMPADcmWZPkQgYTjH+6GDFOp6q+XVXnVNUFzWf1IHBp8xmYf593NCnwEQYz0d8DPr2YExQtYv37DP6Z8y3g0ebnIwzGo/cBB5rlxsWO9ST/DR8AHmxev5PBL/RTwJeANYsd3wwxvxfY3/T7/wA2LIc+B34T+C7wGPDfgDVLtc+BuxiM9R9tEslNM/Uxg3/+/27zmf02g2/yLKW4n2Iw3nz8M/r5Kft/uon7SeDqpdbnJ7z/DD+dFJ13n3vpvyT1hFeKSlJPmNAlqSdM6JLUEyZ0SeoJE7ok9YQJXZJ6woQuST3x/wG/sW2uCVLRzwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# sanity check\n",
    "plt.pcolormesh(area) "
   ]
  },
  {
   "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.6.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}

Great, very useful! It would be nice if this could be merged into master ...

Any idea to merge it into master?

Very nice.. I have been using this development branch and recently had to renovate all my anaconda distribution. Just curious if this is still a branch of if its already in master?