Pythonで世界地図(matplotlib basemap toolkit翻訳)-4

##Plotting data on a map (Example Gallery)
https://matplotlib.org/basemap/users/examples.html

注意）
データのURLにアクセスしてもなかったり、データが新しくなりファイル名が変わっている可能性があります。

次に、Basemapインスタンスメソッドを使用してデータをマップにプロットする方法を示す一連の例を示します。より多くの例は、ベースマップソース配布のexamplesディレクトリに含まれています。データをプロットするためのいくつかのBasemapインスタンスメソッドがあります：


   ・ contour(): 等高線を描きます。
   ・ contourf(): 塗りつぶされた輪郭を描く
   ・ imshow(): 画像を描画する。
   ・ pcolor(): 疑似カラープロットを描きます。
   ・ pcolormesh(): pseudocolorプロット（通常のメッシュの場合はより高速版）を描画します。
   ・ plot(): 線および/またはマーカーを描画する。
   ・ scatter(): マーカーで点を描く。
   ・ quiver(): 描画ベクトル。
   ・ barbs(): 風の棘を引く。（風向・風速の記号）
   ・ drawgreatcircle(): 大きな円を描きなさい。（大圏コース）

これらのインスタンスメソッドの多くは、対応するmatplotlib Axesインスタンスメソッドに単純に転送し、いくつかの追加のプリ/ポスト処理と引数のチェックを行います。また、ベースマップに関連付けられたAxesインスタンスを使用して、matplotlib pyplotインターフェイスまたはOO APIを使用してマップに直接プロットすることもできます。

Basemapインスタンスメソッドの使い方の詳細については、「Matplotlib Basemap Toolkit API」を参照してください。

以下に例を示します（その多くは、httpを介してデータセットを取得するためにnetcdf4-pythonモジュールを利用しています）。

    ベースマップ上に等高線を描画する

```python:
from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np
# 50N、100Wを見下ろす衛星の視点を用いて正射投影図投影を設定します。
# 低解像度の海岸線を使用する。

map = Basemap(projection='ortho',lat_0=45,lon_0=-100,resolution='l')
# 海岸線を描く、国境を塗る、大陸を塗る
map.drawcoastlines(linewidth=0.25)
map.drawcountries(linewidth=0.25)
map.fillcontinents(color='coral',lake_color='aqua')
# 地図投影領域（投影肢）の端を描画し、
map.drawmapboundary(fill_color='aqua')
# 緯度/経度グリッド線を30度ごとに描画します。
map.drawmeridians(np.arange(0,360,30))
map.drawparallels(np.arange(-90,90,30))
# 定期的な緯度/経度グリッドで何らかのデータを作ります。
nlats = 73; nlons = 145; delta = 2.*np.pi/(nlons-1)
lats = (0.5*np.pi-delta*np.indices((nlats,nlons))[0,:,:])
lons = (delta*np.indices((nlats,nlons))[1,:,:])
wave = 0.75*(np.sin(2.*lats)**8*np.cos(4.*lons))
mean = 0.5*np.cos(2.*lats)*((np.sin(2.*lats))**2 + 2.)
# lat / lonグリッドのネイティブマップ投影座標を計算します。
x, y = map(lons*180./np.pi, lats*180./np.pi)
# 地図上の等高線データ。
cs = map.contour(x,y,wave+mean,15,linewidths=1.5)
plt.title('contour lines over filled continent background')
plt.show()
```

![image.png](https://qiita-image-store.s3.amazonaws.com/0/238800/6dba7abe-1539-9010-e069-7770a023a58e.png)

塗りつぶした等高線で降水量をプロット

```python:
from mpl_toolkits.basemap import Basemap, cm
# requires netcdf4-python (netcdf4-python.googlecode.com)
from netCDF4 import Dataset as NetCDFFile
import numpy as np
import matplotlib.pyplot as plt

# NWSから使用される特別な降水量カラーマップを使用して
# NWSから降雨をプロットし、ベースマップに含めます。

nc = NetCDFFile('../../../examples/nws_precip_conus_20061222.nc')
# data from http://water.weather.gov/precip/
prcpvar = nc.variables['amountofprecip']
data = 0.01*prcpvar[:]
latcorners = nc.variables['lat'][:]
loncorners = -nc.variables['lon'][:]
lon_0 = -nc.variables['true_lon'].getValue()
lat_0 = nc.variables['true_lat'].getValue()
# FigureインスタンスとAxesインスタンスを作成する
fig = plt.figure(figsize=(8,8))
ax = fig.add_axes([0.1,0.1,0.8,0.8])
# polar stereographic Basemapインスタンスを作成します。
m = Basemap(projection='stere',lon_0=lon_0,lat_0=90.,lat_ts=lat_0,\
            llcrnrlat=latcorners[0],urcrnrlat=latcorners[2],\
            llcrnrlon=loncorners[0],urcrnrlon=loncorners[2],\
            rsphere=6371200.,resolution='l',area_thresh=10000)
# 海岸線、国境、国境、地図の端を描画します。
m.drawcoastlines()
m.drawstates()
m.drawcountries()
# 緯度線を描く。
parallels = np.arange(0.,90,10.)
m.drawparallels(parallels,labels=[1,0,0,0],fontsize=10)
# 経線を描く
meridians = np.arange(180.,360.,10.)
m.drawmeridians(meridians,labels=[0,0,0,1],fontsize=10)
ny = data.shape[0]; nx = data.shape[1]
lons, lats = m.makegrid(nx, ny) # nxの緯度/経度をnxの等間隔の空間グリッドで取得します。
x, y = m(lons, lats) # map proj座標を計算します。
# 塗りつぶされた輪郭を描く
clevs = [0,1,2.5,5,7.5,10,15,20,30,40,50,70,100,150,200,250,300,400,500,600,750]
cs = m.contourf(x,y,data,clevs,cmap=cm.s3pcpn)
# add colorbar.
cbar = m.colorbar(cs,location='bottom',pad="5%")
cbar.set_label('mm')
# add title
plt.title(prcpvar.long_name+' for period ending '+prcpvar.dateofdata)
plt.show()
```

![image.png](https://qiita-image-store.s3.amazonaws.com/0/238800/b319e075-7450-3401-cc40-2cede840eee5.png)

標高の高い高低の海面気圧の天気図をプロットする

```python:
"""
海面の圧力マップ上にHとLをプロットする
(uses scipy.ndimage.filters and netcdf4-python)
"""
import numpy as np
import matplotlib.pyplot as plt
from datetime import datetime
from mpl_toolkits.basemap import Basemap, addcyclic
from scipy.ndimage.filters import minimum_filter, maximum_filter
from netCDF4 import Dataset

def extrema(mat,mode='wrap',window=10):
    """find the indices of local extrema (min and max)
    in the input array."""
    mn = minimum_filter(mat, size=window, mode=mode)
    mx = maximum_filter(mat, size=window, mode=mode)
    # (mat == mx) true if pixel is equal to the local max
    # (mat == mn) true if pixel is equal to the local in
    # Return the indices of the maxima, minima
    return np.nonzero(mat == mn), np.nonzero(mat == mx)

# plot 00 UTC today.
date = datetime.now().strftime('%Y%m%d')+'00'

# OpenDAPデータセットを開きます。
#data=Dataset("http://nomads.ncep.noaa.gov:9090/dods/gfs/gfs/%s/gfs_%sz_anl" %\
#        (date[0:8],date[8:10]))
data=Dataset("http://nomads.ncep.noaa.gov:9090/dods/gfs_hd/gfs_hd%s/gfs_hd_%sz"%\
        (date[0:8],date[8:10]))



# read lats,lons.
lats = data.variables['lat'][:]
lons1 = data.variables['lon'][:]
nlats = len(lats)
nlons = len(lons1)
# read prmsl, convert to hPa (mb).
prmsl = 0.01*data.variables['prmslmsl'][0]
# ウィンドウパラメータは、検出された高および低の数を制御する。
# (高い価値、少ない高値と低い値)
local_min, local_max = extrema(prmsl, mode='wrap', window=50)
# Basemapインスタンスを作成します。
m =\
Basemap(llcrnrlon=0,llcrnrlat=-80,urcrnrlon=360,urcrnrlat=80,projection='mill')
# 経度のラップアラウンドポイントを追加します。
prmsl, lons = addcyclic(prmsl, lons1)
# 輪郭レベル
clevs = np.arange(900,1100.,5.)
# 地図投影グリッドのx、yを見つける。
lons, lats = np.meshgrid(lons, lats)
x, y = m(lons, lats)
# 図を作成します。
fig=plt.figure(figsize=(8,4.5))
ax = fig.add_axes([0.05,0.05,0.9,0.85])
cs = m.contour(x,y,prmsl,clevs,colors='k',linewidths=1.)
m.drawcoastlines(linewidth=1.25)
m.fillcontinents(color='0.8')
m.drawparallels(np.arange(-80,81,20),labels=[1,1,0,0])
m.drawmeridians(np.arange(0,360,60),labels=[0,0,0,1])
xlows = x[local_min]; xhighs = x[local_max]
ylows = y[local_min]; yhighs = y[local_max]
lowvals = prmsl[local_min]; highvals = prmsl[local_max]
# 下の最小圧力値で、青Lとしてプロットを表示します。
xyplotted = []
# すでにdminメーター内にLまたはHがある場合はプロットしないでください。
yoffset = 0.022*(m.ymax-m.ymin)
dmin = yoffset
for x,y,p in zip(xlows, ylows, lowvals):
    if x < m.xmax and x > m.xmin and y < m.ymax and y > m.ymin:
        dist = [np.sqrt((x-x0)**2+(y-y0)**2) for x0,y0 in xyplotted]
        if not dist or min(dist) > dmin:
            plt.text(x,y,'L',fontsize=14,fontweight='bold',
                    ha='center',va='center',color='b')
            plt.text(x,y-yoffset,repr(int(p)),fontsize=9,
                    ha='center',va='top',color='b',
                    bbox = dict(boxstyle="square",ec='None',fc=(1,1,1,0.5)))
            xyplotted.append((x,y))
# 下の最大圧力値を持つ、赤色のHとしてプロットの高さ。
xyplotted = []
for x,y,p in zip(xhighs, yhighs, highvals):
    if x < m.xmax and x > m.xmin and y < m.ymax and y > m.ymin:
        dist = [np.sqrt((x-x0)**2+(y-y0)**2) for x0,y0 in xyplotted]
        if not dist or min(dist) > dmin:
            plt.text(x,y,'H',fontsize=14,fontweight='bold',
                    ha='center',va='center',color='r')
            plt.text(x,y-yoffset,repr(int(p)),fontsize=9,
                    ha='center',va='top',color='r',
                    bbox = dict(boxstyle="square",ec='None',fc=(1,1,1,0.5)))
            xyplotted.append((x,y))
plt.title('Mean Sea-Level Pressure (with Highs and Lows) %s' % date)
#          平均海面圧力（高低）
plt.show()
```

シェイプファイルからハリケーントラックをプロットする

```python:
"""
Cat 4または5に到達した嵐のためにAtlantic Hurricane Trackを描画します。
嵐がcat 4または5であるトラックの一部は赤く表示されます。
ESRI shapefile data from http://nationalatlas.gov/mld/huralll.html
"""
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
# ランバートコンフォーマルコニックマップ。
m = Basemap(llcrnrlon=-100.,llcrnrlat=0.,urcrnrlon=-20.,urcrnrlat=57.,
            projection='lcc',lat_1=20.,lat_2=40.,lon_0=-60.,
            resolution ='l',area_thresh=1000.)
# shapefileを読んでください。
shp_info = m.readshapefile('../../../examples/huralll020','hurrtracks',drawbounds=False)
# Cat 4に到達した嵐の名前を見つけます。
names = []
for shapedict in m.hurrtracks_info:
    cat = shapedict['CATEGORY']
    name = shapedict['NAME']
    if cat in ['H4','H5'] and name not in names:
        # 名前のついた嵐だけを使う。
        if name != 'NOT NAMED':  names.append(name)
# それらの嵐のプロットトラック。
for shapedict,shape in zip(m.hurrtracks_info,m.hurrtracks):
    name = shapedict['NAME']
    cat = shapedict['CATEGORY']
    if name in names:
        xx,yy = zip(*shape)
        # storm> Cat 4のトラックの一部を濃い赤色で表示します。
        if cat in ['H4','H5']:
            m.plot(xx,yy,linewidth=1.5,color='r')
        elif cat in ['H1','H2','H3']:
            m.plot(xx,yy,color='k')
# 海岸線、子午線および平行線を描画します。
m.drawcoastlines()
m.drawcountries()
m.drawmapboundary(fill_color='#99ffff')
m.fillcontinents(color='#cc9966',lake_color='#99ffff')
m.drawparallels(np.arange(10,70,20),labels=[1,1,0,0])
m.drawmeridians(np.arange(-100,0,20),labels=[0,0,0,1])
plt.title('Atlantic Hurricane Tracks (Storms Reaching Category 4, 1851-2004)')
plt.show()
```

![image.png](https://qiita-image-store.s3.amazonaws.com/0/238800/de698f60-9578-7d90-cf66-5c16dd92e03c.png)

エトポ5トポグラフィ/地形データを画像としてプロットする（指定された光源からのシェーディングの有無を問わない）。

```python:
from mpl_toolkits.basemap import Basemap, shiftgrid, cm
import numpy as np
import matplotlib.pyplot as plt
from netCDF4 import Dataset

# etopo5地形/地形を読んでください。
url = 'http://ferret.pmel.noaa.gov/thredds/catalog/data/PMEL/etopo5.nc'
etopodata = Dataset(url)

topoin = etopodata.variables['ROSE'][:]
lons = etopodata.variables['ETOPO05_X'][:]
lats = etopodata.variables['ETOPO05_Y'][:]
# lonが20から380の代わりに-180から180になるようにデータをシフトします。
topoin,lons = shiftgrid(180.,topoin,lons,start=False)

# 地形/地形を画像としてプロットする。

# FigureインスタンスとAxesインスタンスを作成します。
fig = plt.figure()
ax = fig.add_axes([0.1,0.1,0.8,0.8])
# ベースマップの設定（ 'lcc' = lambert等角錐）。
# WGS84楕円体の大小の球の半径を使用します。
m = Basemap(llcrnrlon=-145.5,llcrnrlat=1.,urcrnrlon=-2.566,urcrnrlat=46.352,\
            rsphere=(6378137.00,6356752.3142),\
            resolution='l',area_thresh=1000.,projection='lcc',\
            lat_1=50.,lon_0=-107.,ax=ax)
# nx x ny規則的に間隔を置いた5kmのネイティブ投影格子に変換する
nx = int((m.xmax-m.xmin)/5000.)+1; ny = int((m.ymax-m.ymin)/5000.)+1
topodat = m.transform_scalar(topoin,lons,lats,nx,ny)
# imshowを使って地図上にプロット画像を作成する。
im = m.imshow(topodat,cm.GMT_haxby)
# 海岸線と政治的境界を描く。
m.drawcoastlines()
m.drawcountries()
m.drawstates()
# 緯度と経度を描きます。
# 地図の左と下のラベル。
parallels = np.arange(0.,80,20.)
m.drawparallels(parallels,labels=[1,0,0,1])
meridians = np.arange(10.,360.,30.)
m.drawmeridians(meridians,labels=[1,0,0,1])
# add colorbar
cb = m.colorbar(im,"right", size="5%", pad='2%')
ax.set_title('ETOPO5 Topography - Lambert Conformal Conic')
plt.show()

# 影付きの救済プロットを作成します。

# 新しいFigure、axesインスタンスを作成します。
fig = plt.figure()
ax = fig.add_axes([0.1,0.1,0.8,0.8])
# 既存のBasemapインスタンスに新しいAxesイメージを添付します。
m.ax = ax
# 光源オブジェクトを作成します。
from matplotlib.colors import LightSource
ls = LightSource(azdeg = 90, altdeg = 20)
# 光源からの陰影を含むデータをRGB配列に変換します。
# (カラーマップを指定する必要があります)
rgb = ls.shade(topodat, cm.GMT_haxby)
im = m.imshow(rgb)
# 海岸線と政治的境界を描く。
m.drawcoastlines()
m.drawcountries()
m.drawstates()
ax.set_title('Shaded ETOPO5 Topography - Lambert Conformal Conic')
plt.show()
```


    ARGOの位置にマーカーをプロットします。


```python:
from netCDF4 import Dataset, num2date
import time, calendar, datetime, numpy
from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import urllib, os
# フォームからダウンロードしたデータ
# http://coastwatch.pfeg.noaa.gov/erddap/tabledap/apdrcArgoAll.html
filename, headers = urllib.urlretrieve('http://coastwatch.pfeg.noaa.gov/erddap/tabledap/apdrcArgoAll.nc?longitude,latitude,time&longitude>=0&longitude<=360&latitude>=-90&latitude<=90&time>=2010-01-01&time<=2010-01-08&distinct()')
dset = Dataset(filename)
lats = dset.variables['latitude'][:]
lons = dset.variables['longitude'][:]
time = dset.variables['time']
times = time[:]
t1 = times.min(); t2 = times.max()
date1 = num2date(t1, units=time.units)
date2 = num2date(t2, units=time.units)
dset.close()
os.remove(filename)
# 浮動小数点数の位置を示すマーカー付きマップ描画
m = Basemap(projection='hammer',lon_0=180)
x, y = m(lons,lats)
m.drawmapboundary(fill_color='#99ffff')
m.fillcontinents(color='#cc9966',lake_color='#99ffff')
m.scatter(x,y,3,marker='o',color='k')
plt.title('Locations of %s ARGO floats active between %s and %s' %\
        (len(lats),date1,date2),fontsize=12)
plt.show()
```

![image.png](https://qiita-image-store.s3.amazonaws.com/0/238800/1ead3ffc-bf0b-86e7-ee55-a4a6012dcf23.png)

SSTの擬似カラープロットと海氷解析。

```python:
from mpl_toolkits.basemap import Basemap
from netCDF4 import Dataset, date2index
import numpy as np
import matplotlib.pyplot as plt
from datetime import datetime
date = datetime(2007,12,15,0) # date to plot.
# open dataset.
dataset = \
Dataset('http://www.ncdc.noaa.gov/thredds/dodsC/OISST-V2-AVHRR_agg')
timevar = dataset.variables['time']
timeindex = date2index(date,timevar) # 希望の日付の時間インデックスを見つける。
# sstを読んでください。 missing_value変数属性を使用してマスクされた配列を自動的に作成します。 
# シングルトン次元を「絞り出す」。
sst = dataset.variables['sst'][timeindex,:].squeeze()
# 氷を読む。
ice = dataset.variables['ice'][timeindex,:].squeeze()
# （グリッドボックスの中心を表す）latとlonを読む。
lats = dataset.variables['lat'][:]
lons = dataset.variables['lon'][:]
lons, lats = np.meshgrid(lons,lats)
# figure、axesインスタンスを作成します。
fig = plt.figure()
ax = fig.add_axes([0.05,0.05,0.9,0.9])
# Basemapインスタンスを作成します。
# 海岸線は使用されていないため、解像度をNoneに設定して大陸の処理をスキップします
#（これにより処理速度が向上します）
m = Basemap(projection='kav7',lon_0=0,resolution=None)
# 地図投影肢の周りに線を引く。
# 地図投影領域の色の背景。
# 土地上の欠損値がこの色を表示します。
m.drawmapboundary(fill_color='0.3')
# プロットsst、次にpcolorの氷
im1 = m.pcolormesh(lons,lats,sst,shading='flat',cmap=plt.cm.jet,latlon=True)
im2 = m.pcolormesh(lons,lats,ice,shading='flat',cmap=plt.cm.gist_gray,latlon=True)
# 平行線と子午線を描きますが、それらにラベルを付けても構いません。
m.drawparallels(np.arange(-90.,99.,30.))
m.drawmeridians(np.arange(-180.,180.,60.))
# add colorbar
cb = m.colorbar(im1,"bottom", size="5%", pad="2%")
# add a title.
ax.set_title('SST and ICE analysis for %s'%date)
plt.show()
```
風ベクトルと風のバーブをプロットする。

```python:
import numpy as np
import matplotlib.pyplot as plt
import datetime
from mpl_toolkits.basemap import Basemap, shiftgrid
from netCDF4 import Dataset
# specify date to plot.
yyyy=1993; mm=03; dd=14; hh=00
date = datetime.datetime(yyyy,mm,dd,hh)
# set OpenDAP server URL.
URLbase="http://nomads.ncdc.noaa.gov/thredds/dodsC/modeldata/cmd_pgbh/"
URL=URLbase+"%04i/%04i%02i/%04i%02i%02i/pgbh00.gdas.%04i%02i%02i%02i.grb2" %\
             (yyyy,yyyy,mm,yyyy,mm,dd,yyyy,mm,dd,hh)
data = Dataset(URL)
# read lats,lons
# 逆の緯度で南から北へ行く。
latitudes = data.variables['lat'][::-1]
longitudes = data.variables['lon'][:].tolist()
# 海面の圧力と10 mの風のデータを得る。
# hpの単位で0.01単位で多項式を入力します。
slpin = 0.01*data.variables['Pressure_msl'][:].squeeze()
uin = data.variables['U-component_of_wind_height_above_ground'][:].squeeze()
vin = data.variables['V-component_of_wind_height_above_ground'][:].squeeze()
# 周期的な点を手動で追加する（余弦関数を使うことができる）
slp = np.zeros((slpin.shape[0],slpin.shape[1]+1),np.float)
slp[:,0:-1] = slpin[::-1]; slp[:,-1] = slpin[::-1,0]
u = np.zeros((uin.shape[0],uin.shape[1]+1),np.float64)
u[:,0:-1] = uin[::-1]; u[:,-1] = uin[::-1,0]
v = np.zeros((vin.shape[0],vin.shape[1]+1),np.float64)
v[:,0:-1] = vin[::-1]; v[:,-1] = vin[::-1,0]
longitudes.append(360.); longitudes = np.array(longitudes)
# make 2-d grid of lons, lats
lons, lats = np.meshgrid(longitudes,latitudes)
# 正射図のベースマップを作る。
m = Basemap(resolution='c',projection='ortho',lat_0=60.,lon_0=-60.)
# 図形の作成、軸の追加
fig1 = plt.figure(figsize=(8,10))
ax = fig1.add_axes([0.1,0.1,0.8,0.8])
# 所望の輪郭レベルを設定する。
clevs = np.arange(960,1061,5)
# グリッドのネイティブx、y座標を計算します。
x, y = m(lons, lats)
# 描画する平行線と子午線を定義します。
parallels = np.arange(-80.,90,20.)
meridians = np.arange(0.,360.,20.)
# SLP輪郭をプロットする。
CS1 = m.contour(x,y,slp,clevs,linewidths=0.5,colors='k',animated=True)
CS2 = m.contourf(x,y,slp,clevs,cmap=plt.cm.RdBu_r,animated=True)
# 投影グリッド上に風ベクトルをプロットする。
# 最初にシフトグリッドを-180から180に変更します
#（経度0から360ではなく）。 さもなければ、補間は乱される。
ugrid,newlons = shiftgrid(180.,u,longitudes,start=False)
vgrid,newlons = shiftgrid(180.,v,longitudes,start=False)
# ベクトルを投影格子に変換する。
uproj,vproj,xx,yy = \
m.transform_vector(ugrid,vgrid,newlons,latitudes,31,31,returnxy=True,masked=True)
# now plot.
Q = m.quiver(xx,yy,uproj,vproj,scale=700)
# キーを作る。
qk = plt.quiverkey(Q, 0.1, 0.1, 20, '20 m/s', labelpos='W')
# 海岸線、平行線、子午線を描画します。
m.drawcoastlines(linewidth=1.5)
m.drawparallels(parallels)
m.drawmeridians(meridians)
# add colorbar
cb = m.colorbar(CS2,"bottom", size="5%", pad="2%")
cb.set_label('hPa')
# set plot title
ax.set_title('SLP and Wind Vectors '+str(date))
plt.show()

# create 2nd figure, add axes
fig2 = plt.figure(figsize=(8,10))
ax = fig2.add_axes([0.1,0.1,0.8,0.8])
# plot SLP contours
CS1 = m.contour(x,y,slp,clevs,linewidths=0.5,colors='k',animated=True)
CS2 = m.contourf(x,y,slp,clevs,cmap=plt.cm.RdBu_r,animated=True)
# マップ上に風の棘をプロットする。
barbs = m.barbs(xx,yy,uproj,vproj,length=5,barbcolor='k',flagcolor='r',linewidth=0.5)
# 海岸線、平行線、子午線を描画します。
m.drawcoastlines(linewidth=1.5)
m.drawparallels(parallels)
m.drawmeridians(meridians)
# add colorbar
cb = m.colorbar(CS2,"bottom", size="5%", pad="2%")
cb.set_label('hPa')
# set plot title.
ax.set_title('SLP and Wind Barbs '+str(date))
plt.show()
```

NYとロンドンの間に大きな円（大圏コース）を描く。

```python:
from mpl_toolkits.basemap import Basemap
import numpy as np
import matplotlib.pyplot as plt
# 新しいFigure、Axesのインスタンスを作成します。
fig=plt.figure()
ax=fig.add_axes([0.1,0.1,0.8,0.8])
# セットアップメルケータマップ投影。
m = Basemap(llcrnrlon=-100.,llcrnrlat=20.,urcrnrlon=20.,urcrnrlat=60.,\
            rsphere=(6378137.00,6356752.3142),\
            resolution='l',projection='merc',\
            lat_0=40.,lon_0=-20.,lat_ts=20.)
# nylat, nylon are lat/lon of New York
nylat = 40.78; nylon = -73.98
# lonlat, lonlon are lat/lon of London.
lonlat = 51.53; lonlon = 0.08
# NYとロンドンの間に大きな円(大圏コース)ルートを描く
m.drawgreatcircle(nylon,nylat,lonlon,lonlat,linewidth=2,color='b')
m.drawcoastlines()
m.fillcontinents()
# draw parallels
m.drawparallels(np.arange(10,90,20),labels=[1,1,0,1])
# draw meridians
m.drawmeridians(np.arange(-180,180,30),labels=[1,1,0,1])
ax.set_title('Great Circle from New York to London')
plt.show()
```

![image.png](https://qiita-image-store.s3.amazonaws.com/0/238800/ad6f4f59-5cc3-2635-126c-1948c297f6e8.png)

昼夜のターミネーターを地図上に描く。

```python:
import numpy as np
from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
from datetime import datetime
# ミラー投影
map = Basemap(projection='mill',lon_0=180)
# 海岸をプロットし、ラベルの子午線と平行を描きます。
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])
# 大陸の 'サンゴ'（zorder = 0）、カラーウェットエリア 'アクア'
map.drawmapboundary(fill_color='aqua')
map.fillcontinents(color='coral',lake_color='aqua')
# 夜間領域を陰影付けし、地図が表示されるようにアルファ透明度を設定します。 
# UTCで現在の時刻を使用します。
date = datetime.utcnow()
CS=map.nightshade(date)
plt.title('Day/Night Map for %s (UTC)' % date.strftime("%d %b %Y %H:%M:%S"))
plt.show()
```

![image.png](https://qiita-image-store.s3.amazonaws.com/0/238800/e735ff7b-6b12-ac73-4d10-2d3985764edc.png)