Skip to content
Snippets Groups Projects
Commit d5400290 authored by Aki Tsuruta's avatar Aki Tsuruta
Browse files

update on methane, mainly changes due to Python2->Python3

parent 89aab8ca
No related branches found
No related tags found
1 merge request!20update on methane, mainly changes due to Python2->Python3
"""CarbonTracker Data Assimilation Shell (CTDAS) Copyright (C) 2017 Wouter Peters. """CarbonTracker Data Assimilation Shell (CTDAS) Copyright (C) 2017 Wouter Peters.
Users are recommended to contact the developers (wouter.peters@wur.nl) to receive Users are recommended to contact the developers (wouter.peters@wur.nl) to receive
updates of the code. See also: http://www.carbontracker.eu. updates of the code. See also: http://www.carbontracker.eu.
This program is free software: you can redistribute it and/or modify it under the This program is free software: you can redistribute it and/or modify it under the
terms of the GNU General Public License as published by the Free Software Foundation, terms of the GNU General Public License as published by the Free Software Foundation,
version 3. This program is distributed in the hope that it will be useful, but version 3. This program is distributed in the hope that it will be useful, but
WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
...@@ -24,7 +24,7 @@ import logging ...@@ -24,7 +24,7 @@ import logging
import os import os
import sys import sys
import numpy as np import numpy as np
from string import join #from string import join
import da.tools.rc as rc import da.tools.rc as rc
from da.tools.initexit import * from da.tools.initexit import *
...@@ -40,35 +40,40 @@ needed_da_items = [ ...@@ -40,35 +40,40 @@ needed_da_items = [
'da.obsoperator', 'da.obsoperator',
'da.obsoperator.rc', 'da.obsoperator.rc',
'da.optimizer.nmembers', 'da.optimizer.nmembers',
'da.suffixname'] 'da.suffixname',
'da.resources.ntime',
'da.resources.npes',
#def validate_rc2(self): 'da.projectaccount',
# """ 'da.tm5.jobfilename'
# Validate the contents of the rc-file given a dictionary of required keys. ]
# Currently required keys are :attr:`~da.tools.initexit.needed_da_items`
# """
# def validate_rc2(self):
# for k, v in self.iteritems(): """
# if v in ['True', 'true', 't', 'T', 'y', 'yes']: Validate the contents of the rc-file given a dictionary of required keys.
# self[k] = True Currently required keys are :attr:`~da.tools.initexit.needed_da_items`
# if v in ['False', 'false', 'f', 'F', 'n', 'no']: """
# self[k] = False
# if 'date' in k : for k, v in self.items():
# self[k] = to_datetime(v) if v in ['True', 'true', 't', 'T', 'y', 'yes']:
# if k in ['time.start', 'time.end', 'time.finish', 'da.restart.tstamp']: self[k] = True
# self[k] = to_datetime(v) if v in ['False', 'false', 'f', 'F', 'n', 'no']:
# for key in needed_da_items: self[k] = False
# if not self.has_key(key): if 'date' in k :
# msg = 'Missing a required value in rc-file : %s' % key self[k] = to_datetime(v)
# logging.error(msg) if k in ['time.start', 'time.end', 'time.finish', 'da.restart.tstamp']:
# logging.error('!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ') self[k] = to_datetime(v)
# logging.error('Please note the update on Dec 02 2011 where rc-file names for DaSystem and ') for key in needed_da_items:
# logging.error('are from now on specified in the main rc-file (see da/rc/da.rc for example)') if key not in self:
# logging.error('!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ') msg = 'Missing a required value in rc-file : %s' % key
# raise IOError, msg logging.error(msg)
# logging.debug('DA Cycle settings have been validated succesfully') logging.error('!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ')
#CycleControl.validate_rc = validate_rc2 logging.error('Please note the update on Dec 02 2011 where rc-file names for DaSystem and ')
logging.error('are from now on specified in the main rc-file (see da/rc/da.rc for example)')
logging.error('!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ')
raise IOError(msg)
logging.debug('DA Cycle settings have been validated succesfully')
CycleControl.validate_rc = validate_rc2
def submit_next_cycle_fmi(self): def submit_next_cycle_fmi(self):
""" """
...@@ -88,21 +93,22 @@ def submit_next_cycle_fmi(self): ...@@ -88,21 +93,22 @@ def submit_next_cycle_fmi(self):
jobid = self['time.end'].strftime('%Y%m%d') jobid = self['time.end'].strftime('%Y%m%d')
targetdir = os.path.join(self['dir.exec']) targetdir = os.path.join(self['dir.exec'])
jobfile = os.path.join(targetdir, 'jb.%s.jb' % jobid) jobfile = os.path.join(targetdir, 'jb.%s.jb' % jobid)
logfile = os.path.join(targetdir, 'jb.%s.log' % jobid)
# Template and commands for job # Template and commands for job
jobparams = {'jobname':"ctdas", 'jobtime':'05:00:00', 'logfile': logfile, 'errfile': logfile, jobtime = self['da.resources.ntime']
'jobpes':'120','jobnodes':'20'} jobprojectaccount = self['da.projectaccount']
jobparams = {'jobname':"ctdas", 'jobtime':jobtime, 'jobnpes':'1', 'projectaccount':jobprojectaccount}
template = self.daplatform.get_job_template(jobparams) template = self.daplatform.get_job_template(jobparams)
execcommand = os.path.join(self['dir.da_submit'], sys.argv[0]) execcommand = os.path.join(self['dir.da_submit'], sys.argv[0])
if '-t' in self.opts: if '-t' in self.opts:
(self.opts).remove('-t') (self.opts).remove('-t')
#template += 'python %s rc=%s %s' % (execcommand, self['da.restart.fname'], join(self.opts, '')) template += 'rcfile=%s \n'%self['da.restart.fname']
template += 'cd %s \n' %self['dir.da_submit'] template += 'cd %s \n' %self['dir.da_submit']
template += '%s rc=%s %s' % (sys.argv[0], self['da.restart.fname'], join(self.opts, '')) template += 'python %s rc=$rcfile %s \n' % (sys.argv[0], ''.join(self.opts))
template += './save2allas/ctdas.sh $rcfile'
# write and submit # write and submit
self.daplatform.write_job(jobfile, template, jobid) self.daplatform.write_job(jobfile, template, jobid)
jobid = self.daplatform.submit_job(jobfile, joblog=logfile) jobid = self.daplatform.submit_job(jobfile)
job_file = 'ctdas.o%s' %jobid[0:-4] job_file = 'ctdas.o%s' %jobid #[0:-4]
logging.info('Jobfile %s' %job_file ) logging.info('Jobfile %s' %job_file )
logging.info('DA cycle has been submitted.') logging.info('DA cycle has been submitted.')
...@@ -122,7 +128,6 @@ def setup_file_structure_fmi(self): ...@@ -122,7 +128,6 @@ def setup_file_structure_fmi(self):
self['dir.input'] = os.path.join(self['dir.da_run'], 'input_%s' %suf ) self['dir.input'] = os.path.join(self['dir.da_run'], 'input_%s' %suf )
self['dir.output'] = os.path.join(self['dir.da_run'], 'output_%s' %suf ) self['dir.output'] = os.path.join(self['dir.da_run'], 'output_%s' %suf )
self['dir.analysis'] = os.path.join(self['dir.da_run'], 'analysis_%s' %suf ) self['dir.analysis'] = os.path.join(self['dir.da_run'], 'analysis_%s' %suf )
#self['dir.jobs'] = os.path.join(self['dir.da_run'], 'jobs')
self['dir.restart'] = os.path.join(self['dir.da_run'], 'restart_%s' %suf ) self['dir.restart'] = os.path.join(self['dir.da_run'], 'restart_%s' %suf )
logging.info("setup_file_structure %s" %self['dir.output']) logging.info("setup_file_structure %s" %self['dir.output'])
...@@ -131,7 +136,6 @@ def setup_file_structure_fmi(self): ...@@ -131,7 +136,6 @@ def setup_file_structure_fmi(self):
create_dirs(os.path.join(self['dir.input'])) create_dirs(os.path.join(self['dir.input']))
create_dirs(os.path.join(self['dir.output'])) create_dirs(os.path.join(self['dir.output']))
create_dirs(os.path.join(self['dir.analysis'])) create_dirs(os.path.join(self['dir.analysis']))
#create_dirs(os.path.join(self['dir.jobs']))
create_dirs(os.path.join(self['dir.restart'])) create_dirs(os.path.join(self['dir.restart']))
logging.info('Succesfully created the file structure for the assimilation job') logging.info('Succesfully created the file structure for the assimilation job')
......
"""CarbonTracker Data Assimilation Shell (CTDAS) Copyright (C) 2017 Wouter Peters. """CarbonTracker Data Assimilation Shell (CTDAS) Copyright (C) 2017 Wouter Peters.
Users are recommended to contact the developers (wouter.peters@wur.nl) to receive Users are recommended to contact the developers (wouter.peters@wur.nl) to receive
updates of the code. See also: http://www.carbontracker.eu. updates of the code. See also: http://www.carbontracker.eu.
This program is free software: you can redistribute it and/or modify it under the This program is free software: you can redistribute it and/or modify it under the
terms of the GNU General Public License as published by the Free Software Foundation, terms of the GNU General Public License as published by the Free Software Foundation,
version 3. This program is distributed in the hope that it will be useful, but version 3. This program is distributed in the hope that it will be useful, but
WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
...@@ -20,18 +20,16 @@ Revision History: ...@@ -20,18 +20,16 @@ Revision History:
File created on Apr 2016. File created on Apr 2016.
""" """
from . import standardvariables
import datetime as dt import datetime as dt
from numpy import array, arange from numpy import array, arange
import os import os
import logging import logging
import sys import sys
sys.path.append('/stornext/store/carbon/CarbonTracker/netcdf4/python_packages_install_dir/lib/python2.7/site-packages/')
import netCDF4 import netCDF4
sys.path.append('/stornext/store/carbon/CarbonTracker/pyhdf/pyhdf-0.8.3/lib/python2.7/site-packages/')
import pyhdf.SD as hdf import pyhdf.SD as hdf
sys.path.append('../') sys.path.append('../')
from da.tools.io4 import * from da.tools.io4 import *
import da.methane.standardvariables
disclaimer = "This data belongs to the CarbonTracker project" disclaimer = "This data belongs to the CarbonTracker project"
...@@ -56,9 +54,9 @@ def add_tc_header2(self): ...@@ -56,9 +54,9 @@ def add_tc_header2(self):
def standard_var2(self,varname): def standard_var2(self,varname):
""" return properties of standard variables """ """ return properties of standard variables """
from . import standardvariables import da.methane.standardvariables as standardvariables
if varname in list(standardvariables.standard_variables.keys()): if varname in standardvariables.standard_variables.keys():
return standardvariables.standard_variables[varname] return standardvariables.standard_variables[varname]
else: else:
return standardvariables.standard_variables['unknown'] return standardvariables.standard_variables['unknown']
......
"""CarbonTracker Data Assimilation Shell (CTDAS) Copyright (C) 2017 Wouter Peters. """CarbonTracker Data Assimilation Shell (CTDAS) Copyright (C) 2017 Wouter Peters.
Users are recommended to contact the developers (wouter.peters@wur.nl) to receive Users are recommended to contact the developers (wouter.peters@wur.nl) to receive
updates of the code. See also: http://www.carbontracker.eu. updates of the code. See also: http://www.carbontracker.eu.
This program is free software: you can redistribute it and/or modify it under the This program is free software: you can redistribute it and/or modify it under the
terms of the GNU General Public License as published by the Free Software Foundation, terms of the GNU General Public License as published by the Free Software Foundation,
version 3. This program is distributed in the hope that it will be useful, but version 3. This program is distributed in the hope that it will be useful, but
WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
...@@ -14,18 +14,16 @@ program. If not, see <http://www.gnu.org/licenses/>.""" ...@@ -14,18 +14,16 @@ program. If not, see <http://www.gnu.org/licenses/>."""
# obs.py # obs.py
""" """
Author : aki Author : peters
Revision History: Revision History:
File revised on Feb 2017. File created on 28 Jul 2010.
""" """
import os import os
import sys import sys
import logging import logging
#from da.baseclasses.statevector import filename
import datetime as dtm import datetime as dtm
from string import strip
from numpy import array, logical_and from numpy import array, logical_and
sys.path.append(os.getcwd()) sys.path.append(os.getcwd())
...@@ -86,8 +84,8 @@ class MethaneObservations(Observations): ...@@ -86,8 +84,8 @@ class MethaneObservations(Observations):
#evn = [s.tostring().lower() for s in evn] #evn = [s.tostring().lower() for s in evn]
#evn = map(strip, evn) #evn = map(strip, evn)
sites = ncf.get_variable('obs_id').take(subselect, axis=0) sites = ncf.get_variable('obs_id').take(subselect, axis=0)
sites = [s.tostring().lower() for s in sites] sites = [s.tostring().lower().decode("utf-8") for s in sites]
sites = list(map(strip, sites)) sites = list(map(str.strip, sites))
lats = ncf.get_variable('latitude').take(subselect, axis=0) lats = ncf.get_variable('latitude').take(subselect, axis=0)
lons = ncf.get_variable('longitude').take(subselect, axis=0) lons = ncf.get_variable('longitude').take(subselect, axis=0)
alts = ncf.get_variable('altitude').take(subselect, axis=0) alts = ncf.get_variable('altitude').take(subselect, axis=0)
...@@ -107,7 +105,6 @@ class MethaneObservations(Observations): ...@@ -107,7 +105,6 @@ class MethaneObservations(Observations):
for n in range(len(dates)): for n in range(len(dates)):
obs[n] = obs[n] obs[n] = obs[n]
#self.datalist.append(MoleFractionSample(ids[n], dates[n], sites[n], obs[n], 0.0, 0.0, 0.0, 0.0, flags[n], alts[n], lats[n], lons[n], evn[n], species[n], strategy[n], 0.0, self.obs_filename))
self.datalist.append( MoleFractionSample(ids[n], dates[n], sites[n], obs[n], 0.0, 0.0, 0.0, 0.0, 0.0, alts[n], lats[n], lons[n], '0000', 'ch4', 1.0, 0.0, self.obs_filename) ) self.datalist.append( MoleFractionSample(ids[n], dates[n], sites[n], obs[n], 0.0, 0.0, 0.0, 0.0, 0.0, alts[n], lats[n], lons[n], '0000', 'ch4', 1.0, 0.0, self.obs_filename) )
logging.debug("Added %d observations to the Data list" % len(dates)) logging.debug("Added %d observations to the Data list" % len(dates))
...@@ -125,8 +122,7 @@ class MethaneObservations(Observations): ...@@ -125,8 +122,7 @@ class MethaneObservations(Observations):
ncf = io.ct_read(filename, method='read') ncf = io.ct_read(filename, method='read')
ids = ncf.get_variable('obs_num') ids = ncf.get_variable('obs_num')
simulated = ncf.get_variable('flask') simulated = ncf.get_variable('flask')
#for i in xrange(simulated.shape[0]):
# print simulated[i,:]
ncf.close() ncf.close()
logging.info("Successfully read data from model sample file (%s)" % filename) logging.info("Successfully read data from model sample file (%s)" % filename)
...@@ -157,7 +153,7 @@ class MethaneObservations(Observations): ...@@ -157,7 +153,7 @@ class MethaneObservations(Observations):
""" """
f = io.CT_CDF(obsinputfile, method='create') f = io.CT_CDF(obsinputfile, method='create')
logging.debug('Creating new observations file for ObservationOperator (%s)' % obsinputfile) logging.info('Creating new observations file for ObservationOperator (%s)' % obsinputfile)
dimid = f.add_dim('obs', len(self.datalist)) dimid = f.add_dim('obs', len(self.datalist))
dim200char = f.add_dim('string_of200chars', 200) dim200char = f.add_dim('string_of200chars', 200)
...@@ -168,7 +164,7 @@ class MethaneObservations(Observations): ...@@ -168,7 +164,7 @@ class MethaneObservations(Observations):
#return obsinputfile #return obsinputfile
data = self.getvalues('id') data = self.getvalues('id')
savedict = io.std_savedict.copy() savedict = io.std_savedict.copy()
savedict['name'] = "obs_num" savedict['name'] = "obs_num"
savedict['dtype'] = "int" savedict['dtype'] = "int"
...@@ -284,7 +280,7 @@ class MethaneObservations(Observations): ...@@ -284,7 +280,7 @@ class MethaneObservations(Observations):
logging.debug('Model-data mismatch active sites : %d ' % self.n_sites_active) logging.debug('Model-data mismatch active sites : %d ' % self.n_sites_active)
logging.debug('Model-data mismatch moved sites : %d ' % self.n_sites_moved) logging.debug('Model-data mismatch moved sites : %d ' % self.n_sites_moved)
cats = [k for k in list(sites_weights.keys()) if 'site.category' in k] cats = [k for k in sites_weights.keys() if 'site.category' in k]
SiteCategories = {} SiteCategories = {}
for key in cats: for key in cats:
...@@ -296,7 +292,7 @@ class MethaneObservations(Observations): ...@@ -296,7 +292,7 @@ class MethaneObservations(Observations):
SiteCategories[name] = {'error':error, 'may_localize':may_localize, 'may_reject':may_reject} SiteCategories[name] = {'error':error, 'may_localize':may_localize, 'may_reject':may_reject}
#print name,SiteCategories[name] #print name,SiteCategories[name]
active = [k for k in list(sites_weights.keys()) if 'site.active' in k] active = [k for k in sites_weights.keys() if 'site.active' in k]
site_info = {} site_info = {}
for key in active: for key in active:
...@@ -308,7 +304,7 @@ class MethaneObservations(Observations): ...@@ -308,7 +304,7 @@ class MethaneObservations(Observations):
for obs in self.datalist: for obs in self.datalist:
obs.mdm = 1000.0 # default is very high model-data-mismatch, until explicitly set by script obs.mdm = 1000.0 # default is very high model-data-mismatch, until explicitly set by script
if obs.code in site_info: if obs.code in site_info:
logging.debug("Observation found (%s)" % obs.code) logging.debug("Observation found (%s)" % obs.code)
obs.mdm = site_info[obs.code]['error'] * self.global_R_scaling obs.mdm = site_info[obs.code]['error'] * self.global_R_scaling
obs.may_localize = site_info[obs.code]['may_localize'] obs.may_localize = site_info[obs.code]['may_localize']
......
"""CarbonTracker Data Assimilation Shell (CTDAS) Copyright (C) 2017 Wouter Peters. """CarbonTracker Data Assimilation Shell (CTDAS) Copyright (C) 2017 Wouter Peters.
Users are recommended to contact the developers (wouter.peters@wur.nl) to receive Users are recommended to contact the developers (wouter.peters@wur.nl) to receive
updates of the code. See also: http://www.carbontracker.eu. updates of the code. See also: http://www.carbontracker.eu.
This program is free software: you can redistribute it and/or modify it under the This program is free software: you can redistribute it and/or modify it under the
terms of the GNU General Public License as published by the Free Software Foundation, terms of the GNU General Public License as published by the Free Software Foundation,
version 3. This program is distributed in the hope that it will be useful, but version 3. This program is distributed in the hope that it will be useful, but
WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
...@@ -43,7 +43,7 @@ class MethaneStateVector(StateVector): ...@@ -43,7 +43,7 @@ class MethaneStateVector(StateVector):
self.speciesdict = {'ch4': np.ones(self.nparams)} self.speciesdict = {'ch4': np.ones(self.nparams)}
logging.debug("A species mask was created, only the following species are recognized in this system:") logging.debug("A species mask was created, only the following species are recognized in this system:")
for k in list(self.speciesdict.keys()): for k in self.speciesdict.keys():
logging.debug(" -> %s" % k) logging.debug(" -> %s" % k)
def get_covariance(self, date, dacycle): def get_covariance(self, date, dacycle):
......
"""CarbonTracker Data Assimilation Shell (CTDAS) Copyright (C) 2017 Wouter Peters. """CarbonTracker Data Assimilation Shell (CTDAS) Copyright (C) 2017 Wouter Peters.
Users are recommended to contact the developers (wouter.peters@wur.nl) to receive Users are recommended to contact the developers (wouter.peters@wur.nl) to receive
updates of the code. See also: http://www.carbontracker.eu. updates of the code. See also: http://www.carbontracker.eu.
This program is free software: you can redistribute it and/or modify it under the This program is free software: you can redistribute it and/or modify it under the
terms of the GNU General Public License as published by the Free Software Foundation, terms of the GNU General Public License as published by the Free Software Foundation,
version 3. This program is distributed in the hope that it will be useful, but version 3. This program is distributed in the hope that it will be useful, but
WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
...@@ -61,7 +61,7 @@ class MethaneStateVector(StateVector): ...@@ -61,7 +61,7 @@ class MethaneStateVector(StateVector):
# These list objects hold the data for each time step of lag in the system. Note that the ensembles for each time step consist # These list objects hold the data for each time step of lag in the system. Note that the ensembles for each time step consist
# of lists of EnsembleMember objects, we define member 0 as the mean of the distribution and n=1,...,nmembers as the spread. # of lists of EnsembleMember objects, we define member 0 as the mean of the distribution and n=1,...,nmembers as the spread.
self.ensemble_members = list(range(self.nlag)) self.ensemble_members = range(self.nlag)
for n in range(self.nlag): for n in range(self.nlag):
self.ensemble_members[n] = [] self.ensemble_members[n] = []
...@@ -119,7 +119,7 @@ class MethaneStateVector(StateVector): ...@@ -119,7 +119,7 @@ class MethaneStateVector(StateVector):
self.speciesdict = {'ch4': np.ones(self.nparams)} self.speciesdict = {'ch4': np.ones(self.nparams)}
logging.debug("A species mask was created, only the following species are recognized in this system:") logging.debug("A species mask was created, only the following species are recognized in this system:")
for k in list(self.speciesdict.keys()): for k in self.speciesdict.keys():
logging.debug(" -> %s" % k) logging.debug(" -> %s" % k)
def get_covariance(self, date, dacycle): def get_covariance(self, date, dacycle):
...@@ -135,7 +135,7 @@ class MethaneStateVector(StateVector): ...@@ -135,7 +135,7 @@ class MethaneStateVector(StateVector):
#----Arbitray covariance matrix ----# #----Arbitray covariance matrix ----#
#fullcov = np.zeros((self.nparams, self.nparams), float) #fullcov = np.zeros((self.nparams, self.nparams), float)
#for i in xrange(self.nparams): #for i in range(self.nparams):
# #fullcov[i,i] = 1. # Identity matrix # #fullcov[i,i] = 1. # Identity matrix
# fullcov[i,i] = 0.08 # fullcov[i,i] = 0.08
#fullcov[self.nparams-1,self.nparams-1] = 1.e-10 #fullcov[self.nparams-1,self.nparams-1] = 1.e-10
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment