1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
|
// Gmsh - Copyright (C) 1997-2021 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// issues on https://gitlab.onelab.info/gmsh/gmsh/issues.
#include "GmshConfig.h"
#include "NearestNeighbor.h"
#if defined(HAVE_ANN)
#include "ANN/ANN.h"
#endif
StringXNumber NearestNeighborOptions_Number[] = {
{GMSH_FULLRC, "View", nullptr, -1.},
};
extern "C" {
GMSH_Plugin *GMSH_RegisterNearestNeighborPlugin()
{
return new GMSH_NearestNeighborPlugin();
}
}
std::string GMSH_NearestNeighborPlugin::getHelp() const
{
return "Plugin(NearestNeighbor) computes the distance from each "
"point in `View' to its nearest neighbor.\n\n"
"If `View' < 0, the plugin is run on the current view.\n\n"
"Plugin(NearestNeighbor) is executed in-place.";
}
int GMSH_NearestNeighborPlugin::getNbOptions() const
{
return sizeof(NearestNeighborOptions_Number) / sizeof(StringXNumber);
}
StringXNumber *GMSH_NearestNeighborPlugin::getOption(int iopt)
{
return &NearestNeighborOptions_Number[iopt];
}
PView *GMSH_NearestNeighborPlugin::execute(PView *v)
{
int iView = (int)NearestNeighborOptions_Number[0].def;
PView *v1 = getView(iView, v);
if(!v1) return v;
PViewData *data1 = v1->getData();
int totpoints = data1->getNumPoints();
if(!totpoints) {
Msg::Error("View[%d] contains no points", iView);
return nullptr;
}
#if defined(HAVE_ANN)
ANNpointArray zeronodes = annAllocPts(totpoints, 3);
int k = 0, step = 0;
for(int ent = 0; ent < data1->getNumEntities(step); ent++) {
for(int ele = 0; ele < data1->getNumElements(step, ent); ele++) {
if(data1->skipElement(step, ent, ele)) continue;
int numNodes = data1->getNumNodes(step, ent, ele);
if(numNodes != 1) continue;
data1->getNode(step, ent, ele, 0, zeronodes[k][0], zeronodes[k][1],
zeronodes[k][2]);
k++;
}
}
ANNkd_tree *kdtree = new ANNkd_tree(zeronodes, totpoints, 3);
ANNidxArray index = new ANNidx[2];
ANNdistArray dist = new ANNdist[2];
v1->setChanged(true);
for(int ent = 0; ent < data1->getNumEntities(step); ent++) {
for(int ele = 0; ele < data1->getNumElements(step, ent); ele++) {
if(data1->skipElement(step, ent, ele)) continue;
int numNodes = data1->getNumNodes(step, ent, ele);
if(numNodes != 1) continue;
double xyz[3];
data1->getNode(step, ent, ele, 0, xyz[0], xyz[1], xyz[2]);
kdtree->annkSearch(xyz, 2, index, dist);
data1->setValue(step, ent, ele, 0, 0, sqrt(dist[1]));
}
}
delete kdtree;
annDeallocPts(zeronodes);
delete[] index;
delete[] dist;
#else
Msg::Error("Nearest neighbor computation requires ANN");
#endif
data1->setName(v1->getData()->getName() + "_NearestNeighbor");
data1->finalize();
return v1;
}
|
