Minor improvements and fixes for internal variables

diff --git a/src/daComposant/daAlgorithms/3DVAR.py b/src/daComposant/daAlgorithms/3DVAR.py
index 1939d6c91529bacfc4c04dc34136b2ef2bd36e57..fa5e47c88fa019a0e084ef7003164fcedccb81a2 100644 (file)
--- a/src/daComposant/daAlgorithms/3DVAR.py
+++ b/src/daComposant/daAlgorithms/3DVAR.py
@@ -129,6 +129,7 @@ class ElementaryAlgorithm(BasicObjects.Algorithm):
                 "MahalanobisConsistency",
                 "OMA",
                 "OMB",
+                "SampledStateForQuantiles",
                 "SigmaObs2",
                 "SimulatedObservationAtBackground",
                 "SimulatedObservationAtCurrentOptimum",

diff --git a/src/daComposant/daCore/NumericObjects.py b/src/daComposant/daCore/NumericObjects.py
index e9e274d79e61f07a9ccb31d7d3c5d9cefb3bf567..fdf94021800e9ca7ef0364e37e2b3191bb1a140c 100644 (file)
--- a/src/daComposant/daCore/NumericObjects.py
+++ b/src/daComposant/daCore/NumericObjects.py
@@ -1583,19 +1583,23 @@ def incr3dvar(selfA, Xb, Y, U, HO, EM, CM, R, B, Q):
     if selfA._toStore("SimulationQuantiles"):
         nech = selfA._parameters["NumberOfSamplesForQuantiles"]
         HXa  = numpy.matrix(numpy.ravel( HXa )).T
+        EXr  = None
         YfQ  = None
         for i in range(nech):
             if selfA._parameters["SimulationForQuantiles"] == "Linear":
                 dXr = numpy.matrix(numpy.random.multivariate_normal(Xa.A1,A) - Xa.A1).T
                 dYr = numpy.matrix(numpy.ravel( HtM * dXr )).T
                 Yr = HXa + dYr
+                if selfA._toStore("SampledStateForQuantiles"): Xr = Xa+dXr
             elif selfA._parameters["SimulationForQuantiles"] == "NonLinear":
                 Xr = numpy.matrix(numpy.random.multivariate_normal(Xa.A1,A)).T
                 Yr = numpy.matrix(numpy.ravel( Hm( Xr ) )).T
             if YfQ is None:
                 YfQ = Yr
+                if selfA._toStore("SampledStateForQuantiles"): EXr = numpy.ravel(Xr)
             else:
                 YfQ = numpy.hstack((YfQ,Yr))
+                if selfA._toStore("SampledStateForQuantiles"): EXr = numpy.vstack((EXr,numpy.ravel(Xr)))
         YfQ.sort(axis=-1)
         YQ = None
         for quantile in selfA._parameters["Quantiles"]:
@@ -1604,6 +1608,8 @@ def incr3dvar(selfA, Xb, Y, U, HO, EM, CM, R, B, Q):
             if YQ is None: YQ = YfQ[:,indice]
             else:          YQ = numpy.hstack((YQ,YfQ[:,indice]))
         selfA.StoredVariables["SimulationQuantiles"].store( YQ )
+        if selfA._toStore("SampledStateForQuantiles"):
+            selfA.StoredVariables["SampledStateForQuantiles"].store( EXr.T )
     if selfA._toStore("SimulatedObservationAtBackground"):
         selfA.StoredVariables["SimulatedObservationAtBackground"].store( numpy.ravel(HXb) )
     if selfA._toStore("SimulatedObservationAtOptimum"):
@@ -2227,19 +2233,23 @@ def psas3dvar(selfA, Xb, Y, U, HO, EM, CM, R, B, Q):
     if selfA._toStore("SimulationQuantiles"):
         nech = selfA._parameters["NumberOfSamplesForQuantiles"]
         HXa  = numpy.matrix(numpy.ravel( HXa )).T
+        EXr  = None
         YfQ  = None
         for i in range(nech):
             if selfA._parameters["SimulationForQuantiles"] == "Linear":
                 dXr = numpy.matrix(numpy.random.multivariate_normal(Xa.A1,A) - Xa.A1).T
                 dYr = numpy.matrix(numpy.ravel( HtM * dXr )).T
                 Yr = HXa + dYr
+                if selfA._toStore("SampledStateForQuantiles"): Xr = Xa+dXr
             elif selfA._parameters["SimulationForQuantiles"] == "NonLinear":
                 Xr = numpy.matrix(numpy.random.multivariate_normal(Xa.A1,A)).T
                 Yr = numpy.matrix(numpy.ravel( Hm( Xr ) )).T
             if YfQ is None:
                 YfQ = Yr
+                if selfA._toStore("SampledStateForQuantiles"): EXr = numpy.ravel(Xr)
             else:
                 YfQ = numpy.hstack((YfQ,Yr))
+                if selfA._toStore("SampledStateForQuantiles"): EXr = numpy.vstack((EXr,numpy.ravel(Xr)))
         YfQ.sort(axis=-1)
         YQ = None
         for quantile in selfA._parameters["Quantiles"]:
@@ -2248,6 +2258,8 @@ def psas3dvar(selfA, Xb, Y, U, HO, EM, CM, R, B, Q):
             if YQ is None: YQ = YfQ[:,indice]
             else:          YQ = numpy.hstack((YQ,YfQ[:,indice]))
         selfA.StoredVariables["SimulationQuantiles"].store( YQ )
+        if selfA._toStore("SampledStateForQuantiles"):
+            selfA.StoredVariables["SampledStateForQuantiles"].store( EXr.T )
     if selfA._toStore("SimulatedObservationAtBackground"):
         selfA.StoredVariables["SimulatedObservationAtBackground"].store( numpy.ravel(HXb) )
     if selfA._toStore("SimulatedObservationAtOptimum"):
@@ -2264,7 +2276,6 @@ def senkf(selfA, Xb, Y, U, HO, EM, CM, R, B, Q, VariantM="KalmanFilterFormula"):
         selfA._parameters["StoreInternalVariables"] = True
     #
     # Opérateurs
-    # ----------
     H = HO["Direct"].appliedControledFormTo
     #
     if selfA._parameters["EstimationOf"] == "State":
@@ -2275,8 +2286,7 @@ def senkf(selfA, Xb, Y, U, HO, EM, CM, R, B, Q, VariantM="KalmanFilterFormula"):
     else:
         Cm = None
     #
-    # Nombre de pas identique au nombre de pas d'observations
-    # -------------------------------------------------------
+    # Durée d'observation et tailles
     if hasattr(Y,"stepnumber"):
         duration = Y.stepnumber()
         __p = numpy.cumprod(Y.shape())[-1]
@@ -2285,7 +2295,6 @@ def senkf(selfA, Xb, Y, U, HO, EM, CM, R, B, Q, VariantM="KalmanFilterFormula"):
         __p = numpy.array(Y).size
     #
     # Précalcul des inversions de B et R
-    # ----------------------------------
     if selfA._parameters["StoreInternalVariables"] \
         or selfA._toStore("CostFunctionJ") \
         or selfA._toStore("CostFunctionJb") \
@@ -2295,10 +2304,9 @@ def senkf(selfA, Xb, Y, U, HO, EM, CM, R, B, Q, VariantM="KalmanFilterFormula"):
         BI = B.getI()
         RI = R.getI()
     #
-    # Initialisation
-    # --------------
     __n = Xb.size
     __m = selfA._parameters["NumberOfMembers"]
+    #
     if hasattr(B,"asfullmatrix"): Pn = B.asfullmatrix(__n)
     else:                         Pn = B
     if hasattr(R,"asfullmatrix"): Rn = R.asfullmatrix(__p)
@@ -2738,19 +2746,23 @@ def std3dvar(selfA, Xb, Y, U, HO, EM, CM, R, B, Q):
     if selfA._toStore("SimulationQuantiles"):
         nech = selfA._parameters["NumberOfSamplesForQuantiles"]
         HXa  = numpy.matrix(numpy.ravel( HXa )).T
+        EXr  = None
         YfQ  = None
         for i in range(nech):
             if selfA._parameters["SimulationForQuantiles"] == "Linear":
                 dXr = numpy.matrix(numpy.random.multivariate_normal(Xa.A1,A) - Xa.A1).T
                 dYr = numpy.matrix(numpy.ravel( HtM * dXr )).T
                 Yr = HXa + dYr
+                if selfA._toStore("SampledStateForQuantiles"): Xr = Xa+dXr
             elif selfA._parameters["SimulationForQuantiles"] == "NonLinear":
                 Xr = numpy.matrix(numpy.random.multivariate_normal(Xa.A1,A)).T
                 Yr = numpy.matrix(numpy.ravel( Hm( Xr ) )).T
             if YfQ is None:
                 YfQ = Yr
+                if selfA._toStore("SampledStateForQuantiles"): EXr = numpy.ravel(Xr)
             else:
                 YfQ = numpy.hstack((YfQ,Yr))
+                if selfA._toStore("SampledStateForQuantiles"): EXr = numpy.vstack((EXr,numpy.ravel(Xr)))
         YfQ.sort(axis=-1)
         YQ = None
         for quantile in selfA._parameters["Quantiles"]:
@@ -2759,6 +2771,8 @@ def std3dvar(selfA, Xb, Y, U, HO, EM, CM, R, B, Q):
             if YQ is None: YQ = YfQ[:,indice]
             else:          YQ = numpy.hstack((YQ,YfQ[:,indice]))
         selfA.StoredVariables["SimulationQuantiles"].store( YQ )
+        if selfA._toStore("SampledStateForQuantiles"):
+            selfA.StoredVariables["SampledStateForQuantiles"].store( EXr.T )
     if selfA._toStore("SimulatedObservationAtBackground"):
         selfA.StoredVariables["SimulatedObservationAtBackground"].store( numpy.ravel(HXb) )
     if selfA._toStore("SimulatedObservationAtOptimum"):
@@ -3242,19 +3256,23 @@ def van3dvar(selfA, Xb, Y, U, HO, EM, CM, R, B, Q):
     if selfA._toStore("SimulationQuantiles"):
         nech = selfA._parameters["NumberOfSamplesForQuantiles"]
         HXa  = numpy.matrix(numpy.ravel( HXa )).T
+        EXr  = None
         YfQ  = None
         for i in range(nech):
             if selfA._parameters["SimulationForQuantiles"] == "Linear":
                 dXr = numpy.matrix(numpy.random.multivariate_normal(Xa.A1,A) - Xa.A1).T
                 dYr = numpy.matrix(numpy.ravel( HtM * dXr )).T
                 Yr = HXa + dYr
+                if selfA._toStore("SampledStateForQuantiles"): Xr = Xa+dXr
             elif selfA._parameters["SimulationForQuantiles"] == "NonLinear":
                 Xr = numpy.matrix(numpy.random.multivariate_normal(Xa.A1,A)).T
                 Yr = numpy.matrix(numpy.ravel( Hm( Xr ) )).T
             if YfQ is None:
                 YfQ = Yr
+                if selfA._toStore("SampledStateForQuantiles"): EXr = numpy.ravel(Xr)
             else:
                 YfQ = numpy.hstack((YfQ,Yr))
+                if selfA._toStore("SampledStateForQuantiles"): EXr = numpy.vstack((EXr,numpy.ravel(Xr)))
         YfQ.sort(axis=-1)
         YQ = None
         for quantile in selfA._parameters["Quantiles"]:
@@ -3263,6 +3281,8 @@ def van3dvar(selfA, Xb, Y, U, HO, EM, CM, R, B, Q):
             if YQ is None: YQ = YfQ[:,indice]
             else:          YQ = numpy.hstack((YQ,YfQ[:,indice]))
         selfA.StoredVariables["SimulationQuantiles"].store( YQ )
+        if selfA._toStore("SampledStateForQuantiles"):
+            selfA.StoredVariables["SampledStateForQuantiles"].store( EXr.T )
     if selfA._toStore("SimulatedObservationAtBackground"):
         selfA.StoredVariables["SimulatedObservationAtBackground"].store( numpy.ravel(HXb) )
     if selfA._toStore("SimulatedObservationAtOptimum"):