Блог

У меня есть 4-мерные данные (время, глубина, y и x), но широта и долгота являются 2d-массивами. y и x — это просто индексы, поэтому просто целые числа, начинающиеся с 0, 1 … end и т.д. Очень похоже на примерный набор данных, предоставленный MetPy:

https://unidata.github.io/MetPy/latest/examples/cross_section.html

К сожалению, это не самый воспроизводимый, потому что он очень специфичен для самих данных. Но у меня возникли проблемы с частью поперечного сечения. Я могу проанализировать данные в соответствии с metpy, но затем я получаю ошибку при получении поперечного сечения:

 data = ds.metpy.parse_cf().squeeze()

print(data)

<xarray.Dataset>
Dimensions:               (y: 1347, x: 1379, deptht: 46, axis_nbounds: 2, time_counter: 12)
Coordinates:
    nav_lat               (y, x) float32 20.92 20.92 20.92 ... 68.49 68.49 68.48
    nav_lon               (y, x) float32 -78.95 -78.9 -78.85 ... -3.614 -3.546
  * deptht                (deptht) float32 3.047 9.454 ... 5.625e 03 5.875e 03
    time_centered         (time_counter) datetime64[ns] 1993-01-16T12:00:00 ....
  * time_counter          (time_counter) datetime64[ns] 1993-01-16T12:00:00 ....
Dimensions without coordinates: y, x, axis_nbounds
Data variables: (12/17)
    deptht_bounds         (deptht, axis_nbounds, y, x) float32 0.0 0.0 ... nan
    time_centered_bounds  (time_counter, axis_nbounds, y, x) datetime64[ns] 1...
    time_counter_bounds   (time_counter, axis_nbounds, y, x) datetime64[ns] 1...
    votemper              (time_counter, deptht, y, x) float32 26.63 ... nan
    vosaline              (time_counter, deptht, y, x) float32 35.88 ... nan
    sosstsst              (time_counter, y, x) float32 26.63 26.58 ... nan nan
    ...                    ...
    sohefldo              (time_counter, y, x) float32 -50.02 -45.85 ... nan nan
    somixhgt              (time_counter, y, x) float32 19.84 19.76 ... nan nan
    sowindsp              (time_counter, y, x) float32 5.591 5.48 ... nan nan
    sohefldp              (time_counter, y, x) float32 nan nan nan ... nan nan
    sowafldp              (time_counter, y, x) float32 3.94e-06 ... nan
    sobowlin              (time_counter, y, x) float32 20.04 20.04 ... nan nan
Attributes:
    name:                      1_VIKING20X.L46-KKG36107B_1d_19930101_19930704...
    description:               ocean T grid variables
    title:                     ocean T grid variables
    Conventions:               CF-1.6
    timeStamp:                 2019-Sep-11 21:02:06 GMT
    uuid:                      a091f081-2943-4c66-aa77-0917f654b802
    history:                   Thu Sep 12 22:15:40 2019: ncrcat -O -F /gfs1/w...
    NCO:                       4.4.8
    nco_openmp_thread_number:  1
 

Затем я пробую часть поперечного сечения (просто чтобы отметить, что 2D широта и долгота помечены как ‘nav_lat’ и ‘nav_lon’, вместо широты и длины):

  start = (40, -40)
 end = (50, -30)
    
 cross = cross_section(data, start, end).set_coords(('nav_lon', 'nav_lat'))
 print(cross)

---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
~/miniconda3/envs/py3_std_maps/lib/python3.9/site-packages/metpy/interpolate/slices.py in cross_section(data, start, end, steps, interp_type)
    165         try:
--> 166             crs_data = data.metpy.pyproj_crs
    167             x = data.metpy.x

~/miniconda3/envs/py3_std_maps/lib/python3.9/site-packages/metpy/xarray.py in pyproj_crs(self)
    252         """Return the coordinate reference system (CRS) as a pyproj object."""
--> 253         return self.crs.to_pyproj()
    254 

~/miniconda3/envs/py3_std_maps/lib/python3.9/site-packages/metpy/xarray.py in crs(self)
    232             return self._data_array.coords['metpy_crs'].item()
--> 233         raise AttributeError('crs attribute is not available.')
    234 

AttributeError: crs attribute is not available.

During handling of the above exception, another exception occurred:

ValueError                                Traceback (most recent call last)
/tmp/ipykernel_2676704/3597335100.py in <module>
----> 1 cross = cross_section(data, start, end).set_coords(('nav_lon', 'nav_lat'))
      2 print(cross)

~/miniconda3/envs/py3_std_maps/lib/python3.9/site-packages/metpy/interpolate/slices.py in cross_section(data, start, end, steps, interp_type)
    154     if isinstance(data, xr.Dataset):
    155         # Recursively apply to dataset
--> 156         return data.map(cross_section, True, (start, end), steps=steps,
    157                         interp_type=interp_type)
    158     elif data.ndim == 0:

~/miniconda3/envs/py3_std_maps/lib/python3.9/site-packages/xarray/core/dataset.py in map(self, func, keep_attrs, args, **kwargs)
   5106         if keep_attrs is None:
   5107             keep_attrs = _get_keep_attrs(default=False)
-> 5108         variables = {
   5109             k: maybe_wrap_array(v, func(v, *args, **kwargs))
   5110             for k, v in self.data_vars.items()

~/miniconda3/envs/py3_std_maps/lib/python3.9/site-packages/xarray/core/dataset.py in <dictcomp>(.0)
   5107             keep_attrs = _get_keep_attrs(default=False)
   5108         variables = {
-> 5109             k: maybe_wrap_array(v, func(v, *args, **kwargs))
   5110             for k, v in self.data_vars.items()
   5111         }

~/miniconda3/envs/py3_std_maps/lib/python3.9/site-packages/metpy/interpolate/slices.py in cross_section(data, start, end, steps, interp_type)
    167             x = data.metpy.x
    168         except AttributeError:
--> 169             raise ValueError('Data missing required coordinate information. Verify that '
    170                              'your data have been parsed by MetPy with proper x and y '
    171                              'dimension coordinates and added crs coordinate of the '

ValueError: Data missing required coordinate information. Verify that your data have been parsed by MetPy with proper x and y dimension coordinates and added crs coordinate of the correct projection for each variable.
 

Я попробовал в качестве альтернативы, потому что, возможно, parse_cf вызывает проблемы:

 data = ds.metpy.assign_latitude_longitude(force=True).squeeze()
 

но тогда я все равно получаю ту же ошибку при применении поперечного сечения.

Есть идеи о том, как это исправить? Еще раз извините за отсутствие воспроизводимости, но любые идеи были бы очень полезны 🙂 Возможно, это связано с проекцией?

Это примерное изображение того, как данные выглядят в одном экземпляре времени и глубины (с учетом температуры):

metpy.interpolate.cross_section требуется, чтобы ваши данные включали координаты измерений x и y и добавленную metpy_crs координату (из любого parse_cf или assign_crs ). В этой ситуации, когда эти координаты измерений x и y отсутствуют, но у вас есть координаты 2D широты и долготы, эти координаты измерений можно вычислить и добавить .metpy.assign_y_x() (вместо assign_latitude_longitude того, что вы заявили, что пытались, что делает обратное — добавление вспомогательных координат широты / долготы из координат измерений y / x).

Итак, если ваш набор данных имеет действительное отображение сетки CF, соответствующее вашей проекции данных, у вас будет:

 data = ds.metpy.parse_cf()
data = data.metpy.assign_y_x()

start = (40, -40)
end = (50, -30)
    
cross = cross_section(data, start, end).set_coords(('nav_lon', 'nav_lat'))
 

Если это не так:

 data = ds.metpy.assign_crs({
    "grid_mapping_name": "lambert_conformal_conic",
    "standard_parallel": 25.0,
    "longitude_of_central_meridian": 265.0,
    "latitude_of_projection_origin": 25.0,
})
data = data.metpy.assign_y_x()

start = (40, -40)
end = (50, -30)
    
cross = cross_section(data, start, end).set_coords(('nav_lon', 'nav_lat'))
 

(изменение атрибутов проекции на те, которые необходимы для вашей заданной проекции)