src/daComposant/daCore/AssimilationStudy.py

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