GDAL_30

gdal_30.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'3000000'"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import gdal\n",
    "import osr\n",
    "\n",
    "gdal.UseExceptions()\n",
    "gdal.VersionInfo()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Define projections"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "srs_4326 = osr.SpatialReference()\n",
    "srs_4326.ImportFromEPSG(4326)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "srs_3857 = osr.SpatialReference()\n",
    "srs_3857.ImportFromEPSG(3857)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Define transformations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "ct_3857_to_4326 = osr.CoordinateTransformation(srs_3857, srs_4326)\n",
    "ct_4326_to_3857 = osr.CoordinateTransformation(srs_4326, srs_3857)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Forward transform (x, y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(inf, inf)"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lat = 50.0\n",
    "lon = -120.0\n",
    "\n",
    "# Incorrect output, and different compared to GDAL 2.4\n",
    "mapx, mapy, z = ct_4326_to_3857.TransformPoint(lon, lat)\n",
    "mapx, mapy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Forward transform (y, x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-13358338.895192828, 6446275.841017158)"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lat = 50.0\n",
    "lon = -120.0\n",
    "\n",
    "# this output is equivalent to the cmd (gdaltransform.exe) (but (x,y) -> (y,x))\n",
    "mapx, mapy, z = ct_4326_to_3857.TransformPoint(lat, lon)\n",
    "mapx, mapy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Backward transform (x, y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(49.99999999999999, -119.99999999999999, 0.0)"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mapx = -13358338.895192828\n",
    "mapy = 6446275.841017158\n",
    "\n",
    "# Wrong output order (?), but correct values # (x,y) -> (y,x)\n",
    "ct_3857_to_4326.TransformPoint(mapx, mapy)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Backward transform (y, x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-75.95934455387825, 57.90788113636135, 0.0)"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mapx = -13358338.895192828\n",
    "mapy = 6446275.841017158\n",
    "\n",
    "# Incorrect results\n",
    "ct_3857_to_4326.TransformPoint(mapy, mapx)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### CMD output reference"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "# (py37) C:\\Users\\Rutger>gdaltransform --version\n",
    "# GDAL 3.0.0, released 2019/05/05\n",
    "\n",
    "# (py37) C:\\Users\\Rutger>gdaltransform -s_srs EPSG:4326 -t_srs EPSG:3857\n",
    "# -120 50\n",
    "# -13358338.8951928 6446275.84101716 0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### CoordinateTransform docstring"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "\u001b[1;31mSignature:\u001b[0m \u001b[0mct_4326_to_3857\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mTransformPoint\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;33m*\u001b[0m\u001b[0margs\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
       "\u001b[1;31mDocstring:\u001b[0m\n",
       "TransformPoint(CoordinateTransformation self, double [3] inout)\n",
       "TransformPoint(CoordinateTransformation self, double [4] inout)\n",
       "TransformPoint(CoordinateTransformation self, double x, double y, double z=0.0)\n",
       "TransformPoint(CoordinateTransformation self, double x, double y, double z, double t)\n",
       "\u001b[1;31mFile:\u001b[0m      c:\\miniconda3\\envs\\py37\\lib\\site-packages\\osgeo\\osr.py\n",
       "\u001b[1;31mType:\u001b[0m      method\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# docstring suggests (x,y,z)\n",
    "ct_4326_to_3857.TransformPoint?"
   ]
  },
  {
   "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.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}