9
$\begingroup$

I have some bytes, say bytes=ByteArray@{189,178,61,188}, read from a file with Import[...,"Byte"]. The big endian bits for these are

Reverse@Mod[NestList[BitShiftRight,#,7],2]&/@Normal@bytes
(* {{1,0,1,1,1,1,0,1},{1,0,1,1,0,0,1,0},{0,0,1,1,1,1,0,1},{1,0,1,1,1,1,0,0}} *)

I happen to know that those four bytes comprise a float. So, according to https://en.wikipedia.org/wiki/Single-precision_floating-point_format, I'll take the first bit as sign, next 8 as exponent, and last 23 as mantissa:

(-1)^#[[1]]2^(PowerRange[2^7,1,1/2].#[[2;;9]]-127)(1+PowerRange[1/2,2^-23,1/2].#[[10;;]])&[Join@@%]

to obtain -2920303/33554432 or -0.0870318.

I haven't handled the special cases when the exponent is all 0's or 1's. Is there a built-in function to convert from a sequence of 32 bits or a list of four bytes as integers to the corresponding float number? Same for double. Also, is there a way to set the accuracy/precision of the result of this function? I would guess Real[123,Accuracy->13], but that doesn't seem to correspond to 123`13.

$\endgroup$
5
  • 4
    $\begingroup$ Not sure -0.0502283 is the correct decoding here. From what I can see, the big endian decoding is -0.08703181 and the little endian decoding is -0.01157826. On my copy of Mathematica, ImportByteArray[ByteArray[{189,178,61,188},"Real32"] gives the little-endian decoding of -0.01157826. $\endgroup$
    – eyorble
    Jul 23, 2022 at 2:47
  • $\begingroup$ @eyorble Yes, my bad with the endianness. Feel free to post as an answer. $\endgroup$
    – Adam
    Jul 23, 2022 at 4:16
  • 4
    $\begingroup$ Does this answer your question? Convert ByteArray to integer and real values $\endgroup$
    – Roman
    Jul 23, 2022 at 9:01
  • $\begingroup$ @Roman yes; your answer is general enough to handle my question also. It may be useful though to have a question that specifically mentions IEEE 754 float in the title (I didn't find that question when searching). $\endgroup$
    – Adam
    Jul 23, 2022 at 16:59
  • 2
    $\begingroup$ I would not close this. It can be handled by existing code, but not very efficiently. $\endgroup$ Jul 23, 2022 at 21:29

3 Answers 3

16
$\begingroup$
ImportByteArray[ByteArray[{189,178,61,188}, "Real32", ByteOrdering -> -1]
ImportByteArray[ByteArray[{189,178,61,188}, "Real32", ByteOrdering -> 1]

-0.01157826

-0.0870318

The -1 ordering is little-endian, the +1 ordering is big-endian. See the documentation for ByteOrdering for details.

Note that double precision floats can be treated similarly using "Real64", and that quad-precision is also supported with "Real128". Other formats can be found in the guide entry guide/ListingOfAllFormats.

I would also additionally recommend seeing ImportString and Import, and note that arrays can be handled gracefully by including an appropriate multiple of the primitive size in bytes or characters.

E.g.

ImportByteArray[ByteArray[{1, 2, 3, 4, 5, 6, 7, 8}], "Integer16",  ByteOrdering -> -1]
FromDigits[Reverse[#], 256] & /@ {{1, 2}, {3, 4}, {5, 6}, {7, 8}}

{513, 1027, 1541, 2055}

{513, 1027, 1541, 2055}

$\endgroup$
9
+50
$\begingroup$

As @DanielLichtblau points out, ImportByteArray is surprisingly slow. Here's a LibraryLink code that uses a bare-metal byte-to-float type cast. Its advantages are:

  • All aspects of IEEE 754 are automatically taken care of because there is no actual conversion happening, only a low-level type cast. Infinities, NaNs, subnormals, etc. should be fine (up to Mathematica's interpretation of these, of course; see below).
  • The conversion is about 6000 times faster than ImportByteArray.

For comparison, we define the slow version as in eyorble's answer:

makefloatE[s_] := First@ImportByteArray[s, "Real32", ByteOrdering -> 1]

and Daniel's bit-fiddling answer

elen = 8;
makefloatD = (bits = Join @@ IntegerDigits[#, 2, 8];
              sign = (-1)^First[bits];
              {exponbits, mantissabits} = TakeDrop[Rest[bits], elen];
              mantissabits = Prepend[mantissabits, 1];
              expon = FromDigits[exponbits, 2] - 2^(elen - 1) + 2;
              mantissa = FromDigits[mantissabits, 2]*2^(-Length[mantissabits]);
              N[sign*2^expon*mantissa]) &;

For the fast version, we use the following LibraryLink C-code:

code = "
#include \"WolframLibrary.h\"
#include \"WolframNumericArrayLibrary.h\"
#include <stdint.h>

// comment out this line if the endianness does not match
#define REVERSE_BYTE_ORDER

DLLEXPORT int makefloat(WolframLibraryData libData, mint Argc,
                        MArgument *Args, MArgument Res) {
  WolframNumericArrayLibrary_Functions naFuns =
      libData->numericarrayLibraryFunctions;
  MNumericArray bytes = MArgument_getMNumericArray(Args[0]);
  mint const * const dimensions = naFuns->MNumericArray_getDimensions(bytes);
  if (dimensions[0] != 4)
    return LIBRARY_DIMENSION_ERROR;
  uint8_t const * const data = naFuns->MNumericArray_getData(bytes);
#ifdef REVERSE_BYTE_ORDER
  uint8_t reversed_data[4] = {data[3], data[2], data[1], data[0]};
  float *x = (float *)reversed_data;
#else
  float *x = (float *)data;
#endif
  MArgument_setReal(Res, *x);
  return LIBRARY_NO_ERROR;
}";

Compile the code with the C compiler driver:

Needs["CCompilerDriver`"]
lib = CreateLibrary[code, "byteconverter", "CompileOptions" -> "-O3", 
   "ShellOutputFunction" -> Print, "ShellCommandFunction" -> Print];
makefloatR = LibraryFunctionLoad[lib, 
   "makefloat", {{LibraryDataType[ByteArray], "Constant"}}, Real];

Now we can compare speeds:

B = ByteArray@{189, 178, 61, 188};

RepeatedTiming[makefloatR@B]
(*    {2.31452*10^-7, -0.0870318}    *)

RepeatedTiming[makefloatE@B]
(*    {0.00135563, -0.0870318}    *)

RepeatedTiming[makefloatD@Normal@B]
(*    {6.48921*10^-6, -0.0870318}    *)

The compiled code is thus about 6000 times faster than ImportByteArray, and 30 times faster than the Mathematica-level bit-fiddling code.

It looks like special cases of IEEE 754 are covered well:

makefloatR@ByteArray@{127, 192, 0, 0}    (* NAN *)
(*    Indeterminate    *)

makefloatR@ByteArray@{127, 128, 0, 0}    (* INFINITY *)
(*    ∞    *)

vectorized version

It may be desirable to convert a long list of numbers in one go. The ImportByteArray answer becomes

makefloatE[s_] := ImportByteArray[s, "Real32", ByteOrdering -> 1]

and the C code becomes

code = "
#include \"WolframLibrary.h\"
#include \"WolframNumericArrayLibrary.h\"
#include <stdint.h>

// comment out this line if the endianness does not match
#define REVERSE_BYTE_ORDER

DLLEXPORT int makefloat_vec(WolframLibraryData libData, mint Argc,
                            MArgument *Args, MArgument Res) {
  WolframNumericArrayLibrary_Functions naFuns = libData->numericarrayLibraryFunctions;
  MNumericArray bytes = MArgument_getMNumericArray(Args[0]);
  mint const * const dimensions = naFuns->MNumericArray_getDimensions(bytes);
  if (dimensions[0] & 3)
    return LIBRARY_DIMENSION_ERROR;
  const mint length = dimensions[0]/4;
  uint8_t const * const data = naFuns->MNumericArray_getData(bytes);
  MTensor result;
  int err = libData->MTensor_new(MType_Real, 1, &length, &result);
  if (err != LIBRARY_NO_ERROR)
    return err;
  double * const data_out = libData->MTensor_getRealData(result);
  for (mint i=0; i<length; i++) {
#ifdef REVERSE_BYTE_ORDER
    uint8_t reversed_data[4] = {data[4*i+3], data[4*i+2], data[4*i+1], data[4*i]};
    float *x = (float *)reversed_data;
#else
    float *x = (float *)(data+4*i);
#endif
    data_out[i] = *x;
  }
  MArgument_setMTensor(Res, result);
  return LIBRARY_NO_ERROR;
}";
Needs["CCompilerDriver`"]
lib = CreateLibrary[code, "byteconverter", "CompileOptions" -> "-O3", 
   "ShellOutputFunction" -> Print, "ShellCommandFunction" -> Print];
makefloatR = LibraryFunctionLoad[lib, 
   "makefloat_vec", {{LibraryDataType[ByteArray], "Constant"}}, {Real, 1}];

Compare execution speed:

SeedRandom[1234];
B = ByteArray@RandomInteger[{0, 127}, 4*10^8];

RepeatedTiming[fE = makefloatE@B;]
(*    {4.58079, Null}    *)

RepeatedTiming[fR = makefloatR@B;]
(*    {0.116399, Null}    *)

fE == fR
(*    True    *)

The C code is still 40 times faster than ImportByteArray but no longer supports IEEE 754 edge cases like NaN and Inf:

makefloatR@ByteArray@{127, 192, 0, 0}    (* NAN *)
(* Numeric data containing a floating point exception
   (NaN or Inf) encountered. *)
(*    $Failed    *)

makefloatR@ByteArray@{127, 128, 0, 0}    (* INFINITY *)
(* Numeric data containing a floating point exception
   (NaN or Inf) encountered. *)
(*    $Failed    *)
$\endgroup$
13
  • 1
    $\begingroup$ @eyorble I've added a vectorized version. $\endgroup$
    – Roman
    Jul 26, 2022 at 7:47
  • 4
    $\begingroup$ Nice method. One thing of possible concern might be alignment of the byte array(s). That said, this casting approach is quite fast. $\endgroup$ Jul 26, 2022 at 13:44
  • 3
    $\begingroup$ @DanielLichtblau I'm concerned about alignments too. But C compilers/linkers these days tend to be aware of programmers' sloppiness and align everything on 4- or 8-byte boundaries, just to prevent some of these very common mistakes. __attribute__ ((aligned (4))) could be used to align reversed_data but we don't have access to the alignments of the data and data_out fields as provided by MNumericArray and MTensor. Keeping my fingers crossed that CPU+kernel will trap bad alignment & fix it (with speed penalty). $\endgroup$
    – Roman
    Jul 26, 2022 at 14:26
  • 1
    $\begingroup$ In our code environment mint = "machine integer". $\endgroup$ Jul 26, 2022 at 18:24
  • 1
    $\begingroup$ @Adam I think mint and int are the same thing. For some reason LibraryLink insists on renaming int but is fine with leaving double alone. As for the const * const part, please see cdecl.org: for example, int const * const data translates as "declare data as const pointer to const int". Thus making sure that the compiler understands that the pointer is constant and also the data pointed to is constant. With these compiler hints I'm trying to speed up the C code by telling the compiler which quantities can be assumed to be constant. $\endgroup$
    – Roman
    Jul 26, 2022 at 20:17
7
$\begingroup$

As shown already, this can be done by ImportByteArray. If you have long lists, this can be slow. Here is a test.

SeedRandom[1234];
Table[
 bytes = Partition[RandomInteger[{0, 255}, 2^n], 4];
 First[Timing[
   Map[ImportByteArray[ByteArray[#], "Real32", ByteOrdering -> -1] &, 
    bytes]
   ]],
 {n, 2, 12}
 ]

(* Out[302]= {0.004421, 0.007021, 0.013251, 0.025017, 0.045063, \
0.084426, 0.172539, 0.333143, 0.673713, 1.36853, 2.72413} *)

One can instead do the conversion the pedestrian way, using the bit fields of a single precision 32-bit real. First is the sign bit, next eight give the exponent between $2^{-127}$ and $2^{127}-1$. There is an implied first bit set to 1 in the mantissa, and the remaining mantissa bits are filled in by our byte array bits.

elen = 8;
SeedRandom[1234];
Table[
 bytes = Partition[RandomInteger[{0, 255}, 2^n], 4];
 Timing[Length[Map[(
      bits = Join @@ IntegerDigits[Reverse@#, 2, 8];
      sign = (-1)^First[bits];
      {exponbits, mantissabits} = TakeDrop[Rest[bits], elen];
      mantissabits = Prepend[mantissabits, 1];
      expon = FromDigits[exponbits, 2] - 2^(elen - 1) + 2;
      mantissa = FromDigits[mantissabits, 2]*2^(-Length[mantissabits]);
      N[sign*2^expon*mantissa]) &, bytes]]
  ],
 {n, 2, 12}
 ]

(* Out[329]= {{0.000095, 1}, {0.000037, 2}, {0.000069, 4}, {0.000142, 
  8}, {0.000234, 16}, {0.000558, 32}, {0.001268, 64}, {0.003829, 
  128}, {0.005531, 256}, {0.008505, 512}, {0.016775, 1024}} *)

That's two orders of magnitude. Oh, and they do give the same results.

In[332]:= SeedRandom[1234];
Table[
 bytes = Partition[RandomInteger[{0, 255}, 2^n], 4];
 Map[ImportByteArray[ByteArray[#], "Real32", 
     ByteOrdering -> -1][[1]] &, bytes],
 {n, 2, 4}
 ]

(* Out[333]= {{1.76353*10^-34}, {3.52751*10^15, -4230.}, {-0.035666, 
  4.6259*10^23, 0.000014599, -6.77262*10^35}} *)

In[334]:= SeedRandom[1234];
Table[
 bytes = Partition[RandomInteger[{0, 255}, 2^n], 4];
 Map[(
    bits = Join @@ IntegerDigits[Reverse@#, 2, 8];
    sign = (-1)^First[bits];
    {exponbits, mantissabits} = TakeDrop[Rest[bits], elen];
    mantissabits = Prepend[mantissabits, 1];
    expon = FromDigits[exponbits, 2] - 2^(elen - 1) + 2;
    mantissa = FromDigits[mantissabits, 2]*2^(-Length[mantissabits]);
    N[sign*2^expon*mantissa]) &, bytes],
 {n, 2, 4}
 ]

(* Out[335]= {{1.76353*10^-34}, {3.52751*10^15, -4230.}, {-0.035666, 
  4.6259*10^23, 0.000014599, -6.77262*10^35}} *)

--- edit ---

As note in a response by @Roman, there are faster ways, for example using a link to C or similar and casting as a 32 bit machine real. I do not have that level of casting ability (I'd need to rise above assistant assistant producer). But the bit fiddling approach can be sped up. First the reference implementation.

toReal32a[bytes_, elen_ : 8] := Module[
  {bits, sign, exponbits, mantissabits, expon, mantissa},
  bits = Join @@ IntegerDigits[Reverse@bytes, 2, 8];
  sign = (-1)^First[bits];
  {exponbits, mantissabits} = TakeDrop[Rest[bits], elen];
  mantissabits = Prepend[mantissabits, 1];
  expon = FromDigits[exponbits, 2] - 2^(elen - 1) + 2;
  mantissa = FromDigits[mantissabits, 2]*2^(-Length[mantissabits]);
  N[sign*2^expon*mantissa]
  ]

First redo using bitwise operations. This will not give a speed gain.

toReal32b[bytes_, elen_ : 8] := Module[
  {uint, signmask = 2^31, offset = 31 - elen, emask, mmask, sign, 
   expon, mantissa},
  emask = 2^31 - 2^offset;
  mmask = 2^offset - 1;
  uint = 2^Range[0, 24, 8] . bytes;
  sign = If[BitAnd[uint, signmask] == 0, 1, -1];
  expon = 
   BitShiftRight[BitAnd[uint, emask], offset] - 2^(elen - 1) + 2;
  mantissa = BitOr[2^offset, BitAnd[uint, mmask]];
  N[sign*2^expon*2^(-offset - 1)*mantissa]
  ]

This version is easy to run through Compile. (Maybe so is the first version, I just doubt it will be quite as fast.)

toReal32c = Compile[{{bytes, _Integer, 1}, {elen, _Integer}}, Module[
    {uint, signmask = 2^31, offset = 31 - elen, emask, mmask, sign, 
     expon, mantissa},
    emask = 2^31 - 2^offset;
    mmask = 2^offset - 1;
    uint = 2^Range[0, 24, 8] . bytes;
    sign = If[BitAnd[uint, signmask] == 0, 1, -1];
    expon = 
     BitShiftRight[BitAnd[uint, emask], offset] - 2^(elen - 1) + 2;
    mantissa = BitOr[2^offset, BitAnd[uint, mmask]];
    sign*2.^(expon - offset - 1)*mantissa
    ], CompilationTarget -> "C", RuntimeOptions -> "Speed"];

Check that they agree.

In[36]:= SeedRandom[1234];
bytes = RandomInteger[{0, 255}, 4]
Map[#[bytes, 8] &, {toReal32a, toReal32b, toReal32c}]

(* Out[37]= {224, 105, 106, 7}

Out[38]= {1.76353*10^-34, 1.76353*10^-34, 1.76353*10^-34} *)

This last delivers far better timings than the reference.

In[39]:= SeedRandom[1234];
Table[
 bytes = Partition[RandomInteger[{0, 255}, 2^n], 4];
 Timing[Length[Map[toReal32a, bytes]]], {n, 2, 12}
 ]

Out[40]= {{0.000112, 1}, {0.000067, 2}, {0.000144, 4}, {0.000226, 
  8}, {0.000448, 16}, {0.001035, 32}, {0.001595, 64}, {0.003396, 
  128}, {0.006998, 256}, {0.014437, 512}, {0.028745, 1024}}

In[41]:= SeedRandom[1234];
Table[
 bytes = Partition[RandomInteger[{0, 255}, 2^n], 4];
 Timing[Length[Map[toReal32b, bytes]]], {n, 2, 12}
 ]

Out[42]= {{0.000149, 1}, {0.000098, 2}, {0.000202, 4}, {0.000424, 
  8}, {0.000737, 16}, {0.001521, 32}, {0.00296, 64}, {0.005906, 
  128}, {0.012203, 256}, {0.020636, 512}, {0.02759, 1024}}

In[43]:= SeedRandom[1234];
Table[
 bytes = Partition[RandomInteger[{0, 255}, 2^n], 4];
 Timing[Length[Map[toReal32c[#, 8] &, bytes]]], {n, 2, 12}
 ]

Out[44]= {{0.000034, 1}, {9.*10^-6, 2}, {0.00001, 4}, {0.000016, 
  8}, {0.000027, 16}, {0.000054, 32}, {0.0001, 64}, {0.000226, 
  128}, {0.00013, 256}, {0.000254, 512}, {0.000502, 1024}}

--- end edit ---

$\endgroup$
7
  • $\begingroup$ Your code does not handle edge cases of IEEE-754 (signed zeros, subnormals, infinities, NaNs with payload, etc.). Maybe they aren't necessary here; but sooner or later someone will apply this code to a non-standard number and run into trouble. $\endgroup$
    – Roman
    Jul 24, 2022 at 5:26
  • 1
    $\begingroup$ @Roman (1) I'd hope one would recognize the need to modify the code first. Like everything else on MSE, it comes as is, and not necessarily suitable for all possible usage. (2) I doubt modifications to handle these cases will amount to even a factor of two in terms of slowdown. $\endgroup$ Jul 24, 2022 at 15:33
  • 2
    $\begingroup$ (1) Actually elen is defined recursively as elen=NoteToSelf["Do not forget to include definition of",elen]. (I added it in an edit.) $\endgroup$ Jul 25, 2022 at 14:07
  • 2
    $\begingroup$ (2) As @Roman tried to indicate*, your request included mention of edge cases (denormalized, infinite, maybe signed zeros). So these will need extra code. (*) I didn't read far enough. I borrowed my attention span from a teenager. $\endgroup$ Jul 25, 2022 at 14:09
  • 1
    $\begingroup$ Ha! I'm not too interested in the edge cases; my data didn't have any. I like the bit fiddling approach, but I'll keep the accept with eyorble since I was looking for a built in or otherwise succinct way. $\endgroup$
    – Adam
    Jul 26, 2022 at 1:16

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.