{ "cells": [ { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import folium\n", "import matplotlib\n", "from matplotlib.pyplot import imread\n", "from matplotlib.colors import Normalize\n", "from folium import raster_layers" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "image = '/notebooks/resources/gpm/gpm_1d.20190531.tif'" ] }, { "cell_type": "code", "execution_count": 55, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-180.0 10.0\n", "-180.0 -60.0\n", "-30.0 -60.0\n", "-30.0 10.0\n", "([-60.0, -30.0], [10.0, -180.0])\n" ] } ], "source": [ "from osgeo import gdal,ogr,osr\n", "\n", "def GetExtent(gt,cols,rows):\n", " ''' Return list of corner coordinates from a geotransform\n", "\n", " @type gt: C{tuple/list}\n", " @param gt: geotransform\n", " @type cols: C{int}\n", " @param cols: number of columns in the dataset\n", " @type rows: C{int}\n", " @param rows: number of rows in the dataset\n", " @rtype: C{[float,...,float]}\n", " @return: coordinates of each corner\n", " '''\n", " ext=[]\n", " xarr=[0,cols]\n", " yarr=[0,rows]\n", "\n", " for px in xarr:\n", " for py in yarr:\n", " x=gt[0]+(px*gt[1])+(py*gt[2])\n", " y=gt[3]+(px*gt[4])+(py*gt[5])\n", " ext.append([x,y])\n", " print (x,y)\n", " yarr.reverse()\n", " return ext\n", "\n", "def ReprojectCoords(coords,src_srs,tgt_srs):\n", " ''' Reproject a list of x,y coordinates.\n", "\n", " @type geom: C{tuple/list}\n", " @param geom: List of [[x,y],...[x,y]] coordinates\n", " @type src_srs: C{osr.SpatialReference}\n", " @param src_srs: OSR SpatialReference object\n", " @type tgt_srs: C{osr.SpatialReference}\n", " @param tgt_srs: OSR SpatialReference object\n", " @rtype: C{tuple/list}\n", " @return: List of transformed [[x,y],...[x,y]] coordinates\n", " '''\n", " trans_coords=[]\n", " transform = osr.CoordinateTransformation( src_srs, tgt_srs)\n", " for x,y in coords:\n", " x,y,z = transform.TransformPoint(x,y)\n", " trans_coords.append([x,y])\n", " return trans_coords\n", "\n", "raster=image\n", "ds=gdal.Open(raster)\n", "\n", "gt=ds.GetGeoTransform()\n", "cols = ds.RasterXSize\n", "rows = ds.RasterYSize\n", "ext=GetExtent(gt,cols,rows)\n", "\n", "src_srs=osr.SpatialReference()\n", "src_srs.ImportFromWkt(ds.GetProjection())\n", "# tgt_srs=osr.SpatialReference()\n", "# tgt_srs.ImportFromEPSG(4326)\n", "tgt_srs = src_srs.CloneGeogCS()\n", "\n", "geo_ext=ReprojectCoords(ext,src_srs,tgt_srs)\n", "\n", "geo_fix = geo_ext\n", "geo_fix[0] = geo_fix[0][::-1]\n", "geo_fix[2] = geo_fix[2][::-1]\n", "extent = geo_fix[2], geo_fix[0]\n", "print(extent)\n" ] }, { "cell_type": "code", "execution_count": 58, "metadata": {}, "outputs": [ { "ename": "AttributeError", "evalue": "'tuple' object has no attribute 'reverse'", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mAttributeError\u001b[0m Traceback (most recent call last)", "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mextent\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mreverse\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m", "\u001b[0;31mAttributeError\u001b[0m: 'tuple' object has no attribute 'reverse'" ] } ], "source": [] }, { "cell_type": "code", "execution_count": 56, "metadata": { "scrolled": false }, "outputs": [ { "data": { "text/html": [ "
" ], "text/plain": [ "" ] }, "execution_count": 56, "metadata": {}, "output_type": "execute_result" } ], "source": [ "m = folium.Map(\n", " control_scale=True \n", ")\n", "\n", "data = matplotlib.pyplot.imread(image)\n", "\n", "# Image bounds on the map in the form\n", "# [[lat_min, lon_min], [lat_max, lon_max]]\n", "m.add_child(raster_layers.ImageOverlay(\n", " data\n", "# , bounds = [[-60.0000000000000000, -30.0000000000000000], [10.0000000000000000, -180.0000000000000000]]\n", " , bounds = extent\n", " , mercator_project = True\n", ")\n", " )" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "print(data)\n", "print(data.shape)\n", "# data = data.transpose()\n", "# print(data)\n", "# print(data.shape)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "!gdalinfo '/notebooks/resources/gpm/gpm_1d.20190531.tif'" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# First: read the geotiff image with GDAL.\n", "from osgeo import gdal, osr\n", "\n", "gdal.UseExceptions()\n", "\n", "ds = gdal.Open(image)\n", "data = ds.ReadAsArray()\n", "gt = ds.GetGeoTransform()\n", "proj = ds.GetProjection()\n", "\n", "inproj = osr.SpatialReference()\n", "inproj.ImportFromWkt(proj)\n", "\n", "print('\\n\\n## ds ##:\\n\\n' + str(ds))\n", "print('\\n\\n## data ##:\\n\\n' + str(data))\n", "print('\\n\\n## gt ##:\\n\\n' + str(gt))\n", "print('\\n\\n## proj ##:\\n\\n' + str(proj))\n", "print('\\n\\n## inproj ##:\\n\\n' + str(inproj))" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "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.5.2" } }, "nbformat": 4, "nbformat_minor": 2 }