GRASS 8 Programmer's Manual 8.6.0dev(2026)-f6f2c534ea
Loading...
Searching...
No Matches
func2d.c
Go to the documentation of this file.
1/*!
2 * \file func2d.c
3 *
4 * \author
5 * Lubos Mitas (original program and various modifications)
6 *
7 * \author
8 * H. Mitasova,
9 * I. Kosinovsky, D. Gerdes,
10 * D. McCauley
11 * (GRASS4.1 version of the program and GRASS4.2 modifications)
12 *
13 * \author
14 * L. Mitas ,
15 * H. Mitasova ,
16 * I. Kosinovsky,
17 * D.Gerdes
18 * D. McCauley (1993, 1995)
19 *
20 * \author modified by McCauley in August 1995
21 * \author modified by Mitasova in August 1995, Nov. 1996
22 *
23 * \copyright
24 * (C) 1993-1999 by Lubos Mitas and the GRASS Development Team
25 *
26 * \copyright
27 * This program is free software under the
28 * GNU General Public License (>=v2).
29 * Read the file COPYING that comes with GRASS
30 * for details.
31 */
32
33#include <stdio.h>
34#include <math.h>
35#include <grass/gis.h>
36#include <grass/interpf.h>
37
38/**
39 * @brief Radial basis function
40 *
41 * Radial basis function - completely regularized spline with tension (d=2)
42 *
43 * parameter description from DESCRIPTION.INTERP
44 *
45 * @param r distance squared
46 * @param fi tension
47 * @return
48 */
49double IL_crst(double r, double fi)
50{
51 double rfsta2 = fi * fi * r / 4.;
52
53 static double c[4] = {8.5733287401, 18.0590169730, 8.6347608925,
54 0.2677737343};
55 static double b[4] = {9.5733223454, 25.6329561486, 21.0996530827,
56 3.9584969228};
57 double ce = 0.57721566;
58
59 static double u[10] = {
60 1.e+00,
61 -.25e+00,
62 .055555555555556e+00,
63 -.010416666666667e+00, /*fixed bug 415.. repl. by 416.. */
64 .166666666666667e-02,
65 -2.31481481481482e-04,
66 2.83446712018141e-05,
67 -3.10019841269841e-06,
68 3.06192435822065e-07,
69 -2.75573192239859e-08};
70 double x = rfsta2;
71 double res;
72
73 double e1, ea, eb;
74
75 if (x < 1.e+00) {
76 res = x *
77 (u[0] +
78 x * (u[1] +
79 x * (u[2] +
80 x * (u[3] +
81 x * (u[4] +
82 x * (u[5] +
83 x * (u[6] +
84 x * (u[7] +
85 x * (u[8] + x * u[9])))))))));
86 return (res);
87 }
88
89 if (x > 25.e+00)
90 e1 = 0.00;
91 else {
92 ea = c[3] + x * (c[2] + x * (c[1] + x * (c[0] + x)));
93 eb = b[3] + x * (b[2] + x * (b[1] + x * (b[0] + x)));
94 e1 = (ea / eb) / (x * exp(x));
95 }
96 res = e1 + ce + log(x);
97 return (res);
98}
99
100/**
101 * @brief Function for calculating derivatives (d=2)
102 *
103 * Derivatives of radial basis function - regularized spline with tension(d=2)
104 *
105 * @param r distance squared
106 * @param fi tension
107 * @param gd1 G1(r)
108 * @param gd2 G2(r)
109 * @return
110 */
111int IL_crstg(double r, double fi, double *gd1, double *gd2)
112{
113 double r2 = r;
114 double rfsta2 = fi * fi * r / 4.;
115 double x, exm, oneme, hold;
116 double fsta2 = fi * fi / 2.;
117
118 x = rfsta2;
119 if (x < 0.001) {
120 *gd1 = 1. - x / 2. + x * x / 6. - x * x * x / 24.;
121 *gd2 = fsta2 * (-.5 + x / 3. - x * x / 8. + x * x * x / 30.);
122 }
123 else {
124 if (x < 35.e+00) {
125 exm = exp(-x);
126 oneme = 1. - exm;
127 *gd1 = oneme / x;
128 hold = x * exm - oneme;
129 *gd2 = (hold + hold) / (r2 * x);
130 }
131 else {
132 *gd1 = 1. / x;
133 *gd2 = -2. / (x * r2);
134 }
135 }
136 return 1;
137}
double IL_crst(double r, double fi)
Radial basis function.
Definition func2d.c:49
int IL_crstg(double r, double fi, double *gd1, double *gd2)
Function for calculating derivatives (d=2)
Definition func2d.c:111
double b
Definition r_raster.c:39
double r
Definition r_raster.c:39
#define x