From 7abed404b478c91d9a201ae805039dca3d7ea1b7 Mon Sep 17 00:00:00 2001
From: Wouter Peters <wouter.peters@wur.nl>
Date: Fri, 6 Aug 2010 09:28:32 +0000
Subject: [PATCH] changed the CtMember and StateVector object so that there is
 no longer a separate 'mean' in the statevector, it is simply member 0

---
 da/ct/statevector.py | 31 +++++++++++++------------------
 1 file changed, 13 insertions(+), 18 deletions(-)

diff --git a/da/ct/statevector.py b/da/ct/statevector.py
index 16ce6348..373df692 100755
--- a/da/ct/statevector.py
+++ b/da/ct/statevector.py
@@ -26,7 +26,7 @@ class CtMember(object):
 
     def __init__(self, membernumber):
         self.membernumber           = membernumber   # the member number
-        self.MeanDeviations         = None           # Parameter values of this member
+        self.ParameterValues        = None           # Parameter values of this member
         self.ModelSample            = None           # Model Sampled Parameter values of this member
         self.ParameterOutputFile    = None           # File to which the member info should be written
 
@@ -40,9 +40,9 @@ class CtMember(object):
 
         filename        = os.path.join(outdir,'parameters.%03d.nc'% self.membernumber)
         f               = CT_CDF(filename,'create')
-        dimparams       = f.AddParamsDim(len(self.MeanDeviations))
+        dimparams       = f.AddParamsDim(len(self.ParameterValues))
 
-        data =  self.MeanDeviations
+        data =  self.ParameterValues
 
         savedict                = std_savedict.copy()
         savedict['name']        = "parametervalues"
@@ -75,9 +75,8 @@ class CtStateVector(object):
         self.isOptimized        = False
 
         # 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 CtMember objects, whereas the MeanValues will simply be a list of array values.
+        # of lists of CtMember objects, we define member 0 as the mean of the distribution and n=1,...,nmembers as the spread.
 
-        self.MeanValues         = range(self.nlag)
         self.EnsembleMembers    = range(self.nlag)
 
         for n in range(self.nlag):
@@ -158,15 +157,14 @@ class CtStateVector(object):
             plt.close('all')
 
 
-        # Create mean values for the new time step (default = 1.0) and add to the MeanValues list
+        # Create mean values for the new time step (default = 1.0) 
 
         NewMean                     = np.ones(self.nparams,float)
-        self.MeanValues[lag-1]      = NewMean
 
         # Create the first ensemble member with a deviation of 0.0 and add to list
 
         NewMember                   = CtMember(0)
-        NewMember.MeanDeviations    = np.zeros((self.nparams),float)  
+        NewMember.ParameterValues   = np.zeros((self.nparams),float)  + NewMean
         dummy                       = self.EnsembleMembers[lag-1].append(NewMember)
 
         # Create members 1:nmembers and add to EnsembleMembers list
@@ -176,7 +174,7 @@ class CtStateVector(object):
             rands                       = np.random.randn(self.nparams)
 
             NewMember                   = CtMember(member)
-            NewMember.MeanDeviations    = np.dot(C,rands)
+            NewMember.ParameterValues   = np.dot(C,rands) + NewMean
             dummy                       = self.EnsembleMembers[lag-1].append(NewMember)
 
         msg =   '%d new ensemble members were added to the state vector # %d'%(self.nmembers,lag)  ; logging.debug(msg)
@@ -193,10 +191,8 @@ class CtStateVector(object):
         # 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 
 
-        dummy = self.MeanValues.pop(0)
         dummy = self.EnsembleMembers.pop(0)
 
-        dummy = self.MeanValues.append([])
         dummy = self.EnsembleMembers.append([])
 
         # And now create a new week of mean + members for n=nlag
@@ -227,17 +223,18 @@ class CtStateVector(object):
 
         for n in range(self.nlag):
 
-            data =  self.MeanValues[n]
+            members = self.EnsembleMembers[n]
+            MeanState =  members[0].ParameterValues
 
             savedict                = f.StandardVar(varname='meanstate')
             savedict['dims']        = dimlag+dimparams 
-            savedict['values']      = data
+            savedict['values']      = MeanState
             savedict['count']       = n
             savedict['comment']     = 'this represents the mean of the ensemble'
             dummy                   = f.AddData(savedict)
 
             members =  self.EnsembleMembers[n]
-            data = np.array([m.MeanDeviations for m in members])
+            data = np.array([m.ParameterValues for m in members])- MeanState
 
             savedict                = f.StandardVar(varname='ensemblestate')
             savedict['dims']        = dimlag+dimmembers+dimparams 
@@ -261,11 +258,9 @@ class CtStateVector(object):
         dummy           = f.close()
 
         for n in range(self.nlag):
-
-            self.MeanValues[n] = MeanState[n,:]
             for m in range(self.nmembers):
                 NewMember                   = CtMember(m)
-                NewMember.MeanDeviations    = EnsembleMembers[n,m,:]
+                NewMember.ParameterValues   = EnsembleMembers[n,m,:] + MeanState  # add the mean to the deviations to hold the full parameter values
                 dummy                       = self.EnsembleMembers[n].append(NewMember)
 
         msg     = 'Successfully read the State Vector from file (%s) ' % (filename,) ; logging.info(msg)
@@ -361,7 +356,7 @@ if __name__ == "__main__":
 
     members[0].WriteToFile(DaCycle.da_settings['dir.input'])
 
-    data = np.array([m.MeanDeviations for m in members])
+    data = np.array([m.ParameterValues for m in members])
 
     StateVector.WriteToFile()
 
-- 
GitLab