Python Interpolation with matplotlib/basemap

Question

I am rather new to programming and am having a very hard time understanding interpolation. Every single source I can find that attempts to explain it is extremely cryptic (especially the package specific sites for basemap/matplotlib). I am mapping using matplotlib's basemap however the nature of my data is that it comes in 5 degree by 5 degree blocks (lat lon blocks). I want to smooth out the map by interpolation.

So first here is my code.

from netCDF4 import Dataset
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap, addcyclic

#load the netcdf file into a variable
mar120="C:/Users/WillEvo/Desktop/sec_giptie_cpl_mar_120.nc"

#grab the data into a new variable
fh=Dataset(mar120,mode="r")

#assign model variable contents to python variables
lons=fh.variables['lon'][:]
lats=fh.variables['lat'][:]
test=fh.variables['NE'][:]

#specifying which time and elevation to map
ionst=test[12,0]

#close the netCDF file
fh.close()

# get rid of white stripe on map
ionst, lons=addcyclic(ionst, lons)

#map settings
m=Basemap(llcrnrlon=-180, llcrnrlat=-87.5, urcrnrlon=180, urcrnrlat=87.5,rsphere=6467997, resolution='i', projection='cyl',area_thresh=10000, lat_0=0, lon_0=0)

#Creating 2d array of latitude and longitude
lon, lat=np.meshgrid(lons, lats)
xi, yi=m(lon, lat)

#setting plot type and which variable to plot
cs=m.pcolormesh(xi,yi,np.squeeze(ionst))

#drawing grid lines
m.drawparallels(np.arange(-90.,90.,30.),labels=[1,0,0,0],fontsize=10)
m.drawmeridians(np.arange(-180.,181.,30.), labels=[0,0,0,1],fontsize=10)

#drawing coast lines
m.drawcoastlines()

#color bar
cbar=m.colorbar(cs, location='bottom', pad="10%")
cbar.set_label("Elecron Density cm-3")

#showing the plot
plt.show()

So now, how can I easily interpolate my data to smooth it out? I have tried to call Basemap.interp however I get an error saying that basemap has no attribute interp.

I am really impartial to what I use to interpolate the data, I just really need someone to explain this to me like I am dumb.

Also note that I am in the process of learning to map so details like labels and such I am not too worried about as of yet. Below is an example map that the code above outputs.

enter image description here

If you don't need to use Python, a good (and easy) way to regrid is through NCL's functions: ncl.ucar.edu/Applications/regrid.shtml — N1B4
– N1B4, Commented May 28, 2015 at 16:09
NOTICE: the line cs=m.pcolormesh(xi, yi, np.squeeze(smooth)) Should have ionst in place of smooth....that smooth was an attempt to interpolate and I forgot to remove the variable — Will.Evo
– Will.Evo, Commented May 28, 2015 at 16:16
Where can I download this file? sec_giptie_cpl_mar_120.nc Thanks. — tommy.carstensen
– tommy.carstensen, Commented Jan 21, 2016 at 1:15

Anthony Lethuillier · Accepted Answer · 2015-05-28 15:38:50Z

4

To smooth things out I would use imshow rather than pcolormesh

For example :

from pylab import *

data = random((3,3))
figure(1)
imshow(data, interpolation='none')

plt.show()

gives : enter image description here

and

imshow(data, interpolation='bicubic')

gives :

enter image description here

The help page gives a list of all possible interpolations : http://matplotlib.org/api/pyplot_api.html#matplotlib.pyplot.imshow

answered May 28, 2015 at 15:38

Anthony Lethuillier

1,5394 gold badges17 silver badges38 bronze badges

Sign up to request clarification or add additional context in comments.

4 Comments

Will.Evo Over a year ago

You are a genius! Holy cow the quality is insanely better without even applying any interpolation. Thanks!!

cfh Over a year ago

@Will.Evo: It actually applies the interpolation internally. Note the interpolation='bicubic' argument.

Will.Evo Over a year ago

Well I didn't explicitly write the interpolation argument, I just changed to imshow and it smoothed out.

Will.Evo Over a year ago

@tommy.carstensen I added an answer to this question with what I think is the relevant code (this was two years ago so I can't be 100% this is exactly that)

Will.Evo · Accepted Answer · 2017-09-19 18:28:46Z

This contains some extra code but this is what I think I ended up with. This was several years ago so I am not 100% certain this is the exact code from above which solved my answer.

from netCDF4 import Dataset
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap, addcyclic
import matplotlib.animation as animation

plt.rcParams['animation.ffmpeg_path'] = 'C:/FFMPEG/bin/ffmpeg'

file_loc = "C:/Users/Will Evonosky/Dropbox/SOARS/SOARS 2015/Data"
#load the netcdf file into a variable
mar120=file_loc+"/My Datasets To Share/SA120_Iono_Acc_WE.nc"

#grab the data into a new variable
fh=Dataset(mar120,mode="r")

#assign model variable contents to python variables
lons=fh.variables['lon'][:]
lats=fh.variables['lat'][:]
var1=fh.variables['GT'][:]

#specifying which time and elevation to map
ionst=var1[0,18,:,:]
details='(Z=6)'

#close the netCDF file
fh.close()



# get rid of white stripe on map
ionst, lons=addcyclic(ionst, lons)

#Setting figure attributes
fig=plt.figure(figsize=(15,15),frameon=False)

#map settings
m=Basemap(llcrnrlon=-180, llcrnrlat=-87.5, urcrnrlon=180, urcrnrlat=87.5,rsphere=6467997, resolution='l', projection='cyl',area_thresh=10000, lat_0=0, lon_0=0)

#Creating 2d array of latitude and longitude
lon, lat=np.meshgrid(lons, lats)
xi, yi=m(lon, lat)

#plotting data onto basemap
cs=m.imshow(ionst, interpolation=None, alpha=.8)
vert=plt.axvline(x=-75, color='black', linewidth=5)
#drawing grid lines
m.drawparallels(np.arange(-90.,90.,30.),labels=[1,0,0,0],fontsize=15)
m.drawmeridians(np.arange(-180.,181.,30.), labels=[0,0,0,1],fontsize=15)

#drawing coast lines
m.drawcoastlines()

#color bar
cbar=m.colorbar(cs, location='bottom', pad="10%")
cbar.set_label(r"Ion Drag $(cm/s^2)$", size=15)

#Title Preferences
plt.title('Ion Drag at '+details, size=25)


#Function to update the plots data
def updateax1(j):
    cs.set_array(var1[j,18,:,:])
    return cs,

#Animate the plot
ani1=animation.FuncAnimation(fig, updateax1, frames=range(24), interval=250, blit=True)
ani1.save('Iondrag_Map.mp4')

#showing the plot
plt.show()

Collectives™ on Stack Overflow

Python Interpolation with matplotlib/basemap

2 Answers 2

4 Comments

Comments

Your Answer

Linked

Hot Network Questions

Collectives™ on Stack Overflow

2 Answers 2

4 Comments

Comments

Your Answer

Sign up or log in

Post as a guest

Linked

Related