49 void nn_quit(
const char* format, ... );
66 nnai * nn = malloc(
sizeof (
nnai ) );
73 nn_quit(
"nnai_create(): n = %d\n", n );
77 nn->
x = malloc( (
size_t) n *
sizeof (
double ) );
78 memcpy( nn->
x, x, (
size_t) n * sizeof (
double ) );
79 nn->
y = malloc( (
size_t) n *
sizeof (
double ) );
80 memcpy( nn->
y, y, (
size_t) n * sizeof (
double ) );
83 for ( i = 0; i < n; ++i )
103 memcpy( w->
weights, weights, (
size_t) ( w->
nvertices ) * sizeof (
double ) );
119 for ( i = 0; i < nn->
n; ++i )
145 for ( i = 0; i < nn->
n; ++i )
155 if ( weight < nn->wmin )
178 #if defined ( NNAI_TEST )
180 #include <sys/time.h>
182 #define NPOINTSIN 10000
187 #define SQ( x ) ( ( x ) * ( x ) )
189 static double franke(
double x,
double y )
193 return 0.75 * exp( ( -SQ( x - 2.0 ) - SQ( y - 2.0 ) ) / 4.0 )
194 + 0.75 * exp( -SQ( x - 2.0 ) / 49.0 - ( y - 2.0 ) / 10.0 )
195 + 0.5 * exp( ( -SQ( x - 7.0 ) - SQ( y - 3.0 ) ) / 4.0 )
196 - 0.2 * exp( -SQ( x - 4.0 ) - SQ( y - 7.0 ) );
202 "Usage: nn_test [-v|-V] [-n <nin> <nxout>]\n"
204 " -a -- use non-Sibsonian interpolation rule\n"
205 " -n <nin> <nout>:\n"
206 " <nin> -- number of input points (default = 10000)\n"
207 " <nout> -- number of output points per side (default = 64)\n"
209 " -V -- very verbose\n"
223 double * xout = NULL;
224 double * yout = NULL;
225 double * zout = NULL;
227 struct timeval tv0, tv1, tv2;
234 switch ( argv[i][1] )
243 nn_quit(
"no number of data points found after -i\n" );
244 nin = atoi( argv[i] );
247 nn_quit(
"no number of ouput points per side found after -i\n" );
248 nx = atoi( argv[i] );
270 printf(
"\nTest of Natural Neighbours array interpolator:\n\n" );
271 printf(
" %d data points\n", nin );
272 printf(
" %d output points\n", nx * nx );
277 printf(
" generating data:\n" );
279 pin = malloc( nin *
sizeof (
point ) );
280 zin = malloc( nin *
sizeof (
double ) );
281 for ( i = 0; i < nin; ++i )
285 p->
x = (double) random() / RAND_MAX;
286 p->
y = (double) random() / RAND_MAX;
287 p->
z = franke( p->
x, p->
y );
290 printf(
" (%f, %f, %f)\n", p->
x, p->
y, p->
z );
296 printf(
" triangulating:\n" );
304 xout = malloc( nout *
sizeof (
double ) );
305 yout = malloc( nout *
sizeof (
double ) );
306 zout = malloc( nout *
sizeof (
double ) );
307 for ( i = 0; i < nout; ++i )
315 cpi = ( nx / 2 ) * ( nx + 1 );
317 gettimeofday( &tv0, &tz );
322 printf(
" creating interpolator:\n" );
327 gettimeofday( &tv1, &tz );
329 long dt = 1000000 * ( tv1.tv_sec - tv0.tv_sec ) + tv1.tv_usec - tv0.tv_usec;
331 printf(
" interpolator creation time = %ld us (%.2f us / point)\n", dt, (
double) dt / nout );
337 printf(
" interpolating:\n" );
341 for ( i = 0; i < nout; ++i )
342 printf(
" (%f, %f, %f)\n", xout[i], yout[i], zout[i] );
345 gettimeofday( &tv2, &tz );
347 long dt = 1000000.0 * ( tv2.tv_sec - tv1.tv_sec ) + tv2.tv_usec - tv1.tv_usec;
349 printf(
" interpolation time = %ld us (%.2f us / point)\n", dt, (
double) dt / nout );
353 printf(
" control point: (%f, %f, %f) (expected z = %f)\n", xout[cpi], yout[cpi], zout[cpi], franke( xout[cpi], yout[cpi] ) );
355 printf(
" interpolating one more time:\n" );
359 for ( i = 0; i < nout; ++i )
360 printf(
" (%f, %f, %f)\n", xout[i], yout[i], zout[i] );
363 gettimeofday( &tv0, &tz );
365 long dt = 1000000.0 * ( tv0.tv_sec - tv2.tv_sec ) + tv0.tv_usec - tv2.tv_usec;
367 printf(
" interpolation time = %ld us (%.2f us / point)\n", dt, (
double) dt / nout );
371 printf(
" control point: (%f, %f, %f) (expected z = %f)\n", xout[cpi], yout[cpi], zout[cpi], franke( xout[cpi], yout[cpi] ) );
373 printf(
" entering new data:\n" );
375 for ( i = 0; i < nin; ++i )
379 p->
z = p->
x * p->
x - p->
y * p->
y;
382 printf(
" (%f, %f, %f)\n", p->
x, p->
y, p->
z );
385 printf(
" interpolating:\n" );
389 for ( i = 0; i < nout; ++i )
390 printf(
" (%f, %f, %f)\n", xout[i], yout[i], zout[i] );
393 printf(
" control point: (%f, %f, %f) (expected z = %f)\n", xout[cpi], yout[cpi], zout[cpi], xout[cpi] * xout[cpi] - yout[cpi] * yout[cpi] );
395 printf(
" restoring data:\n" );
397 for ( i = 0; i < nin; ++i )
401 p->
z = franke( p->
x, p->
y );
404 printf(
" (%f, %f, %f)\n", p->
x, p->
y, p->
z );
407 printf(
" interpolating:\n" );
411 for ( i = 0; i < nout; ++i )
412 printf(
" (%f, %f, %f)\n", xout[i], yout[i], zout[i] );
415 printf(
" control point: (%f, %f, %f) (expected z = %f)\n", xout[cpi], yout[cpi], zout[cpi], franke( xout[cpi], yout[cpi] ) );
static PLCHAR_VECTOR usage
void nnai_interpolate(nnai *nn, double *zin, double *zout)
void nnpi_destroy(nnpi *nn)
nnai * nnai_build(delaunay *d, int n, double *x, double *y)
int * nnpi_get_vertices(nnpi *nn)
void delaunay_destroy(delaunay *d)
int nnpi_get_nvertices(nnpi *nn)
delaunay * delaunay_build(int np, point points[], int ns, int segments[], int nh, double holes[])
void nnai_destroy(nnai *nn)
void nnpi_set_point(nnpi *nn, point *p)
void nnpi_reset(nnpi *nn)
nnpi * nnpi_create(delaunay *d)
void nn_quit(const char *format,...)
void nnpi_normalize_weights(nnpi *nn)
void nnpi_calculate_weights(nnpi *nn)
void nnai_setwmin(nnai *nn, double wmin)
double * nnpi_get_weights(nnpi *nn)
void points_generate2(double xmin, double xmax, double ymin, double ymax, int nx, int ny, int *nout, point **pout)