1 | /*************************************************************************
|
---|
2 | Copyright (c) Sergey Bochkanov (ALGLIB project).
|
---|
3 |
|
---|
4 | >>> SOURCE LICENSE >>>
|
---|
5 | This program is free software; you can redistribute it and/or modify
|
---|
6 | it under the terms of the GNU General Public License as published by
|
---|
7 | the Free Software Foundation (www.fsf.org); either version 2 of the
|
---|
8 | License, or (at your option) any later version.
|
---|
9 |
|
---|
10 | This program is distributed in the hope that it will be useful,
|
---|
11 | but WITHOUT ANY WARRANTY; without even the implied warranty of
|
---|
12 | MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
---|
13 | GNU General Public License for more details.
|
---|
14 |
|
---|
15 | A copy of the GNU General Public License is available at
|
---|
16 | http://www.fsf.org/licensing/licenses
|
---|
17 | >>> END OF LICENSE >>>
|
---|
18 | *************************************************************************/
|
---|
19 | #pragma warning disable 162
|
---|
20 | #pragma warning disable 219
|
---|
21 | using System;
|
---|
22 |
|
---|
23 | public partial class alglib
|
---|
24 | {
|
---|
25 |
|
---|
26 |
|
---|
27 | /*************************************************************************
|
---|
28 | Portable high quality random number generator state.
|
---|
29 | Initialized with HQRNDRandomize() or HQRNDSeed().
|
---|
30 |
|
---|
31 | Fields:
|
---|
32 | S1, S2 - seed values
|
---|
33 | V - precomputed value
|
---|
34 | MagicV - 'magic' value used to determine whether State structure
|
---|
35 | was correctly initialized.
|
---|
36 | *************************************************************************/
|
---|
37 | public class hqrndstate
|
---|
38 | {
|
---|
39 | //
|
---|
40 | // Public declarations
|
---|
41 | //
|
---|
42 |
|
---|
43 | public hqrndstate()
|
---|
44 | {
|
---|
45 | _innerobj = new hqrnd.hqrndstate();
|
---|
46 | }
|
---|
47 |
|
---|
48 | //
|
---|
49 | // Although some of declarations below are public, you should not use them
|
---|
50 | // They are intended for internal use only
|
---|
51 | //
|
---|
52 | private hqrnd.hqrndstate _innerobj;
|
---|
53 | public hqrnd.hqrndstate innerobj { get { return _innerobj; } }
|
---|
54 | public hqrndstate(hqrnd.hqrndstate obj)
|
---|
55 | {
|
---|
56 | _innerobj = obj;
|
---|
57 | }
|
---|
58 | }
|
---|
59 |
|
---|
60 | /*************************************************************************
|
---|
61 | HQRNDState initialization with random values which come from standard
|
---|
62 | RNG.
|
---|
63 |
|
---|
64 | -- ALGLIB --
|
---|
65 | Copyright 02.12.2009 by Bochkanov Sergey
|
---|
66 | *************************************************************************/
|
---|
67 | public static void hqrndrandomize(out hqrndstate state)
|
---|
68 | {
|
---|
69 | state = new hqrndstate();
|
---|
70 | hqrnd.hqrndrandomize(state.innerobj);
|
---|
71 | return;
|
---|
72 | }
|
---|
73 |
|
---|
74 | /*************************************************************************
|
---|
75 | HQRNDState initialization with seed values
|
---|
76 |
|
---|
77 | -- ALGLIB --
|
---|
78 | Copyright 02.12.2009 by Bochkanov Sergey
|
---|
79 | *************************************************************************/
|
---|
80 | public static void hqrndseed(int s1, int s2, out hqrndstate state)
|
---|
81 | {
|
---|
82 | state = new hqrndstate();
|
---|
83 | hqrnd.hqrndseed(s1, s2, state.innerobj);
|
---|
84 | return;
|
---|
85 | }
|
---|
86 |
|
---|
87 | /*************************************************************************
|
---|
88 | This function generates random real number in (0,1),
|
---|
89 | not including interval boundaries
|
---|
90 |
|
---|
91 | State structure must be initialized with HQRNDRandomize() or HQRNDSeed().
|
---|
92 |
|
---|
93 | -- ALGLIB --
|
---|
94 | Copyright 02.12.2009 by Bochkanov Sergey
|
---|
95 | *************************************************************************/
|
---|
96 | public static double hqrnduniformr(hqrndstate state)
|
---|
97 | {
|
---|
98 |
|
---|
99 | double result = hqrnd.hqrnduniformr(state.innerobj);
|
---|
100 | return result;
|
---|
101 | }
|
---|
102 |
|
---|
103 | /*************************************************************************
|
---|
104 | This function generates random integer number in [0, N)
|
---|
105 |
|
---|
106 | 1. N must be less than HQRNDMax-1.
|
---|
107 | 2. State structure must be initialized with HQRNDRandomize() or HQRNDSeed()
|
---|
108 |
|
---|
109 | -- ALGLIB --
|
---|
110 | Copyright 02.12.2009 by Bochkanov Sergey
|
---|
111 | *************************************************************************/
|
---|
112 | public static int hqrnduniformi(hqrndstate state, int n)
|
---|
113 | {
|
---|
114 |
|
---|
115 | int result = hqrnd.hqrnduniformi(state.innerobj, n);
|
---|
116 | return result;
|
---|
117 | }
|
---|
118 |
|
---|
119 | /*************************************************************************
|
---|
120 | Random number generator: normal numbers
|
---|
121 |
|
---|
122 | This function generates one random number from normal distribution.
|
---|
123 | Its performance is equal to that of HQRNDNormal2()
|
---|
124 |
|
---|
125 | State structure must be initialized with HQRNDRandomize() or HQRNDSeed().
|
---|
126 |
|
---|
127 | -- ALGLIB --
|
---|
128 | Copyright 02.12.2009 by Bochkanov Sergey
|
---|
129 | *************************************************************************/
|
---|
130 | public static double hqrndnormal(hqrndstate state)
|
---|
131 | {
|
---|
132 |
|
---|
133 | double result = hqrnd.hqrndnormal(state.innerobj);
|
---|
134 | return result;
|
---|
135 | }
|
---|
136 |
|
---|
137 | /*************************************************************************
|
---|
138 | Random number generator: random X and Y such that X^2+Y^2=1
|
---|
139 |
|
---|
140 | State structure must be initialized with HQRNDRandomize() or HQRNDSeed().
|
---|
141 |
|
---|
142 | -- ALGLIB --
|
---|
143 | Copyright 02.12.2009 by Bochkanov Sergey
|
---|
144 | *************************************************************************/
|
---|
145 | public static void hqrndunit2(hqrndstate state, out double x, out double y)
|
---|
146 | {
|
---|
147 | x = 0;
|
---|
148 | y = 0;
|
---|
149 | hqrnd.hqrndunit2(state.innerobj, ref x, ref y);
|
---|
150 | return;
|
---|
151 | }
|
---|
152 |
|
---|
153 | /*************************************************************************
|
---|
154 | Random number generator: normal numbers
|
---|
155 |
|
---|
156 | This function generates two independent random numbers from normal
|
---|
157 | distribution. Its performance is equal to that of HQRNDNormal()
|
---|
158 |
|
---|
159 | State structure must be initialized with HQRNDRandomize() or HQRNDSeed().
|
---|
160 |
|
---|
161 | -- ALGLIB --
|
---|
162 | Copyright 02.12.2009 by Bochkanov Sergey
|
---|
163 | *************************************************************************/
|
---|
164 | public static void hqrndnormal2(hqrndstate state, out double x1, out double x2)
|
---|
165 | {
|
---|
166 | x1 = 0;
|
---|
167 | x2 = 0;
|
---|
168 | hqrnd.hqrndnormal2(state.innerobj, ref x1, ref x2);
|
---|
169 | return;
|
---|
170 | }
|
---|
171 |
|
---|
172 | /*************************************************************************
|
---|
173 | Random number generator: exponential distribution
|
---|
174 |
|
---|
175 | State structure must be initialized with HQRNDRandomize() or HQRNDSeed().
|
---|
176 |
|
---|
177 | -- ALGLIB --
|
---|
178 | Copyright 11.08.2007 by Bochkanov Sergey
|
---|
179 | *************************************************************************/
|
---|
180 | public static double hqrndexponential(hqrndstate state, double lambdav)
|
---|
181 | {
|
---|
182 |
|
---|
183 | double result = hqrnd.hqrndexponential(state.innerobj, lambdav);
|
---|
184 | return result;
|
---|
185 | }
|
---|
186 |
|
---|
187 | }
|
---|
188 | public partial class alglib
|
---|
189 | {
|
---|
190 |
|
---|
191 |
|
---|
192 | /*************************************************************************
|
---|
193 |
|
---|
194 | *************************************************************************/
|
---|
195 | public class kdtree
|
---|
196 | {
|
---|
197 | //
|
---|
198 | // Public declarations
|
---|
199 | //
|
---|
200 |
|
---|
201 | public kdtree()
|
---|
202 | {
|
---|
203 | _innerobj = new nearestneighbor.kdtree();
|
---|
204 | }
|
---|
205 |
|
---|
206 | //
|
---|
207 | // Although some of declarations below are public, you should not use them
|
---|
208 | // They are intended for internal use only
|
---|
209 | //
|
---|
210 | private nearestneighbor.kdtree _innerobj;
|
---|
211 | public nearestneighbor.kdtree innerobj { get { return _innerobj; } }
|
---|
212 | public kdtree(nearestneighbor.kdtree obj)
|
---|
213 | {
|
---|
214 | _innerobj = obj;
|
---|
215 | }
|
---|
216 | }
|
---|
217 |
|
---|
218 |
|
---|
219 | /*************************************************************************
|
---|
220 | This function serializes data structure to string.
|
---|
221 |
|
---|
222 | Important properties of s_out:
|
---|
223 | * it contains alphanumeric characters, dots, underscores, minus signs
|
---|
224 | * these symbols are grouped into words, which are separated by spaces
|
---|
225 | and Windows-style (CR+LF) newlines
|
---|
226 | * although serializer uses spaces and CR+LF as separators, you can
|
---|
227 | replace any separator character by arbitrary combination of spaces,
|
---|
228 | tabs, Windows or Unix newlines. It allows flexible reformatting of
|
---|
229 | the string in case you want to include it into text or XML file.
|
---|
230 | But you should not insert separators into the middle of the "words"
|
---|
231 | nor you should change case of letters.
|
---|
232 | * s_out can be freely moved between 32-bit and 64-bit systems, little
|
---|
233 | and big endian machines, and so on. You can serialize structure on
|
---|
234 | 32-bit machine and unserialize it on 64-bit one (or vice versa), or
|
---|
235 | serialize it on SPARC and unserialize on x86. You can also
|
---|
236 | serialize it in C# version of ALGLIB and unserialize in C++ one,
|
---|
237 | and vice versa.
|
---|
238 | *************************************************************************/
|
---|
239 | public static void kdtreeserialize(kdtree obj, out string s_out)
|
---|
240 | {
|
---|
241 | alglib.serializer s = new alglib.serializer();
|
---|
242 | s.alloc_start();
|
---|
243 | nearestneighbor.kdtreealloc(s, obj.innerobj);
|
---|
244 | s.sstart_str();
|
---|
245 | nearestneighbor.kdtreeserialize(s, obj.innerobj);
|
---|
246 | s.stop();
|
---|
247 | s_out = s.get_string();
|
---|
248 | }
|
---|
249 |
|
---|
250 |
|
---|
251 | /*************************************************************************
|
---|
252 | This function unserializes data structure from string.
|
---|
253 | *************************************************************************/
|
---|
254 | public static void kdtreeunserialize(string s_in, out kdtree obj)
|
---|
255 | {
|
---|
256 | alglib.serializer s = new alglib.serializer();
|
---|
257 | obj = new kdtree();
|
---|
258 | s.ustart_str(s_in);
|
---|
259 | nearestneighbor.kdtreeunserialize(s, obj.innerobj);
|
---|
260 | s.stop();
|
---|
261 | }
|
---|
262 |
|
---|
263 | /*************************************************************************
|
---|
264 | KD-tree creation
|
---|
265 |
|
---|
266 | This subroutine creates KD-tree from set of X-values and optional Y-values
|
---|
267 |
|
---|
268 | INPUT PARAMETERS
|
---|
269 | XY - dataset, array[0..N-1,0..NX+NY-1].
|
---|
270 | one row corresponds to one point.
|
---|
271 | first NX columns contain X-values, next NY (NY may be zero)
|
---|
272 | columns may contain associated Y-values
|
---|
273 | N - number of points, N>=1
|
---|
274 | NX - space dimension, NX>=1.
|
---|
275 | NY - number of optional Y-values, NY>=0.
|
---|
276 | NormType- norm type:
|
---|
277 | * 0 denotes infinity-norm
|
---|
278 | * 1 denotes 1-norm
|
---|
279 | * 2 denotes 2-norm (Euclidean norm)
|
---|
280 |
|
---|
281 | OUTPUT PARAMETERS
|
---|
282 | KDT - KD-tree
|
---|
283 |
|
---|
284 |
|
---|
285 | NOTES
|
---|
286 |
|
---|
287 | 1. KD-tree creation have O(N*logN) complexity and O(N*(2*NX+NY)) memory
|
---|
288 | requirements.
|
---|
289 | 2. Although KD-trees may be used with any combination of N and NX, they
|
---|
290 | are more efficient than brute-force search only when N >> 4^NX. So they
|
---|
291 | are most useful in low-dimensional tasks (NX=2, NX=3). NX=1 is another
|
---|
292 | inefficient case, because simple binary search (without additional
|
---|
293 | structures) is much more efficient in such tasks than KD-trees.
|
---|
294 |
|
---|
295 | -- ALGLIB --
|
---|
296 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
297 | *************************************************************************/
|
---|
298 | public static void kdtreebuild(double[,] xy, int n, int nx, int ny, int normtype, out kdtree kdt)
|
---|
299 | {
|
---|
300 | kdt = new kdtree();
|
---|
301 | nearestneighbor.kdtreebuild(xy, n, nx, ny, normtype, kdt.innerobj);
|
---|
302 | return;
|
---|
303 | }
|
---|
304 | public static void kdtreebuild(double[,] xy, int nx, int ny, int normtype, out kdtree kdt)
|
---|
305 | {
|
---|
306 | int n;
|
---|
307 |
|
---|
308 | kdt = new kdtree();
|
---|
309 | n = ap.rows(xy);
|
---|
310 | nearestneighbor.kdtreebuild(xy, n, nx, ny, normtype, kdt.innerobj);
|
---|
311 |
|
---|
312 | return;
|
---|
313 | }
|
---|
314 |
|
---|
315 | /*************************************************************************
|
---|
316 | KD-tree creation
|
---|
317 |
|
---|
318 | This subroutine creates KD-tree from set of X-values, integer tags and
|
---|
319 | optional Y-values
|
---|
320 |
|
---|
321 | INPUT PARAMETERS
|
---|
322 | XY - dataset, array[0..N-1,0..NX+NY-1].
|
---|
323 | one row corresponds to one point.
|
---|
324 | first NX columns contain X-values, next NY (NY may be zero)
|
---|
325 | columns may contain associated Y-values
|
---|
326 | Tags - tags, array[0..N-1], contains integer tags associated
|
---|
327 | with points.
|
---|
328 | N - number of points, N>=1
|
---|
329 | NX - space dimension, NX>=1.
|
---|
330 | NY - number of optional Y-values, NY>=0.
|
---|
331 | NormType- norm type:
|
---|
332 | * 0 denotes infinity-norm
|
---|
333 | * 1 denotes 1-norm
|
---|
334 | * 2 denotes 2-norm (Euclidean norm)
|
---|
335 |
|
---|
336 | OUTPUT PARAMETERS
|
---|
337 | KDT - KD-tree
|
---|
338 |
|
---|
339 | NOTES
|
---|
340 |
|
---|
341 | 1. KD-tree creation have O(N*logN) complexity and O(N*(2*NX+NY)) memory
|
---|
342 | requirements.
|
---|
343 | 2. Although KD-trees may be used with any combination of N and NX, they
|
---|
344 | are more efficient than brute-force search only when N >> 4^NX. So they
|
---|
345 | are most useful in low-dimensional tasks (NX=2, NX=3). NX=1 is another
|
---|
346 | inefficient case, because simple binary search (without additional
|
---|
347 | structures) is much more efficient in such tasks than KD-trees.
|
---|
348 |
|
---|
349 | -- ALGLIB --
|
---|
350 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
351 | *************************************************************************/
|
---|
352 | public static void kdtreebuildtagged(double[,] xy, int[] tags, int n, int nx, int ny, int normtype, out kdtree kdt)
|
---|
353 | {
|
---|
354 | kdt = new kdtree();
|
---|
355 | nearestneighbor.kdtreebuildtagged(xy, tags, n, nx, ny, normtype, kdt.innerobj);
|
---|
356 | return;
|
---|
357 | }
|
---|
358 | public static void kdtreebuildtagged(double[,] xy, int[] tags, int nx, int ny, int normtype, out kdtree kdt)
|
---|
359 | {
|
---|
360 | int n;
|
---|
361 | if( (ap.rows(xy)!=ap.len(tags)))
|
---|
362 | throw new alglibexception("Error while calling 'kdtreebuildtagged': looks like one of arguments has wrong size");
|
---|
363 | kdt = new kdtree();
|
---|
364 | n = ap.rows(xy);
|
---|
365 | nearestneighbor.kdtreebuildtagged(xy, tags, n, nx, ny, normtype, kdt.innerobj);
|
---|
366 |
|
---|
367 | return;
|
---|
368 | }
|
---|
369 |
|
---|
370 | /*************************************************************************
|
---|
371 | K-NN query: K nearest neighbors
|
---|
372 |
|
---|
373 | INPUT PARAMETERS
|
---|
374 | KDT - KD-tree
|
---|
375 | X - point, array[0..NX-1].
|
---|
376 | K - number of neighbors to return, K>=1
|
---|
377 | SelfMatch - whether self-matches are allowed:
|
---|
378 | * if True, nearest neighbor may be the point itself
|
---|
379 | (if it exists in original dataset)
|
---|
380 | * if False, then only points with non-zero distance
|
---|
381 | are returned
|
---|
382 | * if not given, considered True
|
---|
383 |
|
---|
384 | RESULT
|
---|
385 | number of actual neighbors found (either K or N, if K>N).
|
---|
386 |
|
---|
387 | This subroutine performs query and stores its result in the internal
|
---|
388 | structures of the KD-tree. You can use following subroutines to obtain
|
---|
389 | these results:
|
---|
390 | * KDTreeQueryResultsX() to get X-values
|
---|
391 | * KDTreeQueryResultsXY() to get X- and Y-values
|
---|
392 | * KDTreeQueryResultsTags() to get tag values
|
---|
393 | * KDTreeQueryResultsDistances() to get distances
|
---|
394 |
|
---|
395 | -- ALGLIB --
|
---|
396 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
397 | *************************************************************************/
|
---|
398 | public static int kdtreequeryknn(kdtree kdt, double[] x, int k, bool selfmatch)
|
---|
399 | {
|
---|
400 |
|
---|
401 | int result = nearestneighbor.kdtreequeryknn(kdt.innerobj, x, k, selfmatch);
|
---|
402 | return result;
|
---|
403 | }
|
---|
404 | public static int kdtreequeryknn(kdtree kdt, double[] x, int k)
|
---|
405 | {
|
---|
406 | bool selfmatch;
|
---|
407 |
|
---|
408 |
|
---|
409 | selfmatch = true;
|
---|
410 | int result = nearestneighbor.kdtreequeryknn(kdt.innerobj, x, k, selfmatch);
|
---|
411 |
|
---|
412 | return result;
|
---|
413 | }
|
---|
414 |
|
---|
415 | /*************************************************************************
|
---|
416 | R-NN query: all points within R-sphere centered at X
|
---|
417 |
|
---|
418 | INPUT PARAMETERS
|
---|
419 | KDT - KD-tree
|
---|
420 | X - point, array[0..NX-1].
|
---|
421 | R - radius of sphere (in corresponding norm), R>0
|
---|
422 | SelfMatch - whether self-matches are allowed:
|
---|
423 | * if True, nearest neighbor may be the point itself
|
---|
424 | (if it exists in original dataset)
|
---|
425 | * if False, then only points with non-zero distance
|
---|
426 | are returned
|
---|
427 | * if not given, considered True
|
---|
428 |
|
---|
429 | RESULT
|
---|
430 | number of neighbors found, >=0
|
---|
431 |
|
---|
432 | This subroutine performs query and stores its result in the internal
|
---|
433 | structures of the KD-tree. You can use following subroutines to obtain
|
---|
434 | actual results:
|
---|
435 | * KDTreeQueryResultsX() to get X-values
|
---|
436 | * KDTreeQueryResultsXY() to get X- and Y-values
|
---|
437 | * KDTreeQueryResultsTags() to get tag values
|
---|
438 | * KDTreeQueryResultsDistances() to get distances
|
---|
439 |
|
---|
440 | -- ALGLIB --
|
---|
441 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
442 | *************************************************************************/
|
---|
443 | public static int kdtreequeryrnn(kdtree kdt, double[] x, double r, bool selfmatch)
|
---|
444 | {
|
---|
445 |
|
---|
446 | int result = nearestneighbor.kdtreequeryrnn(kdt.innerobj, x, r, selfmatch);
|
---|
447 | return result;
|
---|
448 | }
|
---|
449 | public static int kdtreequeryrnn(kdtree kdt, double[] x, double r)
|
---|
450 | {
|
---|
451 | bool selfmatch;
|
---|
452 |
|
---|
453 |
|
---|
454 | selfmatch = true;
|
---|
455 | int result = nearestneighbor.kdtreequeryrnn(kdt.innerobj, x, r, selfmatch);
|
---|
456 |
|
---|
457 | return result;
|
---|
458 | }
|
---|
459 |
|
---|
460 | /*************************************************************************
|
---|
461 | K-NN query: approximate K nearest neighbors
|
---|
462 |
|
---|
463 | INPUT PARAMETERS
|
---|
464 | KDT - KD-tree
|
---|
465 | X - point, array[0..NX-1].
|
---|
466 | K - number of neighbors to return, K>=1
|
---|
467 | SelfMatch - whether self-matches are allowed:
|
---|
468 | * if True, nearest neighbor may be the point itself
|
---|
469 | (if it exists in original dataset)
|
---|
470 | * if False, then only points with non-zero distance
|
---|
471 | are returned
|
---|
472 | * if not given, considered True
|
---|
473 | Eps - approximation factor, Eps>=0. eps-approximate nearest
|
---|
474 | neighbor is a neighbor whose distance from X is at
|
---|
475 | most (1+eps) times distance of true nearest neighbor.
|
---|
476 |
|
---|
477 | RESULT
|
---|
478 | number of actual neighbors found (either K or N, if K>N).
|
---|
479 |
|
---|
480 | NOTES
|
---|
481 | significant performance gain may be achieved only when Eps is is on
|
---|
482 | the order of magnitude of 1 or larger.
|
---|
483 |
|
---|
484 | This subroutine performs query and stores its result in the internal
|
---|
485 | structures of the KD-tree. You can use following subroutines to obtain
|
---|
486 | these results:
|
---|
487 | * KDTreeQueryResultsX() to get X-values
|
---|
488 | * KDTreeQueryResultsXY() to get X- and Y-values
|
---|
489 | * KDTreeQueryResultsTags() to get tag values
|
---|
490 | * KDTreeQueryResultsDistances() to get distances
|
---|
491 |
|
---|
492 | -- ALGLIB --
|
---|
493 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
494 | *************************************************************************/
|
---|
495 | public static int kdtreequeryaknn(kdtree kdt, double[] x, int k, bool selfmatch, double eps)
|
---|
496 | {
|
---|
497 |
|
---|
498 | int result = nearestneighbor.kdtreequeryaknn(kdt.innerobj, x, k, selfmatch, eps);
|
---|
499 | return result;
|
---|
500 | }
|
---|
501 | public static int kdtreequeryaknn(kdtree kdt, double[] x, int k, double eps)
|
---|
502 | {
|
---|
503 | bool selfmatch;
|
---|
504 |
|
---|
505 |
|
---|
506 | selfmatch = true;
|
---|
507 | int result = nearestneighbor.kdtreequeryaknn(kdt.innerobj, x, k, selfmatch, eps);
|
---|
508 |
|
---|
509 | return result;
|
---|
510 | }
|
---|
511 |
|
---|
512 | /*************************************************************************
|
---|
513 | X-values from last query
|
---|
514 |
|
---|
515 | INPUT PARAMETERS
|
---|
516 | KDT - KD-tree
|
---|
517 | X - possibly pre-allocated buffer. If X is too small to store
|
---|
518 | result, it is resized. If size(X) is enough to store
|
---|
519 | result, it is left unchanged.
|
---|
520 |
|
---|
521 | OUTPUT PARAMETERS
|
---|
522 | X - rows are filled with X-values
|
---|
523 |
|
---|
524 | NOTES
|
---|
525 | 1. points are ordered by distance from the query point (first = closest)
|
---|
526 | 2. if XY is larger than required to store result, only leading part will
|
---|
527 | be overwritten; trailing part will be left unchanged. So if on input
|
---|
528 | XY = [[A,B],[C,D]], and result is [1,2], then on exit we will get
|
---|
529 | XY = [[1,2],[C,D]]. This is done purposely to increase performance; if
|
---|
530 | you want function to resize array according to result size, use
|
---|
531 | function with same name and suffix 'I'.
|
---|
532 |
|
---|
533 | SEE ALSO
|
---|
534 | * KDTreeQueryResultsXY() X- and Y-values
|
---|
535 | * KDTreeQueryResultsTags() tag values
|
---|
536 | * KDTreeQueryResultsDistances() distances
|
---|
537 |
|
---|
538 | -- ALGLIB --
|
---|
539 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
540 | *************************************************************************/
|
---|
541 | public static void kdtreequeryresultsx(kdtree kdt, ref double[,] x)
|
---|
542 | {
|
---|
543 |
|
---|
544 | nearestneighbor.kdtreequeryresultsx(kdt.innerobj, ref x);
|
---|
545 | return;
|
---|
546 | }
|
---|
547 |
|
---|
548 | /*************************************************************************
|
---|
549 | X- and Y-values from last query
|
---|
550 |
|
---|
551 | INPUT PARAMETERS
|
---|
552 | KDT - KD-tree
|
---|
553 | XY - possibly pre-allocated buffer. If XY is too small to store
|
---|
554 | result, it is resized. If size(XY) is enough to store
|
---|
555 | result, it is left unchanged.
|
---|
556 |
|
---|
557 | OUTPUT PARAMETERS
|
---|
558 | XY - rows are filled with points: first NX columns with
|
---|
559 | X-values, next NY columns - with Y-values.
|
---|
560 |
|
---|
561 | NOTES
|
---|
562 | 1. points are ordered by distance from the query point (first = closest)
|
---|
563 | 2. if XY is larger than required to store result, only leading part will
|
---|
564 | be overwritten; trailing part will be left unchanged. So if on input
|
---|
565 | XY = [[A,B],[C,D]], and result is [1,2], then on exit we will get
|
---|
566 | XY = [[1,2],[C,D]]. This is done purposely to increase performance; if
|
---|
567 | you want function to resize array according to result size, use
|
---|
568 | function with same name and suffix 'I'.
|
---|
569 |
|
---|
570 | SEE ALSO
|
---|
571 | * KDTreeQueryResultsX() X-values
|
---|
572 | * KDTreeQueryResultsTags() tag values
|
---|
573 | * KDTreeQueryResultsDistances() distances
|
---|
574 |
|
---|
575 | -- ALGLIB --
|
---|
576 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
577 | *************************************************************************/
|
---|
578 | public static void kdtreequeryresultsxy(kdtree kdt, ref double[,] xy)
|
---|
579 | {
|
---|
580 |
|
---|
581 | nearestneighbor.kdtreequeryresultsxy(kdt.innerobj, ref xy);
|
---|
582 | return;
|
---|
583 | }
|
---|
584 |
|
---|
585 | /*************************************************************************
|
---|
586 | Tags from last query
|
---|
587 |
|
---|
588 | INPUT PARAMETERS
|
---|
589 | KDT - KD-tree
|
---|
590 | Tags - possibly pre-allocated buffer. If X is too small to store
|
---|
591 | result, it is resized. If size(X) is enough to store
|
---|
592 | result, it is left unchanged.
|
---|
593 |
|
---|
594 | OUTPUT PARAMETERS
|
---|
595 | Tags - filled with tags associated with points,
|
---|
596 | or, when no tags were supplied, with zeros
|
---|
597 |
|
---|
598 | NOTES
|
---|
599 | 1. points are ordered by distance from the query point (first = closest)
|
---|
600 | 2. if XY is larger than required to store result, only leading part will
|
---|
601 | be overwritten; trailing part will be left unchanged. So if on input
|
---|
602 | XY = [[A,B],[C,D]], and result is [1,2], then on exit we will get
|
---|
603 | XY = [[1,2],[C,D]]. This is done purposely to increase performance; if
|
---|
604 | you want function to resize array according to result size, use
|
---|
605 | function with same name and suffix 'I'.
|
---|
606 |
|
---|
607 | SEE ALSO
|
---|
608 | * KDTreeQueryResultsX() X-values
|
---|
609 | * KDTreeQueryResultsXY() X- and Y-values
|
---|
610 | * KDTreeQueryResultsDistances() distances
|
---|
611 |
|
---|
612 | -- ALGLIB --
|
---|
613 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
614 | *************************************************************************/
|
---|
615 | public static void kdtreequeryresultstags(kdtree kdt, ref int[] tags)
|
---|
616 | {
|
---|
617 |
|
---|
618 | nearestneighbor.kdtreequeryresultstags(kdt.innerobj, ref tags);
|
---|
619 | return;
|
---|
620 | }
|
---|
621 |
|
---|
622 | /*************************************************************************
|
---|
623 | Distances from last query
|
---|
624 |
|
---|
625 | INPUT PARAMETERS
|
---|
626 | KDT - KD-tree
|
---|
627 | R - possibly pre-allocated buffer. If X is too small to store
|
---|
628 | result, it is resized. If size(X) is enough to store
|
---|
629 | result, it is left unchanged.
|
---|
630 |
|
---|
631 | OUTPUT PARAMETERS
|
---|
632 | R - filled with distances (in corresponding norm)
|
---|
633 |
|
---|
634 | NOTES
|
---|
635 | 1. points are ordered by distance from the query point (first = closest)
|
---|
636 | 2. if XY is larger than required to store result, only leading part will
|
---|
637 | be overwritten; trailing part will be left unchanged. So if on input
|
---|
638 | XY = [[A,B],[C,D]], and result is [1,2], then on exit we will get
|
---|
639 | XY = [[1,2],[C,D]]. This is done purposely to increase performance; if
|
---|
640 | you want function to resize array according to result size, use
|
---|
641 | function with same name and suffix 'I'.
|
---|
642 |
|
---|
643 | SEE ALSO
|
---|
644 | * KDTreeQueryResultsX() X-values
|
---|
645 | * KDTreeQueryResultsXY() X- and Y-values
|
---|
646 | * KDTreeQueryResultsTags() tag values
|
---|
647 |
|
---|
648 | -- ALGLIB --
|
---|
649 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
650 | *************************************************************************/
|
---|
651 | public static void kdtreequeryresultsdistances(kdtree kdt, ref double[] r)
|
---|
652 | {
|
---|
653 |
|
---|
654 | nearestneighbor.kdtreequeryresultsdistances(kdt.innerobj, ref r);
|
---|
655 | return;
|
---|
656 | }
|
---|
657 |
|
---|
658 | /*************************************************************************
|
---|
659 | X-values from last query; 'interactive' variant for languages like Python
|
---|
660 | which support constructs like "X = KDTreeQueryResultsXI(KDT)" and
|
---|
661 | interactive mode of interpreter.
|
---|
662 |
|
---|
663 | This function allocates new array on each call, so it is significantly
|
---|
664 | slower than its 'non-interactive' counterpart, but it is more convenient
|
---|
665 | when you call it from command line.
|
---|
666 |
|
---|
667 | -- ALGLIB --
|
---|
668 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
669 | *************************************************************************/
|
---|
670 | public static void kdtreequeryresultsxi(kdtree kdt, out double[,] x)
|
---|
671 | {
|
---|
672 | x = new double[0,0];
|
---|
673 | nearestneighbor.kdtreequeryresultsxi(kdt.innerobj, ref x);
|
---|
674 | return;
|
---|
675 | }
|
---|
676 |
|
---|
677 | /*************************************************************************
|
---|
678 | XY-values from last query; 'interactive' variant for languages like Python
|
---|
679 | which support constructs like "XY = KDTreeQueryResultsXYI(KDT)" and
|
---|
680 | interactive mode of interpreter.
|
---|
681 |
|
---|
682 | This function allocates new array on each call, so it is significantly
|
---|
683 | slower than its 'non-interactive' counterpart, but it is more convenient
|
---|
684 | when you call it from command line.
|
---|
685 |
|
---|
686 | -- ALGLIB --
|
---|
687 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
688 | *************************************************************************/
|
---|
689 | public static void kdtreequeryresultsxyi(kdtree kdt, out double[,] xy)
|
---|
690 | {
|
---|
691 | xy = new double[0,0];
|
---|
692 | nearestneighbor.kdtreequeryresultsxyi(kdt.innerobj, ref xy);
|
---|
693 | return;
|
---|
694 | }
|
---|
695 |
|
---|
696 | /*************************************************************************
|
---|
697 | Tags from last query; 'interactive' variant for languages like Python
|
---|
698 | which support constructs like "Tags = KDTreeQueryResultsTagsI(KDT)" and
|
---|
699 | interactive mode of interpreter.
|
---|
700 |
|
---|
701 | This function allocates new array on each call, so it is significantly
|
---|
702 | slower than its 'non-interactive' counterpart, but it is more convenient
|
---|
703 | when you call it from command line.
|
---|
704 |
|
---|
705 | -- ALGLIB --
|
---|
706 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
707 | *************************************************************************/
|
---|
708 | public static void kdtreequeryresultstagsi(kdtree kdt, out int[] tags)
|
---|
709 | {
|
---|
710 | tags = new int[0];
|
---|
711 | nearestneighbor.kdtreequeryresultstagsi(kdt.innerobj, ref tags);
|
---|
712 | return;
|
---|
713 | }
|
---|
714 |
|
---|
715 | /*************************************************************************
|
---|
716 | Distances from last query; 'interactive' variant for languages like Python
|
---|
717 | which support constructs like "R = KDTreeQueryResultsDistancesI(KDT)"
|
---|
718 | and interactive mode of interpreter.
|
---|
719 |
|
---|
720 | This function allocates new array on each call, so it is significantly
|
---|
721 | slower than its 'non-interactive' counterpart, but it is more convenient
|
---|
722 | when you call it from command line.
|
---|
723 |
|
---|
724 | -- ALGLIB --
|
---|
725 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
726 | *************************************************************************/
|
---|
727 | public static void kdtreequeryresultsdistancesi(kdtree kdt, out double[] r)
|
---|
728 | {
|
---|
729 | r = new double[0];
|
---|
730 | nearestneighbor.kdtreequeryresultsdistancesi(kdt.innerobj, ref r);
|
---|
731 | return;
|
---|
732 | }
|
---|
733 |
|
---|
734 | }
|
---|
735 | public partial class alglib
|
---|
736 | {
|
---|
737 | public class hqrnd
|
---|
738 | {
|
---|
739 | /*************************************************************************
|
---|
740 | Portable high quality random number generator state.
|
---|
741 | Initialized with HQRNDRandomize() or HQRNDSeed().
|
---|
742 |
|
---|
743 | Fields:
|
---|
744 | S1, S2 - seed values
|
---|
745 | V - precomputed value
|
---|
746 | MagicV - 'magic' value used to determine whether State structure
|
---|
747 | was correctly initialized.
|
---|
748 | *************************************************************************/
|
---|
749 | public class hqrndstate
|
---|
750 | {
|
---|
751 | public int s1;
|
---|
752 | public int s2;
|
---|
753 | public double v;
|
---|
754 | public int magicv;
|
---|
755 | };
|
---|
756 |
|
---|
757 |
|
---|
758 |
|
---|
759 |
|
---|
760 | public const int hqrndmax = 2147483563;
|
---|
761 | public const int hqrndm1 = 2147483563;
|
---|
762 | public const int hqrndm2 = 2147483399;
|
---|
763 | public const int hqrndmagic = 1634357784;
|
---|
764 |
|
---|
765 |
|
---|
766 | /*************************************************************************
|
---|
767 | HQRNDState initialization with random values which come from standard
|
---|
768 | RNG.
|
---|
769 |
|
---|
770 | -- ALGLIB --
|
---|
771 | Copyright 02.12.2009 by Bochkanov Sergey
|
---|
772 | *************************************************************************/
|
---|
773 | public static void hqrndrandomize(hqrndstate state)
|
---|
774 | {
|
---|
775 | hqrndseed(math.randominteger(hqrndm1), math.randominteger(hqrndm2), state);
|
---|
776 | }
|
---|
777 |
|
---|
778 |
|
---|
779 | /*************************************************************************
|
---|
780 | HQRNDState initialization with seed values
|
---|
781 |
|
---|
782 | -- ALGLIB --
|
---|
783 | Copyright 02.12.2009 by Bochkanov Sergey
|
---|
784 | *************************************************************************/
|
---|
785 | public static void hqrndseed(int s1,
|
---|
786 | int s2,
|
---|
787 | hqrndstate state)
|
---|
788 | {
|
---|
789 | state.s1 = s1%(hqrndm1-1)+1;
|
---|
790 | state.s2 = s2%(hqrndm2-1)+1;
|
---|
791 | state.v = (double)1/(double)hqrndmax;
|
---|
792 | state.magicv = hqrndmagic;
|
---|
793 | }
|
---|
794 |
|
---|
795 |
|
---|
796 | /*************************************************************************
|
---|
797 | This function generates random real number in (0,1),
|
---|
798 | not including interval boundaries
|
---|
799 |
|
---|
800 | State structure must be initialized with HQRNDRandomize() or HQRNDSeed().
|
---|
801 |
|
---|
802 | -- ALGLIB --
|
---|
803 | Copyright 02.12.2009 by Bochkanov Sergey
|
---|
804 | *************************************************************************/
|
---|
805 | public static double hqrnduniformr(hqrndstate state)
|
---|
806 | {
|
---|
807 | double result = 0;
|
---|
808 |
|
---|
809 | result = state.v*hqrndintegerbase(state);
|
---|
810 | return result;
|
---|
811 | }
|
---|
812 |
|
---|
813 |
|
---|
814 | /*************************************************************************
|
---|
815 | This function generates random integer number in [0, N)
|
---|
816 |
|
---|
817 | 1. N must be less than HQRNDMax-1.
|
---|
818 | 2. State structure must be initialized with HQRNDRandomize() or HQRNDSeed()
|
---|
819 |
|
---|
820 | -- ALGLIB --
|
---|
821 | Copyright 02.12.2009 by Bochkanov Sergey
|
---|
822 | *************************************************************************/
|
---|
823 | public static int hqrnduniformi(hqrndstate state,
|
---|
824 | int n)
|
---|
825 | {
|
---|
826 | int result = 0;
|
---|
827 | int mx = 0;
|
---|
828 |
|
---|
829 |
|
---|
830 | //
|
---|
831 | // Correct handling of N's close to RNDBaseMax
|
---|
832 | // (avoiding skewed distributions for RNDBaseMax<>K*N)
|
---|
833 | //
|
---|
834 | ap.assert(n>0, "HQRNDUniformI: N<=0!");
|
---|
835 | ap.assert(n<hqrndmax-1, "HQRNDUniformI: N>=RNDBaseMax-1!");
|
---|
836 | mx = hqrndmax-1-(hqrndmax-1)%n;
|
---|
837 | do
|
---|
838 | {
|
---|
839 | result = hqrndintegerbase(state)-1;
|
---|
840 | }
|
---|
841 | while( result>=mx );
|
---|
842 | result = result%n;
|
---|
843 | return result;
|
---|
844 | }
|
---|
845 |
|
---|
846 |
|
---|
847 | /*************************************************************************
|
---|
848 | Random number generator: normal numbers
|
---|
849 |
|
---|
850 | This function generates one random number from normal distribution.
|
---|
851 | Its performance is equal to that of HQRNDNormal2()
|
---|
852 |
|
---|
853 | State structure must be initialized with HQRNDRandomize() or HQRNDSeed().
|
---|
854 |
|
---|
855 | -- ALGLIB --
|
---|
856 | Copyright 02.12.2009 by Bochkanov Sergey
|
---|
857 | *************************************************************************/
|
---|
858 | public static double hqrndnormal(hqrndstate state)
|
---|
859 | {
|
---|
860 | double result = 0;
|
---|
861 | double v1 = 0;
|
---|
862 | double v2 = 0;
|
---|
863 |
|
---|
864 | hqrndnormal2(state, ref v1, ref v2);
|
---|
865 | result = v1;
|
---|
866 | return result;
|
---|
867 | }
|
---|
868 |
|
---|
869 |
|
---|
870 | /*************************************************************************
|
---|
871 | Random number generator: random X and Y such that X^2+Y^2=1
|
---|
872 |
|
---|
873 | State structure must be initialized with HQRNDRandomize() or HQRNDSeed().
|
---|
874 |
|
---|
875 | -- ALGLIB --
|
---|
876 | Copyright 02.12.2009 by Bochkanov Sergey
|
---|
877 | *************************************************************************/
|
---|
878 | public static void hqrndunit2(hqrndstate state,
|
---|
879 | ref double x,
|
---|
880 | ref double y)
|
---|
881 | {
|
---|
882 | double v = 0;
|
---|
883 | double mx = 0;
|
---|
884 | double mn = 0;
|
---|
885 |
|
---|
886 | x = 0;
|
---|
887 | y = 0;
|
---|
888 |
|
---|
889 | do
|
---|
890 | {
|
---|
891 | hqrndnormal2(state, ref x, ref y);
|
---|
892 | }
|
---|
893 | while( !((double)(x)!=(double)(0) | (double)(y)!=(double)(0)) );
|
---|
894 | mx = Math.Max(Math.Abs(x), Math.Abs(y));
|
---|
895 | mn = Math.Min(Math.Abs(x), Math.Abs(y));
|
---|
896 | v = mx*Math.Sqrt(1+math.sqr(mn/mx));
|
---|
897 | x = x/v;
|
---|
898 | y = y/v;
|
---|
899 | }
|
---|
900 |
|
---|
901 |
|
---|
902 | /*************************************************************************
|
---|
903 | Random number generator: normal numbers
|
---|
904 |
|
---|
905 | This function generates two independent random numbers from normal
|
---|
906 | distribution. Its performance is equal to that of HQRNDNormal()
|
---|
907 |
|
---|
908 | State structure must be initialized with HQRNDRandomize() or HQRNDSeed().
|
---|
909 |
|
---|
910 | -- ALGLIB --
|
---|
911 | Copyright 02.12.2009 by Bochkanov Sergey
|
---|
912 | *************************************************************************/
|
---|
913 | public static void hqrndnormal2(hqrndstate state,
|
---|
914 | ref double x1,
|
---|
915 | ref double x2)
|
---|
916 | {
|
---|
917 | double u = 0;
|
---|
918 | double v = 0;
|
---|
919 | double s = 0;
|
---|
920 |
|
---|
921 | x1 = 0;
|
---|
922 | x2 = 0;
|
---|
923 |
|
---|
924 | while( true )
|
---|
925 | {
|
---|
926 | u = 2*hqrnduniformr(state)-1;
|
---|
927 | v = 2*hqrnduniformr(state)-1;
|
---|
928 | s = math.sqr(u)+math.sqr(v);
|
---|
929 | if( (double)(s)>(double)(0) & (double)(s)<(double)(1) )
|
---|
930 | {
|
---|
931 |
|
---|
932 | //
|
---|
933 | // two Sqrt's instead of one to
|
---|
934 | // avoid overflow when S is too small
|
---|
935 | //
|
---|
936 | s = Math.Sqrt(-(2*Math.Log(s)))/Math.Sqrt(s);
|
---|
937 | x1 = u*s;
|
---|
938 | x2 = v*s;
|
---|
939 | return;
|
---|
940 | }
|
---|
941 | }
|
---|
942 | }
|
---|
943 |
|
---|
944 |
|
---|
945 | /*************************************************************************
|
---|
946 | Random number generator: exponential distribution
|
---|
947 |
|
---|
948 | State structure must be initialized with HQRNDRandomize() or HQRNDSeed().
|
---|
949 |
|
---|
950 | -- ALGLIB --
|
---|
951 | Copyright 11.08.2007 by Bochkanov Sergey
|
---|
952 | *************************************************************************/
|
---|
953 | public static double hqrndexponential(hqrndstate state,
|
---|
954 | double lambdav)
|
---|
955 | {
|
---|
956 | double result = 0;
|
---|
957 |
|
---|
958 | ap.assert((double)(lambdav)>(double)(0), "HQRNDExponential: LambdaV<=0!");
|
---|
959 | result = -(Math.Log(hqrnduniformr(state))/lambdav);
|
---|
960 | return result;
|
---|
961 | }
|
---|
962 |
|
---|
963 |
|
---|
964 | /*************************************************************************
|
---|
965 |
|
---|
966 | L'Ecuyer, Efficient and portable combined random number generators
|
---|
967 | *************************************************************************/
|
---|
968 | private static int hqrndintegerbase(hqrndstate state)
|
---|
969 | {
|
---|
970 | int result = 0;
|
---|
971 | int k = 0;
|
---|
972 |
|
---|
973 | ap.assert(state.magicv==hqrndmagic, "HQRNDIntegerBase: State is not correctly initialized!");
|
---|
974 | k = state.s1/53668;
|
---|
975 | state.s1 = 40014*(state.s1-k*53668)-k*12211;
|
---|
976 | if( state.s1<0 )
|
---|
977 | {
|
---|
978 | state.s1 = state.s1+2147483563;
|
---|
979 | }
|
---|
980 | k = state.s2/52774;
|
---|
981 | state.s2 = 40692*(state.s2-k*52774)-k*3791;
|
---|
982 | if( state.s2<0 )
|
---|
983 | {
|
---|
984 | state.s2 = state.s2+2147483399;
|
---|
985 | }
|
---|
986 |
|
---|
987 | //
|
---|
988 | // Result
|
---|
989 | //
|
---|
990 | result = state.s1-state.s2;
|
---|
991 | if( result<1 )
|
---|
992 | {
|
---|
993 | result = result+2147483562;
|
---|
994 | }
|
---|
995 | return result;
|
---|
996 | }
|
---|
997 |
|
---|
998 |
|
---|
999 | }
|
---|
1000 | public class nearestneighbor
|
---|
1001 | {
|
---|
1002 | public class kdtree
|
---|
1003 | {
|
---|
1004 | public int n;
|
---|
1005 | public int nx;
|
---|
1006 | public int ny;
|
---|
1007 | public int normtype;
|
---|
1008 | public double[,] xy;
|
---|
1009 | public int[] tags;
|
---|
1010 | public double[] boxmin;
|
---|
1011 | public double[] boxmax;
|
---|
1012 | public int[] nodes;
|
---|
1013 | public double[] splits;
|
---|
1014 | public double[] x;
|
---|
1015 | public int kneeded;
|
---|
1016 | public double rneeded;
|
---|
1017 | public bool selfmatch;
|
---|
1018 | public double approxf;
|
---|
1019 | public int kcur;
|
---|
1020 | public int[] idx;
|
---|
1021 | public double[] r;
|
---|
1022 | public double[] buf;
|
---|
1023 | public double[] curboxmin;
|
---|
1024 | public double[] curboxmax;
|
---|
1025 | public double curdist;
|
---|
1026 | public int debugcounter;
|
---|
1027 | public kdtree()
|
---|
1028 | {
|
---|
1029 | xy = new double[0,0];
|
---|
1030 | tags = new int[0];
|
---|
1031 | boxmin = new double[0];
|
---|
1032 | boxmax = new double[0];
|
---|
1033 | nodes = new int[0];
|
---|
1034 | splits = new double[0];
|
---|
1035 | x = new double[0];
|
---|
1036 | idx = new int[0];
|
---|
1037 | r = new double[0];
|
---|
1038 | buf = new double[0];
|
---|
1039 | curboxmin = new double[0];
|
---|
1040 | curboxmax = new double[0];
|
---|
1041 | }
|
---|
1042 | };
|
---|
1043 |
|
---|
1044 |
|
---|
1045 |
|
---|
1046 |
|
---|
1047 | public const int splitnodesize = 6;
|
---|
1048 | public const int kdtreefirstversion = 0;
|
---|
1049 |
|
---|
1050 |
|
---|
1051 | /*************************************************************************
|
---|
1052 | KD-tree creation
|
---|
1053 |
|
---|
1054 | This subroutine creates KD-tree from set of X-values and optional Y-values
|
---|
1055 |
|
---|
1056 | INPUT PARAMETERS
|
---|
1057 | XY - dataset, array[0..N-1,0..NX+NY-1].
|
---|
1058 | one row corresponds to one point.
|
---|
1059 | first NX columns contain X-values, next NY (NY may be zero)
|
---|
1060 | columns may contain associated Y-values
|
---|
1061 | N - number of points, N>=1
|
---|
1062 | NX - space dimension, NX>=1.
|
---|
1063 | NY - number of optional Y-values, NY>=0.
|
---|
1064 | NormType- norm type:
|
---|
1065 | * 0 denotes infinity-norm
|
---|
1066 | * 1 denotes 1-norm
|
---|
1067 | * 2 denotes 2-norm (Euclidean norm)
|
---|
1068 |
|
---|
1069 | OUTPUT PARAMETERS
|
---|
1070 | KDT - KD-tree
|
---|
1071 |
|
---|
1072 |
|
---|
1073 | NOTES
|
---|
1074 |
|
---|
1075 | 1. KD-tree creation have O(N*logN) complexity and O(N*(2*NX+NY)) memory
|
---|
1076 | requirements.
|
---|
1077 | 2. Although KD-trees may be used with any combination of N and NX, they
|
---|
1078 | are more efficient than brute-force search only when N >> 4^NX. So they
|
---|
1079 | are most useful in low-dimensional tasks (NX=2, NX=3). NX=1 is another
|
---|
1080 | inefficient case, because simple binary search (without additional
|
---|
1081 | structures) is much more efficient in such tasks than KD-trees.
|
---|
1082 |
|
---|
1083 | -- ALGLIB --
|
---|
1084 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
1085 | *************************************************************************/
|
---|
1086 | public static void kdtreebuild(double[,] xy,
|
---|
1087 | int n,
|
---|
1088 | int nx,
|
---|
1089 | int ny,
|
---|
1090 | int normtype,
|
---|
1091 | kdtree kdt)
|
---|
1092 | {
|
---|
1093 | int[] tags = new int[0];
|
---|
1094 | int i = 0;
|
---|
1095 |
|
---|
1096 | ap.assert(n>=1, "KDTreeBuild: N<1!");
|
---|
1097 | ap.assert(nx>=1, "KDTreeBuild: NX<1!");
|
---|
1098 | ap.assert(ny>=0, "KDTreeBuild: NY<0!");
|
---|
1099 | ap.assert(normtype>=0 & normtype<=2, "KDTreeBuild: incorrect NormType!");
|
---|
1100 | ap.assert(ap.rows(xy)>=n, "KDTreeBuild: rows(X)<N!");
|
---|
1101 | ap.assert(ap.cols(xy)>=nx+ny, "KDTreeBuild: cols(X)<NX+NY!");
|
---|
1102 | ap.assert(apserv.apservisfinitematrix(xy, n, nx+ny), "KDTreeBuild: X contains infinite or NaN values!");
|
---|
1103 | tags = new int[n];
|
---|
1104 | for(i=0; i<=n-1; i++)
|
---|
1105 | {
|
---|
1106 | tags[i] = 0;
|
---|
1107 | }
|
---|
1108 | kdtreebuildtagged(xy, tags, n, nx, ny, normtype, kdt);
|
---|
1109 | }
|
---|
1110 |
|
---|
1111 |
|
---|
1112 | /*************************************************************************
|
---|
1113 | KD-tree creation
|
---|
1114 |
|
---|
1115 | This subroutine creates KD-tree from set of X-values, integer tags and
|
---|
1116 | optional Y-values
|
---|
1117 |
|
---|
1118 | INPUT PARAMETERS
|
---|
1119 | XY - dataset, array[0..N-1,0..NX+NY-1].
|
---|
1120 | one row corresponds to one point.
|
---|
1121 | first NX columns contain X-values, next NY (NY may be zero)
|
---|
1122 | columns may contain associated Y-values
|
---|
1123 | Tags - tags, array[0..N-1], contains integer tags associated
|
---|
1124 | with points.
|
---|
1125 | N - number of points, N>=1
|
---|
1126 | NX - space dimension, NX>=1.
|
---|
1127 | NY - number of optional Y-values, NY>=0.
|
---|
1128 | NormType- norm type:
|
---|
1129 | * 0 denotes infinity-norm
|
---|
1130 | * 1 denotes 1-norm
|
---|
1131 | * 2 denotes 2-norm (Euclidean norm)
|
---|
1132 |
|
---|
1133 | OUTPUT PARAMETERS
|
---|
1134 | KDT - KD-tree
|
---|
1135 |
|
---|
1136 | NOTES
|
---|
1137 |
|
---|
1138 | 1. KD-tree creation have O(N*logN) complexity and O(N*(2*NX+NY)) memory
|
---|
1139 | requirements.
|
---|
1140 | 2. Although KD-trees may be used with any combination of N and NX, they
|
---|
1141 | are more efficient than brute-force search only when N >> 4^NX. So they
|
---|
1142 | are most useful in low-dimensional tasks (NX=2, NX=3). NX=1 is another
|
---|
1143 | inefficient case, because simple binary search (without additional
|
---|
1144 | structures) is much more efficient in such tasks than KD-trees.
|
---|
1145 |
|
---|
1146 | -- ALGLIB --
|
---|
1147 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
1148 | *************************************************************************/
|
---|
1149 | public static void kdtreebuildtagged(double[,] xy,
|
---|
1150 | int[] tags,
|
---|
1151 | int n,
|
---|
1152 | int nx,
|
---|
1153 | int ny,
|
---|
1154 | int normtype,
|
---|
1155 | kdtree kdt)
|
---|
1156 | {
|
---|
1157 | int i = 0;
|
---|
1158 | int j = 0;
|
---|
1159 | int maxnodes = 0;
|
---|
1160 | int nodesoffs = 0;
|
---|
1161 | int splitsoffs = 0;
|
---|
1162 | int i_ = 0;
|
---|
1163 | int i1_ = 0;
|
---|
1164 |
|
---|
1165 | ap.assert(n>=1, "KDTreeBuildTagged: N<1!");
|
---|
1166 | ap.assert(nx>=1, "KDTreeBuildTagged: NX<1!");
|
---|
1167 | ap.assert(ny>=0, "KDTreeBuildTagged: NY<0!");
|
---|
1168 | ap.assert(normtype>=0 & normtype<=2, "KDTreeBuildTagged: incorrect NormType!");
|
---|
1169 | ap.assert(ap.rows(xy)>=n, "KDTreeBuildTagged: rows(X)<N!");
|
---|
1170 | ap.assert(ap.cols(xy)>=nx+ny, "KDTreeBuildTagged: cols(X)<NX+NY!");
|
---|
1171 | ap.assert(apserv.apservisfinitematrix(xy, n, nx+ny), "KDTreeBuildTagged: X contains infinite or NaN values!");
|
---|
1172 |
|
---|
1173 | //
|
---|
1174 | // initialize
|
---|
1175 | //
|
---|
1176 | kdt.n = n;
|
---|
1177 | kdt.nx = nx;
|
---|
1178 | kdt.ny = ny;
|
---|
1179 | kdt.normtype = normtype;
|
---|
1180 |
|
---|
1181 | //
|
---|
1182 | // Allocate
|
---|
1183 | //
|
---|
1184 | kdtreeallocdatasetindependent(kdt, nx, ny);
|
---|
1185 | kdtreeallocdatasetdependent(kdt, n, nx, ny);
|
---|
1186 |
|
---|
1187 | //
|
---|
1188 | // Initial fill
|
---|
1189 | //
|
---|
1190 | for(i=0; i<=n-1; i++)
|
---|
1191 | {
|
---|
1192 | for(i_=0; i_<=nx-1;i_++)
|
---|
1193 | {
|
---|
1194 | kdt.xy[i,i_] = xy[i,i_];
|
---|
1195 | }
|
---|
1196 | i1_ = (0) - (nx);
|
---|
1197 | for(i_=nx; i_<=2*nx+ny-1;i_++)
|
---|
1198 | {
|
---|
1199 | kdt.xy[i,i_] = xy[i,i_+i1_];
|
---|
1200 | }
|
---|
1201 | kdt.tags[i] = tags[i];
|
---|
1202 | }
|
---|
1203 |
|
---|
1204 | //
|
---|
1205 | // Determine bounding box
|
---|
1206 | //
|
---|
1207 | for(i_=0; i_<=nx-1;i_++)
|
---|
1208 | {
|
---|
1209 | kdt.boxmin[i_] = kdt.xy[0,i_];
|
---|
1210 | }
|
---|
1211 | for(i_=0; i_<=nx-1;i_++)
|
---|
1212 | {
|
---|
1213 | kdt.boxmax[i_] = kdt.xy[0,i_];
|
---|
1214 | }
|
---|
1215 | for(i=1; i<=n-1; i++)
|
---|
1216 | {
|
---|
1217 | for(j=0; j<=nx-1; j++)
|
---|
1218 | {
|
---|
1219 | kdt.boxmin[j] = Math.Min(kdt.boxmin[j], kdt.xy[i,j]);
|
---|
1220 | kdt.boxmax[j] = Math.Max(kdt.boxmax[j], kdt.xy[i,j]);
|
---|
1221 | }
|
---|
1222 | }
|
---|
1223 |
|
---|
1224 | //
|
---|
1225 | // prepare tree structure
|
---|
1226 | // * MaxNodes=N because we guarantee no trivial splits, i.e.
|
---|
1227 | // every split will generate two non-empty boxes
|
---|
1228 | //
|
---|
1229 | maxnodes = n;
|
---|
1230 | kdt.nodes = new int[splitnodesize*2*maxnodes];
|
---|
1231 | kdt.splits = new double[2*maxnodes];
|
---|
1232 | nodesoffs = 0;
|
---|
1233 | splitsoffs = 0;
|
---|
1234 | for(i_=0; i_<=nx-1;i_++)
|
---|
1235 | {
|
---|
1236 | kdt.curboxmin[i_] = kdt.boxmin[i_];
|
---|
1237 | }
|
---|
1238 | for(i_=0; i_<=nx-1;i_++)
|
---|
1239 | {
|
---|
1240 | kdt.curboxmax[i_] = kdt.boxmax[i_];
|
---|
1241 | }
|
---|
1242 | kdtreegeneratetreerec(kdt, ref nodesoffs, ref splitsoffs, 0, n, 8);
|
---|
1243 |
|
---|
1244 | //
|
---|
1245 | // Set current query size to 0
|
---|
1246 | //
|
---|
1247 | kdt.kcur = 0;
|
---|
1248 | }
|
---|
1249 |
|
---|
1250 |
|
---|
1251 | /*************************************************************************
|
---|
1252 | K-NN query: K nearest neighbors
|
---|
1253 |
|
---|
1254 | INPUT PARAMETERS
|
---|
1255 | KDT - KD-tree
|
---|
1256 | X - point, array[0..NX-1].
|
---|
1257 | K - number of neighbors to return, K>=1
|
---|
1258 | SelfMatch - whether self-matches are allowed:
|
---|
1259 | * if True, nearest neighbor may be the point itself
|
---|
1260 | (if it exists in original dataset)
|
---|
1261 | * if False, then only points with non-zero distance
|
---|
1262 | are returned
|
---|
1263 | * if not given, considered True
|
---|
1264 |
|
---|
1265 | RESULT
|
---|
1266 | number of actual neighbors found (either K or N, if K>N).
|
---|
1267 |
|
---|
1268 | This subroutine performs query and stores its result in the internal
|
---|
1269 | structures of the KD-tree. You can use following subroutines to obtain
|
---|
1270 | these results:
|
---|
1271 | * KDTreeQueryResultsX() to get X-values
|
---|
1272 | * KDTreeQueryResultsXY() to get X- and Y-values
|
---|
1273 | * KDTreeQueryResultsTags() to get tag values
|
---|
1274 | * KDTreeQueryResultsDistances() to get distances
|
---|
1275 |
|
---|
1276 | -- ALGLIB --
|
---|
1277 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
1278 | *************************************************************************/
|
---|
1279 | public static int kdtreequeryknn(kdtree kdt,
|
---|
1280 | double[] x,
|
---|
1281 | int k,
|
---|
1282 | bool selfmatch)
|
---|
1283 | {
|
---|
1284 | int result = 0;
|
---|
1285 |
|
---|
1286 | ap.assert(k>=1, "KDTreeQueryKNN: K<1!");
|
---|
1287 | ap.assert(ap.len(x)>=kdt.nx, "KDTreeQueryKNN: Length(X)<NX!");
|
---|
1288 | ap.assert(apserv.isfinitevector(x, kdt.nx), "KDTreeQueryKNN: X contains infinite or NaN values!");
|
---|
1289 | result = kdtreequeryaknn(kdt, x, k, selfmatch, 0.0);
|
---|
1290 | return result;
|
---|
1291 | }
|
---|
1292 |
|
---|
1293 |
|
---|
1294 | /*************************************************************************
|
---|
1295 | R-NN query: all points within R-sphere centered at X
|
---|
1296 |
|
---|
1297 | INPUT PARAMETERS
|
---|
1298 | KDT - KD-tree
|
---|
1299 | X - point, array[0..NX-1].
|
---|
1300 | R - radius of sphere (in corresponding norm), R>0
|
---|
1301 | SelfMatch - whether self-matches are allowed:
|
---|
1302 | * if True, nearest neighbor may be the point itself
|
---|
1303 | (if it exists in original dataset)
|
---|
1304 | * if False, then only points with non-zero distance
|
---|
1305 | are returned
|
---|
1306 | * if not given, considered True
|
---|
1307 |
|
---|
1308 | RESULT
|
---|
1309 | number of neighbors found, >=0
|
---|
1310 |
|
---|
1311 | This subroutine performs query and stores its result in the internal
|
---|
1312 | structures of the KD-tree. You can use following subroutines to obtain
|
---|
1313 | actual results:
|
---|
1314 | * KDTreeQueryResultsX() to get X-values
|
---|
1315 | * KDTreeQueryResultsXY() to get X- and Y-values
|
---|
1316 | * KDTreeQueryResultsTags() to get tag values
|
---|
1317 | * KDTreeQueryResultsDistances() to get distances
|
---|
1318 |
|
---|
1319 | -- ALGLIB --
|
---|
1320 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
1321 | *************************************************************************/
|
---|
1322 | public static int kdtreequeryrnn(kdtree kdt,
|
---|
1323 | double[] x,
|
---|
1324 | double r,
|
---|
1325 | bool selfmatch)
|
---|
1326 | {
|
---|
1327 | int result = 0;
|
---|
1328 | int i = 0;
|
---|
1329 | int j = 0;
|
---|
1330 |
|
---|
1331 | ap.assert((double)(r)>(double)(0), "KDTreeQueryRNN: incorrect R!");
|
---|
1332 | ap.assert(ap.len(x)>=kdt.nx, "KDTreeQueryRNN: Length(X)<NX!");
|
---|
1333 | ap.assert(apserv.isfinitevector(x, kdt.nx), "KDTreeQueryRNN: X contains infinite or NaN values!");
|
---|
1334 |
|
---|
1335 | //
|
---|
1336 | // Prepare parameters
|
---|
1337 | //
|
---|
1338 | kdt.kneeded = 0;
|
---|
1339 | if( kdt.normtype!=2 )
|
---|
1340 | {
|
---|
1341 | kdt.rneeded = r;
|
---|
1342 | }
|
---|
1343 | else
|
---|
1344 | {
|
---|
1345 | kdt.rneeded = math.sqr(r);
|
---|
1346 | }
|
---|
1347 | kdt.selfmatch = selfmatch;
|
---|
1348 | kdt.approxf = 1;
|
---|
1349 | kdt.kcur = 0;
|
---|
1350 |
|
---|
1351 | //
|
---|
1352 | // calculate distance from point to current bounding box
|
---|
1353 | //
|
---|
1354 | kdtreeinitbox(kdt, x);
|
---|
1355 |
|
---|
1356 | //
|
---|
1357 | // call recursive search
|
---|
1358 | // results are returned as heap
|
---|
1359 | //
|
---|
1360 | kdtreequerynnrec(kdt, 0);
|
---|
1361 |
|
---|
1362 | //
|
---|
1363 | // pop from heap to generate ordered representation
|
---|
1364 | //
|
---|
1365 | // last element is not pop'ed because it is already in
|
---|
1366 | // its place
|
---|
1367 | //
|
---|
1368 | result = kdt.kcur;
|
---|
1369 | j = kdt.kcur;
|
---|
1370 | for(i=kdt.kcur; i>=2; i--)
|
---|
1371 | {
|
---|
1372 | tsort.tagheappopi(ref kdt.r, ref kdt.idx, ref j);
|
---|
1373 | }
|
---|
1374 | return result;
|
---|
1375 | }
|
---|
1376 |
|
---|
1377 |
|
---|
1378 | /*************************************************************************
|
---|
1379 | K-NN query: approximate K nearest neighbors
|
---|
1380 |
|
---|
1381 | INPUT PARAMETERS
|
---|
1382 | KDT - KD-tree
|
---|
1383 | X - point, array[0..NX-1].
|
---|
1384 | K - number of neighbors to return, K>=1
|
---|
1385 | SelfMatch - whether self-matches are allowed:
|
---|
1386 | * if True, nearest neighbor may be the point itself
|
---|
1387 | (if it exists in original dataset)
|
---|
1388 | * if False, then only points with non-zero distance
|
---|
1389 | are returned
|
---|
1390 | * if not given, considered True
|
---|
1391 | Eps - approximation factor, Eps>=0. eps-approximate nearest
|
---|
1392 | neighbor is a neighbor whose distance from X is at
|
---|
1393 | most (1+eps) times distance of true nearest neighbor.
|
---|
1394 |
|
---|
1395 | RESULT
|
---|
1396 | number of actual neighbors found (either K or N, if K>N).
|
---|
1397 |
|
---|
1398 | NOTES
|
---|
1399 | significant performance gain may be achieved only when Eps is is on
|
---|
1400 | the order of magnitude of 1 or larger.
|
---|
1401 |
|
---|
1402 | This subroutine performs query and stores its result in the internal
|
---|
1403 | structures of the KD-tree. You can use following subroutines to obtain
|
---|
1404 | these results:
|
---|
1405 | * KDTreeQueryResultsX() to get X-values
|
---|
1406 | * KDTreeQueryResultsXY() to get X- and Y-values
|
---|
1407 | * KDTreeQueryResultsTags() to get tag values
|
---|
1408 | * KDTreeQueryResultsDistances() to get distances
|
---|
1409 |
|
---|
1410 | -- ALGLIB --
|
---|
1411 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
1412 | *************************************************************************/
|
---|
1413 | public static int kdtreequeryaknn(kdtree kdt,
|
---|
1414 | double[] x,
|
---|
1415 | int k,
|
---|
1416 | bool selfmatch,
|
---|
1417 | double eps)
|
---|
1418 | {
|
---|
1419 | int result = 0;
|
---|
1420 | int i = 0;
|
---|
1421 | int j = 0;
|
---|
1422 |
|
---|
1423 | ap.assert(k>0, "KDTreeQueryAKNN: incorrect K!");
|
---|
1424 | ap.assert((double)(eps)>=(double)(0), "KDTreeQueryAKNN: incorrect Eps!");
|
---|
1425 | ap.assert(ap.len(x)>=kdt.nx, "KDTreeQueryAKNN: Length(X)<NX!");
|
---|
1426 | ap.assert(apserv.isfinitevector(x, kdt.nx), "KDTreeQueryAKNN: X contains infinite or NaN values!");
|
---|
1427 |
|
---|
1428 | //
|
---|
1429 | // Prepare parameters
|
---|
1430 | //
|
---|
1431 | k = Math.Min(k, kdt.n);
|
---|
1432 | kdt.kneeded = k;
|
---|
1433 | kdt.rneeded = 0;
|
---|
1434 | kdt.selfmatch = selfmatch;
|
---|
1435 | if( kdt.normtype==2 )
|
---|
1436 | {
|
---|
1437 | kdt.approxf = 1/math.sqr(1+eps);
|
---|
1438 | }
|
---|
1439 | else
|
---|
1440 | {
|
---|
1441 | kdt.approxf = 1/(1+eps);
|
---|
1442 | }
|
---|
1443 | kdt.kcur = 0;
|
---|
1444 |
|
---|
1445 | //
|
---|
1446 | // calculate distance from point to current bounding box
|
---|
1447 | //
|
---|
1448 | kdtreeinitbox(kdt, x);
|
---|
1449 |
|
---|
1450 | //
|
---|
1451 | // call recursive search
|
---|
1452 | // results are returned as heap
|
---|
1453 | //
|
---|
1454 | kdtreequerynnrec(kdt, 0);
|
---|
1455 |
|
---|
1456 | //
|
---|
1457 | // pop from heap to generate ordered representation
|
---|
1458 | //
|
---|
1459 | // last element is non pop'ed because it is already in
|
---|
1460 | // its place
|
---|
1461 | //
|
---|
1462 | result = kdt.kcur;
|
---|
1463 | j = kdt.kcur;
|
---|
1464 | for(i=kdt.kcur; i>=2; i--)
|
---|
1465 | {
|
---|
1466 | tsort.tagheappopi(ref kdt.r, ref kdt.idx, ref j);
|
---|
1467 | }
|
---|
1468 | return result;
|
---|
1469 | }
|
---|
1470 |
|
---|
1471 |
|
---|
1472 | /*************************************************************************
|
---|
1473 | X-values from last query
|
---|
1474 |
|
---|
1475 | INPUT PARAMETERS
|
---|
1476 | KDT - KD-tree
|
---|
1477 | X - possibly pre-allocated buffer. If X is too small to store
|
---|
1478 | result, it is resized. If size(X) is enough to store
|
---|
1479 | result, it is left unchanged.
|
---|
1480 |
|
---|
1481 | OUTPUT PARAMETERS
|
---|
1482 | X - rows are filled with X-values
|
---|
1483 |
|
---|
1484 | NOTES
|
---|
1485 | 1. points are ordered by distance from the query point (first = closest)
|
---|
1486 | 2. if XY is larger than required to store result, only leading part will
|
---|
1487 | be overwritten; trailing part will be left unchanged. So if on input
|
---|
1488 | XY = [[A,B],[C,D]], and result is [1,2], then on exit we will get
|
---|
1489 | XY = [[1,2],[C,D]]. This is done purposely to increase performance; if
|
---|
1490 | you want function to resize array according to result size, use
|
---|
1491 | function with same name and suffix 'I'.
|
---|
1492 |
|
---|
1493 | SEE ALSO
|
---|
1494 | * KDTreeQueryResultsXY() X- and Y-values
|
---|
1495 | * KDTreeQueryResultsTags() tag values
|
---|
1496 | * KDTreeQueryResultsDistances() distances
|
---|
1497 |
|
---|
1498 | -- ALGLIB --
|
---|
1499 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
1500 | *************************************************************************/
|
---|
1501 | public static void kdtreequeryresultsx(kdtree kdt,
|
---|
1502 | ref double[,] x)
|
---|
1503 | {
|
---|
1504 | int i = 0;
|
---|
1505 | int k = 0;
|
---|
1506 | int i_ = 0;
|
---|
1507 | int i1_ = 0;
|
---|
1508 |
|
---|
1509 | if( kdt.kcur==0 )
|
---|
1510 | {
|
---|
1511 | return;
|
---|
1512 | }
|
---|
1513 | if( ap.rows(x)<kdt.kcur | ap.cols(x)<kdt.nx )
|
---|
1514 | {
|
---|
1515 | x = new double[kdt.kcur, kdt.nx];
|
---|
1516 | }
|
---|
1517 | k = kdt.kcur;
|
---|
1518 | for(i=0; i<=k-1; i++)
|
---|
1519 | {
|
---|
1520 | i1_ = (kdt.nx) - (0);
|
---|
1521 | for(i_=0; i_<=kdt.nx-1;i_++)
|
---|
1522 | {
|
---|
1523 | x[i,i_] = kdt.xy[kdt.idx[i],i_+i1_];
|
---|
1524 | }
|
---|
1525 | }
|
---|
1526 | }
|
---|
1527 |
|
---|
1528 |
|
---|
1529 | /*************************************************************************
|
---|
1530 | X- and Y-values from last query
|
---|
1531 |
|
---|
1532 | INPUT PARAMETERS
|
---|
1533 | KDT - KD-tree
|
---|
1534 | XY - possibly pre-allocated buffer. If XY is too small to store
|
---|
1535 | result, it is resized. If size(XY) is enough to store
|
---|
1536 | result, it is left unchanged.
|
---|
1537 |
|
---|
1538 | OUTPUT PARAMETERS
|
---|
1539 | XY - rows are filled with points: first NX columns with
|
---|
1540 | X-values, next NY columns - with Y-values.
|
---|
1541 |
|
---|
1542 | NOTES
|
---|
1543 | 1. points are ordered by distance from the query point (first = closest)
|
---|
1544 | 2. if XY is larger than required to store result, only leading part will
|
---|
1545 | be overwritten; trailing part will be left unchanged. So if on input
|
---|
1546 | XY = [[A,B],[C,D]], and result is [1,2], then on exit we will get
|
---|
1547 | XY = [[1,2],[C,D]]. This is done purposely to increase performance; if
|
---|
1548 | you want function to resize array according to result size, use
|
---|
1549 | function with same name and suffix 'I'.
|
---|
1550 |
|
---|
1551 | SEE ALSO
|
---|
1552 | * KDTreeQueryResultsX() X-values
|
---|
1553 | * KDTreeQueryResultsTags() tag values
|
---|
1554 | * KDTreeQueryResultsDistances() distances
|
---|
1555 |
|
---|
1556 | -- ALGLIB --
|
---|
1557 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
1558 | *************************************************************************/
|
---|
1559 | public static void kdtreequeryresultsxy(kdtree kdt,
|
---|
1560 | ref double[,] xy)
|
---|
1561 | {
|
---|
1562 | int i = 0;
|
---|
1563 | int k = 0;
|
---|
1564 | int i_ = 0;
|
---|
1565 | int i1_ = 0;
|
---|
1566 |
|
---|
1567 | if( kdt.kcur==0 )
|
---|
1568 | {
|
---|
1569 | return;
|
---|
1570 | }
|
---|
1571 | if( ap.rows(xy)<kdt.kcur | ap.cols(xy)<kdt.nx+kdt.ny )
|
---|
1572 | {
|
---|
1573 | xy = new double[kdt.kcur, kdt.nx+kdt.ny];
|
---|
1574 | }
|
---|
1575 | k = kdt.kcur;
|
---|
1576 | for(i=0; i<=k-1; i++)
|
---|
1577 | {
|
---|
1578 | i1_ = (kdt.nx) - (0);
|
---|
1579 | for(i_=0; i_<=kdt.nx+kdt.ny-1;i_++)
|
---|
1580 | {
|
---|
1581 | xy[i,i_] = kdt.xy[kdt.idx[i],i_+i1_];
|
---|
1582 | }
|
---|
1583 | }
|
---|
1584 | }
|
---|
1585 |
|
---|
1586 |
|
---|
1587 | /*************************************************************************
|
---|
1588 | Tags from last query
|
---|
1589 |
|
---|
1590 | INPUT PARAMETERS
|
---|
1591 | KDT - KD-tree
|
---|
1592 | Tags - possibly pre-allocated buffer. If X is too small to store
|
---|
1593 | result, it is resized. If size(X) is enough to store
|
---|
1594 | result, it is left unchanged.
|
---|
1595 |
|
---|
1596 | OUTPUT PARAMETERS
|
---|
1597 | Tags - filled with tags associated with points,
|
---|
1598 | or, when no tags were supplied, with zeros
|
---|
1599 |
|
---|
1600 | NOTES
|
---|
1601 | 1. points are ordered by distance from the query point (first = closest)
|
---|
1602 | 2. if XY is larger than required to store result, only leading part will
|
---|
1603 | be overwritten; trailing part will be left unchanged. So if on input
|
---|
1604 | XY = [[A,B],[C,D]], and result is [1,2], then on exit we will get
|
---|
1605 | XY = [[1,2],[C,D]]. This is done purposely to increase performance; if
|
---|
1606 | you want function to resize array according to result size, use
|
---|
1607 | function with same name and suffix 'I'.
|
---|
1608 |
|
---|
1609 | SEE ALSO
|
---|
1610 | * KDTreeQueryResultsX() X-values
|
---|
1611 | * KDTreeQueryResultsXY() X- and Y-values
|
---|
1612 | * KDTreeQueryResultsDistances() distances
|
---|
1613 |
|
---|
1614 | -- ALGLIB --
|
---|
1615 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
1616 | *************************************************************************/
|
---|
1617 | public static void kdtreequeryresultstags(kdtree kdt,
|
---|
1618 | ref int[] tags)
|
---|
1619 | {
|
---|
1620 | int i = 0;
|
---|
1621 | int k = 0;
|
---|
1622 |
|
---|
1623 | if( kdt.kcur==0 )
|
---|
1624 | {
|
---|
1625 | return;
|
---|
1626 | }
|
---|
1627 | if( ap.len(tags)<kdt.kcur )
|
---|
1628 | {
|
---|
1629 | tags = new int[kdt.kcur];
|
---|
1630 | }
|
---|
1631 | k = kdt.kcur;
|
---|
1632 | for(i=0; i<=k-1; i++)
|
---|
1633 | {
|
---|
1634 | tags[i] = kdt.tags[kdt.idx[i]];
|
---|
1635 | }
|
---|
1636 | }
|
---|
1637 |
|
---|
1638 |
|
---|
1639 | /*************************************************************************
|
---|
1640 | Distances from last query
|
---|
1641 |
|
---|
1642 | INPUT PARAMETERS
|
---|
1643 | KDT - KD-tree
|
---|
1644 | R - possibly pre-allocated buffer. If X is too small to store
|
---|
1645 | result, it is resized. If size(X) is enough to store
|
---|
1646 | result, it is left unchanged.
|
---|
1647 |
|
---|
1648 | OUTPUT PARAMETERS
|
---|
1649 | R - filled with distances (in corresponding norm)
|
---|
1650 |
|
---|
1651 | NOTES
|
---|
1652 | 1. points are ordered by distance from the query point (first = closest)
|
---|
1653 | 2. if XY is larger than required to store result, only leading part will
|
---|
1654 | be overwritten; trailing part will be left unchanged. So if on input
|
---|
1655 | XY = [[A,B],[C,D]], and result is [1,2], then on exit we will get
|
---|
1656 | XY = [[1,2],[C,D]]. This is done purposely to increase performance; if
|
---|
1657 | you want function to resize array according to result size, use
|
---|
1658 | function with same name and suffix 'I'.
|
---|
1659 |
|
---|
1660 | SEE ALSO
|
---|
1661 | * KDTreeQueryResultsX() X-values
|
---|
1662 | * KDTreeQueryResultsXY() X- and Y-values
|
---|
1663 | * KDTreeQueryResultsTags() tag values
|
---|
1664 |
|
---|
1665 | -- ALGLIB --
|
---|
1666 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
1667 | *************************************************************************/
|
---|
1668 | public static void kdtreequeryresultsdistances(kdtree kdt,
|
---|
1669 | ref double[] r)
|
---|
1670 | {
|
---|
1671 | int i = 0;
|
---|
1672 | int k = 0;
|
---|
1673 |
|
---|
1674 | if( kdt.kcur==0 )
|
---|
1675 | {
|
---|
1676 | return;
|
---|
1677 | }
|
---|
1678 | if( ap.len(r)<kdt.kcur )
|
---|
1679 | {
|
---|
1680 | r = new double[kdt.kcur];
|
---|
1681 | }
|
---|
1682 | k = kdt.kcur;
|
---|
1683 |
|
---|
1684 | //
|
---|
1685 | // unload norms
|
---|
1686 | //
|
---|
1687 | // Abs() call is used to handle cases with negative norms
|
---|
1688 | // (generated during KFN requests)
|
---|
1689 | //
|
---|
1690 | if( kdt.normtype==0 )
|
---|
1691 | {
|
---|
1692 | for(i=0; i<=k-1; i++)
|
---|
1693 | {
|
---|
1694 | r[i] = Math.Abs(kdt.r[i]);
|
---|
1695 | }
|
---|
1696 | }
|
---|
1697 | if( kdt.normtype==1 )
|
---|
1698 | {
|
---|
1699 | for(i=0; i<=k-1; i++)
|
---|
1700 | {
|
---|
1701 | r[i] = Math.Abs(kdt.r[i]);
|
---|
1702 | }
|
---|
1703 | }
|
---|
1704 | if( kdt.normtype==2 )
|
---|
1705 | {
|
---|
1706 | for(i=0; i<=k-1; i++)
|
---|
1707 | {
|
---|
1708 | r[i] = Math.Sqrt(Math.Abs(kdt.r[i]));
|
---|
1709 | }
|
---|
1710 | }
|
---|
1711 | }
|
---|
1712 |
|
---|
1713 |
|
---|
1714 | /*************************************************************************
|
---|
1715 | X-values from last query; 'interactive' variant for languages like Python
|
---|
1716 | which support constructs like "X = KDTreeQueryResultsXI(KDT)" and
|
---|
1717 | interactive mode of interpreter.
|
---|
1718 |
|
---|
1719 | This function allocates new array on each call, so it is significantly
|
---|
1720 | slower than its 'non-interactive' counterpart, but it is more convenient
|
---|
1721 | when you call it from command line.
|
---|
1722 |
|
---|
1723 | -- ALGLIB --
|
---|
1724 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
1725 | *************************************************************************/
|
---|
1726 | public static void kdtreequeryresultsxi(kdtree kdt,
|
---|
1727 | ref double[,] x)
|
---|
1728 | {
|
---|
1729 | x = new double[0,0];
|
---|
1730 |
|
---|
1731 | kdtreequeryresultsx(kdt, ref x);
|
---|
1732 | }
|
---|
1733 |
|
---|
1734 |
|
---|
1735 | /*************************************************************************
|
---|
1736 | XY-values from last query; 'interactive' variant for languages like Python
|
---|
1737 | which support constructs like "XY = KDTreeQueryResultsXYI(KDT)" and
|
---|
1738 | interactive mode of interpreter.
|
---|
1739 |
|
---|
1740 | This function allocates new array on each call, so it is significantly
|
---|
1741 | slower than its 'non-interactive' counterpart, but it is more convenient
|
---|
1742 | when you call it from command line.
|
---|
1743 |
|
---|
1744 | -- ALGLIB --
|
---|
1745 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
1746 | *************************************************************************/
|
---|
1747 | public static void kdtreequeryresultsxyi(kdtree kdt,
|
---|
1748 | ref double[,] xy)
|
---|
1749 | {
|
---|
1750 | xy = new double[0,0];
|
---|
1751 |
|
---|
1752 | kdtreequeryresultsxy(kdt, ref xy);
|
---|
1753 | }
|
---|
1754 |
|
---|
1755 |
|
---|
1756 | /*************************************************************************
|
---|
1757 | Tags from last query; 'interactive' variant for languages like Python
|
---|
1758 | which support constructs like "Tags = KDTreeQueryResultsTagsI(KDT)" and
|
---|
1759 | interactive mode of interpreter.
|
---|
1760 |
|
---|
1761 | This function allocates new array on each call, so it is significantly
|
---|
1762 | slower than its 'non-interactive' counterpart, but it is more convenient
|
---|
1763 | when you call it from command line.
|
---|
1764 |
|
---|
1765 | -- ALGLIB --
|
---|
1766 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
1767 | *************************************************************************/
|
---|
1768 | public static void kdtreequeryresultstagsi(kdtree kdt,
|
---|
1769 | ref int[] tags)
|
---|
1770 | {
|
---|
1771 | tags = new int[0];
|
---|
1772 |
|
---|
1773 | kdtreequeryresultstags(kdt, ref tags);
|
---|
1774 | }
|
---|
1775 |
|
---|
1776 |
|
---|
1777 | /*************************************************************************
|
---|
1778 | Distances from last query; 'interactive' variant for languages like Python
|
---|
1779 | which support constructs like "R = KDTreeQueryResultsDistancesI(KDT)"
|
---|
1780 | and interactive mode of interpreter.
|
---|
1781 |
|
---|
1782 | This function allocates new array on each call, so it is significantly
|
---|
1783 | slower than its 'non-interactive' counterpart, but it is more convenient
|
---|
1784 | when you call it from command line.
|
---|
1785 |
|
---|
1786 | -- ALGLIB --
|
---|
1787 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
1788 | *************************************************************************/
|
---|
1789 | public static void kdtreequeryresultsdistancesi(kdtree kdt,
|
---|
1790 | ref double[] r)
|
---|
1791 | {
|
---|
1792 | r = new double[0];
|
---|
1793 |
|
---|
1794 | kdtreequeryresultsdistances(kdt, ref r);
|
---|
1795 | }
|
---|
1796 |
|
---|
1797 |
|
---|
1798 | /*************************************************************************
|
---|
1799 | Serializer: allocation
|
---|
1800 |
|
---|
1801 | -- ALGLIB --
|
---|
1802 | Copyright 14.03.2011 by Bochkanov Sergey
|
---|
1803 | *************************************************************************/
|
---|
1804 | public static void kdtreealloc(alglib.serializer s,
|
---|
1805 | kdtree tree)
|
---|
1806 | {
|
---|
1807 |
|
---|
1808 | //
|
---|
1809 | // Header
|
---|
1810 | //
|
---|
1811 | s.alloc_entry();
|
---|
1812 | s.alloc_entry();
|
---|
1813 |
|
---|
1814 | //
|
---|
1815 | // Data
|
---|
1816 | //
|
---|
1817 | s.alloc_entry();
|
---|
1818 | s.alloc_entry();
|
---|
1819 | s.alloc_entry();
|
---|
1820 | s.alloc_entry();
|
---|
1821 | apserv.allocrealmatrix(s, tree.xy, -1, -1);
|
---|
1822 | apserv.allocintegerarray(s, tree.tags, -1);
|
---|
1823 | apserv.allocrealarray(s, tree.boxmin, -1);
|
---|
1824 | apserv.allocrealarray(s, tree.boxmax, -1);
|
---|
1825 | apserv.allocintegerarray(s, tree.nodes, -1);
|
---|
1826 | apserv.allocrealarray(s, tree.splits, -1);
|
---|
1827 | }
|
---|
1828 |
|
---|
1829 |
|
---|
1830 | /*************************************************************************
|
---|
1831 | Serializer: serialization
|
---|
1832 |
|
---|
1833 | -- ALGLIB --
|
---|
1834 | Copyright 14.03.2011 by Bochkanov Sergey
|
---|
1835 | *************************************************************************/
|
---|
1836 | public static void kdtreeserialize(alglib.serializer s,
|
---|
1837 | kdtree tree)
|
---|
1838 | {
|
---|
1839 |
|
---|
1840 | //
|
---|
1841 | // Header
|
---|
1842 | //
|
---|
1843 | s.serialize_int(scodes.getkdtreeserializationcode());
|
---|
1844 | s.serialize_int(kdtreefirstversion);
|
---|
1845 |
|
---|
1846 | //
|
---|
1847 | // Data
|
---|
1848 | //
|
---|
1849 | s.serialize_int(tree.n);
|
---|
1850 | s.serialize_int(tree.nx);
|
---|
1851 | s.serialize_int(tree.ny);
|
---|
1852 | s.serialize_int(tree.normtype);
|
---|
1853 | apserv.serializerealmatrix(s, tree.xy, -1, -1);
|
---|
1854 | apserv.serializeintegerarray(s, tree.tags, -1);
|
---|
1855 | apserv.serializerealarray(s, tree.boxmin, -1);
|
---|
1856 | apserv.serializerealarray(s, tree.boxmax, -1);
|
---|
1857 | apserv.serializeintegerarray(s, tree.nodes, -1);
|
---|
1858 | apserv.serializerealarray(s, tree.splits, -1);
|
---|
1859 | }
|
---|
1860 |
|
---|
1861 |
|
---|
1862 | /*************************************************************************
|
---|
1863 | Serializer: unserialization
|
---|
1864 |
|
---|
1865 | -- ALGLIB --
|
---|
1866 | Copyright 14.03.2011 by Bochkanov Sergey
|
---|
1867 | *************************************************************************/
|
---|
1868 | public static void kdtreeunserialize(alglib.serializer s,
|
---|
1869 | kdtree tree)
|
---|
1870 | {
|
---|
1871 | int i0 = 0;
|
---|
1872 | int i1 = 0;
|
---|
1873 |
|
---|
1874 |
|
---|
1875 | //
|
---|
1876 | // check correctness of header
|
---|
1877 | //
|
---|
1878 | i0 = s.unserialize_int();
|
---|
1879 | ap.assert(i0==scodes.getkdtreeserializationcode(), "KDTreeUnserialize: stream header corrupted");
|
---|
1880 | i1 = s.unserialize_int();
|
---|
1881 | ap.assert(i1==kdtreefirstversion, "KDTreeUnserialize: stream header corrupted");
|
---|
1882 |
|
---|
1883 | //
|
---|
1884 | // Unserialize data
|
---|
1885 | //
|
---|
1886 | tree.n = s.unserialize_int();
|
---|
1887 | tree.nx = s.unserialize_int();
|
---|
1888 | tree.ny = s.unserialize_int();
|
---|
1889 | tree.normtype = s.unserialize_int();
|
---|
1890 | apserv.unserializerealmatrix(s, ref tree.xy);
|
---|
1891 | apserv.unserializeintegerarray(s, ref tree.tags);
|
---|
1892 | apserv.unserializerealarray(s, ref tree.boxmin);
|
---|
1893 | apserv.unserializerealarray(s, ref tree.boxmax);
|
---|
1894 | apserv.unserializeintegerarray(s, ref tree.nodes);
|
---|
1895 | apserv.unserializerealarray(s, ref tree.splits);
|
---|
1896 | kdtreealloctemporaries(tree, tree.n, tree.nx, tree.ny);
|
---|
1897 | }
|
---|
1898 |
|
---|
1899 |
|
---|
1900 | /*************************************************************************
|
---|
1901 | Rearranges nodes [I1,I2) using partition in D-th dimension with S as threshold.
|
---|
1902 | Returns split position I3: [I1,I3) and [I3,I2) are created as result.
|
---|
1903 |
|
---|
1904 | This subroutine doesn't create tree structures, just rearranges nodes.
|
---|
1905 | *************************************************************************/
|
---|
1906 | private static void kdtreesplit(kdtree kdt,
|
---|
1907 | int i1,
|
---|
1908 | int i2,
|
---|
1909 | int d,
|
---|
1910 | double s,
|
---|
1911 | ref int i3)
|
---|
1912 | {
|
---|
1913 | int i = 0;
|
---|
1914 | int j = 0;
|
---|
1915 | int ileft = 0;
|
---|
1916 | int iright = 0;
|
---|
1917 | double v = 0;
|
---|
1918 |
|
---|
1919 | i3 = 0;
|
---|
1920 |
|
---|
1921 |
|
---|
1922 | //
|
---|
1923 | // split XY/Tags in two parts:
|
---|
1924 | // * [ILeft,IRight] is non-processed part of XY/Tags
|
---|
1925 | //
|
---|
1926 | // After cycle is done, we have Ileft=IRight. We deal with
|
---|
1927 | // this element separately.
|
---|
1928 | //
|
---|
1929 | // After this, [I1,ILeft) contains left part, and [ILeft,I2)
|
---|
1930 | // contains right part.
|
---|
1931 | //
|
---|
1932 | ileft = i1;
|
---|
1933 | iright = i2-1;
|
---|
1934 | while( ileft<iright )
|
---|
1935 | {
|
---|
1936 | if( (double)(kdt.xy[ileft,d])<=(double)(s) )
|
---|
1937 | {
|
---|
1938 |
|
---|
1939 | //
|
---|
1940 | // XY[ILeft] is on its place.
|
---|
1941 | // Advance ILeft.
|
---|
1942 | //
|
---|
1943 | ileft = ileft+1;
|
---|
1944 | }
|
---|
1945 | else
|
---|
1946 | {
|
---|
1947 |
|
---|
1948 | //
|
---|
1949 | // XY[ILeft,..] must be at IRight.
|
---|
1950 | // Swap and advance IRight.
|
---|
1951 | //
|
---|
1952 | for(i=0; i<=2*kdt.nx+kdt.ny-1; i++)
|
---|
1953 | {
|
---|
1954 | v = kdt.xy[ileft,i];
|
---|
1955 | kdt.xy[ileft,i] = kdt.xy[iright,i];
|
---|
1956 | kdt.xy[iright,i] = v;
|
---|
1957 | }
|
---|
1958 | j = kdt.tags[ileft];
|
---|
1959 | kdt.tags[ileft] = kdt.tags[iright];
|
---|
1960 | kdt.tags[iright] = j;
|
---|
1961 | iright = iright-1;
|
---|
1962 | }
|
---|
1963 | }
|
---|
1964 | if( (double)(kdt.xy[ileft,d])<=(double)(s) )
|
---|
1965 | {
|
---|
1966 | ileft = ileft+1;
|
---|
1967 | }
|
---|
1968 | else
|
---|
1969 | {
|
---|
1970 | iright = iright-1;
|
---|
1971 | }
|
---|
1972 | i3 = ileft;
|
---|
1973 | }
|
---|
1974 |
|
---|
1975 |
|
---|
1976 | /*************************************************************************
|
---|
1977 | Recursive kd-tree generation subroutine.
|
---|
1978 |
|
---|
1979 | PARAMETERS
|
---|
1980 | KDT tree
|
---|
1981 | NodesOffs unused part of Nodes[] which must be filled by tree
|
---|
1982 | SplitsOffs unused part of Splits[]
|
---|
1983 | I1, I2 points from [I1,I2) are processed
|
---|
1984 |
|
---|
1985 | NodesOffs[] and SplitsOffs[] must be large enough.
|
---|
1986 |
|
---|
1987 | -- ALGLIB --
|
---|
1988 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
1989 | *************************************************************************/
|
---|
1990 | private static void kdtreegeneratetreerec(kdtree kdt,
|
---|
1991 | ref int nodesoffs,
|
---|
1992 | ref int splitsoffs,
|
---|
1993 | int i1,
|
---|
1994 | int i2,
|
---|
1995 | int maxleafsize)
|
---|
1996 | {
|
---|
1997 | int n = 0;
|
---|
1998 | int nx = 0;
|
---|
1999 | int ny = 0;
|
---|
2000 | int i = 0;
|
---|
2001 | int j = 0;
|
---|
2002 | int oldoffs = 0;
|
---|
2003 | int i3 = 0;
|
---|
2004 | int cntless = 0;
|
---|
2005 | int cntgreater = 0;
|
---|
2006 | double minv = 0;
|
---|
2007 | double maxv = 0;
|
---|
2008 | int minidx = 0;
|
---|
2009 | int maxidx = 0;
|
---|
2010 | int d = 0;
|
---|
2011 | double ds = 0;
|
---|
2012 | double s = 0;
|
---|
2013 | double v = 0;
|
---|
2014 | int i_ = 0;
|
---|
2015 | int i1_ = 0;
|
---|
2016 |
|
---|
2017 | ap.assert(i2>i1, "KDTreeGenerateTreeRec: internal error");
|
---|
2018 |
|
---|
2019 | //
|
---|
2020 | // Generate leaf if needed
|
---|
2021 | //
|
---|
2022 | if( i2-i1<=maxleafsize )
|
---|
2023 | {
|
---|
2024 | kdt.nodes[nodesoffs+0] = i2-i1;
|
---|
2025 | kdt.nodes[nodesoffs+1] = i1;
|
---|
2026 | nodesoffs = nodesoffs+2;
|
---|
2027 | return;
|
---|
2028 | }
|
---|
2029 |
|
---|
2030 | //
|
---|
2031 | // Load values for easier access
|
---|
2032 | //
|
---|
2033 | nx = kdt.nx;
|
---|
2034 | ny = kdt.ny;
|
---|
2035 |
|
---|
2036 | //
|
---|
2037 | // select dimension to split:
|
---|
2038 | // * D is a dimension number
|
---|
2039 | //
|
---|
2040 | d = 0;
|
---|
2041 | ds = kdt.curboxmax[0]-kdt.curboxmin[0];
|
---|
2042 | for(i=1; i<=nx-1; i++)
|
---|
2043 | {
|
---|
2044 | v = kdt.curboxmax[i]-kdt.curboxmin[i];
|
---|
2045 | if( (double)(v)>(double)(ds) )
|
---|
2046 | {
|
---|
2047 | ds = v;
|
---|
2048 | d = i;
|
---|
2049 | }
|
---|
2050 | }
|
---|
2051 |
|
---|
2052 | //
|
---|
2053 | // Select split position S using sliding midpoint rule,
|
---|
2054 | // rearrange points into [I1,I3) and [I3,I2)
|
---|
2055 | //
|
---|
2056 | s = kdt.curboxmin[d]+0.5*ds;
|
---|
2057 | i1_ = (i1) - (0);
|
---|
2058 | for(i_=0; i_<=i2-i1-1;i_++)
|
---|
2059 | {
|
---|
2060 | kdt.buf[i_] = kdt.xy[i_+i1_,d];
|
---|
2061 | }
|
---|
2062 | n = i2-i1;
|
---|
2063 | cntless = 0;
|
---|
2064 | cntgreater = 0;
|
---|
2065 | minv = kdt.buf[0];
|
---|
2066 | maxv = kdt.buf[0];
|
---|
2067 | minidx = i1;
|
---|
2068 | maxidx = i1;
|
---|
2069 | for(i=0; i<=n-1; i++)
|
---|
2070 | {
|
---|
2071 | v = kdt.buf[i];
|
---|
2072 | if( (double)(v)<(double)(minv) )
|
---|
2073 | {
|
---|
2074 | minv = v;
|
---|
2075 | minidx = i1+i;
|
---|
2076 | }
|
---|
2077 | if( (double)(v)>(double)(maxv) )
|
---|
2078 | {
|
---|
2079 | maxv = v;
|
---|
2080 | maxidx = i1+i;
|
---|
2081 | }
|
---|
2082 | if( (double)(v)<(double)(s) )
|
---|
2083 | {
|
---|
2084 | cntless = cntless+1;
|
---|
2085 | }
|
---|
2086 | if( (double)(v)>(double)(s) )
|
---|
2087 | {
|
---|
2088 | cntgreater = cntgreater+1;
|
---|
2089 | }
|
---|
2090 | }
|
---|
2091 | if( cntless>0 & cntgreater>0 )
|
---|
2092 | {
|
---|
2093 |
|
---|
2094 | //
|
---|
2095 | // normal midpoint split
|
---|
2096 | //
|
---|
2097 | kdtreesplit(kdt, i1, i2, d, s, ref i3);
|
---|
2098 | }
|
---|
2099 | else
|
---|
2100 | {
|
---|
2101 |
|
---|
2102 | //
|
---|
2103 | // sliding midpoint
|
---|
2104 | //
|
---|
2105 | if( cntless==0 )
|
---|
2106 | {
|
---|
2107 |
|
---|
2108 | //
|
---|
2109 | // 1. move split to MinV,
|
---|
2110 | // 2. place one point to the left bin (move to I1),
|
---|
2111 | // others - to the right bin
|
---|
2112 | //
|
---|
2113 | s = minv;
|
---|
2114 | if( minidx!=i1 )
|
---|
2115 | {
|
---|
2116 | for(i=0; i<=2*kdt.nx+kdt.ny-1; i++)
|
---|
2117 | {
|
---|
2118 | v = kdt.xy[minidx,i];
|
---|
2119 | kdt.xy[minidx,i] = kdt.xy[i1,i];
|
---|
2120 | kdt.xy[i1,i] = v;
|
---|
2121 | }
|
---|
2122 | j = kdt.tags[minidx];
|
---|
2123 | kdt.tags[minidx] = kdt.tags[i1];
|
---|
2124 | kdt.tags[i1] = j;
|
---|
2125 | }
|
---|
2126 | i3 = i1+1;
|
---|
2127 | }
|
---|
2128 | else
|
---|
2129 | {
|
---|
2130 |
|
---|
2131 | //
|
---|
2132 | // 1. move split to MaxV,
|
---|
2133 | // 2. place one point to the right bin (move to I2-1),
|
---|
2134 | // others - to the left bin
|
---|
2135 | //
|
---|
2136 | s = maxv;
|
---|
2137 | if( maxidx!=i2-1 )
|
---|
2138 | {
|
---|
2139 | for(i=0; i<=2*kdt.nx+kdt.ny-1; i++)
|
---|
2140 | {
|
---|
2141 | v = kdt.xy[maxidx,i];
|
---|
2142 | kdt.xy[maxidx,i] = kdt.xy[i2-1,i];
|
---|
2143 | kdt.xy[i2-1,i] = v;
|
---|
2144 | }
|
---|
2145 | j = kdt.tags[maxidx];
|
---|
2146 | kdt.tags[maxidx] = kdt.tags[i2-1];
|
---|
2147 | kdt.tags[i2-1] = j;
|
---|
2148 | }
|
---|
2149 | i3 = i2-1;
|
---|
2150 | }
|
---|
2151 | }
|
---|
2152 |
|
---|
2153 | //
|
---|
2154 | // Generate 'split' node
|
---|
2155 | //
|
---|
2156 | kdt.nodes[nodesoffs+0] = 0;
|
---|
2157 | kdt.nodes[nodesoffs+1] = d;
|
---|
2158 | kdt.nodes[nodesoffs+2] = splitsoffs;
|
---|
2159 | kdt.splits[splitsoffs+0] = s;
|
---|
2160 | oldoffs = nodesoffs;
|
---|
2161 | nodesoffs = nodesoffs+splitnodesize;
|
---|
2162 | splitsoffs = splitsoffs+1;
|
---|
2163 |
|
---|
2164 | //
|
---|
2165 | // Recirsive generation:
|
---|
2166 | // * update CurBox
|
---|
2167 | // * call subroutine
|
---|
2168 | // * restore CurBox
|
---|
2169 | //
|
---|
2170 | kdt.nodes[oldoffs+3] = nodesoffs;
|
---|
2171 | v = kdt.curboxmax[d];
|
---|
2172 | kdt.curboxmax[d] = s;
|
---|
2173 | kdtreegeneratetreerec(kdt, ref nodesoffs, ref splitsoffs, i1, i3, maxleafsize);
|
---|
2174 | kdt.curboxmax[d] = v;
|
---|
2175 | kdt.nodes[oldoffs+4] = nodesoffs;
|
---|
2176 | v = kdt.curboxmin[d];
|
---|
2177 | kdt.curboxmin[d] = s;
|
---|
2178 | kdtreegeneratetreerec(kdt, ref nodesoffs, ref splitsoffs, i3, i2, maxleafsize);
|
---|
2179 | kdt.curboxmin[d] = v;
|
---|
2180 | }
|
---|
2181 |
|
---|
2182 |
|
---|
2183 | /*************************************************************************
|
---|
2184 | Recursive subroutine for NN queries.
|
---|
2185 |
|
---|
2186 | -- ALGLIB --
|
---|
2187 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
2188 | *************************************************************************/
|
---|
2189 | private static void kdtreequerynnrec(kdtree kdt,
|
---|
2190 | int offs)
|
---|
2191 | {
|
---|
2192 | double ptdist = 0;
|
---|
2193 | int i = 0;
|
---|
2194 | int j = 0;
|
---|
2195 | int nx = 0;
|
---|
2196 | int i1 = 0;
|
---|
2197 | int i2 = 0;
|
---|
2198 | int d = 0;
|
---|
2199 | double s = 0;
|
---|
2200 | double v = 0;
|
---|
2201 | double t1 = 0;
|
---|
2202 | int childbestoffs = 0;
|
---|
2203 | int childworstoffs = 0;
|
---|
2204 | int childoffs = 0;
|
---|
2205 | double prevdist = 0;
|
---|
2206 | bool todive = new bool();
|
---|
2207 | bool bestisleft = new bool();
|
---|
2208 | bool updatemin = new bool();
|
---|
2209 |
|
---|
2210 |
|
---|
2211 | //
|
---|
2212 | // Leaf node.
|
---|
2213 | // Process points.
|
---|
2214 | //
|
---|
2215 | if( kdt.nodes[offs]>0 )
|
---|
2216 | {
|
---|
2217 | i1 = kdt.nodes[offs+1];
|
---|
2218 | i2 = i1+kdt.nodes[offs];
|
---|
2219 | for(i=i1; i<=i2-1; i++)
|
---|
2220 | {
|
---|
2221 |
|
---|
2222 | //
|
---|
2223 | // Calculate distance
|
---|
2224 | //
|
---|
2225 | ptdist = 0;
|
---|
2226 | nx = kdt.nx;
|
---|
2227 | if( kdt.normtype==0 )
|
---|
2228 | {
|
---|
2229 | for(j=0; j<=nx-1; j++)
|
---|
2230 | {
|
---|
2231 | ptdist = Math.Max(ptdist, Math.Abs(kdt.xy[i,j]-kdt.x[j]));
|
---|
2232 | }
|
---|
2233 | }
|
---|
2234 | if( kdt.normtype==1 )
|
---|
2235 | {
|
---|
2236 | for(j=0; j<=nx-1; j++)
|
---|
2237 | {
|
---|
2238 | ptdist = ptdist+Math.Abs(kdt.xy[i,j]-kdt.x[j]);
|
---|
2239 | }
|
---|
2240 | }
|
---|
2241 | if( kdt.normtype==2 )
|
---|
2242 | {
|
---|
2243 | for(j=0; j<=nx-1; j++)
|
---|
2244 | {
|
---|
2245 | ptdist = ptdist+math.sqr(kdt.xy[i,j]-kdt.x[j]);
|
---|
2246 | }
|
---|
2247 | }
|
---|
2248 |
|
---|
2249 | //
|
---|
2250 | // Skip points with zero distance if self-matches are turned off
|
---|
2251 | //
|
---|
2252 | if( (double)(ptdist)==(double)(0) & !kdt.selfmatch )
|
---|
2253 | {
|
---|
2254 | continue;
|
---|
2255 | }
|
---|
2256 |
|
---|
2257 | //
|
---|
2258 | // We CAN'T process point if R-criterion isn't satisfied,
|
---|
2259 | // i.e. (RNeeded<>0) AND (PtDist>R).
|
---|
2260 | //
|
---|
2261 | if( (double)(kdt.rneeded)==(double)(0) | (double)(ptdist)<=(double)(kdt.rneeded) )
|
---|
2262 | {
|
---|
2263 |
|
---|
2264 | //
|
---|
2265 | // R-criterion is satisfied, we must either:
|
---|
2266 | // * replace worst point, if (KNeeded<>0) AND (KCur=KNeeded)
|
---|
2267 | // (or skip, if worst point is better)
|
---|
2268 | // * add point without replacement otherwise
|
---|
2269 | //
|
---|
2270 | if( kdt.kcur<kdt.kneeded | kdt.kneeded==0 )
|
---|
2271 | {
|
---|
2272 |
|
---|
2273 | //
|
---|
2274 | // add current point to heap without replacement
|
---|
2275 | //
|
---|
2276 | tsort.tagheappushi(ref kdt.r, ref kdt.idx, ref kdt.kcur, ptdist, i);
|
---|
2277 | }
|
---|
2278 | else
|
---|
2279 | {
|
---|
2280 |
|
---|
2281 | //
|
---|
2282 | // New points are added or not, depending on their distance.
|
---|
2283 | // If added, they replace element at the top of the heap
|
---|
2284 | //
|
---|
2285 | if( (double)(ptdist)<(double)(kdt.r[0]) )
|
---|
2286 | {
|
---|
2287 | if( kdt.kneeded==1 )
|
---|
2288 | {
|
---|
2289 | kdt.idx[0] = i;
|
---|
2290 | kdt.r[0] = ptdist;
|
---|
2291 | }
|
---|
2292 | else
|
---|
2293 | {
|
---|
2294 | tsort.tagheapreplacetopi(ref kdt.r, ref kdt.idx, kdt.kneeded, ptdist, i);
|
---|
2295 | }
|
---|
2296 | }
|
---|
2297 | }
|
---|
2298 | }
|
---|
2299 | }
|
---|
2300 | return;
|
---|
2301 | }
|
---|
2302 |
|
---|
2303 | //
|
---|
2304 | // Simple split
|
---|
2305 | //
|
---|
2306 | if( kdt.nodes[offs]==0 )
|
---|
2307 | {
|
---|
2308 |
|
---|
2309 | //
|
---|
2310 | // Load:
|
---|
2311 | // * D dimension to split
|
---|
2312 | // * S split position
|
---|
2313 | //
|
---|
2314 | d = kdt.nodes[offs+1];
|
---|
2315 | s = kdt.splits[kdt.nodes[offs+2]];
|
---|
2316 |
|
---|
2317 | //
|
---|
2318 | // Calculate:
|
---|
2319 | // * ChildBestOffs child box with best chances
|
---|
2320 | // * ChildWorstOffs child box with worst chances
|
---|
2321 | //
|
---|
2322 | if( (double)(kdt.x[d])<=(double)(s) )
|
---|
2323 | {
|
---|
2324 | childbestoffs = kdt.nodes[offs+3];
|
---|
2325 | childworstoffs = kdt.nodes[offs+4];
|
---|
2326 | bestisleft = true;
|
---|
2327 | }
|
---|
2328 | else
|
---|
2329 | {
|
---|
2330 | childbestoffs = kdt.nodes[offs+4];
|
---|
2331 | childworstoffs = kdt.nodes[offs+3];
|
---|
2332 | bestisleft = false;
|
---|
2333 | }
|
---|
2334 |
|
---|
2335 | //
|
---|
2336 | // Navigate through childs
|
---|
2337 | //
|
---|
2338 | for(i=0; i<=1; i++)
|
---|
2339 | {
|
---|
2340 |
|
---|
2341 | //
|
---|
2342 | // Select child to process:
|
---|
2343 | // * ChildOffs current child offset in Nodes[]
|
---|
2344 | // * UpdateMin whether minimum or maximum value
|
---|
2345 | // of bounding box is changed on update
|
---|
2346 | //
|
---|
2347 | if( i==0 )
|
---|
2348 | {
|
---|
2349 | childoffs = childbestoffs;
|
---|
2350 | updatemin = !bestisleft;
|
---|
2351 | }
|
---|
2352 | else
|
---|
2353 | {
|
---|
2354 | updatemin = bestisleft;
|
---|
2355 | childoffs = childworstoffs;
|
---|
2356 | }
|
---|
2357 |
|
---|
2358 | //
|
---|
2359 | // Update bounding box and current distance
|
---|
2360 | //
|
---|
2361 | if( updatemin )
|
---|
2362 | {
|
---|
2363 | prevdist = kdt.curdist;
|
---|
2364 | t1 = kdt.x[d];
|
---|
2365 | v = kdt.curboxmin[d];
|
---|
2366 | if( (double)(t1)<=(double)(s) )
|
---|
2367 | {
|
---|
2368 | if( kdt.normtype==0 )
|
---|
2369 | {
|
---|
2370 | kdt.curdist = Math.Max(kdt.curdist, s-t1);
|
---|
2371 | }
|
---|
2372 | if( kdt.normtype==1 )
|
---|
2373 | {
|
---|
2374 | kdt.curdist = kdt.curdist-Math.Max(v-t1, 0)+s-t1;
|
---|
2375 | }
|
---|
2376 | if( kdt.normtype==2 )
|
---|
2377 | {
|
---|
2378 | kdt.curdist = kdt.curdist-math.sqr(Math.Max(v-t1, 0))+math.sqr(s-t1);
|
---|
2379 | }
|
---|
2380 | }
|
---|
2381 | kdt.curboxmin[d] = s;
|
---|
2382 | }
|
---|
2383 | else
|
---|
2384 | {
|
---|
2385 | prevdist = kdt.curdist;
|
---|
2386 | t1 = kdt.x[d];
|
---|
2387 | v = kdt.curboxmax[d];
|
---|
2388 | if( (double)(t1)>=(double)(s) )
|
---|
2389 | {
|
---|
2390 | if( kdt.normtype==0 )
|
---|
2391 | {
|
---|
2392 | kdt.curdist = Math.Max(kdt.curdist, t1-s);
|
---|
2393 | }
|
---|
2394 | if( kdt.normtype==1 )
|
---|
2395 | {
|
---|
2396 | kdt.curdist = kdt.curdist-Math.Max(t1-v, 0)+t1-s;
|
---|
2397 | }
|
---|
2398 | if( kdt.normtype==2 )
|
---|
2399 | {
|
---|
2400 | kdt.curdist = kdt.curdist-math.sqr(Math.Max(t1-v, 0))+math.sqr(t1-s);
|
---|
2401 | }
|
---|
2402 | }
|
---|
2403 | kdt.curboxmax[d] = s;
|
---|
2404 | }
|
---|
2405 |
|
---|
2406 | //
|
---|
2407 | // Decide: to dive into cell or not to dive
|
---|
2408 | //
|
---|
2409 | if( (double)(kdt.rneeded)!=(double)(0) & (double)(kdt.curdist)>(double)(kdt.rneeded) )
|
---|
2410 | {
|
---|
2411 | todive = false;
|
---|
2412 | }
|
---|
2413 | else
|
---|
2414 | {
|
---|
2415 | if( kdt.kcur<kdt.kneeded | kdt.kneeded==0 )
|
---|
2416 | {
|
---|
2417 |
|
---|
2418 | //
|
---|
2419 | // KCur<KNeeded (i.e. not all points are found)
|
---|
2420 | //
|
---|
2421 | todive = true;
|
---|
2422 | }
|
---|
2423 | else
|
---|
2424 | {
|
---|
2425 |
|
---|
2426 | //
|
---|
2427 | // KCur=KNeeded, decide to dive or not to dive
|
---|
2428 | // using point position relative to bounding box.
|
---|
2429 | //
|
---|
2430 | todive = (double)(kdt.curdist)<=(double)(kdt.r[0]*kdt.approxf);
|
---|
2431 | }
|
---|
2432 | }
|
---|
2433 | if( todive )
|
---|
2434 | {
|
---|
2435 | kdtreequerynnrec(kdt, childoffs);
|
---|
2436 | }
|
---|
2437 |
|
---|
2438 | //
|
---|
2439 | // Restore bounding box and distance
|
---|
2440 | //
|
---|
2441 | if( updatemin )
|
---|
2442 | {
|
---|
2443 | kdt.curboxmin[d] = v;
|
---|
2444 | }
|
---|
2445 | else
|
---|
2446 | {
|
---|
2447 | kdt.curboxmax[d] = v;
|
---|
2448 | }
|
---|
2449 | kdt.curdist = prevdist;
|
---|
2450 | }
|
---|
2451 | return;
|
---|
2452 | }
|
---|
2453 | }
|
---|
2454 |
|
---|
2455 |
|
---|
2456 | /*************************************************************************
|
---|
2457 | Copies X[] to KDT.X[]
|
---|
2458 | Loads distance from X[] to bounding box.
|
---|
2459 | Initializes CurBox[].
|
---|
2460 |
|
---|
2461 | -- ALGLIB --
|
---|
2462 | Copyright 28.02.2010 by Bochkanov Sergey
|
---|
2463 | *************************************************************************/
|
---|
2464 | private static void kdtreeinitbox(kdtree kdt,
|
---|
2465 | double[] x)
|
---|
2466 | {
|
---|
2467 | int i = 0;
|
---|
2468 | double vx = 0;
|
---|
2469 | double vmin = 0;
|
---|
2470 | double vmax = 0;
|
---|
2471 |
|
---|
2472 |
|
---|
2473 | //
|
---|
2474 | // calculate distance from point to current bounding box
|
---|
2475 | //
|
---|
2476 | kdt.curdist = 0;
|
---|
2477 | if( kdt.normtype==0 )
|
---|
2478 | {
|
---|
2479 | for(i=0; i<=kdt.nx-1; i++)
|
---|
2480 | {
|
---|
2481 | vx = x[i];
|
---|
2482 | vmin = kdt.boxmin[i];
|
---|
2483 | vmax = kdt.boxmax[i];
|
---|
2484 | kdt.x[i] = vx;
|
---|
2485 | kdt.curboxmin[i] = vmin;
|
---|
2486 | kdt.curboxmax[i] = vmax;
|
---|
2487 | if( (double)(vx)<(double)(vmin) )
|
---|
2488 | {
|
---|
2489 | kdt.curdist = Math.Max(kdt.curdist, vmin-vx);
|
---|
2490 | }
|
---|
2491 | else
|
---|
2492 | {
|
---|
2493 | if( (double)(vx)>(double)(vmax) )
|
---|
2494 | {
|
---|
2495 | kdt.curdist = Math.Max(kdt.curdist, vx-vmax);
|
---|
2496 | }
|
---|
2497 | }
|
---|
2498 | }
|
---|
2499 | }
|
---|
2500 | if( kdt.normtype==1 )
|
---|
2501 | {
|
---|
2502 | for(i=0; i<=kdt.nx-1; i++)
|
---|
2503 | {
|
---|
2504 | vx = x[i];
|
---|
2505 | vmin = kdt.boxmin[i];
|
---|
2506 | vmax = kdt.boxmax[i];
|
---|
2507 | kdt.x[i] = vx;
|
---|
2508 | kdt.curboxmin[i] = vmin;
|
---|
2509 | kdt.curboxmax[i] = vmax;
|
---|
2510 | if( (double)(vx)<(double)(vmin) )
|
---|
2511 | {
|
---|
2512 | kdt.curdist = kdt.curdist+vmin-vx;
|
---|
2513 | }
|
---|
2514 | else
|
---|
2515 | {
|
---|
2516 | if( (double)(vx)>(double)(vmax) )
|
---|
2517 | {
|
---|
2518 | kdt.curdist = kdt.curdist+vx-vmax;
|
---|
2519 | }
|
---|
2520 | }
|
---|
2521 | }
|
---|
2522 | }
|
---|
2523 | if( kdt.normtype==2 )
|
---|
2524 | {
|
---|
2525 | for(i=0; i<=kdt.nx-1; i++)
|
---|
2526 | {
|
---|
2527 | vx = x[i];
|
---|
2528 | vmin = kdt.boxmin[i];
|
---|
2529 | vmax = kdt.boxmax[i];
|
---|
2530 | kdt.x[i] = vx;
|
---|
2531 | kdt.curboxmin[i] = vmin;
|
---|
2532 | kdt.curboxmax[i] = vmax;
|
---|
2533 | if( (double)(vx)<(double)(vmin) )
|
---|
2534 | {
|
---|
2535 | kdt.curdist = kdt.curdist+math.sqr(vmin-vx);
|
---|
2536 | }
|
---|
2537 | else
|
---|
2538 | {
|
---|
2539 | if( (double)(vx)>(double)(vmax) )
|
---|
2540 | {
|
---|
2541 | kdt.curdist = kdt.curdist+math.sqr(vx-vmax);
|
---|
2542 | }
|
---|
2543 | }
|
---|
2544 | }
|
---|
2545 | }
|
---|
2546 | }
|
---|
2547 |
|
---|
2548 |
|
---|
2549 | /*************************************************************************
|
---|
2550 | This function allocates all dataset-independent array fields of KDTree,
|
---|
2551 | i.e. such array fields that their dimensions do not depend on dataset
|
---|
2552 | size.
|
---|
2553 |
|
---|
2554 | This function do not sets KDT.NX or KDT.NY - it just allocates arrays
|
---|
2555 |
|
---|
2556 | -- ALGLIB --
|
---|
2557 | Copyright 14.03.2011 by Bochkanov Sergey
|
---|
2558 | *************************************************************************/
|
---|
2559 | private static void kdtreeallocdatasetindependent(kdtree kdt,
|
---|
2560 | int nx,
|
---|
2561 | int ny)
|
---|
2562 | {
|
---|
2563 | kdt.x = new double[nx];
|
---|
2564 | kdt.boxmin = new double[nx];
|
---|
2565 | kdt.boxmax = new double[nx];
|
---|
2566 | kdt.curboxmin = new double[nx];
|
---|
2567 | kdt.curboxmax = new double[nx];
|
---|
2568 | }
|
---|
2569 |
|
---|
2570 |
|
---|
2571 | /*************************************************************************
|
---|
2572 | This function allocates all dataset-dependent array fields of KDTree, i.e.
|
---|
2573 | such array fields that their dimensions depend on dataset size.
|
---|
2574 |
|
---|
2575 | This function do not sets KDT.N, KDT.NX or KDT.NY -
|
---|
2576 | it just allocates arrays.
|
---|
2577 |
|
---|
2578 | -- ALGLIB --
|
---|
2579 | Copyright 14.03.2011 by Bochkanov Sergey
|
---|
2580 | *************************************************************************/
|
---|
2581 | private static void kdtreeallocdatasetdependent(kdtree kdt,
|
---|
2582 | int n,
|
---|
2583 | int nx,
|
---|
2584 | int ny)
|
---|
2585 | {
|
---|
2586 | kdt.xy = new double[n, 2*nx+ny];
|
---|
2587 | kdt.tags = new int[n];
|
---|
2588 | kdt.idx = new int[n];
|
---|
2589 | kdt.r = new double[n];
|
---|
2590 | kdt.x = new double[nx];
|
---|
2591 | kdt.buf = new double[Math.Max(n, nx)];
|
---|
2592 | kdt.nodes = new int[splitnodesize*2*n];
|
---|
2593 | kdt.splits = new double[2*n];
|
---|
2594 | }
|
---|
2595 |
|
---|
2596 |
|
---|
2597 | /*************************************************************************
|
---|
2598 | This function allocates temporaries.
|
---|
2599 |
|
---|
2600 | This function do not sets KDT.N, KDT.NX or KDT.NY -
|
---|
2601 | it just allocates arrays.
|
---|
2602 |
|
---|
2603 | -- ALGLIB --
|
---|
2604 | Copyright 14.03.2011 by Bochkanov Sergey
|
---|
2605 | *************************************************************************/
|
---|
2606 | private static void kdtreealloctemporaries(kdtree kdt,
|
---|
2607 | int n,
|
---|
2608 | int nx,
|
---|
2609 | int ny)
|
---|
2610 | {
|
---|
2611 | kdt.x = new double[nx];
|
---|
2612 | kdt.idx = new int[n];
|
---|
2613 | kdt.r = new double[n];
|
---|
2614 | kdt.buf = new double[Math.Max(n, nx)];
|
---|
2615 | kdt.curboxmin = new double[nx];
|
---|
2616 | kdt.curboxmax = new double[nx];
|
---|
2617 | }
|
---|
2618 |
|
---|
2619 |
|
---|
2620 | }
|
---|
2621 | }
|
---|
2622 |
|
---|