← Back to branch summary

~ubuntu-branches/ubuntu/vivid/mricron/vivid

« back to all changes in this revision

Viewing changes to npm_precl/dmath/fitexlin.pas

  • Committer: Package Import Robot
  • Author(s): Michael Hanke
  • Date: 2012-05-15 08:59:27 UTC
  • Revision ID: package-import@ubuntu.com-20120515085927-a2aal84xfmt3ej6c
http://bugs.debian.org/671365
* New upstream code (Closes: #671365).
* Update upstream source URL.
* Bumped Standards-version to 3.9.3, no changes necessary.

Show diffs side-by-side

added added

removed removed

Lines of Context:
npm_precl/dmath/fitexlin.pas
1
 
{ **********************************************************************
2
 
  *                         Unit FITEXLIN.PAS                          *
3
 
  *                            Version 1.1                             *
4
 
  *                    (c) J. Debord, August 2000                      *
5
 
  **********************************************************************
6
 
  This unit fits the "exponential + linear" model:
7
 
 
8
 
                        y = A.[1 - exp(-k.x)] + B.x
9
 
 
10
 
  ********************************************************************** }
11
 
 
12
 
unit FitExLin;
13
 
 
14
 
{$F+}
15
 
 
16
 
interface
17
 
 
18
 
uses
19
 
  FMath, Matrices, Stat, Regress;
20
 
 
21
 
function FuncName : String;
22
 
 
23
 
function FirstParam : Integer;
24
 
 
25
 
function LastParam : Integer;
26
 
 
27
 
function ParamName(I : Integer) : String;
28
 
 
29
 
function RegFunc(X : Float; B : PVector) : Float;
30
 
 
31
 
procedure DerivProc(X : Float; B, D : PVector);
32
 
 
33
 
function FitModel(X, Y : PVector; N : Integer; B : PVector) : Integer;
34
 
 
35
 
 
36
 
implementation
37
 
 
38
 
  function FuncName : String;
39
 
  { --------------------------------------------------------------------
40
 
    Returns the name of the regression function
41
 
    -------------------------------------------------------------------- }
42
 
  begin
43
 
    FuncName := 'y = A[1 - exp(-k.x)] + B.x';
44
 
  end;
45
 
 
46
 
  function FirstParam : Integer;
47
 
  { --------------------------------------------------------------------
48
 
    Returns the index of the first parameter to be fitted
49
 
    -------------------------------------------------------------------- }
50
 
  begin
51
 
    FirstParam := 0;
52
 
  end;
53
 
 
54
 
  function LastParam : Integer;
55
 
  { --------------------------------------------------------------------
56
 
    Returns the index of the last parameter to be fitted
57
 
    -------------------------------------------------------------------- }
58
 
  begin
59
 
    LastParam := 2;
60
 
  end;
61
 
 
62
 
  function ParamName(I : Integer) : String;
63
 
  { --------------------------------------------------------------------
64
 
    Returns the name of the I-th parameter
65
 
    -------------------------------------------------------------------- }
66
 
  begin
67
 
    case I of
68
 
      0 : ParamName := 'A';
69
 
      1 : ParamName := 'k';
70
 
      2 : ParamName := 'B';
71
 
    end;
72
 
  end;
73
 
 
74
 
  function RegFunc(X : Float; B : PVector) : Float;
75
 
  { --------------------------------------------------------------------
76
 
    Computes the regression function at point X
77
 
    B is the vector of parameters, such that :
78
 
 
79
 
    B^[0] = A     B^[1] = k     B^[2] = B
80
 
    -------------------------------------------------------------------- }
81
 
  begin
82
 
    RegFunc := B^[0] * (1.0 - Expo(- B^[1] * X)) + B^[2] * X;
83
 
  end;
84
 
 
85
 
  procedure DerivProc(X : Float; B, D : PVector);
86
 
  { --------------------------------------------------------------------
87
 
    Computes the derivatives of the regression function at point X
88
 
    with respect to the parameters B. The results are returned in D.
89
 
    D^[I] contains the derivative with respect to the I-th parameter.
90
 
    -------------------------------------------------------------------- }
91
 
  var
92
 
    E : Float;
93
 
  begin
94
 
    E := Expo(- B^[1] * X);  { exp(-k.x) }
95
 
    D^[0] := 1.0 - E;        { dy/dA = 1 - exp(-k.x) }
96
 
    D^[1] := B^[0] * X * E;  { dy/dk = A.x.exp(-k.x) }
97
 
    D^[2] := X;              { dy/dB = x }
98
 
  end;
99
 
 
100
 
  function FitModel(X, Y : PVector; N : Integer; B : PVector) : Integer;
101
 
  { --------------------------------------------------------------------
102
 
    Computes initial estimates of the regression parameters
103
 
    --------------------------------------------------------------------
104
 
    Input :  N    = number of points
105
 
             X, Y = point coordinates
106
 
    Output : B    = estimated regression parameters
107
 
    -------------------------------------------------------------------- }
108
 
  var
109
 
    K : Integer;
110
 
    D : Float;
111
 
  begin
112
 
    { B is the slope of the last (linear) part of the curve }
113
 
    K := Round(0.9 * N);
114
 
    if K = N then K := Pred(N);
115
 
    B^[2] := (Y^[N] - Y^[K]) / (X^[N] - X^[K]);
116
 
 
117
 
    { A is the intercept of the linear part }
118
 
    B^[0] := Y^[N] - B^[2] * X^[N];
119
 
 
120
 
    { Slope of the tangent at origin = B + k.A }
121
 
    K := Round(0.1 * N);
122
 
    if K = 1 then K := 2;
123
 
    D := (Y^[K] - Y^[1]) / (X^[K] - X^[1]);
124
 
    B^[1] := (D - B^[1]) / B^[0];
125
 
 
126
 
    FitModel := 0;
127
 
  end;
128
 
 
129
 
  end.