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