Comment désactiver efficacement les bits (Morton inverse)
-
29-10-2019 - |
Question
Cette question: Comment désinterner les bits (non mortel?) A une bonne réponse pour extraire l'une des deux moitiés d'un nombre de Morton (juste les bits impairs), mais j'ai besoin d'une solution qui extrait les deux parties (les bits impairs et les bits pair) dans le moins d'opérations possible.
Pour mon usage, j'aurais besoin de prendre un INT 32 bits et d'extraire deux INT à 16 bits, où l'un est les bits pair et l'autre est les bits impairs décalés à droite par 1 bits, par exemple
input, z: 11101101 01010111 11011011 01101110
output, x: 11100001 10110111 // odd bits shifted right by 1
y: 10111111 11011010 // even bits
Il semble y avoir de nombreuses solutions utilisant des changements et des masques avec des nombres magiques pour générer des numéros Morton (c'est-à-dire des bits entrelacés), par exemple Entrelacer les bits par des numéros de magie binaire, mais je n'ai rien trouvé pour faire l'inverse (c'est-à-dire de détendre).
METTRE À JOUR
Après avoir relu la section de Hacker's Delight on Perfect Shuffles / Unswuffles, j'ai trouvé des exemples utiles que j'ai adaptés comme suit:
// morton1 - extract even bits
uint32_t morton1(uint32_t x)
{
x = x & 0x55555555;
x = (x | (x >> 1)) & 0x33333333;
x = (x | (x >> 2)) & 0x0F0F0F0F;
x = (x | (x >> 4)) & 0x00FF00FF;
x = (x | (x >> 8)) & 0x0000FFFF;
return x;
}
// morton2 - extract odd and even bits
void morton2(uint32_t *x, uint32_t *y, uint32_t z)
{
*x = morton1(z);
*y = morton1(z >> 1);
}
Je pense que cela peut toujours être amélioré, à la fois dans sa forme scalaire actuelle et aussi en profitant de SIMD, donc je suis toujours intéressé par de meilleures solutions (Scalar ou SIMD).
La solution
Si votre processeur gère efficacement les INTS 64 bits, vous pouvez combiner les opérations ...
int64 w = (z &0xAAAAAAAA)<<31 | (z &0x55555555 )
w = (w | (w >> 1)) & 0x3333333333333333;
w = (w | (w >> 2)) & 0x0F0F0F0F0F0F0F0F;
...
Autres conseils
Code pour l'Intel Haswell et les processeurs ultérieurs. Vous pouvez utiliser le jeu d'instructions BMI2 qui contient les instructions PEXT et PDEP. Ceux-ci peuvent (entre autres grandes choses) être utilisés pour construire vos fonctions.
#include <immintrin.h>
#include <stdint.h>
// on GCC, compile with option -mbmi2, requires Haswell or better.
uint64_t xy_to_morton (uint32_t x, uint32_t y)
{
return _pdep_u32(x, 0x55555555) | _pdep_u32(y,0xaaaaaaaa);
}
uint64_t morton_to_xy (uint64_t m, uint32_t *x, uint32_t *y)
{
*x = _pext_u64(m, 0x5555555555555555);
*y = _pext_u64(m, 0xaaaaaaaaaaaaaaaa);
}
Dans le cas où quelqu'un utilise des codes Morton en 3D, il doit donc lire un peu tous les 3, et 64 bits ici est la fonction que j'ai utilisée:
uint64_t morton3(uint64_t x) {
x = x & 0x9249249249249249;
x = (x | (x >> 2)) & 0x30c30c30c30c30c3;
x = (x | (x >> 4)) & 0xf00f00f00f00f00f;
x = (x | (x >> 8)) & 0x00ff0000ff0000ff;
x = (x | (x >> 16)) & 0xffff00000000ffff;
x = (x | (x >> 32)) & 0x00000000ffffffff;
return x;
}
uint64_t bits;
uint64_t x = morton3(bits)
uint64_t y = morton3(bits>>1)
uint64_t z = morton3(bits>>2)
Si vous avez besoin de vitesse, vous pouvez utiliser la conversion de table pour une conversion d'un octet à la fois (le tableau des deux octets est plus rapide mais grand). La procédure est effectuée sous Delphi IDE mais l'assembleur / algorithem est le même.
const
MortonTableLookup : array[byte] of byte = ($00, $01, $10, $11, $12, ... ;
procedure DeinterleaveBits(Input: cardinal);
//In: eax
//Out: dx = EvenBits; ax = OddBits;
asm
movzx ecx, al //Use 0th byte
mov dl, byte ptr[MortonTableLookup + ecx]
//
shr eax, 8
movzx ecx, ah //Use 2th byte
mov dh, byte ptr[MortonTableLookup + ecx]
//
shl edx, 16
movzx ecx, al //Use 1th byte
mov dl, byte ptr[MortonTableLookup + ecx]
//
shr eax, 8
movzx ecx, ah //Use 3th byte
mov dh, byte ptr[MortonTableLookup + ecx]
//
mov ecx, edx
and ecx, $F0F0F0F0
mov eax, ecx
rol eax, 12
or eax, ecx
rol edx, 4
and edx, $F0F0F0F0
mov ecx, edx
rol ecx, 12
or edx, ecx
end;
Je ne voulais pas être limité à un entier de taille fixe et faire des listes de commandes similaires avec des constantes codées en dur, j'ai donc développé une solution C ++ 11 qui utilise la métaprogrammation du modèle pour générer les fonctions et les constantes. Le code d'assembly généré avec -O3
Semble aussi serré que possible sans utiliser BMI:
andl $0x55555555, %eax
movl %eax, %ecx
shrl %ecx
orl %eax, %ecx
andl $0x33333333, %ecx
movl %ecx, %eax
shrl $2, %eax
orl %ecx, %eax
andl $0xF0F0F0F, %eax
movl %eax, %ecx
shrl $4, %ecx
orl %eax, %ecx
movzbl %cl, %esi
shrl $8, %ecx
andl $0xFF00, %ecx
orl %ecx, %esi
Tl; dr repo source et démo en direct.
Mise en œuvre
Fondamentalement à chaque étape du morton1
La fonction fonctionne en déplaçant et en ajoutant à une séquence de constantes qui ressemblent à ceci:
0b0101010101010101
(Alternate 1 et 0)0b0011001100110011
(Alternate 2x 1 et 0)0b0000111100001111
(Alternate 4x 1 et 0)0b0000000011111111
(Alternate 8x 1 et 0)
Si nous devions utiliser D
dimensions, nous aurions un motif avec D-1
zéros et 1
une. Donc, pour les générer, il suffit de générer des consécutifs et d'appliquer un peu de bit ou:
/// @brief Generates 0b1...1 with @tparam n ones
template <class T, unsigned n>
using n_ones = std::integral_constant<T, (~static_cast<T>(0) >> (sizeof(T) * 8 - n))>;
/// @brief Performs `@tparam input | (@tparam input << @tparam width` @tparam repeat times.
template <class T, T input, unsigned width, unsigned repeat>
struct lshift_add :
public lshift_add<T, lshift_add<T, input, width, 1>::value, width, repeat - 1> {
};
/// @brief Specialization for 1 repetition, just does the shift-and-add operation.
template <class T, T input, unsigned width>
struct lshift_add<T, input, width, 1> : public std::integral_constant<T,
(input & n_ones<T, width>::value) | (input << (width < sizeof(T) * 8 ? width : 0))> {
};
Maintenant que nous pouvons générer les constantes au moment de la compilation pour des dimensions arbitraires avec ce qui suit:
template <class T, unsigned step, unsigned dimensions = 2u>
using mask = lshift_add<T, n_ones<T, 1 << step>::value, dimensions * (1 << step), sizeof(T) * 8 / (2 << step)>;
Avec le même type de récursivité, nous pouvons générer des fonctions pour chacune des étapes de l'algorithme x = (x | (x >> K)) & M
:
template <class T, unsigned step, unsigned dimensions>
struct deinterleave {
static T work(T input) {
input = deinterleave<T, step - 1, dimensions>::work(input);
return (input | (input >> ((dimensions - 1) * (1 << (step - 1))))) & mask<T, step, dimensions>::value;
}
};
// Omitted specialization for step 0, where there is just a bitwise and
Il reste à répondre à la question "De combien d'étapes avons-nous besoin?". Cela dépend également du nombre de dimensions. En général, k
étapes calculées 2^k - 1
bits de sortie; Le nombre maximum de bits significatifs pour chaque dimension est donné par z = sizeof(T) * 8 / dimensions
, donc il suffit de prendre 1 + log_2 z
pas. Le problème est maintenant que nous avons besoin de cela comme constexpr
Afin de l'utiliser comme paramètre de modèle. La meilleure façon que j'ai trouvé pour contourner ceci est de définir log2
via la métaprogrammation:
template <unsigned arg>
struct log2 : public std::integral_constant<unsigned, log2<(arg >> 1)>::value + 1> {};
template <>
struct log2<1u> : public std::integral_constant<unsigned, 0u> {};
/// @brief Helper constexpr which returns the number of steps needed to fully interleave a type @tparam T.
template <class T, unsigned dimensions>
using num_steps = std::integral_constant<unsigned, log2<sizeof(T) * 8 / dimensions>::value + 1>;
Et enfin, nous pouvons effectuer un seul appel:
/// @brief Helper function which combines @see deinterleave and @see num_steps into a single call.
template <class T, unsigned dimensions>
T deinterleave_first(T n) {
return deinterleave<T, num_steps<T, dimensions>::value - 1, dimensions>::work(n);
}