{
  "cells": [
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "%matplotlib inline"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "\n# Spherical Grid Conversion\n\n* **ACTM Performer:** AIBEDO Team (PARC, UVic and UW);\n* **Author:** Haruki Hirasawa (hhirasawa@uvic.ca) \n\nWe use ``pygsp`` package to interpolate 2D-gridded data into a spherical grid format. This is done in two steps:\n\n**Step 1:** Generate a text file (\"skeleton\") with the desired Icosahedral level using our python script \nlink [here](https://github.com/kramea/aibedo/blob/preprocess_MS3/preprocessing/gen_icosph_gridfile.py). \n\nThis code block generates a text file that will be used to generate the spherical sample for level 6. \nTo generate a text file for another grid level, please change the ``sphlevel`` in the code. \n\n**Step 2:** Once the text file is generated in step 1, we use the ``cdo`` (Climate Data Operator) \ncommand line tool to generate the interpolated ``netCDF`` file. \nPlease see [here](https://www.isimip.org/protocol/preparing-simulation-files/cdo-help/) for instructions to download ``cdo``. \n\nThe following script is given in command line to generate the interpolated file for model training:\n\n``cdo remapbil,icosphere_6.txt in.nc out.nc``\n\nHere, ``in.nc`` is the 2D-gridded file from ESM model ensembles or Reanalysis datasets, and ``out.nc`` is \nthe name of the interpolated file that will be used for model training.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "import pygsp as pg\nimport numpy as np\n### Generates a text file with centre points of a icosahedral sphere grid\n### Formatted to pass as an argument for 'cdo remap,GRIDFILE in.nc out.nc'\ndef xyz2lonlat(x:float,y:float,zLfloat, radius=6371.0e6):\n    \"\"\"From cartesian geocentric coordinates to 2D geographic coordinates.\"\"\"\n    \"\"\" lon: [-180, 180], lat:[-90, 90] \"\"\"\n    latitude = np.arcsin(z / radius)/np.pi*180\n    longitude = np.arctan2(y, x)/np.pi*180\n    return longitude, latitude \ndef gen_icosphere(level:int,radius:float = 6371.0e6):\n    '''\n    Generates points on a icosphere grid (from Soo Kim)\n    level (int): iteration level\n    radius (float): radius of the earth\n    '''\n    #(1) generate graph of SphereIcosahedral\n    graph = pg.graphs.SphereIcosahedral(2**level)\n    #(2) extract lon, lat coordinate\n    vertices = graph.coords # xyz coordinates\n    #print(\"number of vertices for level \"+str(level)+\" (2**\"+str(level)+\")is \"+str(graph.n_vertices))\n    lon = []\n    lat = []\n    radius = 1\n    # convert cartesian points to spherical lon/lat\n    for tmp_xyz in vertices:\n        tmp_lon, tmp_lat = xyz2lonlat(tmp_xyz[0],tmp_xyz[1],tmp_xyz[2], radius=radius)\n        lon.append(tmp_lon)\n        lat.append(tmp_lat)\n    return lon, lat\nsphlevel = 6\nlon,lat = gen_icosphere(sphlevel,radius = 6371.0e6)\nf = open(\"isosphere_{0}.txt\".format(sphlevel), \"a\")\n# write grid to file\nf.write(\"gridtype = unstructured\\n\")\nf.write(\"gridsize = {0}\\n\".format(len(lon)))\nf.write(\"# Longitudes\\n\")\nf.write(\"xvals = \" + ' '.join(map(str, lon))+\"\\n\")\nf.write(\"# Latitudes\\n\")\nf.write(\"yvals = \" + ' '.join(map(str, lat))+\"\\n\")"
      ]
    }
  ],
  "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": 0
}