I want to plot composite map using netCDF4 and python (using basemap) nc file ca

Question

I want to plot composite map using netCDF4 and python (using basemap)

nc file can be found here. https://www.esrl.noaa.gov/psd/cgi-bin/db_search/DBListFiles.pl?did=59&tid=66777&vid=1498

This nc file has an information of location, time and geopotential height for all location on Earth.

I want to plot a map by composing these information.

For example, suppose I downloaded 'hgt.2010.nc' at the link above, so I have an geopotential height information on a daily basis in 2010 on all globe.

Next, I select the date '2010.01.01', '2010.02.01', '2010.03.01'. Then I want to plot one map by just using average of geopotential heights of three days.

How can I do this? Can you help me by giving some code?

Explanation / Answer

Plotting NCEP/NCAR geopotential height in Python with Basemap

#! /usr/bin/env python from mpl_toolkits.basemap import Basemap, shiftgrid, cm import numpy as np import matplotlib.pyplot as plt from netCDF4 import Dataset import datetime datetime_origin = datetime.datetime(1, 1, 1, 0, 0, 0, 0) # In this case, this is 500 mb. Can check below with d.variables["level"][:] level = 5 time = 0 # 277830 is because this is a ~2.5 x 2.5 degree grid; 111132 meters are in about 111 km appart. This is fairly rough. meters_per_grid = 277830 # Download the dataset from: ftp://ftp.cdc.noaa.gov/Datasets/ncep.reanalysis.derived/pressure/hgt.mon.mean.nc d = Dataset("hgt.mon.mean.nc") hgt = d.variables["hgt"][:][time, level, :, :] # Not sure about the dims, here. lon = d.variables["lon"][:] # need to reverse direction of lat dimension so it's increasing. lat = d.variables["lat"][:][::-1] hgt = hgt[::-1, :] hgt, lon = shiftgrid(180, hgt, lon, start = False) fig = plt.figure() ax = fig.add_axes([0.1, 0.1, 0.8, 0.8]) # Lambert Azimuthal Equal Area # Larger area: #m = Basemap(width = 12000000, height = 12000000, resolution = "l", projection = "laea", lat_ts = 50, lat_0 = 50, lon_0 = -107.0) # Smaller area: #m = Basemap(width = 12000000, height = 8000000, resolution = "l", projection = "laea", lat_ts = 50, lat_0 = 50, lon_0 = -107.0) # For north Polar Stereographic projection m = Basemap(projection='npstere',boundinglat=10,lon_0=270,resolution='l') nx = int((m.xmax - m.xmin)/meters_per_grid); ny = int((m.ymax - m.ymin)/meters_per_grid) hgt = m.transform_scalar(hgt, lon, lat, nx, ny) im = m.imshow(hgt, interpolation = "none") m.drawcoastlines() parallels = np.arange(-90, 90, 30) meridians = np.arange(-180, 180, 60) m.drawparallels(parallels, labels = [1, 0, 0, 1]) m.drawmeridians(meridians, labels = [1, 0, 0, 1]) cb = m.colorbar(im, "right", size = "5%", pad = "2%") ax.set_title("500 mb Geopotential Height") plt.show()

Related Questions

Navigate

• Browse (All)
• Browse I
• Subjects

---

---

Integrity-first tutoring: explanations and feedback only — we do not complete graded work. Learn more.