GRASS Programmer's Manual  6.5.svn(2014)-r66266
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
func2d.c
Go to the documentation of this file.
1 
2 /*-
3  *
4  * Original program and various modifications:
5  * Lubos Mitas
6  *
7  * GRASS4.1 version of the program and GRASS4.2 modifications:
8  * H. Mitasova,
9  * I. Kosinovsky, D. Gerdes
10  * D. McCauley
11  *
12  * Copyright 1993, 1995:
13  * L. Mitas ,
14  * H. Mitasova ,
15  * I. Kosinovsky,
16  * D.Gerdes
17  * D. McCauley
18  *
19  * modified by McCauley in August 1995
20  * modified by Mitasova in August 1995, Nov. 1996
21  *
22  */
23 
24 #include <stdio.h>
25 #include <math.h>
26 #include <grass/gis.h>
27 #include <grass/interpf.h>
28 
29 double IL_crst(double r, double fi)
30 /*
31  * Radial basis function - completely regularized spline with
32  * tension (d=2)
33  */
34 {
35  double rfsta2 = fi * fi * r / 4.;
36 
37  static double c[4] = { 8.5733287401, 18.0590169730, 8.6347608925,
38  0.2677737343
39  };
40  static double b[4] = { 9.5733223454, 25.6329561486, 21.0996530827,
41  3.9584969228
42  };
43  double ce = 0.57721566;
44 
45  static double u[10] = { 1.e+00, -.25e+00,
46  .055555555555556e+00, -.010416666666667e+00, /*fixed bug 415.. repl. by 416.. */
47  .166666666666667e-02, -2.31481481481482e-04,
48  2.83446712018141e-05, -3.10019841269841e-06,
49  3.06192435822065e-07, -2.75573192239859e-08
50  };
51  double x = rfsta2;
52  double res;
53 
54  double e1, ea, eb;
55 
56 
57  if (x < 1.e+00) {
58  res = x * (u[0] + x * (u[1] + x * (u[2] + x * (u[3] + x * (u[4] + x *
59  (u[5] +
60  x *
61  (u[6] +
62  x *
63  (u[7] +
64  x *
65  (u[8] +
66  x *
67  u
68  [9])))))))));
69  return (res);
70  }
71 
72  if (x > 25.e+00)
73  e1 = 0.00;
74  else {
75  ea = c[3] + x * (c[2] + x * (c[1] + x * (c[0] + x)));
76  eb = b[3] + x * (b[2] + x * (b[1] + x * (b[0] + x)));
77  e1 = (ea / eb) / (x * exp(x));
78  }
79  res = e1 + ce + log(x);
80  return (res);
81 }
82 
83 
84 
85 int IL_crstg(double r, double fi, double *gd1, /* G1(r) */
86  double *gd2 /* G2(r) */
87  )
88 
89 
90 /*
91  * Function for calculating derivatives (d=2)
92  */
93 {
94  double r2 = r;
95  double rfsta2 = fi * fi * r / 4.;
96  double x, exm, oneme, hold;
97  double fsta2 = fi * fi / 2.;
98 
99  x = rfsta2;
100  if (x < 0.001) {
101  *gd1 = 1. - x / 2. + x * x / 6. - x * x * x / 24.;
102  *gd2 = fsta2 * (-.5 + x / 3. - x * x / 8. + x * x * x / 30.);
103  }
104  else {
105  if (x < 35.e+00) {
106  exm = exp(-x);
107  oneme = 1. - exm;
108  *gd1 = oneme / x;
109  hold = x * exm - oneme;
110  *gd2 = (hold + hold) / (r2 * x);
111  }
112  else {
113  *gd1 = 1. / x;
114  *gd2 = -2. / (x * r2);
115  }
116  }
117  return 1;
118 }
int IL_crstg(double r, double fi, double *gd1, double *gd2)
Definition: func2d.c:85
float b
Definition: named_colr.c:8
float r
Definition: named_colr.c:8
double IL_crst(double r, double fi)
Definition: func2d.c:29
log
Definition: wxnviz.py:54