Extract surface temperature from Argo observations

   1 import scipy.io as io
   2 import glob
   3 from pylab import *
   4 from mpl_toolkits.basemap import Basemap
   5 
   6 tmp_dir='/scratch/clisap/seaice/TMP/u242023/ARGO/'
   7 
   8 
   9 file_liste=glob.glob(tmp_dir+'*.nc')
  10 D={}# Empty dictionary to store selected profiles
  11 for f in file_liste:# Loop over all data
  12 
  13 
  14     # Open netcdf data file
  15     fid=io.netcdf_file(f,'r')
  16 
  17     # Read content into variables
  18     lat=fid.variables['LATITUDE'][:].copy()
  19     lon=fid.variables['LONGITUDE'][:].copy()
  20     T=fid.variables['TEMP'][:].copy()
  21     P=fid.variables['PRES'][:].copy()
  22     
  23     T[T>=99999]=nan # Set 99999.0 to "Not a Number"
  24     P[P>=99999]=nan
  25 
  26     (nr_profs,Z)=T.shape # Get dimension
  27 
  28     for i in range(nr_profs):
  29         D[(lon[i],lat[i])]=(T[i,:],P[i,:])
  30 
  31     # Close data file
  32     fid.close()
  33 
  34 fid=open('lat_lon_T.tab','w')
  35 for k in D.keys():# write position (lat,lon), surface temperature to file
  36     fid.write(str(k[0])+'\t'+str(k[1])+'\t'+str(D[k][0][0])+'\n')
  37 fid.close()
  38 stop
  39 
  40 
  41 # Draw map of positions
  42 m = Basemap(projection='ortho',lon_0=-45,lat_0=0,resolution='l')
  43 m.bluemarble()
  44 for lon,lat in D.keys():
  45     x,y=m(lon,lat) # Coordinate transfer
  46     m.plot(x,y,'r.')
  47     
  48 savefig('argo_position.png',dpi=150)
  49 
  50 # Plot profiles
  51 figure()
  52 for k in D.keys():
  53 #    print k
  54     plot(D[k][0][:],D[k][1][:])
  55 axis([-2,30,2000,0])
  56 xlabel('T')
  57 ylabel('P')
  58 show()
  59 savefig('Argo_plot.png',dpi=75)