Pythonによる気象・気候データ解析I 2章章末問題

[Pythonによる気象・気候データ解析I: Pythonの基礎・気候値と偏差・回帰相関分析](https://amzn.to/3yvtPCG)

<a href="https://amzn.to/3yvtPCG" target="_blank"><img src="https://m.media-amazon.com/images/I/71u7jHR9KgL._SY522_.jpg" width=150 /></a>

注：しばらくAmazonでは品切れなので楽天等で
<a href="https://hb.afl.rakuten.co.jp/ichiba/3c7b1caa.aa54b4a0.3c7b1cab.f01332e0/?pc=https%3A%2F%2Fitem.rakuten.co.jp%2Fbook%2F17792414%2F&link_type=text&ut=eyJwYWdlIjoiaXRlbSIsInR5cGUiOiJ0ZXh0Iiwic2l6ZSI6IjI0MHgyNDAiLCJuYW0iOjEsIm5hbXAiOiJyaWdodCIsImNvbSI6MSwiY29tcCI6ImRvd24iLCJwcmljZSI6MSwiYm9yIjoxLCJjb2wiOjEsImJidG4iOjEsInByb2QiOjAsImFtcCI6ZmFsc2V9" target="_blank" rel="nofollow sponsored noopener" style="word-wrap:break-word;">Pythonによる気象・気候データ解析1 Pythonの基礎・気候値と偏差・回帰相関分析 [ 神山 翼 ]</a>

の章末問題を一つずつ解いていきます。

# A
## A-1 1997−08

- 本文と同じ年で夏の8月を

```python
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.colors import Normalize
#!pip install japanize_matplotlib > /dev/null
import japanize_matplotlib

draw_year = 1997
draw_month = 8

vim = -5
vmax = 35
vint = 5

cm = plt.get_cmap('seismic')
cs = plt.contourf(lon2, lat2,
                  np.squeeze(sst[:, :, (y==draw_year)*(m==draw_month)]),
                  cmap=cm, norm=Normalize(vmin=vim, vmax=vmax),
                  levels=np.arange(vim, vmax+vint, vint), extend='both')
plt.colorbar(cs)
plt.xlabel('Longitude')
plt.ylabel('Latitude')
plt.xlim(0, 360)
plt.ylim(-90, 90)
plt.title(str(draw_year)+'/'+str(draw_month))
plt.show()
```

- 赤い領域が広がっている

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/c5a8907e-1a28-b201-62b5-a8c1dc23f3c8.png)

## A-2 日本近海

```python
draw_year = 1997
draw_month = 8

vim = -5
vmax = 35
vint = 5

cm = plt.get_cmap('seismic')
cs = plt.contourf(lon2, lat2,
                  np.squeeze(sst[:, :, (y==draw_year)*(m==draw_month)]),
                  cmap=cm, norm=Normalize(vmin=vim, vmax=vmax),
                  levels=np.arange(vim, vmax+vint, vint), extend='both')
plt.colorbar(cs)
plt.xlabel('Longitude')
plt.ylabel('Latitude')

# 日本近海
plt.xlim(120, 160)
plt.ylim(30, 50)
plt.title(str(draw_year)+'/'+str(draw_month) + ' 日本近海')

plt.show()

```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/42ffaa74-275d-4789-86c8-cae20c33ecad.png)

- 同じく draw_month = 12 にしてみると海水温の変化がわかりやすい

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/89d826be-2b1c-9cc0-5c74-90d4b569e536.png)

## A-3

- matplotlibのColorMapのマニュアルを見て寒色暖色系のグラデーションを試してみた

```python
# c.f. https://matplotlib.org/stable/users/explain/colors/colormaps.html

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.colors import Normalize
#!pip install japanize_matplotlib > /dev/null
import japanize_matplotlib

draw_year = 1997
draw_month = 8

vim = -5
vmax = 35
vint = 5

def draw_sea_temp(cm_name):
  cm = plt.get_cmap(cm_name)
  cs = plt.contourf(lon2, lat2,
                    np.squeeze(sst[:, :, (y==draw_year)*(m==draw_month)]),
                    cmap=cm, norm=Normalize(vmin=vim, vmax=vmax),
                    levels=np.arange(vim, vmax+vint, vint), extend='both')
  plt.colorbar(cs)
  plt.xlabel('Longitude')
  plt.ylabel('Latitude')
  plt.xlim(0, 360)
  plt.ylim(-90, 90)
  plt.title(str(draw_year)+'/'+str(draw_month)+' '+cm_name)
  plt.show()

draw_sea_temp('seismic')
draw_sea_temp('bwr')
draw_sea_temp('coolwarm')
draw_sea_temp('rainbow')
draw_sea_temp('jet')
draw_sea_temp('turbo') 
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/e8c77e14-7608-a296-924b-740dd81ea1a0.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/c9a493bb-649c-278d-42af-f33103b30392.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/213388e4-4128-8f57-33f9-df0df9bc512b.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/f7daf3c8-eca6-18f4-93b4-530bd5127b34.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/02526634-8525-22e6-72d2-e754225bf2ba.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/c1e8513f-86ef-1a65-fbf1-a28e5cb354bb.png)

# B
## B-1

```python
import numpy as np

A = np.array([[1, 2, 4], 
              [4, 3, 5],
              [6, 2, 9]])
B = np.array([5, 5, 4])

print(B==4)

A[:, B==4]
```

```python
[False False  True]
array([[4],
       [5],
       [9]])
```

## B-2

(draw_year年)かつ(draw_month月)つまり draw_year年draw_month月 のsstの値を取得する。

# C

- 気圧は950hPaより小さい値もあるがグラデーションがうまく出ないので高めに950hPaとした。
- ３時間毎にプロットして台風の動きを見えるようにした

```python
vim = 950
vmax = 1030
vint = 5

def draw_typhoon(t):
  cm = plt.get_cmap('seismic')
  cs = plt.contourf(lon2, lat2,
                    surface_pressure[:, :, t] / 100,
                    cmap=cm, norm=Normalize(vmin=vim, vmax=vmax),
                    levels=np.arange(vim, vmax+vint, vint), extend='both')
  plt.colorbar(cs)
  plt.xlabel('Longitude')
  plt.ylabel('Latitude')
  plt.xlim(120, 150)
  plt.ylim(23, 45)
  plt.title(str(t) + '時')
  plt.show()

for t in range(0, 24, 3):
  draw_typhoon(t)
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/5ad05b80-7565-09e4-c30f-ecaa8e96c705.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/5261490c-4bfb-e72a-4e01-d58ef96cc089.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/cd0ed875-33ac-edf4-5ea0-7ea0226cbfde.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/0943d295-f10e-2201-3266-b550b7ba7872.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/bab77129-0c55-8702-a4c7-3b7dbfa06434.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/7378ef67-41c1-d4a1-9298-4e3d4333813d.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/2ef1a5ac-b34d-2ad8-245a-22a07e9b67d5.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/5f1e820e-562e-a1fc-9baa-2718d5b1b2f2.png)

- 雲量

```python
cloud = msm_dataset['cloud']

plt.hist(cloud.flatten(), bins=100)
plt.show()
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/0bda5814-0b5b-2109-1085-aead78099a0b.png)

```python
vim = 0
vmax = 100
vint = 5

def draw_cloud(t):
  cm = plt.get_cmap('seismic')
  cs = plt.contourf(lon2, lat2,
                    cloud[:, :, t],
                    cmap=cm, norm=Normalize(vmin=vim, vmax=vmax),
                    levels=np.arange(vim, vmax+vint, vint), extend='both')
  plt.colorbar(cs)
  plt.xlabel('Longitude')
  plt.ylabel('Latitude')
  plt.xlim(120, 150)
  plt.ylim(23, 45)
  plt.title(str(t) + '時')
  plt.show()

for t in range(0, 24, 3):
  draw_cloud(t)
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/35fad0aa-4761-1d5a-f9b5-e7c8d51139a8.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/4a658eb4-65c1-78aa-62df-b6bde916f4c2.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/95597e5f-902e-7372-5a8d-4a6424794d4e.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/b225cd02-c219-0b09-45ab-2ee722796a07.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/dafa17c9-d7b3-5854-9e15-748b99a1d63b.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/f338ac26-2306-592b-cd8e-a99ec5abcf2e.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/c3377f80-f1f9-79f1-38a5-bbd9a82f8aff.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/f9bf65ff-a426-9589-0a46-635e2b0f5556.png)


# D

- 本書発売日 2024-05-01 14:00 の気温を描いてみた
    - 「なんらかの現象」までは思いつかず
    - 27-30.5 とちょっと高めの気がするがこういうものか

[NetCDF化した数値予報GPVデータ](https://database.rish.kyoto-u.ac.jp/arch/jmadata/gpv-netcdf.html)

```python
!pip install netCDF4 > /dev/null

import numpy as np
import netCDF4 as nc4
import datetime
import matplotlib.pyplot as plt
from matplotlib.colors import Normalize
#!pip install japanize_matplotlib > /dev/null
import japanize_matplotlib

loadfile = 'drive/MyDrive/0501.nc'

msm_data = nc4.Dataset(loadfile)

temperature = msm_data.variables['temp'][:]
temperature = temperature.transpose(2, 1, 0)
temperature = temperature.astype(float)

[imt, jmt, tmt] = temperature.shape

lon = msm_data.variables['lon'][:]
lat = msm_data.variables['lat'][:]

[xgrid2, ygrid2] = np.meshgrid(lon, lat)
lon2 = (xgrid2.T).astype(float)
lat2 = (ygrid2.T).astype(float)
del lon, lat, xgrid2, ygrid2

y = np.zeros(tmt)
m = np.zeros(tmt)
d = np.zeros(tmt)
h = np.zeros(tmt)

date_offset = datetime.datetime(2019, 10, 12)
time = msm_data.variables['time'][:]
for tt in range(0, tmt):
  time_data = date_offset + datetime.timedelta(hours=int(time[tt]))
  y[tt] = time_data.year
  m[tt] = time_data.month
  d[tt] = time_data.day
  h[tt] = time_data.hour

plt.hist(temperature.flatten())
plt.show()
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/22b00d23-fdcc-5c79-1a41-82f6e95f4f3c.png)

```python
vim = 27
vmax = 31
vint = 0.5

def draw_temp(t):
  cm = plt.get_cmap('seismic')
  cs = plt.contourf(lon2, lat2,
                    temperature[:, :, t] / 10,
                    cmap=cm, norm=Normalize(vmin=vim, vmax=vmax),
                    levels=np.arange(vim, vmax+vint, vint), extend='both')
  plt.colorbar(cs)
  plt.xlabel('Longitude')
  plt.ylabel('Latitude')
  plt.xlim(120, 150)
  plt.ylim(23, 45)
  plt.title(str(t) + '時')
  plt.show()

draw_temp(14)
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/116049/0a79bb0e-e325-c47c-1b76-08c4646bcbc7.png)

---

- [Pythonによる気象・気候データ解析I 1章章末問題](https://qiita.com/aokikenichi/items/755925f2f16e4f9a022c)
- [Pythonによる気象・気候データ解析I 3章章末問題](https://qiita.com/aokikenichi/items/4be2cb4a8ebb4648fb4d)
- [Pythonによる気象・気候データ解析II 1章章末問題](https://qiita.com/aokikenichi/items/3dd97d989cdd31043e2f)