Индекс Мортона от 2D точки с поплавками
У меня есть 2D точка, которая выглядит примерно так:
class Point
{
float m_x, m_y;
public:
int mortonIndex()
{
// what would go here?
}
};
Я знаю, что делать с целыми числами, но мне нужно использовать числа с плавающей точкой. Я также хочу избежать масштабирования для любого конкретного размера сетки.
Связанная страница в Википедии:
Решение
Есть два способа посмотреть на это:
- «Число с плавающей запятой», являющееся 256-битным числом (или аналогично, «двойное число» как 2100-битное число, реализовано [здесь] [1])
- «Float» — странное 32-разрядное целое число.
Я буду использовать последнее, так как это легче реализовать.
Этот подход использует тот факт, что IEEE floatИзначально они были разработаны для совместимости со старыми целочисленными механизмами баз данных, что позволяло им рассматривать числа с плавающей запятой как целые числа, дополняющие 1.
Точнее, в смысле дополнения 1, порядок значений с плавающей запятой соответствует порядку целых чисел одинаковой ширины (фактически, непосредственное добавление 1 к плавающему плавающему даст вам смежное значение с большей абсолютной величиной **).
class Point
{
float m_x, m_y;
// This assert is not correct when the floating point model
// is not IEEE-compliant, float is not 32-bit, or both.
//
// Such things are hard to find, so we'll just assume
// mostly-sane hardware.
//
static_assert(
(sizeof(int) == sizeof(float)) &&
(sizeof(int)*CHAR_BIT == 32) &&
(sizeof(long long)*CHAR_BIT == 64),
"We need 32-bit ints and floats, and 64-bit long longs!");
public:
// So we don't lose any information, we need 2x the width.
// After all, we're cramming two 32-bit numbers into a single value.
// Lossiness here would mean a container would need to implement
// a binning strategy.
//
// Higher dimensions would require an array, obviously.
//
// Also, we're not going to modify the point, so make this a const routine.
//
long long mortonIndex() const
{
// Pun the x and y coordinates as integers: Just re-interpret the bits.
//
auto ix = reinterpret_cast<const unsigned &>(this->m_x);
auto iy = reinterpret_cast<const unsigned &>(this->m_y);
// Since we're assuming 2s complement arithmetic (99.99% of hardware today),
// we'll need to convert these raw integer-punned floats into
// their corresponding integer "indices".
// Smear their sign bits into these for twiddling below.
//
const auto ixs = static_cast<int>(ix) >> 31;
const auto iys = static_cast<int>(iy) >> 31;
// This is a combination of a fast absolute value and a bias.
//
// We need to adjust the values so -FLT_MAX is close to 0.
//
ix = (((ix & 0x7FFFFFFFL) ^ ixs) - ixs) + 0x7FFFFFFFL;
iy = (((iy & 0x7FFFFFFFL) ^ iys) - iys) + 0x7FFFFFFFL;
// Now we have -FLT_MAX close to 0, and FLT_MAX close to UINT_MAX,
// with everything else in-between.
//
// To make this easy, we'll work with x and y as 64-bit integers.
//
long long xx = ix;
long long yy = iy;
// Dilate and combine as usual...
xx = (xx | (xx << 16)) & 0x0000ffff0000ffffLL;
yy = (yy | (yy << 16)) & 0x0000ffff0000ffffLL;
xx = (xx | (xx << 8)) & 0x00ff00ff00ff00ffLL;
yy = (yy | (yy << 8)) & 0x00ff00ff00ff00ffLL;
xx = (xx | (xx << 4)) & 0x0f0f0f0f0f0f0f0fLL;
yy = (yy | (yy << 4)) & 0x0f0f0f0f0f0f0f0fLL;
xx = (xx | (xx << 2)) & 0x3333333333333333LL;
yy = (yy | (yy << 2)) & 0x3333333333333333LL;
xx = (xx | (xx << 1)) & 0x5555555555555555LL;
yy = (yy | (yy << 1)) & 0x5555555555555555LL;
return xx | (yy << 1);
}
};
Обратите внимание, что вершины результирующей кривой имеют то же распределение, что и позиции в 2D-пространстве с плавающей точкой.
Это может быть проблемой, если вы собираетесь использовать это со структурой на диске, так как кластеризация вблизи координатных осей или начало координат может привести к тому, что запросы диапазона пересекают большое количество полей рядом с ними. В противном случае, IMO, это разумно эффективная альтернатива генерации единых индексов (и без ветвей!).
** Специальная обработка необходима для бесконечностей и NaN, но вы поняли идею.
Другие решения