doc/fr/scripts/simple_NonLinearLeastSquares.py

   1 # -*- coding: utf-8 -*-
   2 #
   3 from numpy import array, ravel
   4 def QuadFunction( coefficients ):
   5     """
   6     Simulation : y = a x^2 + b x + c
   7     """
   8     a, b, c = list(ravel(coefficients))
   9     x_points = (-5, 0, 1, 3, 10)
  10     y_points = []
  11     for x in x_points:
  12         y_points.append( a*x*x + b*x + c )
  13     return array(y_points)
  14 #
  15 Xb   = array([1., 1., 1.])
  16 Yobs = array([57, 2, 3, 17, 192])
  17 #
  18 print("Résolution itérative du problème de calibration")
  19 print("-----------------------.-----------------------")
  20 print("")
  21 from adao import adaoBuilder
  22 case = adaoBuilder.New('')
  23 case.setBackground( Vector = Xb, Stored=True )
  24 case.setObservation( Vector = Yobs, Stored=True )
  25 case.setObservationError( ScalarSparseMatrix = 1. )
  26 case.setObservationOperator( OneFunction = QuadFunction )
  27 case.setAlgorithmParameters(
  28     Algorithm='NonLinearLeastSquares',
  29     Parameters={
  30         'StoreSupplementaryCalculations': [
  31             'CurrentState',
  32             ],
  33         },
  34     )
  35 case.setObserver(
  36     Info="  État intermédiaire en itération courante :",
  37     Template='ValuePrinter',
  38     Variable='CurrentState',
  39     )
  40 case.execute()
  41 print("")
  42 #
  43 #-------------------------------------------------------------------------------
  44 #
  45 print("Calibration de %i coefficients pour une forme quadratique 1D sur %i mesures"%(
  46     len(case.get('Background')),
  47     len(case.get('Observation')),
  48     ))
  49 print("------------------------------------------------------------------------")
  50 print("")
  51 print("Vecteur observation.......................:", ravel(case.get('Observation')))
  52 print("État ébauche a priori.....................:", ravel(case.get('Background')))
  53 print("")
  54 print("Coefficients théoriques attendus..........:", ravel((2,-1,2)))
  55 print("")
  56 print("Coefficients résultants de la calibration.:", ravel(case.get('Analysis')[-1]))