Had to replace a few lines to ensure proper writing of prior fluxes to the...

Had to replace a few lines to ensure proper writing of prior fluxes to the files. Until now, prior fluxes 
were retrieved from a savestate.nc file (and erroneously, the poterior was read by default) while the "old"
method to provide priors was as unscaled fluxes (all multipliers equal to 1.0). This behavior is now back.

da/analysis/expand_fluxes.py

@@ -111,16 +111,20 @@ def SaveWeeklyAvg1x1Data(DaCycle, StateVector):
                     savedir         = DaCycle['dir.output'].replace(startdate.strftime('%Y%m%d'),priordate.strftime('%Y%m%d') )
                     filename        = os.path.join(savedir,'savestate.nc')
                     if os.path.exists(filename):
-                        dummy                   = StateVector.ReadFromFile(filename)
+                        dummy                   = StateVector.ReadFromFile(filename,qual=qual_short)
                         gridmean,gridvariance   = StateVector.StateToGrid(lag=n)
 
+# Replace the mean statevector by all ones (assumed priors)
+
+                        gridmean                = StateVector.VectorToGrid(vectordata = np.ones(StateVector.nparams,))
+
                         msg = 'Read prior dataset from file %s, sds %d: '%(filename,n) ; logging.debug(msg)
                         break
             else:
                 qual_short='opt'
                 savedir                 = DaCycle['dir.output']
                 filename                = os.path.join(savedir,'savestate.nc')
-                dummy                   = StateVector.ReadFromFile(filename)
+                dummy                   = StateVector.ReadFromFile(filename,qual=qual_short)
                 gridmean,gridvariance   = StateVector.StateToGrid(lag=nlag)
 
                 msg = 'Read posterior dataset from file %s, sds %d: '%(filename,nlag) ; logging.debug(msg)
@@ -291,7 +295,12 @@ def SaveWeeklyAvgStateData(DaCycle, StateVector):
                     savedir         = DaCycle['dir.output'].replace(startdate.strftime('%Y%m%d'),priordate.strftime('%Y%m%d') )
                     filename        = os.path.join(savedir,'savestate.nc')
                     if os.path.exists(filename):
-                        dummy                   = StateVector.ReadFromFile(filename)
+                        dummy                   = StateVector.ReadFromFile(filename,qual=qual_short)
+
+# Replace the mean statevector by all ones (assumed priors)
+
+                        statemean               = np.ones((StateVector.nparams,))
+
                         choicelag               = n
 
                         msg = 'Read prior dataset from file %s, sds %d: '%(filename,n) ; logging.debug(msg)
@@ -302,12 +311,13 @@ def SaveWeeklyAvgStateData(DaCycle, StateVector):
                 filename                = os.path.join(savedir,'savestate.nc')
                 dummy                   = StateVector.ReadFromFile(filename)
                 choicelag               = nlag
+                statemean               = StateVector.EnsembleMembers[choicelag-1][0].ParameterValues
 
                 msg = 'Read posterior dataset from file %s, sds %d: '%(filename,nlag) ; logging.debug(msg)
 #
 # if prior, do not multiply fluxes with parameters, otherwise do
 #
-            data                = StateVector.EnsembleMembers[choicelag-1][0].ParameterValues*vectorbio # units of mole region-1 s-1
+            data                = statemean*vectorbio # units of mole region-1 s-1
 
             savedict            = ncf.StandardVar(varname='bio_flux_%s'%qual_short)
             savedict['values']  = data
@@ -337,7 +347,7 @@ def SaveWeeklyAvgStateData(DaCycle, StateVector):
             savedict['count']=next
             ncf.AddData(savedict)
 
-            data                = StateVector.EnsembleMembers[choicelag-1][0].ParameterValues*vectorocn # units of mole region-1 s-1
+            data                = statemean*vectorocn # units of mole region-1 s-1
 
             savedict=ncf.StandardVar(varname='ocn_flux_%s'%qual_short)
             savedict['values']=data
@@ -500,6 +510,7 @@ def SaveWeeklyAvgTCData(DaCycle, StateVector):
 
             savedict            = ncf.StandardVar(varname=vname)
             savedict['name']    = vname
+            savedict['count']   = index
 
             if vname=='latitude': continue
             elif vname=='longitude': continue
@@ -507,8 +518,8 @@ def SaveWeeklyAvgTCData(DaCycle, StateVector):
             elif vname=='idate': continue
             elif 'cov' in vname:
                 
-                #cov                 = data.transpose().dot(data) 
-                cov                 = np.dot(data.transpose(),data) 
+                cov                 = data.transpose().dot(data) 
+                #cov                 = np.dot(data.transpose(),data) 
                 tcdata              = StateVector.VectorToTC(vectordata=cov,cov=True) # vector to TC
 
                 savedict['units']   = '[mol/region/s]**2'
@@ -809,11 +820,11 @@ if __name__ == "__main__":
 
     logging.root.setLevel(logging.DEBUG)
 
-    DaCycle = CycleControl(args={'rc':'../../dajet.rc'})
+    DaCycle = CycleControl(args={'rc':'../../dagriddedjet.rc'})
     DaCycle.Initialize()
     DaCycle.ParseTimes()
 
-    DaSystem    = CtDaSystem('../rc/carbontrackerjet.rc')
+    DaSystem    = CtDaSystem('../rc/carbontrackergridded.rc')
     DaSystem.Initialize()
 
     DaCycle.DaSystem    = DaSystem