I have a xarray dataset called dens and I want to plot it.

This is the data set:

<xarray.dataset>
dimensions:  (time: 641, lat: 30, lon: 30)
coordinates:
  * time     (time) datetime64[ns] 2013-07-01t12:00:00 ... 2013-08-02t12:00:00
  * lon      (lon) float64 32.73 32.83 32.94 33.05 ... 35.53 35.64 35.75 35.85
  * lat      (lat) float64 31.08 31.27 31.47 31.66 ... 36.06 36.25 36.44 36.63
data variables:
    density  (time, lat, lon) float64 2e+03 0.0 0.0 0.0 0.0 ... 0.0 0.0 0.0 0.0

I am using the command

plt.contourf(dens.density.values[-1,:,:]);

plotted it and it was working, but since I wanted the coastline to be plotted on the plot as well, I also tried using

m = basemap(llcrnrlon=data['lon'].min(), llcrnrlat=data['lat'].min(),
                urcrnrlon=data['lon'].max(), urcrnrlat=data['lat'].max(), resolution='i', suppress_ticks=1)
m.drawcoastlines();
m.fillcontinents(color='gray',lake_color='gray');

But when I run all the commands and then run ​​plt.show() the contour plot disappears and all it shows me is the coastline.

How to solve this problem to get a contour map in the same coastline map?

Sorry if this is a stupid question but I'm pretty new to python

thanks for your help,

Yutam

Edit: I just realized now that I was trying to combine two different "toolkits" and it is possible to do all of this using just the basemap toolkit, but just trying to write

m.contourf(dens.density.values[-1,:,:]);

gives me this error:

---------------------------------------------------------------------------
typeerror                                 traceback (most recent call last)
cell in[21], line 1
----> 1 m.contourf(dens.density.values[-1,:,:])

typeerror: basemap.contourf() missing 2 required positional arguments: 'y' and 'data'

Another edit: I keep discovering more and after reading the documentation for basemap I realized the syntax of the command should be like this

m.contourf(dens.lon.values,dens.lat.values,dens.density.values[-1,:,:]);

But now I get this error:

IndexError: too many indices for array: array is 1-dimensional, but 2 were indexed

I'm guessing this is because my density array is 2D, but how do I extract the density values ​​from it? I'm using [-1] in the time dimension because I actually only need the last time step

Thanks again in advance, Yotam

Last edit

This is the final plot, how can I make the surrounding land gray instead of purple? Also, is there a way to describe a larger geographic area, slightly larger, without messing up the data?

Here are the new numbers for my actual data

Correct answer

Plotting using basemaps and xarray has been discussedhere.

m = basemap(llcrnrlon=data['lon'].min(), llcrnrlat=data['lat'].min(),
            urcrnrlon=data['lon'].max(), urcrnrlat=data['lat'].max(), 
  resolution='i', suppress_ticks=1)
m.drawcoastlines();
m.fillcontinents(color='gray',lake_color='gray')
dens.density[-1,:,:].plot.contourf()
plt.show()

import xarray as xr
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.feature as cf

ds = xr.open_dataset('filename.nc')
fig = plt.figure(figsize=(8,8))
crs=ccrs.platecarree()
ax = fig.add_subplot(1,1,1, projection=crs)
gl = ax.gridlines(crs=crs, draw_labels=true,
                linewidth=0.01, color='gray', alpha=0.5, linestyle='-.')

ax.add_feature(cf.coastline.with_scale("50m"), lw=0.5)
ax.add_feature(cf.borders.with_scale("50m"), lw=0.3)

ds.density[-1,:,:].plot.contourf()
plt.show()

from matplotlib.colors import LinearSegmentedColormap
cm = LinearSegmentedColormap.from_list('', ['white', *plt.cm.Blues(np.arange(255))])
ds.density[-1,:,:].plot.contourf(cmap=cm)