src/daComposant/daCore/AssimilationStudy.py

   1 #-*-coding:iso-8859-1-*-
   2 #
   3 #  Copyright (C) 2008-2011  EDF R&D
   4 #
   5 #  This library is free software; you can redistribute it and/or
   6 #  modify it under the terms of the GNU Lesser General Public
   7 #  License as published by the Free Software Foundation; either
   8 #  version 2.1 of the License.
   9 #
  10 #  This library is distributed in the hope that it will be useful,
  11 #  but WITHOUT ANY WARRANTY; without even the implied warranty of
  12 #  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13 #  Lesser General Public License for more details.
  14 #
  15 #  You should have received a copy of the GNU Lesser General Public
  16 #  License along with this library; if not, write to the Free Software
  17 #  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  18 #
  19 #  See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  20 #
  21 #  Author: Jean-Philippe Argaud, jean-philippe.argaud@edf.fr, EDF R&D
  22
  23 __doc__ = """
  24     Définit les outils généraux élémentaires.
  25
  26     Ce module est destiné à etre appelée par AssimilationStudy pour constituer
  27     les objets élémentaires de l'étude.
  28 """
  29 __author__ = "Jean-Philippe ARGAUD"
  30
  31 import os, sys
  32 import numpy
  33 import Logging ; Logging.Logging() # A importer en premier
  34 import Persistence
  35 from BasicObjects import Operator
  36
  37 # ==============================================================================
  38 class AssimilationStudy:
  39     """
  40     Cette classe sert d'interface pour l'utilisation de l'assimilation de
  41     données. Elle contient les méthodes ou accesseurs nécessaires à la
  42     construction d'un calcul d'assimilation.
  43     """
  44     def __init__(self, name=""):
  45         """
  46         Prévoit de conserver l'ensemble des variables nécssaires à un algorithme
  47         élémentaire. Ces variables sont ensuite disponibles pour implémenter un
  48         algorithme élémentaire particulier.
  49
  50         Background............: vecteur Xb
  51         Observation...........: vecteur Y (potentiellement temporel)
  52             d'observations
  53         State.................: vecteur d'état dont une partie est le vecteur de
  54             contrôle. Cette information n'est utile que si l'on veut faire des
  55             calculs sur l'état complet, mais elle n'est pas indispensable pour
  56             l'assimilation.
  57         Control...............: vecteur X contenant toutes les variables de
  58             contrôle, i.e. les paramètres ou l'état dont on veut estimer la
  59             valeur pour obtenir les observations
  60         ObservationOperator...: opérateur d'observation H
  61
  62         Les observations présentent une erreur dont la matrice de covariance est
  63         R. L'ébauche du vecteur de contrôle présente une erreur dont la matrice
  64         de covariance est B.
  65         """
  66         self.__name = str(name)
  67         self.__Xb = None
  68         self.__Y  = None
  69         self.__B  = None
  70         self.__R  = None
  71         self.__Q  = None
  72         self.__H  = {}
  73         self.__M  = {}
  74         #
  75         self.__X  = Persistence.OneVector()
  76         self.__Parameters        = {}
  77         self.__StoredDiagnostics = {}
  78         self.__StoredInputs      = {}
  79         #
  80         self.__B_scalar = None
  81         self.__R_scalar = None
  82         self.__Q_scalar = None
  83         self.__Parameters["B_scalar"] = None
  84         self.__Parameters["R_scalar"] = None
  85         self.__Parameters["Q_scalar"] = None
  86         #
  87         # Variables temporaires
  88         self.__algorithm         = {}
  89         self.__algorithmFile     = None
  90         self.__algorithmName     = None
  91         self.__diagnosticFile    = None
  92         #
  93         # Récupère le chemin du répertoire parent et l'ajoute au path
  94         # (Cela complète l'action de la classe PathManagement dans PlatformInfo,
  95         # qui est activée dans Persistence)
  96         self.__parent = os.path.abspath(os.path.join(os.path.dirname(__file__),".."))
  97         sys.path.insert(0, self.__parent)
  98         sys.path = list(set(sys.path)) # Conserve en unique exemplaire chaque chemin
  99
 100     # ---------------------------------------------------------
 101     def setBackground(self,
 102             asVector           = None,
 103             asPersistentVector = None,
 104             Scheduler          = None,
 105             toBeStored         = False,
 106             ):
 107         """
 108         Permet de définir l'estimation a priori :
 109         - asVector : entrée des données, comme un vecteur compatible avec le
 110           constructeur de numpy.matrix
 111         - asPersistentVector : entrée des données, comme un vecteur de type
 112           persistent contruit avec la classe ad-hoc "Persistence"
 113         - Scheduler est le contrôle temporel des données
 114         - toBeStored : booléen indiquant si la donnée d'entrée est sauvée pour
 115           être rendue disponible au même titre que les variables de calcul
 116         """
 117         if asVector is not None:
 118             if type( asVector ) is type( numpy.matrix([]) ):
 119                 self.__Xb = numpy.matrix( asVector.A1, numpy.float ).T
 120             else:
 121                 self.__Xb = numpy.matrix( asVector,    numpy.float ).T
 122         elif asPersistentVector is not None:
 123             self.__Xb = asPersistentVector
 124         else:
 125             raise ValueError("Error: improperly defined background")
 126         if toBeStored:
 127            self.__StoredInputs["Background"] = self.__Xb
 128         return 0
 129
 130     def setBackgroundError(self,
 131             asCovariance  = None,
 132             asEyeByScalar = None,
 133             toBeStored    = False,
 134             ):
 135         """
 136         Permet de définir la covariance des erreurs d'ébauche :
 137         - asCovariance : entrée des données, comme une matrice compatible avec
 138           le constructeur de numpy.matrix
 139         - asEyeByScalar : entrée des données comme un seul scalaire de variance,
 140           multiplicatif d'une matrice de corrélation identité, aucune matrice
 141           n'étant donc explicitement à donner
 142         - toBeStored : booléen indiquant si la donnée d'entrée est sauvée pour
 143           être rendue disponible au même titre que les variables de calcul
 144         """
 145         if asEyeByScalar is not None:
 146             self.__B_scalar = float(asEyeByScalar)
 147             self.__B        = None
 148         else:
 149             self.__B_scalar = None
 150             self.__B  = numpy.matrix( asCovariance, numpy.float )
 151         #
 152         self.__Parameters["B_scalar"] = self.__B_scalar
 153         if toBeStored:
 154             self.__StoredInputs["BackgroundError"] = self.__B
 155         return 0
 156
 157     # -----------------------------------------------------------
 158     def setObservation(self,
 159             asVector           = None,
 160             asPersistentVector = None,
 161             Scheduler          = None,
 162             toBeStored         = False,
 163             ):
 164         """
 165         Permet de définir les observations :
 166         - asVector : entrée des données, comme un vecteur compatible avec le
 167           constructeur de numpy.matrix
 168         - asPersistentVector : entrée des données, comme un vecteur de type
 169           persistent contruit avec la classe ad-hoc "Persistence"
 170         - Scheduler est le contrôle temporel des données disponibles
 171         - toBeStored : booléen indiquant si la donnée d'entrée est sauvée pour
 172           être rendue disponible au même titre que les variables de calcul
 173         """
 174         if asVector is not None:
 175             if type( asVector ) is type( numpy.matrix([]) ):
 176                 self.__Y = numpy.matrix( asVector.A1, numpy.float ).T
 177             else:
 178                 self.__Y = numpy.matrix( asVector,    numpy.float ).T
 179         elif asPersistentVector is not None:
 180             self.__Y = asPersistentVector
 181         else:
 182             raise ValueError("Error: improperly defined observations")
 183         if toBeStored:
 184             self.__StoredInputs["Observation"] = self.__Y
 185         return 0
 186
 187     def setObservationError(self,
 188             asCovariance  = None,
 189             asEyeByScalar = None,
 190             toBeStored    = False,
 191             ):
 192         """
 193         Permet de définir la covariance des erreurs d'observations :
 194         - asCovariance : entrée des données, comme une matrice compatible avec
 195           le constructeur de numpy.matrix
 196         - asEyeByScalar : entrée des données comme un seul scalaire de variance,
 197           multiplicatif d'une matrice de corrélation identité, aucune matrice
 198           n'étant donc explicitement à donner
 199         - toBeStored : booléen indiquant si la donnée d'entrée est sauvée pour
 200           être rendue disponible au même titre que les variables de calcul
 201         """
 202         if asEyeByScalar is not None:
 203             self.__R_scalar = float(asEyeByScalar)
 204             self.__R        = None
 205         else:
 206             self.__R_scalar = None
 207             self.__R  = numpy.matrix( asCovariance, numpy.float )
 208         #
 209         self.__Parameters["R_scalar"] = self.__R_scalar
 210         if toBeStored:
 211             self.__StoredInputs["ObservationError"] = self.__R
 212         return 0
 213
 214     def setObservationOperator(self,
 215             asFunction = {"Direct":None, "Tangent":None, "Adjoint":None},
 216             asMatrix   = None,
 217             appliedToX = None,
 218             toBeStored = False,
 219             ):
 220         """
 221         Permet de définir un opérateur d'observation H. L'ordre de priorité des
 222         définitions et leur sens sont les suivants :
 223         - si asFunction["Tangent"] et asFunction["Adjoint"] ne sont pas None
 224           alors on définit l'opérateur à l'aide de fonctions. Si la fonction
 225           "Direct" n'est pas définie, on prend la fonction "Tangent".
 226         - si les fonctions ne sont pas disponibles et si asMatrix n'est pas
 227           None, alors on définit l'opérateur "Direct" et "Tangent" à l'aide de
 228           la matrice, et l'opérateur "Adjoint" à l'aide de la transposée. La
 229           matrice fournie doit être sous une forme compatible avec le
 230           constructeur de numpy.matrix.
 231         - si l'argument "appliedToX" n'est pas None, alors on définit, pour des
 232           X divers, l'opérateur par sa valeur appliquée à cet X particulier,
 233           sous la forme d'un dictionnaire appliedToX[NAME] avec NAME un nom.
 234           L'opérateur doit néanmoins déjà avoir été défini comme d'habitude.
 235         - toBeStored : booléen indiquant si la donnée d'entrée est sauvée pour
 236           être rendue disponible au même titre que les variables de calcul
 237         """
 238         if (type(asFunction) is type({})) and \
 239                 asFunction.has_key("Tangent") and asFunction.has_key("Adjoint") and \
 240                 (asFunction["Tangent"] is not None) and (asFunction["Adjoint"] is not None):
 241             if not asFunction.has_key("Direct") or (asFunction["Direct"] is None):
 242                 self.__H["Direct"]  = Operator( fromMethod = asFunction["Tangent"]  )
 243             else:
 244                 self.__H["Direct"] = Operator( fromMethod = asFunction["Direct"]  )
 245             self.__H["Tangent"]    = Operator( fromMethod = asFunction["Tangent"] )
 246             self.__H["Adjoint"]    = Operator( fromMethod = asFunction["Adjoint"] )
 247         elif asMatrix is not None:
 248             mat = numpy.matrix( asMatrix, numpy.float )
 249             self.__H["Direct"]  = Operator( fromMatrix = mat )
 250             self.__H["Tangent"] = Operator( fromMatrix = mat )
 251             self.__H["Adjoint"] = Operator( fromMatrix = mat.T )
 252         else:
 253             raise ValueError("Error: improperly defined observation operator")
 254         #
 255         if appliedToX is not None:
 256             self.__H["AppliedToX"] = {}
 257             if type(appliedToX) is not dict:
 258                 raise ValueError("Error: observation operator defined by \"appliedToX\" need a dictionary as argument.")
 259             for key in appliedToX.keys():
 260                 if type( appliedToX[key] ) is type( numpy.matrix([]) ):
 261                     # Pour le cas où l'on a une vraie matrice
 262                     self.__H["AppliedToX"][key] = numpy.matrix( appliedToX[key].A1, numpy.float ).T
 263                 elif type( appliedToX[key] ) is type( numpy.array([]) ) and len(appliedToX[key].shape) > 1:
 264                     # Pour le cas où l'on a un vecteur représenté en array avec 2 dimensions
 265                     self.__H["AppliedToX"][key] = numpy.matrix( appliedToX[key].reshape(len(appliedToX[key]),), numpy.float ).T
 266                 else:
 267                     self.__H["AppliedToX"][key] = numpy.matrix( appliedToX[key],    numpy.float ).T
 268         else:
 269             self.__H["AppliedToX"] = None
 270         #
 271         if toBeStored:
 272             self.__StoredInputs["ObservationOperator"] = self.__H
 273         return 0
 274
 275     # -----------------------------------------------------------
 276     def setEvolutionModel(self,
 277             asFunction = {"Direct":None, "Tangent":None, "Adjoint":None},
 278             asMatrix   = None,
 279             Scheduler  = None,
 280             toBeStored = False,
 281             ):
 282         """
 283         Permet de définir un opérateur d'évolution M. L'ordre de priorité des
 284         définitions et leur sens sont les suivants :
 285         - si asFunction["Tangent"] et asFunction["Adjoint"] ne sont pas None
 286           alors on définit l'opérateur à l'aide de fonctions. Si la fonction
 287           "Direct" n'est pas définie, on prend la fonction "Tangent".
 288         - si les fonctions ne sont pas disponibles et si asMatrix n'est pas
 289           None, alors on définit l'opérateur "Direct" et "Tangent" à l'aide de
 290           la matrice, et l'opérateur "Adjoint" à l'aide de la transposée. La
 291           matrice fournie doit être sous une forme compatible avec le
 292           constructeur de numpy.matrix.
 293         - toBeStored : booléen indiquant si la donnée d'entrée est sauvée pour
 294           être rendue disponible au même titre que les variables de calcul
 295         """
 296         if (type(asFunction) is type({})) and (asFunction["Tangent"] is not None) and (asFunction["Adjoint"] is not None):
 297             if not asFunction.has_key("Direct") or (asFunction["Direct"] is None):
 298                 self.__M["Direct"] = Operator( fromMethod = asFunction["Tangent"]  )
 299             else:
 300                 self.__M["Direct"] = Operator( fromMethod = asFunction["Direct"]  )
 301             self.__M["Tangent"]    = Operator( fromMethod = asFunction["Tangent"] )
 302             self.__M["Adjoint"]    = Operator( fromMethod = asFunction["Adjoint"] )
 303         elif asMatrix is not None:
 304             matrice = numpy.matrix( asMatrix, numpy.float )
 305             self.__M["Direct"]  = Operator( fromMatrix = matrice )
 306             self.__M["Tangent"] = Operator( fromMatrix = matrice )
 307             self.__M["Adjoint"] = Operator( fromMatrix = matrice.T )
 308         else:
 309             raise ValueError("Error: improperly defined evolution operator")
 310         #
 311         if toBeStored:
 312             self.__StoredInputs["EvolutionModel"] = self.__M
 313         return 0
 314
 315     def setEvolutionError(self,
 316             asCovariance  = None,
 317             asEyeByScalar = None,
 318             toBeStored    = False,
 319             ):
 320         """
 321         Permet de définir la covariance des erreurs de modèle :
 322         - asCovariance : entrée des données, comme une matrice compatible avec
 323           le constructeur de numpy.matrix
 324         - asEyeByScalar : entrée des données comme un seul scalaire de variance,
 325           multiplicatif d'une matrice de corrélation identité, aucune matrice
 326           n'étant donc explicitement à donner
 327         - toBeStored : booléen indiquant si la donnée d'entrée est sauvée pour
 328           être rendue disponible au même titre que les variables de calcul
 329         """
 330         if asEyeByScalar is not None:
 331             self.__Q_scalar = float(asEyeByScalar)
 332             self.__Q        = None
 333         else:
 334             self.__Q_scalar = None
 335             self.__Q  = numpy.matrix( asCovariance, numpy.float )
 336         #
 337         self.__Parameters["Q_scalar"] = self.__Q_scalar
 338         if toBeStored:
 339             self.__StoredInputs["EvolutionError"] = self.__Q
 340         return 0
 341
 342     # -----------------------------------------------------------
 343     def setControls (self,
 344             asVector = None,
 345             toBeStored   = False,
 346             ):
 347         """
 348         Permet de définir la valeur initiale du vecteur X contenant toutes les
 349         variables de contrôle, i.e. les paramètres ou l'état dont on veut
 350         estimer la valeur pour obtenir les observations. C'est utile pour un
 351         algorithme itératif/incrémental.
 352         - asVector : entrée des données, comme un vecteur compatible avec le
 353           constructeur de numpy.matrix.
 354         - toBeStored : booléen indiquant si la donnée d'entrée est sauvée pour
 355           être rendue disponible au même titre que les variables de calcul
 356         """
 357         if asVector is not None:
 358             self.__X.store( asVector )
 359         if toBeStored:
 360             self.__StoredInputs["Controls"] = self.__X
 361         return 0
 362
 363     # -----------------------------------------------------------
 364     def setAlgorithm(self, choice = None ):
 365         """
 366         Permet de sélectionner l'algorithme à utiliser pour mener à bien l'étude
 367         d'assimilation. L'argument est un champ caractère se rapportant au nom
 368         d'un fichier contenu dans "../daAlgorithms" et réalisant l'opération
 369         d'assimilation sur les arguments fixes.
 370         """
 371         if choice is None:
 372             raise ValueError("Error: algorithm choice has to be given")
 373         if self.__algorithmName is not None:
 374             raise ValueError("Error: algorithm choice has already been done as \"%s\", it can't be changed."%self.__algorithmName)
 375         daDirectory = "daAlgorithms"
 376         #
 377         # Recherche explicitement le fichier complet
 378         # ------------------------------------------
 379         module_path = None
 380         for directory in sys.path:
 381             if os.path.isfile(os.path.join(directory, daDirectory, str(choice)+'.py')):
 382                 module_path = os.path.abspath(os.path.join(directory, daDirectory))
 383         if module_path is None:
 384             raise ImportError("No algorithm module named \"%s\" was found in a \"%s\" subdirectory\n             The search path is %s"%(choice, daDirectory, sys.path))
 385         #
 386         # Importe le fichier complet comme un module
 387         # ------------------------------------------
 388         try:
 389             sys_path_tmp = sys.path ; sys.path.insert(0,module_path)
 390             self.__algorithmFile = __import__(str(choice), globals(), locals(), [])
 391             self.__algorithmName = str(choice)
 392             sys.path = sys_path_tmp ; del sys_path_tmp
 393         except ImportError, e:
 394             raise ImportError("The module named \"%s\" was found, but is incorrect at the import stage.\n             The import error message is: %s"%(choice,e))
 395         #
 396         # Instancie un objet du type élémentaire du fichier
 397         # -------------------------------------------------
 398         self.__algorithm = self.__algorithmFile.ElementaryAlgorithm()
 399         self.__StoredInputs["AlgorithmName"] = str(choice)
 400         return 0
 401
 402     def setAlgorithmParameters(self, asDico=None):
 403         """
 404         Permet de définir les paramètres de l'algorithme, sous la forme d'un
 405         dictionnaire.
 406         """
 407         if asDico is not None:
 408             self.__Parameters.update( dict( asDico ) )
 409         #
 410         self.__StoredInputs["AlgorithmParameters"] = self.__Parameters
 411         return 0
 412
 413     # -----------------------------------------------------------
 414     def setDiagnostic(self, choice = None, name = "", unit = "", basetype = None, parameters = {} ):
 415         """
 416         Permet de sélectionner un diagnostic a effectuer.
 417         """
 418         if choice is None:
 419             raise ValueError("Error: diagnostic choice has to be given")
 420         daDirectory = "daDiagnostics"
 421         #
 422         # Recherche explicitement le fichier complet
 423         # ------------------------------------------
 424         module_path = None
 425         for directory in sys.path:
 426             if os.path.isfile(os.path.join(directory, daDirectory, str(choice)+'.py')):
 427                 module_path = os.path.abspath(os.path.join(directory, daDirectory))
 428         if module_path is None:
 429             raise ImportError("No diagnostic module named \"%s\" was found in a \"%s\" subdirectory\n             The search path is %s"%(choice, daDirectory, sys.path))
 430         #
 431         # Importe le fichier complet comme un module
 432         # ------------------------------------------
 433         try:
 434             sys_path_tmp = sys.path ; sys.path.insert(0,module_path)
 435             self.__diagnosticFile = __import__(str(choice), globals(), locals(), [])
 436             sys.path = sys_path_tmp ; del sys_path_tmp
 437         except ImportError, e:
 438             raise ImportError("The module named \"%s\" was found, but is incorrect at the import stage.\n             The import error message is: %s"%(choice,e))
 439         #
 440         # Instancie un objet du type élémentaire du fichier
 441         # -------------------------------------------------
 442         if self.__StoredInputs.has_key(name):
 443             raise ValueError("A default input with the same name \"%s\" already exists."%str(name))
 444         elif self.__StoredDiagnostics.has_key(name):
 445             raise ValueError("A diagnostic with the same name \"%s\" already exists."%str(name))
 446         else:
 447             self.__StoredDiagnostics[name] = self.__diagnosticFile.ElementaryDiagnostic(
 448                 name       = name,
 449                 unit       = unit,
 450                 basetype   = basetype,
 451                 parameters = parameters )
 452         return 0
 453
 454     # -----------------------------------------------------------
 455     def shape_validate(self):
 456         """
 457         Validation de la correspondance correcte des tailles des variables et
 458         des matrices s'il y en a.
 459         """
 460         if self.__Xb is None:                  __Xb_shape = (0,)
 461         elif hasattr(self.__Xb,"shape"):
 462             if type(self.__Xb.shape) is tuple: __Xb_shape = self.__Xb.shape
 463             else:                              __Xb_shape = self.__Xb.shape()
 464         else: raise TypeError("Xb has no attribute of shape: problem !")
 465         #
 466         if self.__Y is None:                  __Y_shape = (0,)
 467         elif hasattr(self.__Y,"shape"):
 468             if type(self.__Y.shape) is tuple: __Y_shape = self.__Y.shape
 469             else:                             __Y_shape = self.__Y.shape()
 470         else: raise TypeError("Y has no attribute of shape: problem !")
 471         #
 472         if self.__B is None and self.__B_scalar is None:
 473             __B_shape = (0,0)
 474         elif self.__B is None and self.__B_scalar is not None:
 475             __B_shape = (max(__Xb_shape),max(__Xb_shape))
 476         elif hasattr(self.__B,"shape"):
 477             if type(self.__B.shape) is tuple: __B_shape = self.__B.shape
 478             else:                             __B_shape = self.__B.shape()
 479         else: raise TypeError("B has no attribute of shape: problem !")
 480         #
 481         if self.__R is None and self.__R_scalar is None:
 482             __R_shape = (0,0)
 483         elif self.__R is None and self.__R_scalar is not None:
 484             __R_shape = (max(__Y_shape),max(__Y_shape))
 485         elif hasattr(self.__R,"shape"):
 486             if type(self.__R.shape) is tuple: __R_shape = self.__R.shape
 487             else:                             __R_shape = self.__R.shape()
 488         else: raise TypeError("R has no attribute of shape: problem !")
 489         #
 490         if self.__Q is None:                  __Q_shape = (0,0)
 491         elif hasattr(self.__Q,"shape"):
 492             if type(self.__Q.shape) is tuple: __Q_shape = self.__Q.shape
 493             else:                             __Q_shape = self.__Q.shape()
 494         else: raise TypeError("Q has no attribute of shape: problem !")
 495         #
 496         if len(self.__H) == 0:                          __H_shape = (0,0)
 497         elif type(self.__H) is type({}):                __H_shape = (0,0)
 498         elif hasattr(self.__H["Direct"],"shape"):
 499             if type(self.__H["Direct"].shape) is tuple: __H_shape = self.__H["Direct"].shape
 500             else:                                       __H_shape = self.__H["Direct"].shape()
 501         else: raise TypeError("H has no attribute of shape: problem !")
 502         #
 503         if len(self.__M) == 0:                          __M_shape = (0,0)
 504         elif type(self.__M) is type({}):                __M_shape = (0,0)
 505         elif hasattr(self.__M["Direct"],"shape"):
 506             if type(self.__M["Direct"].shape) is tuple: __M_shape = self.__M["Direct"].shape
 507             else:                                       __M_shape = self.__M["Direct"].shape()
 508         else: raise TypeError("M has no attribute of shape: problem !")
 509         #
 510         # Vérification des conditions
 511         # ---------------------------
 512         if not( len(__Xb_shape) == 1 or min(__Xb_shape) == 1 ):
 513             raise ValueError("Shape characteristic of Xb is incorrect: \"%s\""%(__Xb_shape,))
 514         if not( len(__Y_shape) == 1 or min(__Y_shape) == 1 ):
 515             raise ValueError("Shape characteristic of Y is incorrect: \"%s\""%(__Y_shape,))
 516         #
 517         if not( min(__B_shape) == max(__B_shape) ):
 518             raise ValueError("Shape characteristic of B is incorrect: \"%s\""%(__B_shape,))
 519         if not( min(__R_shape) == max(__R_shape) ):
 520             raise ValueError("Shape characteristic of R is incorrect: \"%s\""%(__R_shape,))
 521         if not( min(__Q_shape) == max(__Q_shape) ):
 522             raise ValueError("Shape characteristic of Q is incorrect: \"%s\""%(__Q_shape,))
 523         if not( min(__M_shape) == max(__M_shape) ):
 524             raise ValueError("Shape characteristic of M is incorrect: \"%s\""%(__M_shape,))
 525         #
 526         if len(self.__H) > 0 and not(type(self.__H) is type({})) and not( __H_shape[1] == max(__Xb_shape) ):
 527             raise ValueError("Shape characteristic of H \"%s\" and X \"%s\" are incompatible"%(__H_shape,__Xb_shape))
 528         if len(self.__H) > 0 and not(type(self.__H) is type({})) and not( __H_shape[0] == max(__Y_shape) ):
 529             raise ValueError("Shape characteristic of H \"%s\" and Y \"%s\" are incompatible"%(__H_shape,__Y_shape))
 530         if len(self.__H) > 0 and not(type(self.__H) is type({})) and len(self.__B) > 0 and not( __H_shape[1] == __B_shape[0] ):
 531             raise ValueError("Shape characteristic of H \"%s\" and B \"%s\" are incompatible"%(__H_shape,__B_shape))
 532         if len(self.__H) > 0 and not(type(self.__H) is type({})) and len(self.__R) > 0 and not( __H_shape[0] == __R_shape[1] ):
 533             raise ValueError("Shape characteristic of H \"%s\" and R \"%s\" are incompatible"%(__H_shape,__R_shape))
 534         #
 535         if self.__B is not None and len(self.__B) > 0 and not( __B_shape[1] == max(__Xb_shape) ):
 536             raise ValueError("Shape characteristic of B \"%s\" and Xb \"%s\" are incompatible"%(__B_shape,__Xb_shape))
 537         #
 538         if self.__R is not None and len(self.__R) > 0 and not( __R_shape[1] == max(__Y_shape) ):
 539             raise ValueError("Shape characteristic of R \"%s\" and Y \"%s\" are incompatible"%(__R_shape,__Y_shape))
 540         #
 541         if self.__M is not None and len(self.__M) > 0 and not(type(self.__M) is type({})) and not( __M_shape[1] == max(__Xb_shape) ):
 542             raise ValueError("Shape characteristic of M \"%s\" and X \"%s\" are incompatible"%(__M_shape,__Xb_shape))
 543         #
 544         return 1
 545
 546     # -----------------------------------------------------------
 547     def analyze(self):
 548         """
 549         Permet de lancer le calcul d'assimilation.
 550
 551         Le nom de la méthode à activer est toujours "run". Les paramètres en
 552         arguments de la méthode sont fixés. En sortie, on obtient les résultats
 553         dans la variable de type dictionnaire "StoredVariables", qui contient en
 554         particulier des objets de Persistence pour les analyses, OMA...
 555         """
 556         self.shape_validate()
 557         #
 558         self.__algorithm.run(
 559             Xb         = self.__Xb,
 560             Y          = self.__Y,
 561             H          = self.__H,
 562             M          = self.__M,
 563             R          = self.__R,
 564             B          = self.__B,
 565             Q          = self.__Q,
 566             Parameters = self.__Parameters,
 567             )
 568         return 0
 569
 570     # -----------------------------------------------------------
 571     def get(self, key=None):
 572         """
 573         Renvoie les résultats disponibles après l'exécution de la méthode
 574         d'assimilation, ou les diagnostics disponibles. Attention, quand un
 575         diagnostic porte le même nom qu'une variable stockée, c'est la variable
 576         stockée qui est renvoyée, et le diagnostic est inatteignable.
 577         """
 578         if key is not None:
 579             if self.__algorithm.has_key(key):
 580                 return self.__algorithm.get( key )
 581             elif self.__StoredInputs.has_key(key):
 582                 return self.__StoredInputs[key]
 583             elif self.__StoredDiagnostics.has_key(key):
 584                 return self.__StoredDiagnostics[key]
 585             else:
 586                 raise ValueError("The requested key \"%s\" does not exists as an input, a diagnostic or a stored variable."%key)
 587         else:
 588             allvariables = self.__algorithm.get()
 589             allvariables.update( self.__StoredDiagnostics )
 590             allvariables.update( self.__StoredInputs )
 591             return allvariables
 592
 593     def get_available_variables(self):
 594         """
 595         Renvoie les variables potentiellement utilisables pour l'étude,
 596         initialement stockées comme données d'entrées ou dans les algorithmes,
 597         identifiés par les chaînes de caractères. L'algorithme doit avoir été
 598         préalablement choisi sinon la méthode renvoie "None".
 599         """
 600         if len( self.__algorithm.keys()) == 0 and len( self.__StoredInputs.keys() ) == 0:
 601             return None
 602         else:
 603             variables = []
 604             if len( self.__algorithm.keys()) > 0:
 605                 variables.extend( self.__algorithm.get().keys() )
 606             if len( self.__StoredInputs.keys() ) > 0:
 607                 variables.extend( self.__StoredInputs.keys() )
 608             variables.sort()
 609             return variables
 610
 611     def get_available_algorithms(self):
 612         """
 613         Renvoie la liste des algorithmes potentiellement utilisables, identifiés
 614         par les chaînes de caractères.
 615         """
 616         files = []
 617         for directory in sys.path:
 618             if os.path.isdir(os.path.join(directory,"daAlgorithms")):
 619                 for fname in os.listdir(os.path.join(directory,"daAlgorithms")):
 620                     root, ext = os.path.splitext(fname)
 621                     if ext == '.py' and root != '__init__':
 622                         files.append(root)
 623         files.sort()
 624         return files
 625
 626     def get_available_diagnostics(self):
 627         """
 628         Renvoie la liste des diagnostics potentiellement utilisables, identifiés
 629         par les chaînes de caractères.
 630         """
 631         files = []
 632         for directory in sys.path:
 633             if os.path.isdir(os.path.join(directory,"daDiagnostics")):
 634                 for fname in os.listdir(os.path.join(directory,"daDiagnostics")):
 635                     root, ext = os.path.splitext(fname)
 636                     if ext == '.py' and root != '__init__':
 637                         files.append(root)
 638         files.sort()
 639         return files
 640
 641     # -----------------------------------------------------------
 642     def get_algorithms_main_path(self):
 643         """
 644         Renvoie le chemin pour le répertoire principal contenant les algorithmes
 645         dans un sous-répertoire "daAlgorithms"
 646         """
 647         return self.__parent
 648
 649     def add_algorithms_path(self, asPath=None):
 650         """
 651         Ajoute au chemin de recherche des algorithmes un répertoire dans lequel
 652         se trouve un sous-répertoire "daAlgorithms"
 653
 654         Remarque : si le chemin a déjà été ajouté pour les diagnostics, il n'est
 655         pas indispensable de le rajouter ici.
 656         """
 657         if not os.path.isdir(asPath):
 658             raise ValueError("The given "+asPath+" argument must exist as a directory")
 659         if not os.path.isdir(os.path.join(asPath,"daAlgorithms")):
 660             raise ValueError("The given \""+asPath+"\" argument must contain a subdirectory named \"daAlgorithms\"")
 661         if not os.path.isfile(os.path.join(asPath,"daAlgorithms","__init__.py")):
 662             raise ValueError("The given \""+asPath+"/daAlgorithms\" path must contain a file named \"__init__.py\"")
 663         sys.path.insert(0, os.path.abspath(asPath))
 664         sys.path = list(set(sys.path)) # Conserve en unique exemplaire chaque chemin
 665         return 1
 666
 667     def get_diagnostics_main_path(self):
 668         """
 669         Renvoie le chemin pour le répertoire principal contenant les diagnostics
 670         dans un sous-répertoire "daDiagnostics"
 671         """
 672         return self.__parent
 673
 674     def add_diagnostics_path(self, asPath=None):
 675         """
 676         Ajoute au chemin de recherche des algorithmes un répertoire dans lequel
 677         se trouve un sous-répertoire "daDiagnostics"
 678
 679         Remarque : si le chemin a déjà été ajouté pour les algorithmes, il n'est
 680         pas indispensable de le rajouter ici.
 681         """
 682         if not os.path.isdir(asPath):
 683             raise ValueError("The given "+asPath+" argument must exist as a directory")
 684         if not os.path.isdir(os.path.join(asPath,"daDiagnostics")):
 685             raise ValueError("The given \""+asPath+"\" argument must contain a subdirectory named \"daDiagnostics\"")
 686         if not os.path.isfile(os.path.join(asPath,"daDiagnostics","__init__.py")):
 687             raise ValueError("The given \""+asPath+"/daDiagnostics\" path must contain a file named \"__init__.py\"")
 688         sys.path.insert(0, os.path.abspath(asPath))
 689         sys.path = list(set(sys.path)) # Conserve en unique exemplaire chaque chemin
 690         return 1
 691
 692     # -----------------------------------------------------------
 693     def setDataObserver(self,
 694             VariableName   = None,
 695             HookFunction   = None,
 696             HookParameters = None,
 697             Scheduler      = None,
 698             ):
 699         """
 700         Permet d'associer un observer à une ou des variables nommées gérées en
 701         interne, activable selon des règles définies dans le Scheduler. A chaque
 702         pas demandé dans le Scheduler, il effectue la fonction HookFunction avec
 703         les arguments (variable persistante VariableName, paramètres HookParameters).
 704         """
 705         #
 706         if type( self.__algorithm ) is dict:
 707             raise ValueError("No observer can be build before choosing an algorithm.")
 708         #
 709         # Vérification du nom de variable et typage
 710         # -----------------------------------------
 711         if type( VariableName ) is str:
 712             VariableNames = [VariableName,]
 713         elif type( VariableName ) is list:
 714             VariableNames = map( str, VariableName )
 715         else:
 716             raise ValueError("The observer requires a name or a list of names of variables.")
 717         #
 718         # Association interne de l'observer à la variable
 719         # -----------------------------------------------
 720         for n in VariableNames:
 721             if not self.__algorithm.has_key( n ):
 722                 raise ValueError("An observer requires to be set on a variable named %s which does not exist."%n)
 723         else:
 724             self.__algorithm.StoredVariables[ n ].setDataObserver(
 725                 Scheduler      = Scheduler,
 726                 HookFunction   = HookFunction,
 727                 HookParameters = HookParameters,
 728                 )
 729
 730     def removeDataObserver(self,
 731             VariableName   = None,
 732             HookFunction   = None,
 733             ):
 734         """
 735         Permet de retirer un observer à une ou des variable nommée.
 736         """
 737         #
 738         if type( self.__algorithm ) is dict:
 739             raise ValueError("No observer can be removed before choosing an algorithm.")
 740         #
 741         # Vérification du nom de variable et typage
 742         # -----------------------------------------
 743         if type( VariableName ) is str:
 744             VariableNames = [VariableName,]
 745         elif type( VariableName ) is list:
 746             VariableNames = map( str, VariableName )
 747         else:
 748             raise ValueError("The observer requires a name or a list of names of variables.")
 749         #
 750         # Association interne de l'observer à la variable
 751         # -----------------------------------------------
 752         for n in VariableNames:
 753             if not self.__algorithm.has_key( n ):
 754                 raise ValueError("An observer requires to be removed on a variable named %s which does not exist."%n)
 755         else:
 756             self.__algorithm.StoredVariables[ n ].removeDataObserver(
 757                 HookFunction   = HookFunction,
 758                 )
 759
 760     # -----------------------------------------------------------
 761     def setDebug(self, level=10):
 762         """
 763         Utiliser par exemple "import logging ; level = logging.DEBUG" avant cet
 764         appel pour changer le niveau de verbosité, avec :
 765         NOTSET=0 < DEBUG=10 < INFO=20 < WARNING=30 < ERROR=40 < CRITICAL=50
 766         """
 767         import logging
 768         log = logging.getLogger()
 769         log.setLevel( level )
 770
 771     def unsetDebug(self):
 772         """
 773         Remet le logger au niveau par défaut
 774         """
 775         import logging
 776         log = logging.getLogger()
 777         log.setLevel( logging.WARNING )
 778
 779     def prepare_to_pickle(self):
 780         self.__algorithmFile = None
 781         self.__diagnosticFile = None
 782         self.__H  = {}
 783
 784 # ==============================================================================
 785 if __name__ == "__main__":
 786     print '\n AUTODIAGNOSTIC \n'
 787
 788     ADD = AssimilationStudy("Ma premiere etude BLUE")
 789
 790     ADD.setBackground         (asVector     = [0, 1, 2])
 791     ADD.setBackgroundError    (asCovariance = "1 0 0;0 1 0;0 0 1")
 792     ADD.setObservation        (asVector     = [0.5, 1.5, 2.5])
 793     ADD.setObservationError   (asCovariance = "1 0 0;0 1 0;0 0 1")
 794     ADD.setObservationOperator(asMatrix     = "1 0 0;0 1 0;0 0 1")
 795
 796     ADD.setAlgorithm(choice="Blue")
 797
 798     ADD.analyze()
 799
 800     print "Nombre d'analyses  :", ADD.get("Analysis").stepnumber()
 801     print "Ebauche            :", [0, 1, 2]
 802     print "Observation        :", [0.5, 1.5, 2.5]
 803     print "Demi-somme         :", list((numpy.array([0, 1, 2])+numpy.array([0.5, 1.5, 2.5]))/2)
 804     print "  qui doit être identique à :"
 805     print "Analyse résultante :", ADD.get("Analysis").valueserie(0)
 806     print "Innovation         :", ADD.get("Innovation").valueserie(0)
 807     print
 808
 809     print "Algorithmes disponibles.......................:", ADD.get_available_algorithms()
 810     # print " Chemin des algorithmes.....................:", ADD.get_algorithms_main_path()
 811     print "Diagnostics types disponibles.................:", ADD.get_available_diagnostics()
 812     # print " Chemin des diagnostics.....................:", ADD.get_diagnostics_main_path()
 813     print "Variables disponibles.........................:", ADD.get_available_variables()
 814     print
 815
 816     ADD.setDiagnostic("RMS", "Ma RMS")
 817
 818     liste = ADD.get().keys()
 819     liste.sort()
 820     print "Variables et diagnostics nommés disponibles...:", liste
 821
 822     print
 823     print "Exemple de mise en place d'un observeur :"
 824     def obs(var=None,info=None):
 825         print "  ---> Mise en oeuvre de l'observer"
 826         print "       var  =",var.valueserie(-1)
 827         print "       info =",info
 828     ADD.setDataObserver( 'Analysis', HookFunction=obs, Scheduler = [2, 4], HookParameters = "Second observer")
 829     # Attention, il faut décaler le stockage de 1 pour suivre le pas interne
 830     # car le pas 0 correspond à l'analyse ci-dessus.
 831     for i in range(1,6):
 832         print
 833         print "Action sur la variable observée, étape :",i
 834         ADD.get('Analysis').store( [i, i, i] )
 835     print
 836
 837     print "Mise en debug et hors debug"
 838     print "Nombre d'analyses  :", ADD.get("Analysis").stepnumber()
 839     ADD.setDebug()
 840     ADD.analyze()
 841     ADD.unsetDebug()
 842     print "Nombre d'analyses  :", ADD.get("Analysis").stepnumber()
 843     ADD.analyze()
 844     print "Nombre d'analyses  :", ADD.get("Analysis").stepnumber()
 845     print