Path: 07. Jupyter-ostalo / basemap.ipynb

Views: ⁹⁴⁰

Kernel: Python 3 (Anaconda)

Basemap

Basemap je dio Matplotlib biblioteke s funkcijama i metodama za pravljenje geografskih karata. Korištenje Basemapa ilustrirat ćemo na nekim primjerima.

In [19]:

# osnovne biblioteke koje ćemo koristiti
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
%matplotlib inline
import numpy as np

Crtanje najkraće udaljenosti

In [2]:

fig = plt.figure(figsize=(14,10))
ax = fig.add_axes([0.1,0.1,0.9,0.9])
# Mercatorova projekcija
m = Basemap(projection='merc',llcrnrlat=-80,urcrnrlat=80,\
            llcrnrlon=-180,urcrnrlon=180,lat_ts=20,resolution='c')
# geografske koordinate Zagreba
zaglat = 45.48; zaglon = 15.58
# geografske koordinate Pekinga
pelat = 39.55; pelon = 116.25

m.drawgreatcircle(zaglon,zaglat,pelon,pelat,linewidth=2,color='b')
m.drawcoastlines()
m.fillcontinents(color='coral',lake_color='aqua')
m.drawmapboundary(fill_color='aqua')
m.drawparallels(np.arange(-90.,91.,30.),labels=[1,1,0,1])
m.drawmeridians(np.arange(-180.,181.,60.),labels=[1,1,0,1])
zx,zy=m(zaglon,zaglat)
px,py=m(pelon,pelat)
ax.text(zx,zy,'Zagreb',fontsize=14,fontweight='bold',ha='center',va='center',color='b')
ax.text(px,py,'Peking',fontsize=14,fontweight='bold',ha='center',va='center',color='r')
ax.set_title('Zagreb - Peking');

/projects/anaconda3/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:3327: MatplotlibDeprecationWarning: The ishold function was deprecated in version 2.0.
  b = ax.ishold()
/projects/anaconda3/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:3336: MatplotlibDeprecationWarning: axes.hold is deprecated.
    See the API Changes document (http://matplotlib.org/api/api_changes.html)
    for more details.
  ax.hold(b)
/projects/anaconda3/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:1800: MatplotlibDeprecationWarning: The get_axis_bgcolor function was deprecated in version 2.0. Use get_facecolor instead.
  axisbgc = ax.get_axis_bgcolor()

Dan i noć

In [3]:

from datetime import datetime
# Millerova projekcija
fig = plt.figure(figsize=(8,8))
ax = fig.add_axes([0,0,1,1])
map = Basemap(projection='mill',lon_0=90)

map.drawcoastlines()
map.drawparallels(np.arange(-90,90,30),labels=[1,0,0,0])
map.drawmeridians(np.arange(map.lonmin,map.lonmax+30,60),labels=[0,0,0,1])

map.drawmapboundary(fill_color='aqua')
map.fillcontinents(color='coral',lake_color='aqua')

date = datetime.utcnow()
CS=map.nightshade(date)
ax.set_title(u'Mapa dana/noći za %s (UTC)' % date.strftime("%d %b %Y %H:%M:%S"));

/projects/anaconda3/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:1800: MatplotlibDeprecationWarning: The get_axis_bgcolor function was deprecated in version 2.0. Use get_facecolor instead.
  axisbgc = ax.get_axis_bgcolor()
/projects/anaconda3/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:3683: MatplotlibDeprecationWarning: The ishold function was deprecated in version 2.0.
  b = ax.ishold()
/projects/anaconda3/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:3719: VisibleDeprecationWarning: using a non-integer number instead of an integer will result in an error in the future
  xx = x[x.shape[0]/2,:]
/projects/anaconda3/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:3752: MatplotlibDeprecationWarning: axes.hold is deprecated.
    See the API Changes document (http://matplotlib.org/api/api_changes.html)
    for more details.
  ax.hold(b)

Crtanje potresa

In [4]:

# Standardne biblioteke
from dateutil.tz import tzutc, tzlocal

# Ostalo što nam treba za crtanje
from matplotlib.colors import Normalize
from matplotlib.cm import get_cmap

Podaci se dohvaćaju sa stranice earthquake.usgs.gov. U datoteku query.csv spremamo podatke o zemljotresima u zadnjih tjedan dana.

In [85]:

url = "https://earthquake.usgs.gov/fdsnws/event/1/query.csv"
import time
import datetime as dt
today = dt.datetime.now()
now = today - dt.timedelta(hours=1)
endtime= now.isoformat()
starttime = (now - dt.timedelta(days=7)).isoformat()
params = dict(minmagnitude='2.5',orderby='time',starttime=starttime,endtime=endtime)

In [86]:

import requests
r = requests.get(url, params = params)
with open('query.csv','w') as file:
    file.write(r.text)

In [87]:

# Parsiranje
import csv

lat, lon, mag = [], [], []
with open('query.csv', 'r') as csvfile:
    reader = csv.DictReader(csvfile)
    for redak in reader:
        lon.append(float(redak['longitude']))
        lat.append(float(redak['latitude']))
        mag.append(float(redak['mag']))

In [88]:

# Normalizacija magnituda
norm = Normalize()
mag_norms = norm(np.array(mag))
z = (mag_norms * 10.0)**2.0

Crtanje

In [89]:

fig = plt.figure(figsize=(14,10))


m = Basemap(projection='moll', lat_0=0, lon_0=0, resolution='l', area_thresh=1000.0)
m.drawcoastlines()
m.fillcontinents(color='coral',lake_color='aqua')
m.drawmapboundary(fill_color='aqua')
m.drawmeridians(np.arange(0, 360, 30))
m.drawparallels(np.arange(-90, 90, 30))
cmap = get_cmap('jet')
x,y = m(lon, lat)
sc = m.scatter(x, y, s=z, cmap=cmap, c=mag)
plt.colorbar(sc, orientation='horizontal')

plt.title("Potresi");

/projects/anaconda3/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:1656: MatplotlibDeprecationWarning: The get_axis_bgcolor function was deprecated in version 2.0. Use get_facecolor instead.
  fill_color = ax.get_axis_bgcolor()
/projects/anaconda3/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:1800: MatplotlibDeprecationWarning: The get_axis_bgcolor function was deprecated in version 2.0. Use get_facecolor instead.
  axisbgc = ax.get_axis_bgcolor()
/projects/anaconda3/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:3289: MatplotlibDeprecationWarning: The ishold function was deprecated in version 2.0.
  b = ax.ishold()
/projects/anaconda3/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:3298: MatplotlibDeprecationWarning: axes.hold is deprecated.
    See the API Changes document (http://matplotlib.org/api/api_changes.html)
    for more details.
  ax.hold(b)

Crtanje SST (sea surface temperature)

Podaci se preuzimaju sa stranice CISL Research Data Archive (http://rda.ucar.edu/). Podaci su u netCDF formatu, koji je u verziji 4 u biti podskup HDF formata (Hierarchical Data Format). Opis formata možete naći ovdje.

Opis konkretnih podataka je ovdje, a sami podaci se nalaze ovdje.

In [90]:

from netCDF4 import Dataset
#f = Dataset('avhrr-only-v2.20150324.nc')
f = Dataset('avhrr-only-v2.20170330.nc')

In [91]:

# Kako su strukturirani podaci 
print (f)

<class 'netCDF4._netCDF4.Dataset'>
root group (NETCDF3_CLASSIC data model, file format NETCDF3):
    Conventions: CF-1.0
    title: Daily-OI-V2, Final, Data (Ship, Buoy, AVHRR: NOAA19, METOP, NCEP-ice)
    History: Version 2.0
    creation_date: 2017-04-14 11:10
    Description: Reynolds, et al.(2007) Daily High-resolution Blended Analyses. Available at ftp://eclipse.ncdc.noaa.gov/pub/OI-daily/daily-sst.pdf  Climatology is based on 1971-2000 OI.v2 SST, Satellite data: Navy  NOAA19 METOP AVHRR, Ice data: NCEP ice
    Source: NOAA/National Climatic Data Center
    Contact: Dick Reynolds, email: [email protected] & Chunying Liu, email: [email protected]
    dimensions(sizes): time(1), zlev(1), lat(720), lon(1440)
    variables(dimensions): float32 time(time), float32 zlev(zlev), float32 lat(lat), float32 lon(lon), int16 sst(time,zlev,lat,lon), int16 anom(time,zlev,lat,lon), int16 err(time,zlev,lat,lon), int16 ice(time,zlev,lat,lon)
    groups: 

In [92]:

print(f.variables.keys()) # imena varijabli
temp = f.variables['sst']  # temperatura
print(temp)

odict_keys(['time', 'zlev', 'lat', 'lon', 'sst', 'anom', 'err', 'ice'])
<class 'netCDF4._netCDF4.Variable'>
int16 sst(time, zlev, lat, lon)
    long_name: Daily sea surface temperature
    units: degrees C
    _FillValue: -999
    add_offset: 0.0
    scale_factor: 0.01
    valid_min: -300
    valid_max: 4500
unlimited dimensions: 
current shape = (1, 1, 720, 1440)
filling off

In [93]:

temp.dimensions

('time', 'zlev', 'lat', 'lon')

In [94]:

data = temp[0,0]
data.shape

(720, 1440)

In [95]:

lonvals = f.variables['lon'][:]
latvals = f.variables['lat'][:]

In [96]:

plt.rcParams['figure.figsize'] = (12.0, 8.0)
plt.contourf(lonvals, latvals, data, 20, cmap=plt.get_cmap('YlGnBu_r'))
plt.colorbar()
plt.show()

Sad ćemo nacrati iste podatke na zemaljskoj kugli.

In [97]:

m = Basemap(projection='ortho', lon_0=-50, lat_0=40, resolution='l')

In [98]:

X,Y = np.meshgrid(lonvals, latvals)
x,y = m(X,Y)

In [99]:

pc = m.contourf(x, y, data, 30, cmap=plt.get_cmap('YlGnBu_r'))
m.bluemarble()
m.drawmapboundary()
m.drawcoastlines()
plt.title('Ortografska projekcija')
plt.colorbar(pc, orientation='vertical')
plt.show()

/projects/anaconda3/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:3683: MatplotlibDeprecationWarning: The ishold function was deprecated in version 2.0.
  b = ax.ishold()
/projects/anaconda3/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:3752: MatplotlibDeprecationWarning: axes.hold is deprecated.
    See the API Changes document (http://matplotlib.org/api/api_changes.html)
    for more details.
  ax.hold(b)
/projects/anaconda3/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:1656: MatplotlibDeprecationWarning: The get_axis_bgcolor function was deprecated in version 2.0. Use get_facecolor instead.
  fill_color = ax.get_axis_bgcolor()
/projects/anaconda3/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:3363: MatplotlibDeprecationWarning: The ishold function was deprecated in version 2.0.
  b = ax.ishold()
/projects/anaconda3/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py:3372: MatplotlibDeprecationWarning: axes.hold is deprecated.
    See the API Changes document (http://matplotlib.org/api/api_changes.html)
    for more details.
  ax.hold(b)

In [18]:

from verzije import *
from IPython.display import HTML
HTML(print_sysinfo()+info_packages('matplotlib,basemap'))

Python verzija	3.5.3
kompajler	GCC 4.8.2 20140120 (Red Hat 4.8.2-15)
sustav	Linux
broj CPU-a	8
interpreter	64bit

matplotlib verzija	2.0.0
basemap verzija	1.0.8