Skip to content

GitLab

"...lib/td11/git@git.noc.ac.uk:brecinosrivas/mdf_reader.git" did not exist on "515c10ebdc68f38af2c7ec6e04f262ac07337422"
Commit 218eacca authored by thopri's avatar thopri
Browse files

updated Extract Class Init function

parent eecd7693
Hide whitespace changes
Inline Side-by-side
Showing with 136 additions and 118 deletions
inputs/namelist_cmems.bdy
View file @ 218eacca
......@@ -52,7 +52,7 @@
! CMEMS Data Source Configuration
!------------------------------------------------------------------------------
ln_use_cmems = .true.
ln_download_cmems = .true.
ln_download_cmems = .false.
sn_cmems_dir = '/Users/thopri/Projects/PyNEMO/inputs/' ! where to download CMEMS input files (static and variable)
!------------------------------------------------------------------------------
! CMEMS MOTU Configuration (for Boundary Data)
......@@ -115,8 +115,8 @@
!------------------------------------------------------------------------------
nn_year_000 = 2017 ! year start
nn_year_end = 2017 ! year end
nn_month_000 = 01 ! month start (default = 1 is years>1)
nn_month_end = 04 ! month end (default = 12 is years>1)
nn_month_000 = 02 ! month start (default = 1 is years>1)
nn_month_end = 03 ! month end (default = 12 is years>1)
sn_dst_calendar = 'gregorian' ! output calendar format
nn_base_year = 1960 ! base year for time counter
sn_tide_grid = './src_data/tide/grid_tpxo7.2.nc'
......
pynemo/config/motu_config_template.ini
View file @ 218eacca
......@@ -28,7 +28,7 @@ log_level=20
# block_size block used to download file (integer expressing bytes) default=65535
# block_size=65535
socket_timeout=120000
socket_timeout=120000
# Http proxy to connect to Motu server
# proxy_server=proxy.domain.net:8080
......
pynemo/nemo_bdy_dl_cmems.py
View file @ 218eacca
......@@ -118,14 +118,15 @@ def MWD_request_cmems(args,date_min,date_max,F):
for m in range(len(month_end)):
mnth_dl = request_cmems(args, month_start[m], month_end[m])
if mnth_dl == 0:
logger.info('CMEMS month request ' + str((m + 1)) + 'of' + (str(len(month_end) + 1)) + ' successful')
logger.info('CMEMS month request ' + str((m + 1)) + 'of' + (str(len(month_end))) + ' successful')
if type(mnth_dl) == str:
logger.error(
'CMEMS month request ' + str((m + 1)) + 'of' + str((len(month_end) + 1)) + ' unsuccessful: Error Msg below')
'CMEMS month request ' + str((m + 1)) + 'of' + str((len(month_end))) + ' unsuccessful: Error Msg below')
logger.error(mnth_dl)
return mnth_dl
if mnth_dl == 1:
return 1
if F == 'W':
week_start = pd.date_range(date_min, date_max,
freq='W').strftime("%Y-%m-%d").tolist()
......@@ -136,14 +137,15 @@ def MWD_request_cmems(args,date_min,date_max,F):
for w in range(len(week_end)):
wk_dl = request_cmems(args, week_start[w], week_end[w])
if wk_dl == 0:
logger.info('CMEMS week request ' + str((w + 1)) + 'of' + str((len(week_end) + 1)) + ' successful')
logger.info('CMEMS week request ' + str((w + 1)) + 'of' + str((len(week_end))) + ' successful')
if type(wk_dl) == str:
logger.error(
'CMEMS week request ' + str((m + 1)) + 'of' + str((len(week_end)) + 1) + ' unsuccessful: Error Msg below')
'CMEMS week request ' + str((m + 1)) + 'of' + str((len(week_end))) + ' unsuccessful: Error Msg below')
logger.error(wk_dl)
return wk_dl
if wk_dl == 1:
return 1
if F == 'D':
days = pd.date_range(date_min, date_max,
freq='D').strftime("%Y-%m-%d").tolist()
......@@ -156,13 +158,15 @@ def MWD_request_cmems(args,date_min,date_max,F):
return
if type(dy_dl) == str:
logger.error('CMEMS day request ' + str((d + 1)) + 'of' + (
str(len(days) + 1)) + ' unsuccessful: Error Msg below')
str(len(days))) + ' unsuccessful: Error Msg below')
logger.error(dy_dl)
return dy_dl
else:
if (not F == 'M' or not F == 'W' or not F == 'D'):
time_int_err = 'incorrect character used to define time download interval please use M, W or D'
logger.error(time_int_err)
return time_int_err
return 0
'''
......@@ -246,7 +250,11 @@ def request_cmems(args, date_min, date_max):
logger.info('download of request xml file for variable ' + ' '.join(grids[key]) + ' successful')
xml = locs[key]+args['dl_prefix']+'_'+str(date_min)+'_'+str(date_max)+'_'+str(key)+ '.xml'
root = ET.parse(xml).getroot()
try:
root = ET.parse(xml).getroot()
except ET.ParseError:
return 'Parse Error in XML file, This generally occurs when CMEMS service is down and returns an unexpected XML.'
logger.info('size of request ' + root.attrib['size'])
if 'OK' in root.attrib['msg']:
......
pynemo/nemo_bdy_extr_tm3.py
View file @ 218eacca
......@@ -74,6 +74,14 @@ class Extract:
self.g_type = grd
self.settings = setup
self.key_vec = False
self.t_mask = None
self.u_mask = None
self.v_mask = None
self.dist_wei = None
self.dist_fac = None
self.tmp_valid = None
self.data_ind = None
self.nan_ind = None
# TODO: Why are we deepcopying the coordinates???
......@@ -452,7 +460,7 @@ class Extract:
# Get first and last date within range, init to cover entire range
first_date = 0
last_date = len(sc_time.time_counter) - 1
last_date = len(sc_time.time_counter) - 1
rev_seq = range(len(sc_time.time_counter))
rev_seq.reverse()
for date in rev_seq:
......@@ -471,19 +479,19 @@ class Extract:
# TODO: sort generic mask variable name
try:
varid_3 = nc_3['tmask']
t_mask = varid_3[:1, :sc_z_len, j_run, i_run]
self.t_mask = varid_3[:1, :sc_z_len, j_run, i_run]
except:
varid_3 = nc_3['mask']
varid_3 = np.expand_dims(varid_3, axis=0)
t_mask = varid_3[:1, :sc_z_len, np.min(j_run):np.max(j_run) + 1, np.min(i_run):np.max(i_run) + 1]
self.t_mask = varid_3[:1, :sc_z_len, np.min(j_run):np.max(j_run) + 1, np.min(i_run):np.max(i_run) + 1]
# TODO: Sort out issue with j_run and i_run not broadcasting to varid_3
if self.key_vec:
#varid_3 = nc_3['umask']
varid_3 = nc_3['mask']
u_mask = varid_3[:1, :sc_z_len, j_run, extended_i]
self.u_mask = varid_3[:1, :sc_z_len, j_run, extended_i]
#varid_3 = nc_3['vmask']
varid_3 = nc_3['mask']
v_mask = varid_3[:1, :sc_z_len, extended_j, i_run]
self.v_mask = varid_3[:1, :sc_z_len, extended_j, i_run]
nc_3.close()
# Identify missing values and scale factors if defined
......@@ -547,8 +555,8 @@ class Extract:
# Average vector vars onto T-grid
if self.key_vec:
# First make sure land points have a zero val
sc_alt_arr[0] *= u_mask
sc_alt_arr[1] *= v_mask
sc_alt_arr[0] *= self.u_mask
sc_alt_arr[1] *= self.v_mask
# Average from to T-grid assuming C-grid stagger
sc_array[0] = 0.5 * (sc_alt_arr[0][:,:,:,:-1] +
sc_alt_arr[0][:,:,:,1:])
......@@ -560,7 +568,7 @@ class Extract:
# Note using isnan/sum is relatively fast, but less than
# bottleneck external lib
self.logger.info('SC ARRAY MIN MAX : %s %s', np.nanmin(sc_array[0]), np.nanmax(sc_array[0]))
sc_array[0][t_mask == 0] = np.NaN
sc_array[0][self.t_mask == 0] = np.NaN
self.logger.info( 'SC ARRAY MIN MAX : %s %s', np.nanmin(sc_array[0]), np.nanmax(sc_array[0]))
if not np.isnan(np.sum(meta_data[vn]['sf'])):
sc_array[0] *= meta_data[vn]['sf']
......@@ -568,7 +576,7 @@ class Extract:
sc_array[0] += meta_data[vn]['os']
if self.key_vec:
sc_array[1][t_mask == 0] = np.NaN
sc_array[1][self.t_mask == 0] = np.NaN
if not np.isnan(np.sum(meta_data[vn + 1]['sf'])):
sc_array[1] *= meta_data[vn + 1]['sf']
if not np.isnan(np.sum(meta_data[vn + 1]['os'])):
......@@ -611,7 +619,7 @@ class Extract:
# source data to dest bdy pts. Only need do once.
if self.first:
# identify valid pts
data_ind = np.invert(np.isnan(sc_bdy[0,:,:,:]))
self.data_ind = np.invert(np.isnan(sc_bdy[0,:,:,:]))
# dist_tot is currently 2D so extend along depth
# axis to allow single array calc later, also remove
# any invalid pts using our eldritch data_ind
......@@ -619,10 +627,10 @@ class Extract:
self.dist_tot = (np.repeat(self.dist_tot, sc_z_len).reshape(
self.dist_tot.shape[0],
self.dist_tot.shape[1], sc_z_len)).transpose(2,0,1)
self.dist_tot *= data_ind
self.dist_tot *= self.data_ind
self.logger.info('DIST TOT ZEROS %s', np.sum(self.dist_tot == 0))
self.logger.info('DIST IND ZEROS %s', np.sum(data_ind == 0))
self.logger.info('DIST IND ZEROS %s', np.sum(self.data_ind == 0))
# Identify problem pts due to grid discontinuities
# using dists > lat
......@@ -634,22 +642,22 @@ class Extract:
# Calculate guassian weighting with correlation dist
r0 = self.settings['r0']
dist_wei = (1/(r0 * np.power(2 * np.pi, 0.5)))*(np.exp( -0.5 *np.power(self.dist_tot / r0, 2)))
self.dist_wei = (1/(r0 * np.power(2 * np.pi, 0.5)))*(np.exp( -0.5 *np.power(self.dist_tot / r0, 2)))
# Calculate sum of weightings
dist_fac = np.sum(dist_wei * data_ind, 2)
self.dist_fac = np.sum(self.dist_wei * self.data_ind, 2)
# identify loc where all sc pts are land
nan_ind = np.sum(data_ind, 2) == 0
self.logger.info('NAN IND : %s ', np.sum(nan_ind))
self.nan_ind = np.sum(self.data_ind, 2) == 0
self.logger.info('NAN IND : %s ', np.sum(self.nan_ind))
# Calc max zlevel to which data available on sc grid
data_ind = np.sum(nan_ind == 0, 0) - 1
self.data_ind = np.sum(self.nan_ind == 0, 0) - 1
# set land val to level 1 otherwise indexing problems
# may occur- should not affect later results because
# land is masked in weightings array
data_ind[data_ind == -1] = 0
self.data_ind[self.data_ind == -1] = 0
# transform depth levels at each bdy pt to vector
# index that can be used to speed up calcs
data_ind += np.arange(0, sc_z_len * self.num_bdy, sc_z_len)
self.data_ind += np.arange(0, sc_z_len * self.num_bdy, sc_z_len)
# ? Attribute only used on first run so clear.
del self.dist_tot
......@@ -657,8 +665,8 @@ class Extract:
# weighted averaged onto new horizontal grid
for vn in range(self.nvar):
self.logger.info(' sc_bdy %s %s', np.nanmin(sc_bdy), np.nanmax(sc_bdy))
dst_bdy = (np.nansum(sc_bdy[vn,:,:,:] * dist_wei, 2) /
dist_fac)
dst_bdy = (np.nansum(sc_bdy[vn,:,:,:] * self.dist_wei, 2) /
self.dist_fac)
self.logger.info(' dst_bdy %s %s', np.nanmin(dst_bdy), np.nanmax(dst_bdy))
# Quick check to see we have not got bad values
if np.sum(dst_bdy == np.inf) > 0:
......@@ -667,8 +675,8 @@ class Extract:
# weight vector array and rotate onto dest grid
if self.key_vec:
# [:,:,:,vn+1]
dst_bdy_2 = (np.nansum(sc_bdy[vn+1,:,:,:] * dist_wei, 2) /
dist_fac)
dst_bdy_2 = (np.nansum(sc_bdy[vn+1,:,:,:] * self.dist_wei, 2) /
self.dist_fac)
self.logger.info('time to to rot and rep ')
self.logger.info('%s %s', np.nanmin(dst_bdy), np.nanmax(dst_bdy))
self.logger.info( '%s en to %s %s' , self.rot_str,self.rot_dir, dst_bdy.shape)
......@@ -677,18 +685,18 @@ class Extract:
self.logger.info('%s %s', np.nanmin(dst_bdy), np.nanmax(dst_bdy))
# Apply 1-2-1 filter along bdy pts using NN ind self.id_121
if self.first:
tmp_valid = np.invert(np.isnan(
self.tmp_valid = np.invert(np.isnan(
dst_bdy.flatten(1)[self.id_121]))
# Finished first run operations
self.first = False
dst_bdy = (np.nansum(dst_bdy.flatten(1)[self.id_121] *
self.tmp_filt, 2) / np.sum(self.tmp_filt *
tmp_valid, 2))
self.tmp_valid, 2))
# Set land pts to zero
self.logger.info(' pre dst_bdy[nan_ind] %s %s', np.nanmin(dst_bdy), np.nanmax(dst_bdy))
dst_bdy[nan_ind] = 0
self.logger.info(' pre dst_bdy %s %s', np.nanmin(dst_bdy), np.nanmax(dst_bdy))
dst_bdy[self.nan_ind] = 0
self.logger.info(' post dst_bdy %s %s', np.nanmin(dst_bdy), np.nanmax(dst_bdy))
# Remove any data on dst grid that is in land
dst_bdy[:,np.isnan(self.bdy_z)] = 0
......@@ -699,7 +707,7 @@ class Extract:
# If all else fails fill down using deepest pt
dst_bdy = dst_bdy.flatten(1)
dst_bdy += ((dst_bdy == 0) *
dst_bdy[data_ind].repeat(sc_z_len))
dst_bdy[self.data_ind].repeat(sc_z_len))
# Weighted averaged on new vertical grid
dst_bdy = (dst_bdy[self.z_ind[:,0]] * self.z_dist[:,0] +
dst_bdy[self.z_ind[:,1]] * self.z_dist[:,1])
......
pynemo/profile.py
View file @ 218eacca
......@@ -142,8 +142,8 @@ def process_bdy(setup_filepath=0, mask_gui=False):
elif settings['year_end'] - settings['year_000'] == 0:
days_mth = monthrange(settings['year_end'],settings['month_end'])
date_min = str(settings['year_000'])+'-'+str(settings['month_000'])+'-01'
date_max = str(settings['year_end'])+'-'+str(settings['month_end'])+'-'+str(days_mth[1])
date_min = str(settings['year_000'])+'-'+str(settings['month_000']).zfill(2)+'-01'
date_max = str(settings['year_end'])+'-'+str(settings['month_end']).zfill(2)+'-'+str(days_mth[1])
elif settings['year_end'] - settings['year_000'] < 0:
error_msg = 'end date before start date please ammend bdy file'
......@@ -168,6 +168,7 @@ def process_bdy(setup_filepath=0, mask_gui=False):
sys.exit(error_msg)
# download request for CMEMS data, try whole time interval first.
logger.info('starting CMES download now (this can take a while)...')
dl = dl_cmems.request_cmems(settings, date_min, date_max)
if dl == 0:
logger.info('CMES data downloaded successfully')
......
test_scripts/test_cmems.py
View file @ 218eacca
......@@ -8,89 +8,90 @@ into profile.py
@author: thopri
"""
# import nemo_bdy_dl_cmems as dl_cmems
# import logging
# import nemo_bdy_setup as setup
# from calendar import monthrange
#
# logger = logging.getLogger(__name__)
# logging.basicConfig(filename='/Users/thopri/Projects/PyNEMO/test_scripts/cmems_test.log', level=logging.INFO)
#
# Setup = setup.Setup('/Users/thopri/Projects/PyNEMO/inputs/namelist_cmems.bdy') # default settings file
# settings = Setup.settings
#
# if settings['use_cmems'] == True:
#
# if settings['year_end'] - settings['year_000'] > 0:
# date_min = settings['year_000']+'-01-01'
# date_max = settings['year_end']+'-12-31'
#
# days_mth = monthrange(settings['year_end'],settings['month_end'])
#
# date_min = str(settings['year_000'])+'-'+str(settings['month_000'])+'-01'
#
# date_max = str(settings['year_end'])+'-'+str(settings['month_end'])+'-'+str(days_mth[1])
#
# cmes_config= {
# 'ini_config_template' : settings['cmes_config_template'],
# 'user' : settings['cmes_usr'],
# 'pwd' : settings['cmes_pwd'],
# 'motu_server' : settings['motu_server'],
# 'service_id' : settings['cmes_model'],
# 'product_id' : settings['cmes_product'],
# 'date_min' : date_min,
# 'date_max' : date_max,
# 'latitude_min' : settings['latitude_min'],
# 'latitude_max' : settings['latitude_max'],
# 'longitude_min' : settings['longitude_min'],
# 'longitude_max' : settings['longitude_max'],
# 'depth_min' : settings['depth_min'],
# 'depth_max' : settings['depth_max'],
# 'variable' : settings['cmes_variable'],
# 'src_dir' : settings['src_dir'],
# 'out_name' : settings['cmes_output'],
# 'config_out' : settings['cmes_config']
# }
#
# chk = dl_cmems.chk_motu()
#
# if chk == 1:
# logger.error('motuclient not installed, please install by running $ pip install motuclient')
#
# else:
# logger.info('version ' +chk+ ' of motuclient is installed')
# logger.info('requesting CMES download now (this can take a while)...')
# dl = dl_cmems.request_CMEMS(cmes_config)
#
# if dl == 0:
# logger.info('CMES data downloaded successfully')
#
# if type(dl) == str:
# logger.error(dl)
#
from pynemo import nemo_bdy_dl_cmems as dl_cmems
from pynemo import nemo_bdy_setup as setup
from calendar import monthrange
import sys
import requests, json
import time
import threading
Setup = setup.Setup('/Users/thopri/Projects/PyNEMO/inputs/namelist_cmems.bdy') # default settings file
settings = Setup.settings
if settings['use_cmems'] == True:
if settings['year_end'] - settings['year_000'] > 0:
date_min = settings['year_000']+'-01-01'
date_max = settings['year_end']+'-12-31'
days_mth = monthrange(settings['year_end'],settings['month_end'])
date_min = str(settings['year_000'])+'-'+str(settings['month_000'])+'-01'
date_max = str(settings['year_end'])+'-'+str(settings['month_end'])+'-'+str(days_mth[1])
cmes_config= {
'ini_config_template' : settings['cmes_config_template'],
'user' : settings['cmes_usr'],
'pwd' : settings['cmes_pwd'],
'motu_server' : settings['motu_server'],
'service_id' : settings['cmes_model'],
'product_id' : settings['cmes_product'],
'date_min' : date_min,
'date_max' : date_max,
'latitude_min' : settings['latitude_min'],
'latitude_max' : settings['latitude_max'],
'longitude_min' : settings['longitude_min'],
'longitude_max' : settings['longitude_max'],
'depth_min' : settings['depth_min'],
'depth_max' : settings['depth_max'],
'variable' : settings['cmes_variable'],
'src_dir' : settings['src_dir'],
'out_name' : settings['cmes_output'],
'config_out' : settings['cmes_config']
}
chk = dl_cmems.chk_motu()
#my_config = {'verbose': sys.stderr }
cmems_URL = 'http://nrt.cmems-du.eu/motu-web/Motu'
cas_url = 'https://cmems-cas.cls.fr/cas/login'
dl = dl_cmems.request_CMEMS(cmes_config)
#payload = {'user':'tprime','pwd':'*5TrsWI8i&Ds','action':'describeproduct','service':'GLOBAL_ANALYSIS_FORECAST_PHY_001_024-TDS','product':'global-analysis-forecast-phy-001-024'}
class progress_bar_loading(threading.Thread):
#with requests.Session() as session:
# r = session.post(cmems_URL,data=payload)
def run(self):
global stop
global kill
print 'Loading.... ',
sys.stdout.flush()
i = 0
while stop != True:
if (i%4) == 0:
sys.stdout.write('\b/')
elif (i%4) == 1:
sys.stdout.write('\b-')
elif (i%4) == 2:
sys.stdout.write('\b\\')
elif (i%4) == 3:
sys.stdout.write('\b|')
response = requests.get(cmems_URL, verify=False)
sys.stdout.flush()
time.sleep(0.2)
i+=1
cookies = response.cookies
if kill == True:
print '\b\b\b\b ABORT!',
else:
print '\b\b done!',
payload = {'_eventId': 'submit','lt':'e1s1', 'submit': 'LOGIN', 'username': 'tprime', 'password': '*5TrsWI8i&Ds'}
sessionResp = requests.post(cas_url, data=payload,params=cookies,verify=False)
print sessionResp.status_code
print sessionResp.content
kill = False
stop = False
p = progress_bar_loading()
p.start()
#data = { 'user' : 'tprime', 'pwd' : '*5TrsWI8i&Ds' }
#r = requests.post(cmems_URL, data=json.dumps(data), verify=False)
#token = json.loads(r.text)['session']
\ No newline at end of file
try:
#anything you want to run.
time.sleep(1)
stop = True
except KeyboardInterrupt or EOFError:
kill = True
stop = True
\ No newline at end of file
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