我按照代码绘制了这个网站上的横截面图。 (https://unidata.github.io/MetPy/latest/examples/cross_section.html#sphx-glr-examples-cross-section-py)
在示例中，运行以下代码将生成此结果(交叉)。cross = cross_section(data, start, end).set_coords(('lat', 'lon'))在此处输入图像描述
但是，与示例不同的是，运行横截面代码后结果值的 x，y 坐标已更改。并且经纬度是固定的。 在此处输入图像描述 我不明白为什么这个结果会这样。
我的代码
import cartopy.crs as ccrs
import cartopy.feature as cfeature
import matplotlib.pyplot as plt
import numpy as np
import os
import xarray as xr
import metpy.calc as mpcalc
from metpy.cbook import get_test_data
from metpy.interpolate import cross_section
os.chdir ("D:/PPL/CrossSection/20190404/")
ncfile = r'./2019040400pres.nc'
NC = xr.open_dataset(ncfile)
pres_list = ['1000','975','950','925','900','875','850','800','750','700','650','600','550','500','450','400','350','300','250','200','150','100','70','50'];
p_list = np.array(pres_list, dtype = np.float64)
NC = NC.metpy.parse_cf().squeeze()
NC = NC.assign_coords(isobaric = p_list).set_coords('isobaric')
我的代码必须根据 nc 文件除以同量异位压层进行分析，因此我合并了除以每层的值。 喜欢这个
tmp1 = NC['TMP_1000mb']
tmp2 = NC['TMP_975mb']
...
tmp23 = NC['TMP_70mb']
tmp24 = NC['TMP_50mb']
TMP = xr.concat([tmp1,tmp2],'isobaric')
TMP = xr.concat([TMP,tmp3],'isobaric')
...
TMP = xr.concat([TMP,tmp23],'isobaric')
TMP = xr.concat([TMP,tmp24],'isobaric')
TMP['isobaric']=p_list
NC['Temperature'] = TMP
所以，我打印NC。结果是
<xarray.Dataset>
Dimensions: (isobaric: 24, x: 602, y: 781)
Coordinates:
* y (y) float64 0.0 1.5e+03 3e+03 ... 1.168e+06 1.17e+06
* x (x) float64 0.0 1.5e+03 3e+03 ... 9e+05 9.015e+05
latitude (y, x) float64 32.26 32.26 32.26 ... 42.94 42.94 42.93
longitude (y, x) float64 121.8 121.9 121.9 ... 132.5 132.5 132.5
time datetime64[ns] 2019-04-04
metpy_crs object Projection: latitude_longitude
* isobaric (isobaric) float64 1e+03 975.0 950.0 ... 100.0 70.0 50.0
Data variables:
...
...
Temperature (isobaric, y, x) float32 282.87515 282.87515 ... nan nan
Attributes:
Conventions: CF-1.0
History: created by wgrib2
GRIB2_grid_template: 30
但是，当我运行横截面时，x，y更改为开始和结束，并且原始经度和纬度是固定的。
start = (37.5, 126.63)
end = (37.8 ,128.86)
cross = cross_section(NC, start, end)
#I didn't write .set_codes here because unlike the example, coordinates are given
#latitude and longitude.
print(cross)
<xarray.Dataset>
Dimensions: (index: 100, isobaric: 24)
Coordinates:
latitude (index) float64 32.26 32.26 32.26 ... 32.26 32.26 32.26
longitude (index) float64 121.8 121.8 121.8 ... 121.8 121.8 121.8
time datetime64[ns] 2019-04-04
metpy_crs object Projection: latitude_longitude
x (index) float64 126.6 126.7 126.7 ... 128.8 128.8 128.9
y (index) float64 37.5 37.5 37.51 37.51 ... 37.79 37.8 37.8
* index (index) int32 0 1 2 3 4 5 6 7 ... 92 93 94 95 96 97 98 99
* isobaric (isobaric) float64 1e+03 975.0 950.0 ... 100.0 70.0 50.0
Data variables:
...
...
Temperature (isobaric, index) float64 282.9 282.9 ... 209.6 209.6
Attributes:
Conventions: CF-1.0
History: created by wgrib2
GRIB2_grid_template: 30
我希望保持x，y，并且纬度和经度在开始和结束范围内，如示例中所示。 我真的不明白为什么会这样。
感谢您的关注。
作为参考，我附加了数据源和初始数据。 数据来源 ： 来自 KMA 的 LDAPS 数据
完整的NC数据(初始数据)
<xarray.Dataset>
Dimensions: (time: 1, x: 602, y: 781)
Coordinates:
* y (y) float64 0.0 1.5e+03 3e+03 ... 1.168e+06 1.17e+06
* x (x) float64 0.0 1.5e+03 3e+03 ... 9e+05 9.015e+05
latitude (y, x) float64 ...
longitude (y, x) float64 ...
* time (time) datetime64[ns] 2019-04-04
Data variables:
DZDT_1000mb (time, y, x) float32 ...
DZDT_975mb (time, y, x) float32 ...
...
DZDT_70mb (time, y, x) float32 ...
DZDT_50mb (time, y, x) float32 ...
UGRD_1000mb (time, y, x) float32 ...
VGRD_1000mb (time, y, x) float32 ...
UGRD_975mb (time, y, x) float32 ...
VGRD_975mb (time, y, x) float32 ...
...
UGRD_70mb (time, y, x) float32 ...
VGRD_70mb (time, y, x) float32 ...
UGRD_50mb (time, y, x) float32 ...
VGRD_50mb (time, y, x) float32 ...
HGT_1000mb (time, y, x) float32 ...
HGT_975mb (time, y, x) float32 ...
...
HGT_70mb (time, y, x) float32 ...
HGT_50mb (time, y, x) float32 ...
TMP_1000mb (time, y, x) float32 ...
TMP_975mb (time, y, x) float32 ...
...
TMP_70mb (time, y, x) float32 ...
TMP_50mb (time, y, x) float32 ...
var0_1_194_1000mb (time, y, x) float32 ...
var0_1_194_975mb (time, y, x) float32 ...
...
var0_1_194_70mb (time, y, x) float32 ...
var0_1_194_50mb (time, y, x) float32 ...
RH_1000mb (time, y, x) float32 ...
RH_975mb (time, y, x) float32 ...
...
RH_70mb (time, y, x) float32 ...
RH_50mb (time, y, x) float32 ...
Attributes:
Conventions: CF-1.0
History: created by wgrib2
GRIB2_grid_template: 30
请注意，数据集中从.metpy.parse_cf()中标识的metpy_crs为Projection: latitude_longitude，这是不正确的，因为数据集在投影网格空间中具有 2D 纬度/经度坐标和 1D y/x 坐标。这可能是由于数据集的元数据不符合 CF 标准而触发的。现在，正因为如此，MetPy 将您的x坐标视为经度，y视为纬度(请注意这些坐标中的范围)，即使它们不是。这里的结果绝对远非理想(没有错误消息或任何东西)，所以我建议在 MetPy 的问题跟踪器上提交一个问题。
要直接解决问题，您需要手动指定数据的 CRS/投影。为此(在 MetPy v1.0 及更高版本中)，请在 MetPy 的数据集访问器上使用assign_crs方法，而不是使用parse_cf(仅适用于符合 CF 的数据集)。