get_wgs84.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "import json\n",
    "\n",
    "import gdal\n",
    "from osgeo import ogr\n",
    "from osgeo import osr"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Sample code to extract WGS84 bounds from Geotif files"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_wgs84_bounds(fname):\n",
    "    src = gdal.Open(fname)\n",
    "    ulx, xres, xskew, uly, yskew, yres  = src.GetGeoTransform()\n",
    "    lrx = ulx + (src.RasterXSize * xres)\n",
    "    lry = uly + (src.RasterYSize * yres)\n",
    "    \n",
    "    # Setup the source projection - you can also import from epsg, proj4...\n",
    "    source = osr.SpatialReference()\n",
    "    source.ImportFromWkt(src.GetProjection())\n",
    "\n",
    "    # The target projection\n",
    "    target = osr.SpatialReference()\n",
    "    target.ImportFromEPSG(4326)\n",
    "\n",
    "    # Create the transform from source to WGS84\n",
    "    transform = osr.CoordinateTransformation(source, target)\n",
    "\n",
    "    # Transform the points\n",
    "    # top-left corner\n",
    "    lng1, lat1, z1 = transform.TransformPoint(ulx, uly)\n",
    "    # top-right corner\n",
    "    lng2,lat2,z2 = transform.TransformPoint(lrx, uly)\n",
    "    # bot-right corner\n",
    "    lng3, lat3, z2 = transform.TransformPoint(lrx, lry)\n",
    "    #bot-left corner\n",
    "    lng4, lat4, z2 = transform.TransformPoint(ulx, lry)\n",
    "    \n",
    "    \n",
    "    return [lng1,lat1,lng2,lat2,lng3,lat3,lng4,lat4]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Original tif files not supplied, but can be downloaded via link in Readme\n",
    "\n",
    "### For each map, the WGS84 bounds were already extracted and saved in json dir."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "data_dir = 'dataset/TOP50raster'\n",
    "json_dir = os.path.join(data_dir,'json')\n",
    "tif_dir = os.path.join(data_dir,'tif')\n",
    "\n",
    "for fname in os.listdir(tif_dir):\n",
    "    if '.tif' in fname:\n",
    "        map_fname = os.path.join(tif_dir,fname)\n",
    "        \n",
    "        json_fname = os.path.join(json_dir,fname.replace('.tif','.json'))\n",
    "        bounds = get_wgs84_bounds(map_fname)\n",
    "        obj = {'label_coordinates_WGS84': bounds}\n",
    "        json.dump(obj,open(json_fname,'w'))\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.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}