Skip to content

Commit baa0dd0

Browse files
zeuxzeux
authored
Merge pull request zeux#889 from zeux/spatial-points
spatialorder: Introduce experimental point clusterizer
2 parents 756c652 + 6fc5ccb commit baa0dd0

3 files changed

Lines changed: 208 additions & 32 deletions

File tree

‎demo/main.cpp‎

Lines changed: 30 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1237,6 +1237,35 @@ void spatialSortTriangles(const Mesh& mesh)
12371237
(end - start) * 1000);
12381238
}
12391239

1240+
void spatialClusterPoints(const Mesh& mesh, size_t cluster_size)
1241+
{
1242+
typedef PackedVertexOct PV;
1243+
1244+
double start = timestamp();
1245+
1246+
std::vector<unsigned int> index(mesh.vertices.size());
1247+
meshopt_spatialClusterPoints(&index[0], &mesh.vertices[0].px, mesh.vertices.size(), sizeof(Vertex), cluster_size);
1248+
1249+
double end = timestamp();
1250+
1251+
std::vector<PV> pv(mesh.vertices.size());
1252+
packMesh(pv, mesh.vertices);
1253+
1254+
std::vector<PV> pvo(mesh.vertices.size());
1255+
for (size_t i = 0; i < index.size(); ++i)
1256+
pvo[i] = pv[index[i]];
1257+
1258+
std::vector<unsigned char> vbuf(meshopt_encodeVertexBufferBound(mesh.vertices.size(), sizeof(PV)));
1259+
vbuf.resize(meshopt_encodeVertexBuffer(&vbuf[0], vbuf.size(), &pvo[0], mesh.vertices.size(), sizeof(PV)));
1260+
1261+
size_t csize = compress(vbuf);
1262+
1263+
printf("SpatialCP: %.1f bits/vertex (post-deflate %.1f bits/vertex); sort %.2f msec\n",
1264+
double(vbuf.size() * 8) / double(mesh.vertices.size()),
1265+
double(csize * 8) / double(mesh.vertices.size()),
1266+
(end - start) * 1000);
1267+
}
1268+
12401269
void tessellationAdjacency(const Mesh& mesh)
12411270
{
12421271
double start = timestamp();
@@ -1525,6 +1554,7 @@ void process(const char* path)
15251554

15261555
spatialSort(mesh);
15271556
spatialSortTriangles(mesh);
1557+
spatialClusterPoints(mesh, 64);
15281558

15291559
reindexFuzzy(mesh);
15301560
coverage(mesh);

‎src/meshoptimizer.h‎

Lines changed: 10 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -714,6 +714,16 @@ MESHOPTIMIZER_API void meshopt_spatialSortRemap(unsigned int* destination, const
714714
*/
715715
MESHOPTIMIZER_EXPERIMENTAL void meshopt_spatialSortTriangles(unsigned int* destination, const unsigned int* indices, size_t index_count, const float* vertex_positions, size_t vertex_count, size_t vertex_positions_stride);
716716

717+
/**
718+
* Experimental: Spatial clusterizer
719+
* Reorders points into clusters optimized for spatial locality, and generates a new index buffer.
720+
* Ensures the output can be split into cluster_size chunks where each chunk has good positional locality. Only the last chunk will be smaller than cluster_size.
721+
*
722+
* destination must contain enough space for the resulting index buffer (vertex_count elements)
723+
* vertex_positions should have float3 position in the first 12 bytes of each vertex
724+
*/
725+
MESHOPTIMIZER_EXPERIMENTAL void meshopt_spatialClusterPoints(unsigned int* destination, const float* vertex_positions, size_t vertex_count, size_t vertex_positions_stride, size_t cluster_size);
726+
717727
/**
718728
* Quantize a float into half-precision (as defined by IEEE-754 fp16) floating point value
719729
* Generates +-inf for overflow, preserves NaN, flushes denormals to zero, rounds to nearest

‎src/spatialorder.cpp‎

Lines changed: 168 additions & 32 deletions
Original file line numberDiff line numberDiff line change
@@ -22,7 +22,7 @@ inline unsigned long long part1By2(unsigned long long x)
2222
return x;
2323
}
2424

25-
static void computeOrder(unsigned long long* result, const float* vertex_positions_data, size_t vertex_count, size_t vertex_positions_stride)
25+
static void computeOrder(unsigned long long* result, const float* vertex_positions_data, size_t vertex_count, size_t vertex_positions_stride, bool morton)
2626
{
2727
size_t vertex_stride_float = vertex_positions_stride / sizeof(float);
2828

@@ -60,61 +60,158 @@ static void computeOrder(unsigned long long* result, const float* vertex_positio
6060
int y = int((v[1] - minv[1]) * scale + 0.5f);
6161
int z = int((v[2] - minv[2]) * scale + 0.5f);
6262

63-
result[i] = part1By2(x) | (part1By2(y) << 1) | (part1By2(z) << 2);
63+
if (morton)
64+
result[i] = part1By2(x) | (part1By2(y) << 1) | (part1By2(z) << 2);
65+
else
66+
result[i] = ((unsigned long long)x << 0) | ((unsigned long long)y << 20) | ((unsigned long long)z << 40);
6467
}
6568
}
6669

67-
static void computeHistogram(unsigned int (&hist)[1024][5], const unsigned long long* data, size_t count)
70+
static void radixSort10(unsigned int* destination, const unsigned int* source, const unsigned short* keys, size_t count)
6871
{
72+
unsigned int hist[1024];
6973
memset(hist, 0, sizeof(hist));
7074

71-
// compute 5 10-bit histograms in parallel
75+
// compute histogram (assume keys are 10-bit)
76+
for (size_t i = 0; i < count; ++i)
77+
hist[keys[i]]++;
78+
79+
unsigned int sum = 0;
80+
81+
// replace histogram data with prefix histogram sums in-place
82+
for (int i = 0; i < 1024; ++i)
83+
{
84+
unsigned int h = hist[i];
85+
hist[i] = sum;
86+
sum += h;
87+
}
88+
89+
assert(sum == count);
90+
91+
// reorder values
92+
for (size_t i = 0; i < count; ++i)
93+
{
94+
unsigned int id = keys[source[i]];
95+
96+
destination[hist[id]++] = source[i];
97+
}
98+
}
99+
100+
static void computeHistogram(unsigned int (&hist)[256][2], const unsigned short* data, size_t count)
101+
{
102+
memset(hist, 0, sizeof(hist));
103+
104+
// compute 2 8-bit histograms in parallel
72105
for (size_t i = 0; i < count; ++i)
73106
{
74107
unsigned long long id = data[i];
75108

76-
hist[(id >> 0) & 1023][0]++;
77-
hist[(id >> 10) & 1023][1]++;
78-
hist[(id >> 20) & 1023][2]++;
79-
hist[(id >> 30) & 1023][3]++;
80-
hist[(id >> 40) & 1023][4]++;
109+
hist[(id >> 0) & 255][0]++;
110+
hist[(id >> 8) & 255][1]++;
81111
}
82112

83-
unsigned int sum0 = 0, sum1 = 0, sum2 = 0, sum3 = 0, sum4 = 0;
113+
unsigned int sum0 = 0, sum1 = 0;
84114

85115
// replace histogram data with prefix histogram sums in-place
86-
for (int i = 0; i < 1024; ++i)
116+
for (int i = 0; i < 256; ++i)
87117
{
88-
unsigned int h0 = hist[i][0], h1 = hist[i][1], h2 = hist[i][2], h3 = hist[i][3], h4 = hist[i][4];
118+
unsigned int h0 = hist[i][0], h1 = hist[i][1];
89119

90120
hist[i][0] = sum0;
91121
hist[i][1] = sum1;
92-
hist[i][2] = sum2;
93-
hist[i][3] = sum3;
94-
hist[i][4] = sum4;
95122

96123
sum0 += h0;
97124
sum1 += h1;
98-
sum2 += h2;
99-
sum3 += h3;
100-
sum4 += h4;
101125
}
102126

103-
assert(sum0 == count && sum1 == count && sum2 == count && sum3 == count && sum4 == count);
127+
assert(sum0 == count && sum1 == count);
104128
}
105129

106-
static void radixPass(unsigned int* destination, const unsigned int* source, const unsigned long long* keys, size_t count, unsigned int (&hist)[1024][5], int pass)
130+
static void radixPass(unsigned int* destination, const unsigned int* source, const unsigned short* keys, size_t count, unsigned int (&hist)[256][2], int pass)
107131
{
108-
int bitoff = pass * 10;
132+
int bitoff = pass * 8;
109133

110134
for (size_t i = 0; i < count; ++i)
111135
{
112-
unsigned int id = unsigned(keys[source[i]] >> bitoff) & 1023;
136+
unsigned int id = unsigned(keys[source[i]] >> bitoff) & 255;
113137

114138
destination[hist[id][pass]++] = source[i];
115139
}
116140
}
117141

142+
static void partitionPoints(unsigned int* target, const unsigned int* order, const unsigned char* sides, size_t split, size_t count)
143+
{
144+
size_t l = 0, r = split;
145+
146+
for (size_t i = 0; i < count; ++i)
147+
{
148+
unsigned char side = sides[order[i]];
149+
target[side ? r : l] = order[i];
150+
l += 1;
151+
l -= side;
152+
r += side;
153+
}
154+
155+
assert(l == split && r == count);
156+
}
157+
158+
static void splitPoints(unsigned int* destination, unsigned int* orderx, unsigned int* ordery, unsigned int* orderz, const unsigned long long* keys, size_t count, void* scratch, size_t cluster_size)
159+
{
160+
if (count <= cluster_size)
161+
{
162+
memcpy(destination, orderx, count * sizeof(unsigned int));
163+
return;
164+
}
165+
166+
unsigned int* axes[3] = {orderx, ordery, orderz};
167+
168+
int bestk = -1;
169+
unsigned int bestdim = 0;
170+
171+
for (int k = 0; k < 3; ++k)
172+
{
173+
const unsigned int mask = (1 << 20) - 1;
174+
unsigned int dim = (unsigned(keys[axes[k][count - 1]] >> (k * 20)) & mask) - (unsigned(keys[axes[k][0]] >> (k * 20)) & mask);
175+
176+
if (dim >= bestdim)
177+
{
178+
bestk = k;
179+
bestdim = dim;
180+
}
181+
}
182+
183+
assert(bestk >= 0);
184+
185+
// split roughly in half, with the left split always being aligned to cluster size
186+
size_t split = ((count / 2) + cluster_size - 1) / cluster_size * cluster_size;
187+
assert(split > 0 && split < count);
188+
189+
// mark sides of split for partitioning
190+
unsigned char* sides = static_cast<unsigned char*>(scratch) + count * sizeof(unsigned int);
191+
192+
for (size_t i = 0; i < split; ++i)
193+
sides[axes[bestk][i]] = 0;
194+
195+
for (size_t i = split; i < count; ++i)
196+
sides[axes[bestk][i]] = 1;
197+
198+
// partition all axes into two sides, maintaining order
199+
unsigned int* temp = static_cast<unsigned int*>(scratch);
200+
201+
for (int k = 0; k < 3; ++k)
202+
{
203+
if (k == bestk)
204+
continue;
205+
206+
unsigned int* axis = axes[k];
207+
memcpy(temp, axis, sizeof(unsigned int) * count);
208+
partitionPoints(axis, temp, sides, split, count);
209+
}
210+
211+
splitPoints(destination, orderx, ordery, orderz, keys, split, scratch, cluster_size);
212+
splitPoints(destination + split, orderx + split, ordery + split, orderz + split, keys, count - split, scratch, cluster_size);
213+
}
214+
118215
} // namespace meshopt
119216

120217
void meshopt_spatialSortRemap(unsigned int* destination, const float* vertex_positions, size_t vertex_count, size_t vertex_positions_stride)
@@ -127,22 +224,25 @@ void meshopt_spatialSortRemap(unsigned int* destination, const float* vertex_pos
127224
meshopt_Allocator allocator;
128225

129226
unsigned long long* keys = allocator.allocate<unsigned long long>(vertex_count);
130-
computeOrder(keys, vertex_positions, vertex_count, vertex_positions_stride);
227+
computeOrder(keys, vertex_positions, vertex_count, vertex_positions_stride, /* morton= */ true);
131228

132-
unsigned int hist[1024][5];
133-
computeHistogram(hist, keys, vertex_count);
134-
135-
unsigned int* scratch = allocator.allocate<unsigned int>(vertex_count);
229+
unsigned int* scratch = allocator.allocate<unsigned int>(vertex_count * 2); // 4b for order + 2b for keys
230+
unsigned short* keyk = (unsigned short*)(scratch + vertex_count);
136231

137232
for (size_t i = 0; i < vertex_count; ++i)
138233
destination[i] = unsigned(i);
139234

235+
unsigned int* order[] = {scratch, destination};
236+
140237
// 5-pass radix sort computes the resulting order into scratch
141-
radixPass(scratch, destination, keys, vertex_count, hist, 0);
142-
radixPass(destination, scratch, keys, vertex_count, hist, 1);
143-
radixPass(scratch, destination, keys, vertex_count, hist, 2);
144-
radixPass(destination, scratch, keys, vertex_count, hist, 3);
145-
radixPass(scratch, destination, keys, vertex_count, hist, 4);
238+
for (int k = 0; k < 5; ++k)
239+
{
240+
// copy 10-bit key segments into keyk to reduce cache pressure during radix pass
241+
for (size_t i = 0; i < vertex_count; ++i)
242+
keyk[i] = (unsigned short)((keys[i] >> (k * 10)) & 1023);
243+
244+
radixSort10(order[k % 2], order[(k + 1) % 2], keyk, vertex_count);
245+
}
146246

147247
// since our remap table is mapping old=>new, we need to reverse it
148248
for (size_t i = 0; i < vertex_count; ++i)
@@ -202,3 +302,39 @@ void meshopt_spatialSortTriangles(unsigned int* destination, const unsigned int*
202302
destination[r * 3 + 2] = c;
203303
}
204304
}
305+
306+
void meshopt_spatialClusterPoints(unsigned int* destination, const float* vertex_positions, size_t vertex_count, size_t vertex_positions_stride, size_t cluster_size)
307+
{
308+
using namespace meshopt;
309+
310+
assert(vertex_positions_stride >= 12 && vertex_positions_stride <= 256);
311+
assert(vertex_positions_stride % sizeof(float) == 0);
312+
assert(cluster_size > 0);
313+
314+
meshopt_Allocator allocator;
315+
316+
unsigned long long* keys = allocator.allocate<unsigned long long>(vertex_count);
317+
computeOrder(keys, vertex_positions, vertex_count, vertex_positions_stride, /* morton= */ false);
318+
319+
unsigned int* order = allocator.allocate<unsigned int>(vertex_count * 3);
320+
unsigned int* scratch = allocator.allocate<unsigned int>(vertex_count * 2); // 4b for order + 1b for side or 2b for keys
321+
unsigned short* keyk = reinterpret_cast<unsigned short*>(scratch + vertex_count);
322+
323+
for (int k = 0; k < 3; ++k)
324+
{
325+
// copy 16-bit key segments into keyk to reduce cache pressure during radix pass
326+
for (size_t i = 0; i < vertex_count; ++i)
327+
keyk[i] = (unsigned short)(keys[i] >> (k * 20));
328+
329+
unsigned int hist[256][2];
330+
computeHistogram(hist, keyk, vertex_count);
331+
332+
for (size_t i = 0; i < vertex_count; ++i)
333+
order[k * vertex_count + i] = unsigned(i);
334+
335+
radixPass(scratch, order + k * vertex_count, keyk, vertex_count, hist, 0);
336+
radixPass(order + k * vertex_count, scratch, keyk, vertex_count, hist, 1);
337+
}
338+
339+
splitPoints(destination, order, order + vertex_count, order + 2 * vertex_count, keys, vertex_count, scratch, cluster_size);
340+
}

0 commit comments

Comments
 (0)