diff --git a/gridded/da/analysis/expand_fluxes.py b/gridded/da/analysis/expand_fluxes.py index 46b3a87ed763fa41d5494dfb1792be735222647d..c07d461a0deaec073d0d7fe64f8a8015807364b0 100755 --- a/gridded/da/analysis/expand_fluxes.py +++ b/gridded/da/analysis/expand_fluxes.py @@ -1042,8 +1042,8 @@ if __name__ == "__main__": DaCycle.DaSystem = DaSystem - StateVector = CtStateVector(DaCycle) - StateVector.initialize() + StateVector = CtStateVector() + StateVector.initialize(DaCycle) while DaCycle['time.end'] < DaCycle['time.finish']: diff --git a/gridded/da/baseclasses/statevector.py b/gridded/da/baseclasses/statevector.py index 70ea63cc99967e90fe57c7f5940bc7f2b489ad4e..f5c9a9ddeea7a242b61a32606184eb89074b44f9 100755 --- a/gridded/da/baseclasses/statevector.py +++ b/gridded/da/baseclasses/statevector.py @@ -22,9 +22,9 @@ your own baseclass StateVector we refer to :ref:`tut_chapter5`. """ import os -import sys import logging import numpy as np +from datetime import timedelta import da.tools.io4 as io identifier = 'Baseclass Statevector ' @@ -114,20 +114,16 @@ class StateVector(object): """ - def __init__(self, DaCycle=None): + def __init__(self): self.Identifier = identifier self.Version = version # The following code allows the object to be initialized with a DaCycle object already present. Otherwise, it can # be added at a later moment. - if DaCycle != None: - self.DaCycle = DaCycle - else: - self.DaCycle = {} logging.info('Statevector object initialized: %s' % self.Identifier) - def initialize(self): + def initialize(self, DaCycle): """ initialize the object by specifying the dimensions. There are two major requirements for each statvector that you want to build: @@ -138,9 +134,9 @@ class StateVector(object): An example is given below. """ - self.nlag = int(self.DaCycle['time.nlag']) - self.nmembers = int(self.DaCycle['da.optimizer.nmembers']) - self.nparams = int(self.DaCycle.DaSystem['nparameters']) + self.nlag = int(DaCycle['time.nlag']) + self.nmembers = int(DaCycle['da.optimizer.nmembers']) + self.nparams = int(DaCycle.DaSystem['nparameters']) self.nobs = 0 self.ObsToAssimmilate = () # empty containter to hold observations to assimilate later on @@ -157,14 +153,14 @@ class StateVector(object): # This specifies the file to read with the gridded mask at 1x1 degrees. Each gridbox holds a number that specifies the parametermember # that maps onto it. From this map, a dictionary is created that allows a reverse look-up so that we can map parameters to a grid. - mapfile = os.path.join(self.DaCycle.DaSystem['regionsfile']) + mapfile = os.path.join(DaCycle.DaSystem['regionsfile']) ncf = io.CT_Read(mapfile, 'read') self.gridmap = ncf.get_variable('regions') self.tcmap = ncf.get_variable('transcom_regions') ncf.close() - logging.debug("A TransCom map on 1x1 degree was read from file %s" % self.DaCycle.DaSystem['regionsfile']) - logging.debug("A parameter map on 1x1 degree was read from file %s" % self.DaCycle.DaSystem['regionsfile']) + logging.debug("A TransCom map on 1x1 degree was read from file %s" % DaCycle.DaSystem['regionsfile']) + logging.debug("A parameter map on 1x1 degree was read from file %s" % DaCycle.DaSystem['regionsfile']) # Create a dictionary for state <-> gridded map conversions @@ -281,7 +277,7 @@ class StateVector(object): logging.debug('%d new ensemble members were added to the state vector # %d' % (self.nmembers, lag)) - def propagate(self, startdate, cycle): + def propagate(self, DaCycle): """ :rtype: None @@ -292,7 +288,7 @@ class StateVector(object): In the future, this routine can incorporate a formal propagation of the statevector. """ - from datetime import timedelta + # Remove State Vector n=1 by simply "popping" it from the list and appending a new empty list at the front. This empty list will # hold the new ensemble for the new cycle @@ -300,8 +296,8 @@ class StateVector(object): self.EnsembleMembers.append([]) # And now create a new time step of mean + members for n=nlag - date = startdate + timedelta(days=(self.nlag - 0.5) * int(cycle)) - cov = self.get_covariance(date) + date = DaCycle['time.start'] + timedelta(days=(self.nlag - 0.5) * int(DaCycle['time.cycle'])) + cov = self.get_covariance(date, DaCycle) self.make_new_ensemble(self.nlag, cov) logging.info('The state vector has been propagated by one cycle') @@ -590,7 +586,9 @@ class StateVector(object): return (mean, cov) - + def get_covariance(self, date, cycleparams): + pass + ################### End Class StateVector ################### if __name__ == "__main__": diff --git a/gridded/da/ct/statevector.py b/gridded/da/ct/statevector.py index 65163f8272f1f540f9f9735b31237e7628ac6367..16a21324ed8abd9cc3c36dc9fd8164c6132dd173 100755 --- a/gridded/da/ct/statevector.py +++ b/gridded/da/ct/statevector.py @@ -27,7 +27,7 @@ version = '0.0' class CtStateVector(StateVector): """ This is a StateVector object for CarbonTracker. It has a private method to make new ensemble members """ - def get_covariance(self, date): + def get_covariance(self, date, DaCycle): """ Make a new ensemble from specified matrices, the attribute lag refers to the position in the state vector. Note that lag=1 means an index of 0 in python, hence the notation lag-1 in the indexing below. The argument is thus referring to the lagged state vector as [1,2,3,4,5,..., nlag] @@ -39,8 +39,8 @@ class CtStateVector(StateVector): # Get the needed matrices from the specified covariance files - file_ocn_cov = self.DaCycle.DaSystem['ocn.covariance'] - file_bio_cov = self.DaCycle.DaSystem['bio.covariance'] #LU logika tego to powrot do dacycle zeby potem z systemu(CT) pobrac parametr + file_ocn_cov = DaCycle.DaSystem['ocn.covariance'] + file_bio_cov = DaCycle.DaSystem['bio.covariance'] #LU logika tego to powrot do dacycle zeby potem z systemu(CT) pobrac parametr # replace YYYY.MM in the ocean covariance file string diff --git a/gridded/da/ctgridded/statevector.py b/gridded/da/ctgridded/statevector.py index 6e2d35ad70703006cbe8fcd95c4e4764c9e027c0..6f49238b39b54a2ee60d8b96c933dffa1155701d 100755 --- a/gridded/da/ctgridded/statevector.py +++ b/gridded/da/ctgridded/statevector.py @@ -27,7 +27,7 @@ version = '0.0' class CtGriddedStateVector(StateVector): """ This is a StateVector object for CarbonTracker. It has a private method to make new ensemble members """ - def get_covariance(self, date): + def get_covariance(self, date, DaCycle): """ Make a new ensemble from specified matrices, the attribute lag refers to the position in the state vector. Note that lag=1 means an index of 0 in python, hence the notation lag-1 in the indexing below. The argument is thus referring to the lagged state vector as [1,2,3,4,5,..., nlag] @@ -41,11 +41,11 @@ class CtGriddedStateVector(StateVector): # Get the needed matrices from the specified covariance files - file_ocn_cov = self.DaCycle.DaSystem['ocn.covariance'] + file_ocn_cov = DaCycle.DaSystem['ocn.covariance'] - cov_files = os.listdir(self.DaCycle.DaSystem['bio.cov.dir']) + cov_files = os.listdir(DaCycle.DaSystem['bio.cov.dir']) - cov_files = [os.path.join(self.DaCycle.DaSystem['bio.cov.dir'], f) for f in cov_files if self.DaCycle.DaSystem['bio.cov.prefix'] in f] + cov_files = [os.path.join(DaCycle.DaSystem['bio.cov.dir'], f) for f in cov_files if DaCycle.DaSystem['bio.cov.prefix'] in f] logging.debug("Found %d covariances to use for biosphere" % len(cov_files)) diff --git a/gridded/da/tm5/observationoperator.py b/gridded/da/tm5/observationoperator.py index 1b4f04408bbc8878e2d0de0e277c9c0649227700..bcc0903246ad0fe89903e6dfe36802ca8a6c10de 100755 --- a/gridded/da/tm5/observationoperator.py +++ b/gridded/da/tm5/observationoperator.py @@ -85,12 +85,13 @@ class TM5ObservationOperator(ObservationOperator): logging.info('Observation Operator initialized: %s (%s)' % (self.Identifier, self.Version)) - def initialize(self): + def initialize(self, DaCycle): """ Execute all steps needed to prepare the ObsOperator for use inside CTDAS, only done at the very first cycle normally """ - + self.DaCycle = DaCycle + if self.DaCycle['time.restart'] == False : logging.info('First time step, setting up and compiling the TM5 model before proceeding!') # Modify the rc-file to reflect directory structure defined by CTDAS @@ -124,7 +125,7 @@ class TM5ObservationOperator(ObservationOperator): logging.debug('Reloading the da.obsoperator.rc file for this DaCycle') - self.load_rc(self.DaCycle['da.obsoperator.rc']) + self.load_rc(self.DaCycle['da.obsoperator.rc']) #LU czy tutaj jest to kiedykolwiek zapisywane? chyba tak, jak zapisujemy daCycle... logging.debug('Note that the obsoperator is not recompiled if this is a recovery from a crash!!!') diff --git a/gridded/da/tools/pipeline.py b/gridded/da/tools/pipeline.py index 142364fd7bf3cba067dd75209d4fdee72c6ec888..74e15d51468127ee8a18931d18146eaba47eac57 100755 --- a/gridded/da/tools/pipeline.py +++ b/gridded/da/tools/pipeline.py @@ -46,7 +46,7 @@ def forward_pipeline(DaCycle, PlatForm, DaSystem, Samples, StateVector, ObsOpera # Create State Vector - StateVector.initialize() #LU w prepare state inicjalizujemy wektor stanu ktoremu dopiero co przypisalismy wartosci. + StateVector.initialize(DaCycle) #LU w prepare state inicjalizujemy wektor stanu ktoremu dopiero co przypisalismy wartosci. # Read from other simulation and write priors, then read posteriors and propagate @@ -87,10 +87,10 @@ def start_job(DaCycle, DaSystem, DaPlatForm, StateVector, Samples, ObsOperator): DaCycle.DaSystem = DaSystem #LU laczy te listy parametrow ale na zasadzie hierarchii DaCycle.DaPlatForm = DaPlatForm #LU przypisuje cyklowi platforme DaCycle.initialize() #LU setup file structure etc - StateVector.DaCycle = DaCycle # also embed object in StateVector so it can access cycle information for I/O etc #LU cykl zostaje przypisany state vectorowi + #StateVector.DaCycle = DaCycle # also embed object in StateVector so it can access cycle information for I/O etc #LU cykl zostaje przypisany state vectorowi #Samples.DaCycle = DaCycle # also embed object in Samples object so it can access cycle information for I/O etc #LU cykl zostaje przypisany probkom - ObsOperator.DaCycle = DaCycle # also embed object in ObsOperator object so it can access cycle information for I/O etc #LU cykl zostaje przypisany obsoperatorowi - ObsOperator.initialize() # Setup Observation Operator #LU a pote mobsoperator jest inicjalizowany + #ObsOperator.DaCycle = DaCycle # also embed object in ObsOperator object so it can access cycle information for I/O etc #LU cykl zostaje przypisany obsoperatorowi + ObsOperator.initialize(DaCycle) # Setup Observation Operator #LU a pote mobsoperator jest inicjalizowany def prepare_state(DaCycle, StateVector): @@ -103,7 +103,7 @@ def prepare_state(DaCycle, StateVector): logging.info(header + "starting prepare_state" + footer) - StateVector.initialize() #LU w prepare state inicjalizujemy wektor stanu ktoremu dopiero co przypisalismy wartosci. + StateVector.initialize(DaCycle) #LU w prepare state inicjalizujemy wektor stanu ktoremu dopiero co przypisalismy wartosci. #LU to jest zalezne od cyklu, i cykl pojedynczy moze miec time.restart lub moze nie miec. if not DaCycle['time.restart']: #LU jesli jest to pierwszy cykl @@ -113,7 +113,7 @@ def prepare_state(DaCycle, StateVector): nlag = StateVector.nlag #LU dla kazdego od zera do dlugosc(cykl) wyznaczamy date oraz znajdujemy kowariancje i robimy nowa ensemble. for n in range(nlag): date = DaCycle['time.start'] + datetime.timedelta(days=(n + 0.5) * int(DaCycle['time.cycle'])) #LU ta data jest tutaj potrzebna tylko do znalezienia odpowiedniego pliku z kowariancja. - cov = StateVector.get_covariance(date) + cov = StateVector.get_covariance(date, DaCycle) StateVector.make_new_ensemble(n + 1, cov) else: @@ -125,7 +125,7 @@ def prepare_state(DaCycle, StateVector): # Now propagate the ensemble by one cycle to prepare for the current cycle - StateVector.propagate(DaCycle['time.start'], DaCycle['time.cycle']) + StateVector.propagate(DaCycle) # Finally, also write the StateVector to a file so that we can always access the a-priori information