0

I'm looking to show surface wind data from "https://cds.climate.copernicus.eu/cdsapp#!/dataset/reanalysis-era5-land?tab=form"

I'm using the code below to unpack the data in the form of a diagram "https://confluence.ecmwf.int/display/CKB/How+to+plot+GRIB+files+with+Python+and+matplotlib"

However I get the error 'TypeError: Dimensions of C (1801, 3600) are incompatible with X (3600) and/or Y (1801); see help(pcolormesh)', which should work as C is (rows, columns) and X represents columns and Y represents rows, but this error suggests the data does not fit?

The code is below, please any suggestions would be massively appreciated, thanks!

import pygrib
import matplotlib.pyplot as plt
import matplotlib.colors as colors
from mpl_toolkits.basemap import Basemap
from mpl_toolkits.basemap import shiftgrid
import numpy as np
 
plt.figure(figsize=(12,8))
 
grib = 'adaptor.mars.internal-1669570066.0148444-9941-17-702abb2d-e37e-4ef7-a19e-a04fb24e5a20.grib' # Set the file name of your input GRIB file
grbs = pygrib.open(grib)
 
grb = grbs.select()[0]
data = grb.values
 
# need to shift data grid longitudes from (0..360) to (-180..180)
lons = np.linspace(float(grb['longitudeOfFirstGridPointInDegrees']), \
float(grb['longitudeOfLastGridPointInDegrees']), int(grb['Ni']) )
lats = np.linspace(float(grb['latitudeOfFirstGridPointInDegrees']), \
float(grb['latitudeOfLastGridPointInDegrees']), int(grb['Nj']) )
data, lons = shiftgrid(180., data, lons, start=False)
grid_lon, grid_lat = np.meshgrid(lons, lats) #regularly spaced 2D grid
 
m = Basemap(projection='cyl', llcrnrlon=-180, \
    urcrnrlon=180.,llcrnrlat=lats.min(),urcrnrlat=lats.max(), \
    resolution='c')
    
m.drawcoastlines()
m.drawmapboundary()
m.drawparallels(np.arange(-90.,120.,30.),labels=[1,0,0,0])
m.drawmeridians(np.arange(-180.,180.,60.),labels=[0,0,0,1])
    
x, y = m(grid_lon, grid_lat)

#x = np.transpose(x)
#y = np.transpose(y)

cs = m.pcolormesh(x,y,data,shading='flat',cmap=plt.cm.gist_stern_r)
 
plt.colorbar(cs,orientation='vertical', shrink=0.5)
plt.title('CAMS AOD forecast') # Set the name of the variable to plot
plt.savefig(grib+'.png') # Set the output file name

However, even if I transpose the x and y data I get the same error type but: X (1801) and/or Y (3600). What ideally should happen is a map loads with North-South and East-West wind surface speeds.

Numpy pcolormesh: TypeError: Dimensions of C are incompatible with X and/or Y

Above is a similar question. However, the suggested transposing data solution which worked in that case did not work here.

Azhar Khan
  • 3,829
  • 11
  • 26
  • 32
Samuel B
  • 13
  • 3
  • The x and the y are the *boundaries* of your cells. So, there should be 3601 boundaries if there are 3600 cells in the x direction. – JohanC Dec 01 '22 at 10:26
  • @JohanC Thanks for the response, would you by any chance know how to add these extra columns and rows? Many thanks! Sam – Samuel B Dec 02 '22 at 01:06
  • Now you do `lons = np.linspace(..., int(grb['Ni']) )`. It would make sense to add 1: `lons = np.linspace(..., int(grb['Ni']) + 1)`. And similar for `lats`. – JohanC Dec 02 '22 at 06:32
  • @JohanC Hi, thanks very much for the suggestion, it looks like the approach works on Y, changing it to 1802 (from 1801), however modifying the X value this way instead (int(grb['Ni']) + 1) changes the lons array from (-180, ... , 180) to (0 , ... , 360), with grid_lat and grid_lon both staying at (1802, 3600)? – Samuel B Dec 02 '22 at 08:55

0 Answers0