Minor internal clarification and performance improvements

diff --git a/src/daComposant/daAlgorithms/FunctionTest.py b/src/daComposant/daAlgorithms/FunctionTest.py
index c7d622bab0276cfd7ab5d9ec04389560b2b7bbab..8da7096a243a7348826e97e3896d8009eb0ef2ea 100644 (file)
--- a/src/daComposant/daAlgorithms/FunctionTest.py
+++ b/src/daComposant/daAlgorithms/FunctionTest.py
@@ -97,7 +97,7 @@ class ElementaryAlgorithm(BasicObjects.Algorithm):
         msgs += ("     -----------------------------\n")
         msgs += ("     Characteristics of input vector X, internally converted:\n")
         msgs += ("       Type...............: %s\n")%type( Xn )
-        msgs += ("       Lenght of vector...: %i\n")%max(numpy.matrix( Xn ).shape)
+        msgs += ("       Lenght of vector...: %i\n")%max(numpy.ravel( Xn ).shape)
         msgs += ("       Minimum value......: %."+str(_p)+"e\n")%numpy.min( Xn )
         msgs += ("       Maximum value......: %."+str(_p)+"e\n")%numpy.max( Xn )
         msgs += ("       Mean of vector.....: %."+str(_p)+"e\n")%numpy.mean( Xn, dtype=mfp )
@@ -131,7 +131,7 @@ class ElementaryAlgorithm(BasicObjects.Algorithm):
             msgs  = ("===> Information after evaluation:\n")
             msgs += ("\n     Characteristics of simulated output vector Y=H(X), to compare to others:\n")
             msgs += ("       Type...............: %s\n")%type( Yn )
-            msgs += ("       Lenght of vector...: %i\n")%max(numpy.matrix( Yn ).shape)
+            msgs += ("       Lenght of vector...: %i\n")%max(numpy.ravel( Yn ).shape)
             msgs += ("       Minimum value......: %."+str(_p)+"e\n")%numpy.min( Yn )
             msgs += ("       Maximum value......: %."+str(_p)+"e\n")%numpy.max( Yn )
             msgs += ("       Mean of vector.....: %."+str(_p)+"e\n")%numpy.mean( Yn, dtype=mfp )

diff --git a/src/daComposant/daAlgorithms/ObserverTest.py b/src/daComposant/daAlgorithms/ObserverTest.py
index 39446d279573ed18997b1cba9aa6d4bbf2bc3418..323370adfaa60d57ede0346fd463182a7ad05ae0 100644 (file)
--- a/src/daComposant/daAlgorithms/ObserverTest.py
+++ b/src/daComposant/daAlgorithms/ObserverTest.py
@@ -68,7 +68,7 @@ class ElementaryAlgorithm(BasicObjects.Algorithm):
         self.StoredVariables["SigmaObs2"].store( 1. )
         self.StoredVariables["SigmaBck2"].store( 1. )
         self.StoredVariables["MahalanobisConsistency"].store( 1. )
-        self.StoredVariables["SimulationQuantiles"].store( numpy.matrix((__YY,__YY,__YY)) )
+        self.StoredVariables["SimulationQuantiles"].store( numpy.array((__YY,__YY,__YY)) )
         self.StoredVariables["SimulatedObservationAtBackground"].store( __YY )
         self.StoredVariables["SimulatedObservationAtCurrentState"].store( __YY )
         self.StoredVariables["SimulatedObservationAtOptimum"].store( __YY )

diff --git a/src/daComposant/daCore/Aidsm.py b/src/daComposant/daCore/Aidsm.py
index 1ca48cd27210de758b7d770f983a75eaefa65d40..79ddef25c1c5e5507d76078c734340055ada8027 100644 (file)
--- a/src/daComposant/daCore/Aidsm.py
+++ b/src/daComposant/daCore/Aidsm.py
@@ -763,9 +763,9 @@ class Aidsm(object):
         self.__adaoObject["AlgorithmParameters"].executePythonScheme( self.__adaoObject )
         if "UserPostAnalysis" in self.__adaoObject and len(self.__adaoObject["UserPostAnalysis"])>0:
             self.__objname = self.__retrieve_objname()
-            __Upa = [str(val).replace("ADD.","self.").replace(self.__objname+".","self.") for val in self.__adaoObject["UserPostAnalysis"]]
+            __Upa = map(str, self.__adaoObject["UserPostAnalysis"])
             __Upa = eval("\n".join(__Upa))
-            exec(__Upa, {}, {'self':self})
+            exec(__Upa, {}, {'self':self, 'ADD':self, 'case':self, 'adaopy':self, self.__objname:self})
         return 0
 
     def __executeYACSScheme(self, FileName=None):

diff --git a/src/daComposant/daCore/BasicObjects.py b/src/daComposant/daCore/BasicObjects.py
index bb0bec8883fd517cd2cdfca899caa1a3d4742b4d..aba702c66b5c44ffe56df8bc375f22eeb8801bdf 100644 (file)
--- a/src/daComposant/daCore/BasicObjects.py
+++ b/src/daComposant/daCore/BasicObjects.py
@@ -226,7 +226,7 @@ class Operator(object):
                 else:
                     if self.__Matrix is not None:
                         self.__addOneMatrixCall()
-                        _xv = numpy.matrix(numpy.ravel(xv)).T
+                        _xv = numpy.ravel(xv).reshape((-1,1))
                         _hv = self.__Matrix * _xv
                     else:
                         self.__addOneMethodCall()

diff --git a/src/daComposant/daCore/NumericObjects.py b/src/daComposant/daCore/NumericObjects.py
index c5b6306065be21dbc94b1261905156d0cecba92e..85ecff2048846829367362bc70444f4bf66b1b06 100644 (file)
--- a/src/daComposant/daCore/NumericObjects.py
+++ b/src/daComposant/daCore/NumericObjects.py
@@ -890,7 +890,7 @@ def c2ukf(selfA, Xb, Y, U, HO, EM, CM, R, B, Q):
         XEtnnp = []
         for point in range(nbSpts):
             if selfA._parameters["EstimationOf"] == "State":
-                XEtnnpi = numpy.asmatrix(numpy.ravel( Mm( (Xnp[:,point], Un) ) )).T
+                XEtnnpi = numpy.asarray( Mm( (Xnp[:,point], Un) ) ).reshape((-1,1))
                 if Cm is not None and Un is not None: # Attention : si Cm est aussi dans M, doublon !
                     Cm = Cm.reshape(Xn.size,Un.size) # ADAO & check shape
                     XEtnnpi = XEtnnpi + Cm * Un
@@ -899,10 +899,10 @@ def c2ukf(selfA, Xb, Y, U, HO, EM, CM, R, B, Q):
             elif selfA._parameters["EstimationOf"] == "Parameters":
                 # --- > Par principe, M = Id, Q = 0
                 XEtnnpi = Xnp[:,point]
-            XEtnnp.append( XEtnnpi )
-        XEtnnp = numpy.hstack( XEtnnp )
+            XEtnnp.append( numpy.ravel(XEtnnpi).reshape((-1,1)) )
+        XEtnnp = numpy.concatenate( XEtnnp, axis=1 )
         #
-        Xncm = numpy.matrix( XEtnnp.getA()*numpy.array(Wm) ).sum(axis=1)
+        Xncm = ( XEtnnp * Wm ).sum(axis=1)
         #
         if selfA._parameters["Bounds"] is not None and selfA._parameters["ConstrainedBy"] == "EstimateProjection":
             Xncm = ApplyBounds( Xncm, selfA._parameters["Bounds"] )
@@ -910,14 +910,14 @@ def c2ukf(selfA, Xb, Y, U, HO, EM, CM, R, B, Q):
         if selfA._parameters["EstimationOf"] == "State":        Pnm = Q
         elif selfA._parameters["EstimationOf"] == "Parameters": Pnm = 0.
         for point in range(nbSpts):
-            Pnm += Wc[i] * (XEtnnp[:,point]-Xncm) * (XEtnnp[:,point]-Xncm).T
+            Pnm += Wc[i] * ((XEtnnp[:,point]-Xncm).reshape((-1,1)) * (XEtnnp[:,point]-Xncm))
         #
         if selfA._parameters["EstimationOf"] == "Parameters" and selfA._parameters["Bounds"] is not None:
             Pnmdemi = selfA._parameters["Reconditioner"] * numpy.real(scipy.linalg.sqrtm(Pnm))
         else:
             Pnmdemi = numpy.real(scipy.linalg.sqrtm(Pnm))
         #
-        Xnnp = numpy.hstack([Xncm, Xncm+Gamma*Pnmdemi, Xncm-Gamma*Pnmdemi])
+        Xnnp = numpy.hstack([Xncm.reshape((-1,1)), Xncm.reshape((-1,1))+Gamma*Pnmdemi, Xncm.reshape((-1,1))-Gamma*Pnmdemi])
         #
         if selfA._parameters["Bounds"] is not None and selfA._parameters["ConstrainedBy"] == "EstimateProjection":
             for point in range(nbSpts):
@@ -926,27 +926,27 @@ def c2ukf(selfA, Xb, Y, U, HO, EM, CM, R, B, Q):
         Ynnp = []
         for point in range(nbSpts):
             if selfA._parameters["EstimationOf"] == "State":
-                Ynnpi = numpy.asmatrix(numpy.ravel( Hm( (Xnnp[:,point], None) ) )).T
+                Ynnpi = Hm( (Xnnp[:,point], None) )
             elif selfA._parameters["EstimationOf"] == "Parameters":
-                Ynnpi = numpy.asmatrix(numpy.ravel( Hm( (Xnnp[:,point], Un) ) )).T
-            Ynnp.append( Ynnpi )
-        Ynnp = numpy.hstack( Ynnp )
+                Ynnpi = Hm( (Xnnp[:,point], Un) )
+            Ynnp.append( numpy.ravel(Ynnpi).reshape((-1,1)) )
+        Ynnp = numpy.concatenate( Ynnp, axis=1 )
         #
-        Yncm = numpy.matrix( Ynnp.getA()*numpy.array(Wm) ).sum(axis=1)
+        Yncm = ( Ynnp * Wm ).sum(axis=1)
         #
         Pyyn = R
         Pxyn = 0.
         for point in range(nbSpts):
-            Pyyn += Wc[i] * (Ynnp[:,point]-Yncm) * (Ynnp[:,point]-Yncm).T
-            Pxyn += Wc[i] * (Xnnp[:,point]-Xncm) * (Ynnp[:,point]-Yncm).T
+            Pyyn += Wc[i] * ((Ynnp[:,point]-Yncm).reshape((-1,1)) * (Ynnp[:,point]-Yncm))
+            Pxyn += Wc[i] * ((Xnnp[:,point]-Xncm).reshape((-1,1)) * (Ynnp[:,point]-Yncm))
         #
-        _Innovation  = Ynpu - Yncm
+        _Innovation  = Ynpu - Yncm.reshape((-1,1))
         if selfA._parameters["EstimationOf"] == "Parameters":
             if Cm is not None and Un is not None: # Attention : si Cm est aussi dans H, doublon !
                 _Innovation = _Innovation - Cm * Un
         #
         Kn = Pxyn * Pyyn.I
-        Xn = Xncm + Kn * _Innovation
+        Xn = Xncm.reshape((-1,1)) + Kn * _Innovation
         Pn = Pnm - Kn * Pyyn * Kn.T
         #
         if selfA._parameters["Bounds"] is not None and selfA._parameters["ConstrainedBy"] == "EstimateProjection":
@@ -988,7 +988,7 @@ def c2ukf(selfA, Xb, Y, U, HO, EM, CM, R, B, Q):
             or selfA._toStore("CurrentOptimum") \
             or selfA._toStore("APosterioriCovariance"):
             Jb  = float( 0.5 * (Xa - Xb).T * BI * (Xa - Xb) )
-            Jo  = float( 0.5 * _Innovation.T * RI * _Innovation )
+            Jo  = float( 0.5 * _Innovation.T * (RI * _Innovation) )
             J   = Jb + Jo
             selfA.StoredVariables["CostFunctionJb"].store( Jb )
             selfA.StoredVariables["CostFunctionJo"].store( Jo )
@@ -1300,11 +1300,11 @@ def enks(selfA, Xb, Y, U, HO, EM, CM, R, B, Q, VariantM="EnKS16-KalmanFilterForm
         #
         if U is not None:
             if hasattr(U,"store") and len(U)>1:
-                Un = numpy.asmatrix(numpy.ravel( U[step] )).T
+                Un = numpy.ravel( U[step] ).reshape((-1,1))
             elif hasattr(U,"store") and len(U)==1:
-                Un = numpy.asmatrix(numpy.ravel( U[0] )).T
+                Un = numpy.ravel( U[0] ).reshape((-1,1))
             else:
-                Un = numpy.asmatrix(numpy.ravel( U )).T
+                Un = numpy.ravel( U ).reshape((-1,1))
         else:
             Un = None
         #
@@ -1442,11 +1442,11 @@ def etkf(selfA, Xb, Y, U, HO, EM, CM, R, B, Q, VariantM="KalmanFilterFormula"):
         #
         if U is not None:
             if hasattr(U,"store") and len(U)>1:
-                Un = numpy.asmatrix(numpy.ravel( U[step] )).T
+                Un = numpy.ravel( U[step] ).reshape((-1,1))
             elif hasattr(U,"store") and len(U)==1:
-                Un = numpy.asmatrix(numpy.ravel( U[0] )).T
+                Un = numpy.ravel( U[0] ).reshape((-1,1))
             else:
-                Un = numpy.asmatrix(numpy.ravel( U )).T
+                Un = numpy.ravel( U ).reshape((-1,1))
         else:
             Un = None
         #
@@ -4001,51 +4001,51 @@ def uskf(selfA, Xb, Y, U, HO, EM, CM, R, B, Q):
         XEtnnp = []
         for point in range(nbSpts):
             if selfA._parameters["EstimationOf"] == "State":
-                XEtnnpi = numpy.asmatrix(numpy.ravel( Mm( (Xnp[:,point], Un) ) )).T
+                XEtnnpi = numpy.asarray( Mm( (Xnp[:,point], Un) ) ).reshape((-1,1))
                 if Cm is not None and Un is not None: # Attention : si Cm est aussi dans M, doublon !
                     Cm = Cm.reshape(Xn.size,Un.size) # ADAO & check shape
                     XEtnnpi = XEtnnpi + Cm * Un
             elif selfA._parameters["EstimationOf"] == "Parameters":
                 # --- > Par principe, M = Id, Q = 0
                 XEtnnpi = Xnp[:,point]
-            XEtnnp.append( XEtnnpi )
-        XEtnnp = numpy.hstack( XEtnnp )
+            XEtnnp.append( numpy.ravel(XEtnnpi).reshape((-1,1)) )
+        XEtnnp = numpy.concatenate( XEtnnp, axis=1 )
         #
-        Xncm = numpy.matrix( XEtnnp.getA()*numpy.array(Wm) ).sum(axis=1)
+        Xncm = ( XEtnnp * Wm ).sum(axis=1)
         #
         if selfA._parameters["EstimationOf"] == "State":        Pnm = Q
         elif selfA._parameters["EstimationOf"] == "Parameters": Pnm = 0.
         for point in range(nbSpts):
-            Pnm += Wc[i] * (XEtnnp[:,point]-Xncm) * (XEtnnp[:,point]-Xncm).T
+            Pnm += Wc[i] * ((XEtnnp[:,point]-Xncm).reshape((-1,1)) * (XEtnnp[:,point]-Xncm))
         #
         Pnmdemi = numpy.real(scipy.linalg.sqrtm(Pnm))
         #
-        Xnnp = numpy.hstack([Xncm, Xncm+Gamma*Pnmdemi, Xncm-Gamma*Pnmdemi])
+        Xnnp = numpy.hstack([Xncm.reshape((-1,1)), Xncm.reshape((-1,1))+Gamma*Pnmdemi, Xncm.reshape((-1,1))-Gamma*Pnmdemi])
         #
         Ynnp = []
         for point in range(nbSpts):
             if selfA._parameters["EstimationOf"] == "State":
-                Ynnpi = numpy.asmatrix(numpy.ravel( Hm( (Xnnp[:,point], None) ) )).T
+                Ynnpi = Hm( (Xnnp[:,point], None) )
             elif selfA._parameters["EstimationOf"] == "Parameters":
-                Ynnpi = numpy.asmatrix(numpy.ravel( Hm( (Xnnp[:,point], Un) ) )).T
-            Ynnp.append( Ynnpi )
-        Ynnp = numpy.hstack( Ynnp )
+                Ynnpi = Hm( (Xnnp[:,point], Un) )
+            Ynnp.append( numpy.ravel(Ynnpi).reshape((-1,1)) )
+        Ynnp = numpy.concatenate( Ynnp, axis=1 )
         #
-        Yncm = numpy.matrix( Ynnp.getA()*numpy.array(Wm) ).sum(axis=1)
+        Yncm = ( Ynnp * Wm ).sum(axis=1)
         #
         Pyyn = R
         Pxyn = 0.
         for point in range(nbSpts):
-            Pyyn += Wc[i] * (Ynnp[:,point]-Yncm) * (Ynnp[:,point]-Yncm).T
-            Pxyn += Wc[i] * (Xnnp[:,point]-Xncm) * (Ynnp[:,point]-Yncm).T
+            Pyyn += Wc[i] * ((Ynnp[:,point]-Yncm).reshape((-1,1)) * (Ynnp[:,point]-Yncm))
+            Pxyn += Wc[i] * ((Xnnp[:,point]-Xncm).reshape((-1,1)) * (Ynnp[:,point]-Yncm))
         #
-        _Innovation  = Ynpu - Yncm
+        _Innovation  = Ynpu - Yncm.reshape((-1,1))
         if selfA._parameters["EstimationOf"] == "Parameters":
             if Cm is not None and Un is not None: # Attention : si Cm est aussi dans H, doublon !
                 _Innovation = _Innovation - Cm * Un
         #
         Kn = Pxyn * Pyyn.I
-        Xn = Xncm + Kn * _Innovation
+        Xn = Xncm.reshape((-1,1)) + Kn * _Innovation
         Pn = Pnm - Kn * Pyyn * Kn.T
         #
         Xa = Xn # Pointeurs
@@ -4084,7 +4084,7 @@ def uskf(selfA, Xb, Y, U, HO, EM, CM, R, B, Q):
             or selfA._toStore("CurrentOptimum") \
             or selfA._toStore("APosterioriCovariance"):
             Jb  = float( 0.5 * (Xa - Xb).T * BI * (Xa - Xb) )
-            Jo  = float( 0.5 * _Innovation.T * RI * _Innovation )
+            Jo  = float( 0.5 * _Innovation.T * (RI * _Innovation) )
             J   = Jb + Jo
             selfA.StoredVariables["CostFunctionJb"].store( Jb )
             selfA.StoredVariables["CostFunctionJo"].store( Jo )