2d lat/lon 数据为什么我的 pcolor 图没有绘制数据
2d lat/lon data why is my pcolor plot not drawing the data
希望这是我天真地忽略的东西,但我无法弄清楚为什么我的 pcolor 图没有显示我的任何数据。我有 2D lat/lon,其中一个时间步长是 2D 海冰浓度。任何提示将不胜感激这是我在下面尝试的:
from netCDF4 import Dataset
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
# ifile='upper_box.nc' # doesn't produce any plot
ifile='1979_sfc_out.nc'
fh = Dataset(ifile,mode='r')
lons=fh.variables['lon'][:]
lats=fh.variables['lat'][:]
seaice=fh.variables['SEAICE'][:]
seaice_units = fh.variables['SEAICE'].units
fh.close()
# set basemap instance, specifying our desired map and projection settings
lon_0 = lons.mean()
lat_0 = lats.mean()
m = Basemap(width=5000000,height=3500000,
resolution='l',projection='stere',\
lat_ts=40,lat_0=lat_0,lon_0=lon_0)
# Plot Data
# one timestep
seaice_firstTimestep = seaice[0,:,:]
cs = m.pcolor(lons,lats,seaice_firstTimestep)
# Add Grid Lines
#m.drawparallels(np.arange(-80., 81., 10.), labels=[1,0,0,0], fontsize=10)
#m.drawmeridians(np.arange(-180., 181., 10.), labels=[0,0,0,1], fontsize=10)
# Add Coastlines, States, and Country Boundaries
m.drawcoastlines()
m.drawstates()
m.drawcountries()
# Add Colorbar
cbar = m.colorbar(cs, location='bottom', pad="10%")
cbar.set_label(seaice_units)
# Add Title
plt.title('Sea ice upper box')
plt.show()
我的数据在那里,颜色条的范围是正确的:
>>> seaice_firstTimestep[70:80,70:80]
array([[ 0. , 0. , 0. , 0. , 0. ,
0. , 0. , 0. , 0. , 0.90625 ],
[ 0. , 0. , 0. , 0. , 0. ,
0. , 0. , 0.9140625, 0.9140625, 0.9140625],
[ 0. , 0.9375 , 0.9375 , 0.9296875, 0.9296875,
0.921875 , 0.9140625, 0.9140625, 0.9140625, 0.9140625],
[ 0.953125 , 0.9453125, 0.9375 , 0.9375 , 0.921875 ,
0.921875 , 0.9140625, 0.9140625, 0.9140625, 0.9140625],
[ 0.9453125, 0.9453125, 0.9453125, 0.9375 , 0.9296875,
0.921875 , 0.9140625, 0.9140625, 0.9140625, 0.9140625],
[ 0.9453125, 0.9453125, 0.9453125, 0.9453125, 0.9375 ,
0.9296875, 0.921875 , 0.921875 , 0.921875 , 0.921875 ],
[ 0.9453125, 0.953125 , 0.953125 , 0.9453125, 0.9375 ,
0.9375 , 0.9296875, 0.9296875, 0.9296875, 0.9296875],
[ 0.953125 , 0.953125 , 0.953125 , 0.9453125, 0.9453125,
0.9375 , 0.9375 , 0.9375 , 0.9375 , 0.9375 ],
[ 0.9453125, 0.953125 , 0.953125 , 0.9453125, 0.9453125,
0.9453125, 0.9453125, 0.9453125, 0.9453125, 0.9453125],
[ 0.9453125, 0.9453125, 0.953125 , 0.9453125, 0.9453125,
0.9453125, 0.9453125, 0.9453125, 0.9453125, 0.9453125]], dtype=float32)
>>>
我的 ncdump 在这里,也许它不知道在哪里绘制海冰浓度,因为我的尺寸是 x,y 而不是 lat/lon?
netcdf 79_sfc_out {
dimensions:
x = 83 ;
y = 94 ;
time = UNLIMITED ; // (8736 currently)
nv4 = 4 ;
variables:
float time(time) ;
time:axis = "T" ;
time:long_name = "time" ;
time:standard_name = "time" ;
time:units = "hours since 1979-1-2 00:00:00" ;
time:calendar = "standard" ;
float x(x) ;
x:axis = "x" ;
x:long_name = "X-coordinate in Cartesian system" ;
x:standard_name = "projection_x_coordinate" ;
x:units = "meters" ;
float y(y) ;
y:axis = "y" ;
y:long_name = "Y-coordinate in Cartesian system" ;
y:standard_name = "projection_y_coordinate" ;
y:units = "meters" ;
float lon(y, x) ;
lon:units = "degrees_east" ;
lon:valid_range = -180., 180. ;
lon:standard_name = "longitude" ;
lon:bounds = "lon_bnds" ;
float lat(y, x) ;
lat:units = "degrees_north" ;
lat:valid_range = -90., 90. ;
lat:standard_name = "latitude" ;
lat:bounds = "lat_bnds" ;
float lon_bnds(y, x, nv4) ;
lon_bnds:units = "degreesE" ;
float lat_bnds(y, x, nv4) ;
lat_bnds:units = "degreesN" ;
char mapping ;
mapping:false_easting = 0. ;
mapping:false_northing = 0. ;
mapping:grid_mapping_name = "polar_stereographic" ;
mapping:latitude_of_projection_origin = 90. ;
mapping:standard_parallel = 64. ;
mapping:straight_vertical_longitude_from_pole = -152. ;
mapping:semi_major_axis = 6370000. ;
mapping:semi_minor_axis = 6370000. ;
float SEAICE(time, y, x) ;
SEAICE:_FillValue = -9999.f ;
SEAICE:units = "fraction" ;
SEAICE:long_name = "Ice concentration (ice=1;no ice=0)" ;
SEAICE:grid_mapping = "mapping" ;
SEAICE:coordinates = "lon lat" ;
float U10(time, y, x) ;
U10:_FillValue = -9999.f ;
U10:units = "m/s" ;
U10:long_name = "U-component of wind at 10m height" ;
U10:grid_mapping = "mapping" ;
U10:coordinates = "lon lat" ;
float V10(time, y, x) ;
V10:_FillValue = -9999.f ;
V10:units = "m/s" ;
V10:long_name = "V-component of wind at 10m height" ;
V10:grid_mapping = "mapping" ;
V10:coordinates = "lon lat" ;
在绘图之前,您必须将经纬度网格转换为投影经纬度,如下所示:
m = Basemap(width=5000000,height=3500000,
resolution='l',projection='stere',\
lat_ts=40,lat_0=lat_0,lon_0=lon_0)
LON,LAT = m(lons,lats)
然后在绘图时使用新的 LAT、LON 以适应投影:
cs = m.pcolor(LON,LAT,seaice_firstTimestep)
有关详细信息,请阅读本手册:https://matplotlib.org/basemap/users/mapcoords.html
from netCDF4 import Dataset
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
# ifile='upper_box.nc' # doesn't produce any plot
ifile='1979_sfc_out.nc'
fh = Dataset(ifile,mode='r')
lons=fh.variables['lon'][:]
lats=fh.variables['lat'][:]
seaice=fh.variables['SEAICE'][:]
seaice_units = fh.variables['SEAICE'].units
fh.close()
# set basemap instance, specifying our desired map and projection settings
lon_0 = lons.mean()
lat_0 = lats.mean()
m = Basemap(width=5000000,height=3500000,
resolution='l',projection='stere',\
lat_ts=40,lat_0=lat_0,lon_0=lon_0)
# Plot Data
# one timestep
seaice_firstTimestep = seaice[0,:,:]
cs = m.pcolor(lons,lats,seaice_firstTimestep)
# Add Grid Lines
#m.drawparallels(np.arange(-80., 81., 10.), labels=[1,0,0,0], fontsize=10)
#m.drawmeridians(np.arange(-180., 181., 10.), labels=[0,0,0,1], fontsize=10)
# Add Coastlines, States, and Country Boundaries
m.drawcoastlines()
m.drawstates()
m.drawcountries()
# Add Colorbar
cbar = m.colorbar(cs, location='bottom', pad="10%")
cbar.set_label(seaice_units)
# Add Title
plt.title('Sea ice upper box')
plt.show()
我的数据在那里,颜色条的范围是正确的:
>>> seaice_firstTimestep[70:80,70:80]
array([[ 0. , 0. , 0. , 0. , 0. ,
0. , 0. , 0. , 0. , 0.90625 ],
[ 0. , 0. , 0. , 0. , 0. ,
0. , 0. , 0.9140625, 0.9140625, 0.9140625],
[ 0. , 0.9375 , 0.9375 , 0.9296875, 0.9296875,
0.921875 , 0.9140625, 0.9140625, 0.9140625, 0.9140625],
[ 0.953125 , 0.9453125, 0.9375 , 0.9375 , 0.921875 ,
0.921875 , 0.9140625, 0.9140625, 0.9140625, 0.9140625],
[ 0.9453125, 0.9453125, 0.9453125, 0.9375 , 0.9296875,
0.921875 , 0.9140625, 0.9140625, 0.9140625, 0.9140625],
[ 0.9453125, 0.9453125, 0.9453125, 0.9453125, 0.9375 ,
0.9296875, 0.921875 , 0.921875 , 0.921875 , 0.921875 ],
[ 0.9453125, 0.953125 , 0.953125 , 0.9453125, 0.9375 ,
0.9375 , 0.9296875, 0.9296875, 0.9296875, 0.9296875],
[ 0.953125 , 0.953125 , 0.953125 , 0.9453125, 0.9453125,
0.9375 , 0.9375 , 0.9375 , 0.9375 , 0.9375 ],
[ 0.9453125, 0.953125 , 0.953125 , 0.9453125, 0.9453125,
0.9453125, 0.9453125, 0.9453125, 0.9453125, 0.9453125],
[ 0.9453125, 0.9453125, 0.953125 , 0.9453125, 0.9453125,
0.9453125, 0.9453125, 0.9453125, 0.9453125, 0.9453125]], dtype=float32)
>>>
我的 ncdump 在这里,也许它不知道在哪里绘制海冰浓度,因为我的尺寸是 x,y 而不是 lat/lon?
netcdf 79_sfc_out {
dimensions:
x = 83 ;
y = 94 ;
time = UNLIMITED ; // (8736 currently)
nv4 = 4 ;
variables:
float time(time) ;
time:axis = "T" ;
time:long_name = "time" ;
time:standard_name = "time" ;
time:units = "hours since 1979-1-2 00:00:00" ;
time:calendar = "standard" ;
float x(x) ;
x:axis = "x" ;
x:long_name = "X-coordinate in Cartesian system" ;
x:standard_name = "projection_x_coordinate" ;
x:units = "meters" ;
float y(y) ;
y:axis = "y" ;
y:long_name = "Y-coordinate in Cartesian system" ;
y:standard_name = "projection_y_coordinate" ;
y:units = "meters" ;
float lon(y, x) ;
lon:units = "degrees_east" ;
lon:valid_range = -180., 180. ;
lon:standard_name = "longitude" ;
lon:bounds = "lon_bnds" ;
float lat(y, x) ;
lat:units = "degrees_north" ;
lat:valid_range = -90., 90. ;
lat:standard_name = "latitude" ;
lat:bounds = "lat_bnds" ;
float lon_bnds(y, x, nv4) ;
lon_bnds:units = "degreesE" ;
float lat_bnds(y, x, nv4) ;
lat_bnds:units = "degreesN" ;
char mapping ;
mapping:false_easting = 0. ;
mapping:false_northing = 0. ;
mapping:grid_mapping_name = "polar_stereographic" ;
mapping:latitude_of_projection_origin = 90. ;
mapping:standard_parallel = 64. ;
mapping:straight_vertical_longitude_from_pole = -152. ;
mapping:semi_major_axis = 6370000. ;
mapping:semi_minor_axis = 6370000. ;
float SEAICE(time, y, x) ;
SEAICE:_FillValue = -9999.f ;
SEAICE:units = "fraction" ;
SEAICE:long_name = "Ice concentration (ice=1;no ice=0)" ;
SEAICE:grid_mapping = "mapping" ;
SEAICE:coordinates = "lon lat" ;
float U10(time, y, x) ;
U10:_FillValue = -9999.f ;
U10:units = "m/s" ;
U10:long_name = "U-component of wind at 10m height" ;
U10:grid_mapping = "mapping" ;
U10:coordinates = "lon lat" ;
float V10(time, y, x) ;
V10:_FillValue = -9999.f ;
V10:units = "m/s" ;
V10:long_name = "V-component of wind at 10m height" ;
V10:grid_mapping = "mapping" ;
V10:coordinates = "lon lat" ;
在绘图之前,您必须将经纬度网格转换为投影经纬度,如下所示:
m = Basemap(width=5000000,height=3500000,
resolution='l',projection='stere',\
lat_ts=40,lat_0=lat_0,lon_0=lon_0)
LON,LAT = m(lons,lats)
然后在绘图时使用新的 LAT、LON 以适应投影:
cs = m.pcolor(LON,LAT,seaice_firstTimestep)
有关详细信息,请阅读本手册:https://matplotlib.org/basemap/users/mapcoords.html