You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
 

209 lines
13 KiB

{
"cells": [
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"scrolled": false
},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAOcAAADnCAYAAADl9EEgAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAEv9JREFUeJzt3XlQU+ffBfBjCCigSEVxY0DcN9DiUq3Sn7hvdUbr2jJ1LCq4oqJYxY6OtSgMWrEu7bRi615cqxZrVYgLWiu4IFoXxFocheogImLEkLz/VF9FdnPz3Jt7PjNMQpLheyQcn3BJ8lQxmUwgIvnRiA5ARMVjOYlkiuUkkimWk0imWE4imdKWdmW9evVMWVlZlspCpFa3TSZTo6IXlrpyZmVlwWQyCf1ISEhQ9XylZXjy5AnatGkDOzs7hIeHC/keNG7cGK1bt1bM/QDAo7j+8WEtmZWDgwNSUlIQHByMhQsXwsXFBdu2bbPY/OnTpyMjIwMJCQkWmykVlpPMTqPRIDIyEjk5OejTpw/8/f3RuHFjJCYmSjp38+bNWL16NWJiYuDq6irpLEtgOUkyDg4O2L59OzIyMuDm5gZfX198/PHHkszas2cPxo4di9mzZ8Pf31+SGZZW6gEhInNo0KABjh8/jvj4ePTt2xfPnj3Drl27Srz9yZMnodPp0LBhQ4wYMQLVq1cv9esvW7YM8+fPR1BQECIjI80dXxiunGQxPXv2xLFjx7B3717s2LHjtevy8/MxZ84cuLi44IMPPsDKlSsxefJk1KhRAzVr1sTVq1cRFhaG06dPo6CgANnZ2diwYQNatWqFsLAwfP3111i7dq2gf5k0WE6yqG7dumHAgAEICQkBAOh0OnTt2hU1atTA999/j9GjRyM7OxsPHjzA06dPkZGRgbCwMGg0Gnz33Xfo3r07qlatChcXF0yePBm1a9dGamoqgoODBf/LzK/Mh7Xr1q2DjY0NbGxsoNVqX54v7rIX54s7Le18cR8aTfH/bxiNxhI/r+j58tyusLAQDx48MNt8k8kEo9EIo9EIg8EAo9GIwsLC1y4zmUyvXZeXl4ejR4++vM2L2xc9fXHeZDIVe5sXs4u7TXHXvXratm1bzJkz5+Wcoh+v/rtevezVr/3i1GAwICMjAw4ODtDr9Xj33Xexf/9+DBw48I37283NDaGhodDpdC/vB6PRWOLPhzUps5xz584t7u8ypZ4v7fSFop+XJCoqCn5+fuW6rRSioqLQu3dvs3/dKlWqlHrZq+cjIyMxd+7c1y4vz2lZlxU9/+LzFz/4r14+b948bNmy5eV1xd3m1etevbzobV7MqV+/PpKTk+Hs7Fyh750aigmUo5y5ubmWyPEGg8EAjUaD48ePo7CwUNgdotPpyv0fiZQZnj9/LjzD3bt3zfb1HB0dMX/+/AoXU02qlPaDV6VKFZPoP+bm5eWVebTOmudba4bk5GS0atUKDg4OQuZXhlQZ/Pz8YDKZ3nwoVcbTikyiJSQkqHq+NWbIy8szATDp9Xoh8ytLqgz/9eyN/qnjwTvJysWLF6HRaFC1alXRUWSN5SSLS0lJgb29vegYssdyksVdu3aNB4LKgeUki9u7dy+8vLxEx5A9PreWLCouLg63b9/GiRMnREeRPa6cZFFTpkxBr1694ObmJjqK7HHlJIvhqlkxXDnJYrhqVgxXTrIIrpoVx5WTLIKrZsVx5STJcdWsHK6cJDmumpXDlZMkxVWz8rhykqS4alYeV06SzOXLl/H3339z1awkrpwkmfDwcLi7u3PVrCSWkyRz8OBBjBo1SnQMxWI5SRL5+fl4+PAhZs+eLTqKYrGcJIn4+HjY2tpaxZ4lorCcJImbN2+iWrVqomMoGstJkujVqxfy8vLeeBNuKj+WkyTRtm1bAMC5c+cEJ1EulpMk4+zsjLi4ONExFIvlJMk0atQIf/zxh+gYisVykmTatWuHq1evio6hWCwnSaZHjx64d++e6BiKxXKSZAYNGgS9Xo/8/HzRURSJ5STJ1K5dG7a2tjh8+LDoKIrEcpKk6tSpgyNHjoiOoUgsJ0mqefPmSE5OFh1DkVhOklSnTp1w8+ZN0TEUiZvnyny+0jM8fvwY169fR4cOHYTMNydunluE6E1TRc9XeoZnz56ZAJhu374tZL45cfNcsip2dnZwdHTEr7/+KjqK4rCcJLkGDRrwfYQqgeUkybVp0wYpKSmiYygOy0mSe//995GRkSE6huKwnCQ5V1dXPH/+XHQMxWE5SXJpaWnC/wyiRCwnSS45ORn169cXHUNxWE6SXFJSEnr37i06huKwnCSpW7du4f79+wgKChIdRXFYTpLUzJkz4eHhgWbNmomOojjcyIgkYzAYEBcXh+joaNFRFIkrJ0kmKioKNjY2CAwMFB1FkVhOksyqVaswbNgwaDT8MasMftdIEmfOnEFmZiaWL18uOopisZwkidmzZ6Nt27aoV6+e6CiKxQNCZHZ5eXk4deoU9uzZIzqKonHlJLMLCwuDk5MThgwZIjqKorGcZHY//fQTxo0bJzqG4rGcZFZ79uzB48ePsWTJEtFRFI/lJLNasGABfH194eDgIDqK4vGAEJnN4cOH8ddffyE2NlZ0FKvAlZPM5tNPP0W/fv3Qpk0b0VGsAstJZrFixQo8ePAA27ZtEx3FarCc9NYKCgoQFhaGqVOnwtnZWXQcq8Fy0lsLDAyEra0tn6pnZjwgRG/l7t272LhxI2JiYvgEdzPjd5PeSu/evdGoUSOMHTtWdBSrw5WTKq1nz564efMmdxGTCFdOqpTRo0cjMTERZ8+ehZubm+g4VokrJ1WYTqdDbGwsjhw5Am9vb9FxrBb355T5fDlmSE1NhZ2dHZo3by5kvijcn7MI0fsyip4vtwxHjx41ValS5a3326zsfJG4PyfJ2qRJk9C9e3e4u7uLjmL1+DsnlVtiYiJu3LiBuLg40VFUgSsnlVtQUBA6d+6MJk2aiI6iClw5qVyysrJw5coVboJrQVw5qUwHDx7EnTt3sGzZMr4czIJYTipVdHQ0Bg0ahNq1a2POnDmi46gKH9ZSicaPH4+YmBhERETAw8NDdBzVYTnpDQaDAf/73//w559/4sCBAxg4cCB0Op3oWKrDctJrsrOz4e3tjdzcXKSmpqJFixaiI6kWf+ekl1JSUuDu7g5bW1vcuXOHxRSM5SQAwK5du+Dj44POnTvj5s2bcHJyEh1J9VhOQlJSEkaMGIFJkyYhPj6e72ggE7wXCBMnToSPjw+++eYb0VHoFTwgREhJScH+/ftFx6AiuHKqXFJSEoxGI/r16yc6ChXBcqpcbGws6tSpw98zZYj3iMolJCTAy8tLdAwqBsupctevX+dDWpliOVUsPz8fubm5GDNmjOgoVAyWU8V27tyJqlWr8q0tZYrlVLH9+/fD09NTdAwqAcupYklJSejSpYvoGFQCllPF7ty5g6FDh4qOQSVgOVXq8uXLMBgMGDhwoOgoVAKWU6W2b98OFxcXaLV8BqdcsZwqFR8fzzfrkjmWU6WuXLmCPn36iI5BpWA5VaigoAA5OTl88oHMsZwqtG/fPtjZ2fGd22WO5VShQ4cOoX79+qJjUBlYThXKzs7mUVoF4Oa5Mp8vRYZHjx4hLS0N7dq1K3dJRX8fRM+XMgM3z1XofKky1K9f3zRixAihGSpC9HwpM4Cb59KrvvjiC+zduxcFBQWio1AJWE6VCgwMhJ2dHRYtWiQ6CpWA5VQpjUaDgIAArFmzRnQUKgHLqWIRERHQ6/WYNWuW6ChUDJZTxapVq4Y1a9YgOjoaFy5cEB2HimA5VW78+PHw8/NDt27dkJmZKToOvYLlJPz+++9o0KAB2rVrB71eLzoO/YflJGg0Gly8eBEGgwGdOnWC0WgUHYnActJ/HBwccP78eaSlpWHIkCGi4xBYTnqFu7s7dDodfvvtNx7BlQGWk17z3nvvYcuWLVi5ciXWrVsnOo6q8aUJ9IZRo0bhxo0bmDp1Kpo2bcp3TBCE5aRiLViwANevX8fAgQORnJwsOo4q8WEtlWjjxo3o2bMnfHx8kJ2dLTqO6rCcVKpDhw5h+vTpuHXrFlq1aoVDhw6JjqQaLCeVacWKFWjTpg1cXV0xYMAANGzYEOvXrxcdy+qxnFQu1apVw7Fjx5CZmYmuXbsiKCgItra2aN68OSZOnIjExETREa0Oy0kV4urqip07d+Lp06eIjY1F+/btcfDgQfj6+sLBwQHDhg3DjRs3RMe0CiwnVYpWq8XQoUMRGxuLjIwM6PV6LFq0CBcuXEDLli2xfPly0REVj+Uks7Czs0NoaCjS09OxdOlShIaGYtq0aaJjKRrLSWYXGhqKn3/+GWvWrMGZM2dEx1EslpMkMXz4cPj6+mLEiBGioygWy0mS2bVrFzIzMzF48GC+DK0SWE6STO3atZGQkID4+Hg4OTlh5MiROHfunOhYisFykqS6deuGnJwcTJ8+HUlJSejYsSOcnJwwa9YsrqZlYDlJcnZ2dggPD0d6ejpycnIQHByMtWvXwtPTk8/ZLQXLSRbl5OSEL7/8Enfv3kVhYSHat2/PFbQELCcJUatWLaSmpiIzMxNLly4VHUeWWE4SxtnZGcHBwViyZAkMBoPoOLLDcpJQERER0Gq1CAgIEB1FdlhOEkqj0WDDhg3YtGkTOnbsiLNnz4qOJBvcPFfm89WS4cmTJ/jnn3+Qn5+PmjVromnTphadXx7cPLcI0Zumip6vtgwnTpww2dvbm1q3bm16/vy5xeeXhpvnkqp1794d169fx61btzBz5kzRcYRiOUl23NzcEBISgvXr16v6b6AsJ8nSwoULUVhYqOoXbbOcJEtarRZjxoxBRESE6CjCsJwkW6tWrcLDhw+xc+dO0VGEYDlJtpycnNC7d2/MnTtXdBQhWE6StW+//Ra3bt3CkydPREexOJaTZM3T0xNdunRBenq66CgWx3KS7B04cADPnz9HaGio6CgWxXKS7NWqVQsNGzbEypUrodfrRcexGJaTFKFu3bpwdHTE5MmTRUexGJaTFGPx4sXYtGkT8vPzRUexCJaTFGPatGlwcnJCUFCQ6CgWwXKSooSHh2Pr1q3Izc0VHUVyLCcpSmBgIN555x0EBgaKjiI5lpMUJzIyErGxsVb/uyfLSYozbtw42NnZYevWraKjSIrlJEVq2bIlduzYITqGpFhOUqQBAwYgKSlJdAxJsZykSBMmTEB2djby8vJER5EMy0mK5OnpCXt7e/z444+io0iG5STFat26NXbv3i06hmRYTlKsDz/8EOfPnxcdQzIsJynWhAkTkJOTg5ycHNFRJMFykmI1aNAAjo6OiImJER1FEiwnKZqXlxf27t0rOoYkWE5StCFDhuDixYuiY0iC5SRFCwgIQG5uLv7991/RUcyO5SRFc3V1RY0aNbB+/XrRUcyO5STFa9++Pfbt2yc6htlxf06Zz2eGsuffv38fGRkZ8PHxEZbhbXB/ToXOZ4ay5xcWFprs7e1N8+bNE5bhbYD7c5K10mg0mDJlCqKjo1FQUCA6jtmwnGQVli5dCltbW/Tt21d0FLNhOckqaLVaHD9+HCdOnMDixYtFxzELlpOshre3N6Kjo7Fo0SKcPHlSdJy3xnKSVZk6dSr69u2Ljz76SHSUt8ZyktXZvn07Hjx4gNjYWNFR3grLSVbH2dkZvr6++Oqrr0RHeSssJ1mlRYsW4dKlSzAYDKKjVBrLSVapR48esLGxQVxcnOgolcZyktVyd3dX9O+dLCdZrW7duuHUqVOiY1Qay0lWa/To0cjIyBAdo9JYTrJa/fv3R2FhIc6ePSs6SqWwnGS1NBoN6tati02bNomOUiksJ1m1Dh06QPRrkiuL5SSrNnLkSFy7dg1Go1F0lApjOcmq+fv7AwA2bNggOEnFsZxk1TQaDXr06KHIl5GxnGT1Nm/ejHv37iEkJER0lAphOcnqubq6Yt26dVi5cqWi9lVhOUkVAgIC4OLiguDgYNFRyo3lJNWYM2cOtm/frphXqrCcpBohISHQaDSKeZ0ny0mqodFoMHbsWERHR4uOUi4sJ6lKVFQUcnJyoNPpREcpE8tJqlK9enV4eHgoYsNdlpNUx8/PjysnkRxNmDABd+7ckf1RW5aTVKdr167QarXYvXu36CilYjlJlZo1a4YtW7aIjlEqlpNUqX///jh9+rToGKXi5rkyn88M0szX6/W4fPkyOnToICzDC9w8V6HzmUGa+YWFhSYApvT0dGEZXgA3zyX6fxqNBvb29khMTBQdpUQsJ6lWzZo1cenSJdExSsRykmqZTCbY2NiIjlEilpNU69mzZ6hRo4boGCViOUmVfvnlFzx69AgjR44UHaVELCdZnW3btsHLywsuLi5o0qQJlixZ8tpbY8bFxWH48OH47LPP0KRJE4FJS8dyktUwGAzo0KED/P39Ua9ePcycOROdOnXCV199BUdHR3h7e8Pd3R2DBw/GJ598gh9++EF05FJpRQcgMge9Xg8vLy9kZWUhPT0dHh4eL68zGAxYvXo1EhISUKtWLXz++edo0aKFwLTlw3KS4hkMBrRs2RKPHz9GWloaXF1dX7teq9VixowZmDFjhqCElVPm0/csmIVIrW6bTKZGRS8stZxEJA4PCBHJFMtJJFMsJ5FMsZxEMsVyEsnU/wEdn4JyqwRzFgAAAABJRU5ErkJggg==\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"import cartopy.crs as ccrs\n",
"import cartopy.feature as cf\n",
"import matplotlib.pyplot as plt\n",
"import numpy as np\n",
"\n",
"\n",
"# proj = ccrs.LambertConformal(central_latitude = 25, \n",
"# central_longitude = 265, \n",
"# standard_parallels = (25, 25))\n",
"\n",
"# proj = ccrs.PlateCarree()\n",
"proj = ccrs.PlateCarree()\n",
"\n",
"# Data and coordinates (from download link above)\n",
"# with np.load('nam_218_20120414_1200_006.npz') as nam:\n",
"# dat = nam['dpc']\n",
"# lat = nam['lat']\n",
"# lon = nam['lon']\n",
"\n",
"ax = plt.axes(projection = proj)\n",
"# ax.pcolormesh(lon, lat, dat, transform = ccrs.PlateCarree())\n",
"# ax.add_feature(cf.NaturalEarthFeature(\n",
"# category='cultural',\n",
"# name='admin_1_states_provinces_lines',\n",
"# scale='10m',\n",
"# facecolor='none'))\n",
"\n",
"\n",
"ax.coastlines('10m')\n",
"# ax.add_feature(cf.BORDERS)\n",
"ax.add_feature(cf.STATES.with_scale('10m'))\n",
"\n",
"# The extent bounds are specified as an array [[x0, y0], [x1, y1]], \n",
"# where x0 is the left side of the extent, y0 is the top, x1 is the right and y1 is the bottom.\n",
"# extent (x0, x1, y0, y1)\n",
"# extent = [-180,180, -90,90] # world\n",
"# extent = [-90, -30, 20, -60] # south america\n",
"# extent = [-74, -31, 5.5, -33] # brazil\n",
"# extent = [-53.5, -45, -11, -20] # brazil\n",
"extent = [-48, -47, -16, -17] # brazil\n",
"# extent = [-90, -30, 20, -60]\n",
"# extent = [-100, 30, 0, 80]\n",
"# extent = [-87.35, -79.5, 24.1, 30.8]\n",
"\n",
"\n",
"ax.set_extent(extent, crs=ccrs.PlateCarree())\n",
"\n",
"ax.gridlines()\n",
"\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"\n",
"## ds ##:\n",
"\n",
"<osgeo.gdal.Dataset; proxy of <Swig Object of type 'GDALDatasetShadow *' at 0x7f00185e8840> >\n",
"\n",
"\n",
"## data ##:\n",
"\n",
"[[[28 28 28 ... 29 28 28]\n",
" [27 29 29 ... 28 28 29]\n",
" [29 32 31 ... 31 32 35]\n",
" ...\n",
" [47 47 46 ... 35 34 34]\n",
" [46 46 48 ... 33 33 35]\n",
" [46 46 48 ... 33 34 35]]\n",
"\n",
" [[46 47 50 ... 51 51 51]\n",
" [46 47 49 ... 51 52 51]\n",
" [48 50 52 ... 51 52 53]\n",
" ...\n",
" [58 57 57 ... 58 57 57]\n",
" [57 57 57 ... 57 57 58]\n",
" [58 57 57 ... 57 58 58]]\n",
"\n",
" [[53 54 54 ... 55 56 55]\n",
" [54 55 55 ... 56 55 55]\n",
" [55 57 57 ... 55 55 57]\n",
" ...\n",
" [62 61 61 ... 59 59 59]\n",
" [61 61 61 ... 59 58 59]\n",
" [62 61 60 ... 59 58 58]]]\n",
"\n",
"\n",
"## gt ##:\n",
"\n",
"(885110.0, 10.0, 0.0, 8176910.0, 0.0, -10.0)\n",
"\n",
"\n",
"## proj ##:\n",
"\n",
"PROJCS[\"WGS 84 / UTM zone 22S\",GEOGCS[\"WGS 84\",DATUM[\"WGS_1984\",SPHEROID[\"WGS 84\",6378137,298.257223563,AUTHORITY[\"EPSG\",\"7030\"]],AUTHORITY[\"EPSG\",\"6326\"]],PRIMEM[\"Greenwich\",0,AUTHORITY[\"EPSG\",\"8901\"]],UNIT[\"degree\",0.0174532925199433,AUTHORITY[\"EPSG\",\"9122\"]],AUTHORITY[\"EPSG\",\"4326\"]],PROJECTION[\"Transverse_Mercator\"],PARAMETER[\"latitude_of_origin\",0],PARAMETER[\"central_meridian\",-51],PARAMETER[\"scale_factor\",0.9996],PARAMETER[\"false_easting\",500000],PARAMETER[\"false_northing\",10000000],UNIT[\"metre\",1,AUTHORITY[\"EPSG\",\"9001\"]],AXIS[\"Easting\",EAST],AXIS[\"Northing\",NORTH],AUTHORITY[\"EPSG\",\"32722\"]]\n",
"\n",
"\n",
"## inproj ##:\n",
"\n",
"PROJCS[\"WGS 84 / UTM zone 22S\",\n",
" GEOGCS[\"WGS 84\",\n",
" DATUM[\"WGS_1984\",\n",
" SPHEROID[\"WGS 84\",6378137,298.257223563,\n",
" AUTHORITY[\"EPSG\",\"7030\"]],\n",
" AUTHORITY[\"EPSG\",\"6326\"]],\n",
" PRIMEM[\"Greenwich\",0,\n",
" AUTHORITY[\"EPSG\",\"8901\"]],\n",
" UNIT[\"degree\",0.0174532925199433,\n",
" AUTHORITY[\"EPSG\",\"9122\"]],\n",
" AUTHORITY[\"EPSG\",\"4326\"]],\n",
" PROJECTION[\"Transverse_Mercator\"],\n",
" PARAMETER[\"latitude_of_origin\",0],\n",
" PARAMETER[\"central_meridian\",-51],\n",
" PARAMETER[\"scale_factor\",0.9996],\n",
" PARAMETER[\"false_easting\",500000],\n",
" PARAMETER[\"false_northing\",10000000],\n",
" UNIT[\"metre\",1,\n",
" AUTHORITY[\"EPSG\",\"9001\"]],\n",
" AXIS[\"Easting\",EAST],\n",
" AXIS[\"Northing\",NORTH],\n",
" AUTHORITY[\"EPSG\",\"32722\"]]\n"
]
}
],
"source": [
"# First: read the geotiff image with GDAL.\n",
"from osgeo import gdal, osr\n",
"\n",
"gdal.UseExceptions()\n",
"\n",
"\n",
"fname = '/notebooks/resources/T22KHG_20190425T132241_TCI_smaller.tif'\n",
"\n",
"ds = gdal.Open(fname)\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": []
}
],
"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
}