построение данных из netcdf с картографическими данными, не выводя данные на 0 долготе - PullRequest
3 голосов
/ 29 марта 2020

Я начинаю свое путешествие с привязанными к сетке данными, и я пытался построить некоторые температурные данные из файла netcdf с картопией. Я следовал за некоторыми примерами, и я не могу понять, почему у моих графиков есть белая линия в середине. (я уже проверил данные и матрицы заполнены числами, без NaN)

import cartopy.crs as ccrs
import matplotlib.pyplot as plt
import numpy as np
import xarray as xr
from cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatterimport glob


data = xr.open_dataset('aux1.nc')
lat = data.lat
lon = data.lon
time = data.time
Temp = data.air


#Calculo la temperatura media anual
Tanual = Temp.resample(time="y").mean()
#Promedio de todos los meses
Tprom = Temp.mean(dim="time").values


#Grafico
fig = plt.figure(figsize=(10, 4))
ax = fig.add_subplot(1, 1, 1, projection=ccrs.PlateCarree())
ax.coastlines()
ax.set_global()
ct = ax.contourf(lon,lat,Tprom,transform=ccrs.PlateCarree(),cmap="bwr")
ax.gridlines()
cb = plt.colorbar(ct,orientation="vertical",extendrect='True')
cb.set_label("Temperatura [°C]")
ax.set_xticks(np.arange(-180,181,60), crs=ccrs.PlateCarree())
ax.set_yticks(np.arange(-90,91,30), crs=ccrs.PlateCarree())
lon_formatter = LongitudeFormatter(zero_direction_label=True)
lat_formatter = LatitudeFormatter()
ax.xaxis.set_major_formatter(lon_formatter)
ax.yaxis.set_major_formatter(lat_formatter)

https://i.stack.imgur.com/aPHz4.jpg

1 Ответ

4 голосов
/ 29 марта 2020

Хороший вопрос! Проблема в том, что для большинства климатических данных с координатной сеткой координата долготы выглядит примерно так:

array([1.25, 3.75, 6.25, ..., 351.25, 353.75, 356.25, 358.75])

Так что нет явной точки longitude=0, и это часто дает вам тонкую белую линию на графике. Я вижу эту проблему в опубликованных работах (даже «Природа») иногда тоже!

Есть много способов обойти это, но самый простой способ - использовать пакет cartopy, который имеет утилиту под названием add_cyclic_point, которая в основном интерполирует данные по обе стороны от точки longitude=0. (Ссылка: https://scitools.org.uk/cartopy/docs/v0.15/cartopy/util/util.html)

Единственным недостатком этого метода является то, что при использовании xarray это означает, что вам нужно извлечь данные вручную, а затем вы потеряете метаданные, поэтому вот функция, которую я написал, чтобы сделать это красивым и простым в использовании, сохраняя при этом метаданные.

from cartopy.util import add_cyclic_coord
import xarray as xr 

def xr_add_cyclic_point(da):
    """
    Inputs
    da: xr.DataArray with dimensions (time,lat,lon)
    """

    # Use add_cyclic_point to interpolate input data
    lon_idx = da.dims.index('lon')
    wrap_data, wrap_lon = add_cyclic_point(da.values, coord=da.lon, axis=lon_idx)

    # Generate output DataArray with new data but same structure as input
    outp_da = xr.DataArray(data=wrap_data, 
                           coords = {'time': da.time, 'lat': da.lat, 'lon': wrap_lon}, 
                           dims=da.dims, 
                           attrs=da.attrs)

    return outp_da

Пример

Так, например, если мой исходный DataArray выглядит следующим образом :

<xarray.DataArray 'tas' (time: 60, lat: 90, lon: 144)>
[777600 values with dtype=float32]
Coordinates:
  * lat      (lat) float64 -89.49 -87.98 -85.96 -83.93 ... 85.96 87.98 89.49
  * lon      (lon) float64 1.25 3.75 6.25 8.75 11.25 ... 351.3 353.8 356.2 358.8
  * time     (time) object 1901-01-16 12:00:00 ... 1905-12-16 12:00:00
Attributes:
    long_name:         Near-Surface Air Temperature
    units:             K
    valid_range:       [100. 400.]
    cell_methods:      time: mean
    standard_name:     air_temperature
    original_units:    deg_k
    original_name:     t_ref
    cell_measures:     area: areacella
    associated_files:  baseURL: http://cmip-pcmdi.llnl.gov/CMIP5/dataLocation...

И когда я строю среднее время, это дает:

tas.mean(dim='time').plot.contourf()

enter image description here

Теперь я могу используйте мою функцию, чтобы сгенерировать новый интерполированный DataArray, например:

wrapped_tas = xr_add_cyclic_point(tas)
wrapped_tas

<xarray.DataArray (time: 60, lat: 90, lon: 145)>
array([[[251.19466, 251.19469, 251.19472, ..., 251.19226, 251.19073,
         251.19466], ...
        [250.39403, 250.39468, 250.39961, ..., 250.39429, 250.39409,
         250.39403]]], dtype=float32)
Coordinates:
  * time     (time) object 1901-01-16 12:00:00 ... 1905-12-16 12:00:00
  * lat      (lat) float64 -89.49 -87.98 -85.96 -83.93 ... 85.96 87.98 89.49
  * lon      (lon) float64 1.25 3.75 6.25 8.75 11.25 ... 353.8 356.2 358.8 361.2
Attributes:
    long_name:         Near-Surface Air Temperature
    units:             K
    valid_range:       [100. 400.]
    cell_methods:      time: mean
    standard_name:     air_temperature
    original_units:    deg_k
    original_name:     t_ref
    cell_measures:     area: areacella
    associated_files:  baseURL: http://cmip-pcmdi.llnl.gov/CMIP5/dataLocation...

Как вы можете видеть, координата долготы была увеличена на одну точку, до go с 144-> 145 длины , это означает, что теперь он «обтекает» точку longitude=0.

Этот новый массив данных при построении дает график без белой линии:)

wrapped_tas.mean(dim='time').plot.contour()

enter image description here

Надеюсь, это поможет! :)

...