PLplot
5.11.0
Main Page
Namespaces
Classes
Files
File List
File Members
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Macros
nn.h
Go to the documentation of this file.
1
//--------------------------------------------------------------------------
2
//
3
// File: nn.h
4
//
5
// Created: 04/08/2000
6
//
7
// Author: Pavel Sakov
8
// CSIRO Marine Research
9
//
10
// Purpose: Header file for nn library
11
//
12
// Description: None
13
//
14
// Revisions: None
15
//
16
//--------------------------------------------------------------------------
17
18
#if !defined ( _NN_H )
19
#define _NN_H
20
21
#include "
nndll.h
"
22
23
typedef
enum
{
SIBSON
,
NON_SIBSONIAN
}
NN_RULE
;
24
25
#if !defined ( _POINT_STRUCT )
26
#define _POINT_STRUCT
27
typedef
struct
28
{
29
double
x;
30
double
y;
31
double
z;
32
}
point
;
33
#endif
34
35
//* Smoothes the input point array by averaging the input x,y and z values
36
//** for each cell within virtual rectangular nx by ny grid. The corners of the
37
//** grid are created from min and max values of the input array. It also frees
38
//** the original array and returns results and new dimension via original
39
//** data and size pointers.
40
//*
41
//* @param pn Pointer to number of points (input/output)
42
//* @param ppoints Pointer to array of points (input/output) [*pn]
43
//* @param nx Number of x nodes in decimation
44
//* @param ny Number of y nodes in decimation
45
//
46
void
points_thin
(
int
* n,
point
** points,
int
nx,
int
ny );
47
48
//* Generates rectangular grid nx by ny using min and max x and y values from
49
//** the input point array. Allocates space for the output point array, be sure
50
//** to free it when necessary!
51
//*
52
//* @param n Number of points
53
//* @param points Array of points [n]
54
//* @param nx Number of x nodes
55
//* @param ny Number of y nodes
56
//* @param nout Pointer to number of output points
57
//* @param pout Ppointer to array of output points [*nout]
58
//
59
void
points_generate1
(
int
n,
point
points[],
int
nx,
int
ny,
double
zoom
,
int
* nout,
point
** pout );
60
61
//* Generates rectangular grid nx by ny using specified min and max x and y
62
//** values. Allocates space for the output point array, be sure to free it
63
//** when necessary!
64
//*
65
//* @param xmin Min x value
66
//* @param xmax Max x value
67
//* @param ymin Min y value
68
//* @param ymax Max y value
69
//* @param nx Number of x nodes
70
//* @param ny Number of y nodes
71
//* @param zoom Zoom coefficient
72
//* @param nout Pointer to number of output points
73
//* @param pout Pointer to array of output points [*nout]
74
//
75
void
points_generate2
(
double
xmin
,
double
xmax
,
double
ymin
,
double
ymax
,
int
nx,
int
ny,
int
* nout,
point
** pout );
76
77
//* Reads array of points from a columnar file.
78
//
79
// @param fname File name (can be "stdin" dor stndard input)
80
// @param dim Number of dimensions (must be 2 or 3)
81
// @param n Pointer to number of points (output)
82
// @param points Pointer to array of points [*n] (output)
83
//
84
void
points_read
(
char
* fname,
int
dim,
int
* n,
point
** points );
85
86
//* Scales Y coordinate so that the resulting set fits into square:
87
//** xmax - xmin = ymax - ymin
88
//*
89
//* @param n Number of points
90
//* @param points The points to scale
91
//* @return Y axis compression coefficient
92
//
93
double
points_scaletosquare
(
int
n,
point
* points );
94
95
//* Compresses Y domain by a given multiple.
96
//
97
// @param n Number of points
98
// @param points The points to scale
99
// @param Y axis compression coefficient as returned by points_scaletosquare()
100
//
101
void
points_scale
(
int
n,
point
* points,
double
k );
102
103
//* Structure to perform the Delaunay triangulation of a given array of points.
104
//
105
// Contains a deep copy of the input array of points.
106
// Contains triangles, circles and edges resulted from the triangulation.
107
// Contains neighbour triangles for each triangle.
108
// Contains point to triangle map.
109
//
110
struct
delaunay
;
111
typedef
struct
delaunay
delaunay
;
112
113
//* Builds Delaunay triangulation of the given array of points.
114
//
115
// @param np Number of points
116
// @param points Array of points [np] (input)
117
// @param ns Number of forced segments
118
// @param segments Array of (forced) segment endpoint indices [2*ns]
119
// @param nh Number of holes
120
// @param holes Array of hole (x,y) coordinates [2*nh]
121
// @return Delaunay triangulation with triangulation results
122
//
123
delaunay
*
delaunay_build
(
int
np,
point
points
[],
int
ns,
int
segments[],
int
nh,
double
holes[] );
124
125
//* Destroys Delaunay triangulation.
126
//
127
// @param d Structure to be destroyed
128
//
129
void
delaunay_destroy
(
delaunay
* d );
130
131
//* `lpi' -- "linear point interpolator" is a structure for
132
// conducting linear interpolation on a given data on a "point-to-point" basis.
133
// It interpolates linearly within each triangle resulted from the Delaunay
134
// triangluation of input data. `lpi' is much faster than all
135
// Natural Neighbours interpolators below.
136
//
137
struct
lpi
;
138
typedef
struct
lpi
lpi
;
139
140
//* Builds linear interpolator.
141
//
142
// @param d Delaunay triangulation
143
// @return Linear interpolator
144
//
145
lpi
*
lpi_build
(
delaunay
*
d
);
146
147
//* Destroys linear interpolator.
148
//
149
// @param l Structure to be destroyed
150
//
151
void
lpi_destroy
(
lpi
* l );
152
153
//* Finds linearly interpolated value in a point.
154
//
155
// @param l Linear point interpolator
156
// @param p Point to be interpolated (p->x, p->y -- input; p->z -- output)
157
//
158
void
lpi_interpolate_point
(
lpi
* l,
point
* p );
159
160
// Linearly interpolates data from one array of points for another array of
161
// points.
162
//
163
// @param nin Number of input points
164
// @param pin Array of input points [pin]
165
// @param nout Number of ouput points
166
// @param pout Array of output points [nout]
167
//
168
NNDLLIMPEXP
169
void
lpi_interpolate_points
(
int
nin,
point
pin[],
int
nout,
point
pout[] );
170
171
//* `nnpi' -- "Natural Neighbours point interpolator" is a
172
// structure for conducting Natural Neighbours interpolation on a given data on
173
// a "point-to-point" basis. Because it involves weight calculation for each
174
// next output point, it is not particularly suitable for consequitive
175
// interpolations on the same set of observation points -- use
176
// `nnhpi' or `nnai' in these cases.
177
//
178
struct
nnpi
;
179
typedef
struct
nnpi
nnpi
;
180
181
//* Creates Natural Neighbours point interpolator.
182
//
183
// @param d Delaunay triangulation
184
// @return Natural Neighbours interpolation
185
//
186
nnpi
*
nnpi_create
(
delaunay
*
d
);
187
188
//* Destroys Natural Neighbours point interpolation.
189
//
190
// @param nn Structure to be destroyed
191
//
192
void
nnpi_destroy
(
nnpi
* nn );
193
194
//* Finds Natural Neighbours-interpolated value in a point.
195
//
196
// @param nn NN point interpolator
197
// @param p Point to be interpolated (p->x, p->y -- input; p->z -- output)
198
//
199
void
nnpi_interpolate_point
(
nnpi
* nn,
point
*
p
);
200
201
//* Natural Neighbours-interpolates data in one array of points for another
202
//** array of points.
203
//*
204
//* @param nin Number of input points
205
//* @param pin Array of input points [pin]
206
//* @param wmin Minimal allowed weight
207
//* @param nout Number of output points
208
//* @param pout Array of output points [nout]
209
//
210
NNDLLIMPEXP
211
void
nnpi_interpolate_points
(
int
nin,
point
pin[],
double
wmin
,
int
nout,
point
pout[] );
212
213
//* Sets minimal allowed weight for Natural Neighbours interpolation.
214
// @param nn Natural Neighbours point interpolator
215
// @param wmin Minimal allowed weight
216
//
217
void
nnpi_setwmin
(
nnpi
* nn,
double
wmin
);
218
219
//* `nnhpi' is a structure for conducting consequitive
220
// Natural Neighbours interpolations on a given spatial data set in a random
221
// sequence of points from a set of finite size, taking advantage of repeated
222
// interpolations in the same point. It allows to modify Z
223
// coordinate of data between interpolations.
224
//
225
struct
nnhpi
;
226
typedef
struct
nnhpi
nnhpi
;
227
228
//* Creates Natural Neighbours hashing point interpolator.
229
//
230
// @param d Delaunay triangulation
231
// @param size Hash table size (should be of order of number of output points)
232
// @return Natural Neighbours interpolation
233
//
234
nnhpi
*
nnhpi_create
(
delaunay
* d,
int
size );
235
236
//* Destroys Natural Neighbours hashing point interpolation.
237
//
238
// @param nn Structure to be destroyed
239
//
240
void
nnhpi_destroy
(
nnhpi
* nn );
241
242
//* Finds Natural Neighbours-interpolated value in a point.
243
//
244
// @param nnhpi NN hashing point interpolator
245
// @param p Point to be interpolated (p->x, p->y -- input; p->z -- output)
246
//
247
void
nnhpi_interpolate
(
nnhpi
* nn,
point
* p );
248
249
//* Modifies interpolated data.
250
// Finds point* pd in the underlying Delaunay triangulation such that
251
// pd->x = p->x and pd->y = p->y, and copies p->z to pd->z. Exits with error
252
// if the point is not found.
253
//
254
// @param nn Natural Neighbours hashing point interpolator
255
// @param p New data
256
//
257
void
nnhpi_modify_data
(
nnhpi
* nn,
point
* p );
258
259
//* Sets minimal allowed weight for Natural Neighbours interpolation.
260
// @param nn Natural Neighbours point hashing interpolator
261
// @param wmin Minimal allowed weight
262
//
263
void
nnhpi_setwmin
(
nnhpi
* nn,
double
wmin );
264
265
// `nnai' is a tructure for conducting consequitive Natural
266
// Neighbours interpolations on a given spatial data set in a given array of
267
// points. It allows to modify Z coordinate of data between interpolations.
268
// `nnai' is the fastest of the three Natural Neighbours
269
// interpolators here.
270
//
271
struct
nnai
;
272
typedef
struct
nnai
nnai
;
273
274
//* Builds Natural Neighbours array interpolator. This includes calculation of
275
// weights used in nnai_interpolate().
276
//
277
// @param d Delaunay triangulation
278
// @return Natural Neighbours interpolation
279
//
280
nnai
*
nnai_build
(
delaunay
*
d
,
int
n
,
double
*
x
,
double
*
y
);
281
282
//* Destroys Natural Neighbours array interpolator.
283
//
284
// @param nn Structure to be destroyed
285
//
286
void
nnai_destroy
(
nnai
* nn );
287
288
//* Conducts NN interpolation in a fixed array of output points using
289
// data specified for a fixed array of input points. Uses pre-calculated
290
// weights.
291
//
292
// @param nn NN array interpolator
293
// @param zin input data [nn->d->npoints]
294
// @param zout output data [nn->n]. Must be pre-allocated!
295
//
296
void
nnai_interpolate
(
nnai
* nn,
double
* zin,
double
* zout );
297
298
//* Sets minimal allowed weight for Natural Neighbours interpolation.
299
// @param nn Natural Neighbours array interpolator
300
// @param wmin Minimal allowed weight
301
//
302
void
nnai_setwmin
(
nnai
* nn,
double
wmin
);
303
304
// Sets the verbosity level within nn package.
305
// 0 (default) - silent
306
// 1 - verbose
307
// 2 - very verbose
308
//
309
extern
int
nn_verbose
;
310
311
// Switches between weight calculation methods.
312
// SIBSON -- classic Sibson method
313
// NON_SIBSONIAN -- simpler and (I think) more robust method
314
//
315
extern
NNDLLIMPEXP_DATA
(
NN_RULE
)
nn_rule
;
316
317
// Contains version string for the nn package.
318
//
319
extern const
char
*
nn_version
;
320
321
// Limits verbose information to a particular vertex (used mainly for
322
// debugging purposes).
323
//
324
extern
int
nn_test_vertice
;
325
326
#endif // _NN_H
plplot_source
lib
nn
nn.h
Generated on Sun Apr 12 2015 03:08:35 for PLplot by
1.8.1.2