/ Published in: Python
Expand |
Embed | Plain Text
Copy this code and paste it in your HTML
from mpl_toolkits.basemap import Basemap from scipy.io.numpyio import fread import matplotlib.pyplot as plt import matplotlib.numerix.ma as M import matplotlib.colors as C from pylab import * import numpy as np import scipy as sc import os, fnmatch import array import sys import datetime d=str(datetime.date.today()) d=d.replace('-', '') path='c:/' #sys.argv[1] for fn in os.listdir(path): if fnmatch.fnmatch(fn, d + '*.b*'): fileName=path + fn break masking=True treshold=0 savename=d + '_SO2.png' dpi=120 smoothness=1 projection='mill' titre='SO2' legende='Delta Brightness Temperature (K)' centre=[-90,0] print('reading binary data') Lat=[] Lon=[] Val=[] nbreligne=long(os.stat(fileName)[6])/(8*int(fileName[-2:])) rawfile=np.fromfile(open(fileName,'rb'),'d',-1) Lat=rawfile[0:nbreligne] Lon=rawfile[nbreligne:nbreligne*2] Val=rawfile[nbreligne*21:nbreligne*22] if(masking==True): i=0 while(i<nbreligne): if(Val[i]>-1.1): Val[i]=-1.1 i+=1 print('plotting data') map=Basemap(projection='mill',lat_0=-90,lon_0=0,llcrnrlat=-90,urcrnrlat=90,\ urcrnrlon=180,llcrnrlon=-180,resolution='i',area_thresh=30000.) xi=np.linspace(-180,180,360*smoothness) #1440 yi=np.linspace(-90,90,180*smoothness) #720 zi=griddata(Lon,Lat,Val,xi,yi) topodat=map.transform_scalar(zi,xi,yi,720,1440) map.imshow(topodat,cm.winter, interpolation='bilinear',aspect='auto',norm = C.Normalize(vmin = -5., vmax = -1.2, clip = True),extent=[-180,180,-90,90],filternorm=1) cb=plt.colorbar(shrink=0.7) cb.ax.set_ylabel(legende,fontsize=11) for t in cb.ax.get_yticklabels(): t.set_fontsize(7) meridians = arange(-180,180,60) parallels = arange(-90,90,30) map.drawparallels(parallels,labels=[1,0,0,0],fontsize=7,linewidth=0.25) map.drawmeridians(meridians,labels=[0,0,0,1],fontsize=7,linewidth=0.25) title(titre) print('saving map') map.drawcoastlines(0.25,antialiased=1) plt.savefig(savename,dpi=dpi) print('done')