I am getting a very strange error using basemap. No error appears, yet my 3rd plot has no data plotted when data does indeed exist. Below is my code. When run, you will see that both modis and seawifs data is plotted, but viirs is not. I can't figure out why.
import numpy as np
import urllib
import urllib2
import netCDF4
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
from datetime import datetime, date, time, timedelta
import json
import math
def indexsearch(datebroken,year, month, day):
for i in range(0,len(datebroken)):
if (datebroken[i,0] == year and datebroken[i,1] == month and datebroken[i,2] == day):
return i
url = 'http://coastwatch.pfeg.noaa.gov/erddap/griddap/erdMWchlamday.nc?chlorophyll' +\
'[(2002-07-16T12:00:00Z):1:(2015-04-16T00:00:00Z)][(0.0):1:(0.0)][(36):1:(39)][(235):1:(240)]'
file = 'erdMWchlamday.nc'
urllib.urlretrieve(url, file)
ncfilemod = netCDF4.Dataset(file)
ncv1 = ncfilemod.variables
print ncv1.keys()
time1=ncv1['time'][:]
inceptiondate = datetime(1970, 1, 1, 0, 0, 0)
timenew1=[]
for i in time1[:]:
newdate = inceptiondate + timedelta(seconds=i)
timenew1.append(newdate.strftime('%Y%m%d%H'))
datebroken1 = np.zeros((len(timenew1),4),dtype=int)
for i in range(0,len(timenew1)):
datebroken1[i,0] = int(timenew1[i][0:4])
datebroken1[i,1] = int(timenew1[i][4:6])
datebroken1[i,2] = int(timenew1[i][6:8])
datebroken1[i,3] = int(timenew1[i][8:10])
lon1= ncv1['longitude'][:]
lat1 = ncv1['latitude'][:]
lons1, lats1 = np.meshgrid(lon1,lat1)
chla1 = ncv1['chlorophyll'][:,0,:,:]
url = 'http://coastwatch.pfeg.noaa.gov/erddap/griddap/erdSWchlamday.nc?chlorophyll' +\
'[(1997-09-16):1:(2010-12-16T12:00:00Z)][(0.0):1:(0.0)][(36):1:(39)][(235):1:(240)]'
file = 'erdSWchlamday.nc'
urllib.urlretrieve(url, file)
#Ncfile 2
ncfilewif = netCDF4.Dataset(file)
ncv2 = ncfilewif.variables
print ncv2.keys()
time2=ncv2['time'][:]
inceptiondate = datetime(1970, 1, 1, 0, 0, 0)
timenew2=[]
for i in time2[:]:
newdate = inceptiondate + timedelta(seconds=i)
timenew2.append(newdate.strftime('%Y%m%d%H'))
datebroken2 = np.zeros((len(timenew2),4),dtype=int)
for i in range(0,len(timenew2)):
datebroken2[i,0] = int(timenew2[i][0:4])
datebroken2[i,1] = int(timenew2[i][4:6])
datebroken2[i,2] = int(timenew2[i][6:8])
datebroken2[i,3] = int(timenew2[i][8:10])
lon2= ncv2['longitude'][:]
lat2 = ncv2['latitude'][:]
lons2, lats2 = np.meshgrid(lon2,lat2)
chla2 = ncv2['chlorophyll'][:,0,:,:]
url = 'http://coastwatch.pfeg.noaa.gov/erddap/griddap/erdVH2chlamday.nc?chla' +\
'[(2012-01-15):1:(2015-05-15T00:00:00Z)][(39):1:(36)][(-125):1:(-120)]'
file = 'erdVH2chlamday.nc'
urllib.urlretrieve(url, file)
ncfileviir = netCDF4.Dataset(file)
ncv3 = ncfileviir.variables
print ncv3.keys()
time3=ncv3['time'][:]
inceptiondate = datetime(1970, 1, 1, 0, 0, 0)
timenew3=[]
for i in time3[:]:
newdate = inceptiondate + timedelta(seconds=i)
timenew3.append(newdate.strftime('%Y%m%d%H'))
datebroken3 = np.zeros((len(timenew3),4),dtype=int)
for i in range(0,len(timenew3)):
datebroken3[i,0] = int(timenew3[i][0:4])
datebroken3[i,1] = int(timenew3[i][4:6])
datebroken3[i,2] = int(timenew3[i][6:8])
datebroken3[i,3] = int(timenew3[i][8:10])
lon3= ncv3['longitude'][:]
lat3 = ncv3['latitude'][:]
lons3, lats3 = np.meshgrid(lon3,lat3)
chla3 = ncv3['chla'][:,:,:]
i1=indexsearch(datebroken1,2012,6,16)
print i1
i2=indexsearch(datebroken2,2010,6,16)
print i2
i3=indexsearch(datebroken3,2012,6,15)
print i3
chla1plot = chla1[i1,:,:]
chla2plot = chla2[i2,:,:]
chla3plot = chla3[i3,:,:]
ncfileviir.close()
ncfilemod.close()
ncfilewif.close()
Important code is below here. All code above is just pulling the data into python to plot.
minlat = 36
maxlat = 39
minlon = 235
maxlon = 240
# Create map
fig = plt.figure()
#####################################################################################################################
#plot figure 1
ax1 = fig.add_subplot(221)
m = Basemap(projection='merc', llcrnrlat=minlat,urcrnrlat=maxlat,llcrnrlon=minlon, urcrnrlon=maxlon,resolution='h')
cs1 = m.pcolormesh(lons1,lats1,chla1plot,cmap=plt.cm.jet,latlon=True)
m.drawcoastlines()
m.drawmapboundary()
m.fillcontinents()
m.drawcountries()
m.drawstates()
m.drawrivers()
#Sets up parallels and meridians.
parallels = np.arange(36.,39,1.)
# labels = [left,right,top,bottom]
m.drawparallels(parallels,labels=[False,True,True,False])
meridians = np.arange(235.,240.,1.)
m.drawmeridians(meridians,labels=[True,False,False,True])
ax1.set_title('Modis')
#####################################################################################################################
#plot figure 2
ax2 = fig.add_subplot(222)
cs2 = m.pcolormesh(lons2,lats2,chla2plot,cmap=plt.cm.jet,latlon=True)
m.drawcoastlines()
m.drawmapboundary()
m.fillcontinents()
m.drawcountries()
m.drawstates()
m.drawrivers()
#Sets up parallels and meridians.
parallels = np.arange(36.,39,1.)
# labels = [left,right,top,bottom]
m.drawparallels(parallels,labels=[False,True,True,False])
meridians = np.arange(235.,240.,1.)
m.drawmeridians(meridians,labels=[True,False,False,True])
ax2.set_title('SeaWIFS')
#####################################################################################################################
#plot figure 3
ax3 = fig.add_subplot(223)
cs3 = m.pcolormesh(lons3,np.flipud(lats3),np.flipud(chla3plot),cmap=plt.cm.jet,latlon=True)
m.drawcoastlines()
m.drawmapboundary()
m.fillcontinents()
m.drawcountries()
m.drawstates()
m.drawrivers()
#Sets up parallels and meridians.
parallels = np.arange(36.,39,1.)
# labels = [left,right,top,bottom]
m.drawparallels(parallels,labels=[False,True,True,False])
meridians = np.arange(235.,240.,1.)
m.drawmeridians(meridians,labels=[True,False,False,True])
ax3.set_title('VIIRS')
# Save figure (without 'white' borders)
#plt.savefig('SSTtest.png', bbox_inches='tight')
plt.show()
My results are shown here!
![results]: https://i.stack.imgur.com/dRjkU.png
The issue that I found was that I had
The final plot was -125 to -120 which basemap did not automatically handle, but instead placed the plot at an area where I did not have data. I added a new m = basemap statement and changed the meridian numbers for the third graph using -125 to -120 as my longitude and the graph plotted just fine.