Commit 8e732cc8 authored by Richard Hofmeister's avatar Richard Hofmeister

barentssea: added woa S,T interpolation

parent d9b8725b
import netCDF4
from pylab import *
from datetime import datetime
import sys
lonlat = {'barents_sea':(45.0,73.0),'greenland_shelf':(-14.,77.5), 'nansen_basin':(40.0,84.0)}
if len(sys.argv)>3:
try:
stationname = sys.argv[1]
lonlat[stationname] = (float(sys.argv[2]),float(sys.argv[3]))
except:
print('usage: get_st_profiles.py stationname lon lat')
sys.exit()
elif len(sys.argv)>1:
stationname=sys.argv[1]
if stationname not in lonlat:
print('known station names: %s'%str(lonlat.keys()))
sys.exit()
else:
stationname='barents_sea'
nc = netCDF4.Dataset('/work/gg0877/KST/MiMeMo/woa/woa_arctic_1.00.nc')
ncv = nc.variables
def get_ij(nc,lon,lat):
nlon=nc.variables['lon'][:]
nlat=nc.variables['lat'][:]
i = argmin(abs(nlon-lon))
j = argmin(abs(nlat-lat))
return i,j
startdate = datetime(1990,1,15,0,0,0)
enddate = datetime(2018,1,20,0,0,0)
for station in [stationname]:#lonlat:
lon,lat = lonlat[station]
i,j = get_ij(nc,lon,lat)
s = nc.variables['s_mn'][:,:,j,i]
t = nc.variables['t_mn'][:,:,j,i]
d = nc.variables['depth'][:]
sf = open('s_profile_woa_%s.dat'%station,'w')
tf = open('t_profile_woa_%s.dat'%station,'w')
currdate=startdate
while currdate < enddate:
midx = currdate.month-1
snum = min(where(s[midx,:].mask)[0])
tnum = min(where(t[midx,:].mask)[0])
num = min(snum,tnum)
sf.write('%s %d 2\n'%(str(currdate),num))
tf.write('%s %d 2\n'%(str(currdate),num))
for n in range(num):
sf.write(' %0.1f %0.2f\n'%(-d[n],s[midx,n]))
tf.write(' %0.1f %0.2f\n'%(-d[n],t[midx,n]))
if currdate.month==12:
currdate = datetime(currdate.year+1,1,15,0,0,0)
else:
currdate = datetime(currdate.year,currdate.month+1,15,0,0,0)
sf.close()
tf.close()
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment