Как получить месяц максимального стока
Я хочу получить месяц максимального стока для каждого года и для временного ряда в целом. Идея состоит в том, чтобы охарактеризовать глобальную сезонность, взглянув на месяц максимального стока. Затем я хочу попытаться определить, имеет ли каждый пиксель унимодальный или бимодальный режим.
Я хочу создать карту, подобную той, что приведена в Примерах Pangeo здесь .
То, что это показывает, является часом максимального количества осадков. Я хочу показать МЕСЯЦ максимального стока (в виде целого числа).
Получение данных
Здесь я загружаю данные стока GRUN и создаю объект xarray.
ПРИМЕЧАНИЕ: Набор данных здесь> 1 ГБ. Я использую его, чтобы сделать этот пример полностью воспроизводимым.
# get the data
import subprocess
command = """
wget -O grun.nc https://www.research-collection.ethz.ch/bitstream/handle/20.500.11850/324386/GRUN_v1_GSWP3_WGS84_05_1902_2014.nc?sequence=1&isAllowed=y
"""
import os
if not os.path.exists('grun.nc'):
process = subprocess.Popen(command.split(), stdout=subprocess.PIPE)
output, error = process.communicate()
# read the data
import xarray as xr
ds = xr.open_dataset('grun.nc')
# select a subset so we can work with it more quickly
ds = ds.isel(time=slice(-100,-1))
ds
Out[]:
<xarray.Dataset>
Dimensions: (lat: 360, lon: 720, time: 99)
Coordinates:
* lon (lon) float64 -179.8 -179.2 -178.8 -178.2 ... 178.8 179.2 179.8
* lat (lat) float64 -89.75 -89.25 -88.75 -88.25 ... 88.75 89.25 89.75
* time (time) datetime64[ns] 2006-09-01 2006-10-01 ... 2014-11-01
Data variables:
Runoff (time, lat, lon) float32 ...
Attributes:
title: GRUN
version: GRUN 1.0
meteorological_forcing: GSWP3
temporal_resolution: monthly
spatial_resolution: 0.5x0.5
crs: WGS84
proj4: +proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs
EPSG: 4326
references: Ghiggi et al.,2019. GRUN: An observation-based g...
authors: Gionata Ghiggi; Lukas Gudmundsson
contacts: gionata.ghiggi@gmail.com; lukas.gudmundsson@env....
institution: Land-Climate Dynamics, Institute for Atmospheric...
institution_id: IAC ETHZ
Что я пробовал
У меня есть значения nan, поэтому я не могу просто применить argmax() к набору данных. Я использую тот же подход, что и @jhamman здесь в сочетании с примерами Pangeo выше. Я не совсем уверен, что это дает мне, но, кажется, дает мне
# Apply argmax where you have NAN values
def my_func(ds, dim=None):
return ds.isel(**{dim: ds['Runoff'].argmax(dim)})
mask = ds['Runoff'].isel(time=0).notnull() # determine where you have valid data
ds2 = ds.fillna(-9999) # fill nans with a missing flag of some kind
new = ds2.reset_coords(drop=True).groupby('time.month').apply(my_func, dim='time').where(mask) # do the groupby operation/reduction and reapply the mask
new
Out[]:
<xarray.Dataset>
Dimensions: (lat: 360, lon: 720, month: 12)
Coordinates:
* lon (lon) float64 -179.8 -179.2 -178.8 -178.2 ... 178.8 179.2 179.8
* lat (lat) float64 -89.75 -89.25 -88.75 -88.25 ... 88.75 89.25 89.75
* month (month) int64 1 2 3 4 5 6 7 8 9 10 11 12
Data variables:
Runoff (month, lat, lon) float32 nan nan nan nan nan ... nan nan nan nan
Attributes:
title: GRUN
version: GRUN 1.0
meteorological_forcing: GSWP3
temporal_resolution: monthly
spatial_resolution: 0.5x0.5
crs: WGS84
proj4: +proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs
EPSG: 4326
references: Ghiggi et al.,2019. GRUN: An observation-based g...
authors: Gionata Ghiggi; Lukas Gudmundsson
contacts: gionata.ghiggi@gmail.com; lukas.gudmundsson@env....
institution: Land-Climate Dynamics, Institute for Atmospheric...
institution_id: IAC ETHZ
Это дает мне
import matplotlib.pyplot as plt
fig,ax = plt.subplots(figsize=(12,8))
new.Runoff.sel(month=10).plot(ax=ax, cmap='twilight')
Идеальный вывод
Я хочу, чтобы значение каждого пикселя было месяцем максимального стока.
Рад преобразовать в pandas при необходимости.
Таким образом, я получу набор данных xr.Dataset с целым числом для месяца максимального стока. В идеале было бы замечательно иметь месяц максимального стока с течением времени, чтобы я также мог видеть, как изменилась эта сезонность.
<xarray.Dataset>
Dimensions: (lat: 360, lon: 720)
Coordinates:
* lon (lon) float64 -179.8 -179.2 -178.8 -178.2 ... 178.8 179.2 179.8
* lat (lat) float64 -89.75 -89.25 -88.75 -88.25 ... 88.75 89.25 89.75
Data variables:
Month_of_max (lat, lon) int32 ...
# OR EVEN BETTER
<xarray.Dataset>
Dimensions: (lat: 360, lon: 720, Year: 10)
Coordinates:
* lon (lon) float64 -179.8 -179.2 -178.8 -178.2 ... 178.8 179.2 179.8
* lat (lat) float64 -89.75 -89.25 -88.75 -88.25 ... 88.75 89.25 89.75
* year (year) float64 2010 2011 2012 2013 ...
Data variables:
Month_of_max (lat, lon, year) int32 ...