I have some files with floating point values in them that I'm having problems figuring out their exact encoding. I've tried several ways to convert to a standard double value, but haven't had much luck. I know what the values are supposed to convert to, but need a method to extract directly from the files.
These are HP-1000 (HP21xx) series floating point values. Similar to 48 bit Pascal, but not the same.
Other than some old, unreliable documentation (couldn't get a conversion using what it claimed was the format), and a list of formats from quadibloc.com I have found nothing else. Format from quadibloc gives the following:
(Note, big-endian order)
MSb: mantissa sign (says it's 2s complement)
39 bits: mantissa
7 bits: exponent
1 bit: exponent sign
The other documentation isn't clear, but seems to say it's all 2s complement. I've tried the Pascal conversion code I found, but it doesn't even come close. (after moving the exponent sign to high bit of the byte).
Examples, converted using a PPC that emulates an HP-1000:
A7EB851EB90A -22.02
A870A3D70A0A -21.89
A8AE147AE20A -21.83
A8E147AE140A -21.78
A9666666670A -21.65
AC70A3D70A0A -20.89
ACC28F5C290A -20.81
ACCCCCCCCC0A -20.8
AE70A3D70B0A -20.39
AEB851EB850A -20.32
This was from only one file. I have at least 100 of these files to extract from.
Any ideas? Sure would be appreciated.
Edit: I guess I should have said what language I am working in. In this case, it's C#. As for additional data, I have a ton of it. More examples here, these include negative exponents, at least when converted to decimal.
400000000002 1
43851EB85204 2.11
451EB851EB04 2.16
4C7AE147AE04 2.39
4EC4EC4EC5F3 4.8076923077E-03
519999999A04 2.55
5838B6BE9BF3 5.3846153846E-03
5B851EB85202 1.43
5BD70A3D7108 11.48
5C7AE147AE08 11.56
5E51EB851E08 11.79
5FD70A3D7108 11.98
62E147AE1508 12.36
64A3D70A3E08 12.58
666666666702 1.6
67AE147AE202 1.62
6B204B9E54F3 6.5384615385E-03
733333333302 1.8
762762762AF3 7.2115384616E-03
794DFB4619F3 7.4038461538E-03
7C74941627F3 7.5961538465E-03
7E07E07E14F3 7.6923076925E-03
Should have included positive numbers before.
Tried the suggestions, but in C#. Wildly different results. I even split out the bits of the calculations, but I get very high numbers compared to the numbers given in the suggestions. Here's the code with the bits expanded.
Int64[] test_data = new Int64[] {
0xA7EB851EB90A, 0xA870A3D70A0A, 0xA8AE147AE20A, 0xA8E147AE140A,
0xA9666666670A, 0xAC70A3D70A0A, 0xACC28F5C290A, 0xACCCCCCCCC0A,
0xAE70A3D70B0A, 0xAEB851EB850A, 0x400000000002, 0x43851EB85204,
0x451EB851EB04, 0x4C7AE147AE04, 0x4EC4EC4EC5F3, 0x519999999A04,
0x5838B6BE9BF3, 0x5B851EB85202, 0x5BD70A3D7108, 0x5C7AE147AE08,
0x5E51EB851E08, 0x5FD70A3D7108, 0x62E147AE1508, 0x64A3D70A3E08,
0x666666666702, 0x67AE147AE202, 0x6B204B9E54F3, 0x733333333302,
0x762762762AF3, 0x794DFB4619F3, 0x7C74941627F3, 0x7E07E07E14F3
};
private void Checkit()
{
Byte[] td = new Byte[6]; // 6 byte array for reversed data
Byte[] ld = new Byte[8]; // 8 byte conversion array
for (Int32 i = 0; i < test_data.Length; i++) // Loop through data
{
ld = BitConverter.GetBytes(test_data[i]); // Get value as byte array
for (Int32 j = 0; j < 6; j++) // Copy the bytes in reverse
td[5 - j] = ld[j];
// Go test them
Console.WriteLine("{0:X6} --> {1:E}", test_data[i], Real48ToDouble(ref td));
}
}
Double Real48ToDouble(ref Byte[] realValue)
{
// Values are using input value of A7EB851EB90A and 7E07E07E14F3
if (realValue[0] == 0)
return 0.0; // Null exponent = 0
Byte[] b = new Byte[8] { 0, 0, 0, 0, 0, 0, 0, 0 }; // 64 bit byte array
for (Int32 i = 0; i < 5; i++) // Copy over the 48 bit info
b[4 - i] = realValue[i];
Int64 mant = BitConverter.ToInt64(b, 0); // 0x000000A7EB851EB9 - Get mantissa with sign
// 0x0000007E07E07E14
Int32 expo = realValue[5]; // 0x0000000A - Grab the exponent
// 0x000000F3
Int32 mant_sign = (Int32)(mant >> 39); // 0x00000001
// 0x00000000
// sign extend mantissa from 40 to 64 bits, then take absolute value
mant = (Int64)((mant ^ (1L << 39)) - (1L << 39)); // 0xFFFFFFA7EB851EB9 - First calc
// 0x0000007E07E07E14
mant = Math.Abs(mant); // 0x00000058147AE147 - Make absolute
// 0x0000007E07E07E14
// convert mantissa to floating-point
Double fmant = mant * Math.Pow(2, -39.0); // 0.68812499999876309 - Second calc
// 0.98461538463743636
// rotate exponent field right by 1 bit
expo = (expo >> 1) | ((expo & 1) << 7); // 0x00000005
// 0x000000F9
// sign extend exponent from 8 to 32 bits
expo = ((expo ^ (1 << 7)) - (1 << 7)); // 0x00000005
// 0xFFFFFFF9
// compute scale factor from exponent field
Double scale = Math.Pow(2, expo); // 32.0 - Scale
// 0.0078125
// scale mantissa, and apply sign bit for final result
Double num = fmant * scale; // 22.019999999960419 - Make the final abs number
// 0.0076923076924799716
return (mant_sign != 0) ? (-num) : num; // -22.019999999960419 - Return with sign
// 0.0076923076924799716
}
Changed code above The code above is now working correctly.
An article in Hewlett-Packard Journal, October 1978, pg. 14 gives the following information about the floating-point format used here:
Extended-precision numbers have a 39-bit signed mantissa and the same seven-bit signed exponent. All mantissas are normalized, which means they are in the ranges [½, 1) and [-1, -½).
Considered together with the information in the question this means the mantissa is a signed, two's-complement, 40-bit integer represented by the first five bytes, and must be divided by 239 to get its numerical floating-point value.
The binary exponent is stored in the last byte. Based on the description of the 40-bit mantissa field as a 39-bit mantissa in the journal article, while the exponent is described a having seven bits, together with information on the exponent sign bit from the question, it appears that the exponent is a signed, two's-complement 8-bit integer, which has been rotated left by one bit so that its sign bit winds up in bit [0] of the exponent byte. To process it, we need to rotate right one bit.
Based on the above information, I wrote the following C99 program that successfully decodes the test vectors from the question:
#include <stdio.h>
#include <stdlib.h>
#include <stdint.h>
#include <math.h>
int64_t test_data [] =
{
0xA7EB851EB90A, 0xA870A3D70A0A, 0xA8AE147AE20A, 0xA8E147AE140A,
0xA9666666670A, 0xAC70A3D70A0A, 0xACC28F5C290A, 0xACCCCCCCCC0A,
0xAE70A3D70B0A, 0xAEB851EB850A, 0x400000000002, 0x43851EB85204,
0x451EB851EB04, 0x4C7AE147AE04, 0x4EC4EC4EC5F3, 0x519999999A04,
0x5838B6BE9BF3, 0x5B851EB85202, 0x5BD70A3D7108, 0x5C7AE147AE08,
0x5E51EB851E08, 0x5FD70A3D7108, 0x62E147AE1508, 0x64A3D70A3E08,
0x666666666702, 0x67AE147AE202, 0x6B204B9E54F3, 0x733333333302,
0x762762762AF3, 0x794DFB4619F3, 0x7C74941627F3, 0x7E07E07E14F3
};
#define NBR_TEST_VECTORS (sizeof(test_data) / sizeof(int64_t))
int main (void)
{
for (int i = 0; i < NBR_TEST_VECTORS; i++) {
// extract mantissa and exponent bits
int64_t mant = test_data[i] >> 8;
int expo = test_data[i] & 0xFF;
int mant_sign = mant >> 39;
// sign extend mantissa from 40 to 64 bits, then take absolute value
mant = ((mant ^ (1LL << 39)) - (1LL << 39));
mant = llabs (mant);
// convert mantissa to floating-point
double fmant = mant * pow (2.0, -39.0);
// rotate exponent field right by 1 bit
expo = (expo >> 1) | ((expo & 1) << 7);
// sign extend exponent from 8 to 32 bits
expo = ((expo ^ (1 << 7)) - (1 << 7));
// compute scale factor from exponent field
double scale = pow (2.0, expo);
// scale mantissa, and apply sign bit for final result
double num = fmant * scale;
num = mant_sign ? -num : num;
printf ("%012llx --> % 20.11e\n", test_data[i], num);
}
return EXIT_SUCCESS;
}
The output of the above program should look as follows:
a7eb851eb90a --> -2.20200000000e+001
a870a3d70a0a --> -2.18900000000e+001
a8ae147ae20a --> -2.18300000000e+001
a8e147ae140a --> -2.17800000000e+001
a9666666670a --> -2.16500000000e+001
ac70a3d70a0a --> -2.08900000000e+001
acc28f5c290a --> -2.08100000000e+001
accccccccc0a --> -2.08000000000e+001
ae70a3d70b0a --> -2.03900000000e+001
aeb851eb850a --> -2.03200000000e+001
400000000002 --> 1.00000000000e+000
43851eb85204 --> 2.11000000000e+000
451eb851eb04 --> 2.16000000000e+000
4c7ae147ae04 --> 2.39000000000e+000
4ec4ec4ec5f3 --> 4.80769230769e-003
519999999a04 --> 2.55000000000e+000
5838b6be9bf3 --> 5.38461538456e-003
5b851eb85202 --> 1.43000000000e+000
5bd70a3d7108 --> 1.14800000000e+001
5c7ae147ae08 --> 1.15600000000e+001
5e51eb851e08 --> 1.17900000000e+001
5fd70a3d7108 --> 1.19800000000e+001
62e147ae1508 --> 1.23600000000e+001
64a3d70a3e08 --> 1.25800000000e+001
666666666702 --> 1.60000000000e+000
67ae147ae202 --> 1.62000000000e+000
6b204b9e54f3 --> 6.53846153847e-003
733333333302 --> 1.80000000000e+000
762762762af3 --> 7.21153846158e-003
794dfb4619f3 --> 7.40384615382e-003
7c74941627f3 --> 7.59615384651e-003
7e07e07e14f3 --> 7.69230769248e-003
Implementation of njuffa approach in Delphi
type
THexabyte = array[0..5] of Byte;
function ConvertHPToDouble(const Data: THexabyte): Double;
var
iexp: Integer;
pbi: PByteArray;
i64: Int64;
begin
iexp := Data[5] shr 1; //unsigned shift
if (Data[5] and 1) <> 0 then //use exp sign
iexp := - iexp;
pbi := @i64;
pbi[0] := 0;
pbi[1] := 0;
pbi[2] := 0;
pbi[3] := Data[4]; //shuffle bytes into intel order
pbi[4] := Data[3];
pbi[5] := Data[2];
pbi[6] := Data[1];
pbi[7] := Data[0];
i64 := i64 div (1 shl 23); //arithmetic shift saves sign
Result := i64 / (Int64(1) shl (40 - iexp));
end;
var
Data: THexabyte;
i: Integer;
s: string;
begin
s := 'A870A3D70A0A';
for i := 0 to 5 do
Data[i] := StrToInt('$' + Copy(s, i * 2 + 1, 2));
Memo1.Lines.Add(Format('%14.7f', [ConvertHPToDouble(Data)]));
// gives -21.8900000
Converting njuffa's
code into Delphi:
- The mantissa is left-shifted 16 bits to form a 64 bit signed integer.
- The exponent is extracted and rotated right 1 bit to form a signed 8 bit integer.
- The resulting double is scaled by multiplying the mantissa with
IntPower(2,-63+Exponent)
program HP_Float_To_Double;
{$APPTYPE CONSOLE}
uses
System.SysUtils,Math;
var
test_data : TArray<Int64> = [
$A7EB851EB90A, $A870A3D70A0A, $A8AE147AE20A, $A8E147AE140A,
$A9666666670A, $AC70A3D70A0A, $ACC28F5C290A, $ACCCCCCCCC0A,
$AE70A3D70B0A, $AEB851EB850A, $400000000002, $43851EB85204,
$451EB851EB04, $4C7AE147AE04, $4EC4EC4EC5F3, $519999999A04,
$5838B6BE9BF3, $5B851EB85202, $5BD70A3D7108, $5C7AE147AE08,
$5E51EB851E08, $5FD70A3D7108, $62E147AE1508, $64A3D70A3E08,
$666666666702, $67AE147AE202, $6B204B9E54F3, $733333333302,
$762762762AF3, $794DFB4619F3, $7C74941627F3, $7E07E07E14F3];
function Exponent( b: Byte): ShortInt;
begin
Result := (b shr 1) or ((b and 1) shl 7);
end;
function HPFloatToDouble( HP: Int64): Double;
var
Exp: ShortInt;
begin
HP := HP shl 16; // Shift left 16 bits
Exp := Exponent(PByte(@HP)[2]); // Rotate exponent right 1 position into signed 8 bit
HP := (HP and $FFFFFFFFFF000000); // Clear exponent part from mantissa
Result := HP*IntPower(2,-63+Exp); // Scale result
end;
var
I64 : Int64;
begin
for I64 in test_data do
WriteLn(Format('%x -> %20.11e',[I64, HPFloatToDouble(I64)]));
ReadLn;
end.
The output is following:
A7EB851EB90A -> -2,2020000000E+001
A870A3D70A0A -> -2,1890000000E+001
A8AE147AE20A -> -2,1830000000E+001
A8E147AE140A -> -2,1780000000E+001
A9666666670A -> -2,1650000000E+001
AC70A3D70A0A -> -2,0890000000E+001
ACC28F5C290A -> -2,0810000000E+001
ACCCCCCCCC0A -> -2,0800000000E+001
AE70A3D70B0A -> -2,0390000000E+001
AEB851EB850A -> -2,0320000000E+001
400000000002 -> 1,0000000000E+000
43851EB85204 -> 2,1100000000E+000
451EB851EB04 -> 2,1600000000E+000
4C7AE147AE04 -> 2,3900000000E+000
4EC4EC4EC5F3 -> 4,8076923077E-003
519999999A04 -> 2,5500000000E+000
5838B6BE9BF3 -> 5,3846153846E-003
5B851EB85202 -> 1,4300000000E+000
5BD70A3D7108 -> 1,1480000000E+001
5C7AE147AE08 -> 1,1560000000E+001
5E51EB851E08 -> 1,1790000000E+001
5FD70A3D7108 -> 1,1980000000E+001
62E147AE1508 -> 1,2360000000E+001
64A3D70A3E08 -> 1,2580000000E+001
666666666702 -> 1,6000000000E+000
67AE147AE202 -> 1,6200000000E+000
6B204B9E54F3 -> 6,5384615385E-003
733333333302 -> 1,8000000000E+000
762762762AF3 -> 7,2115384616E-003
794DFB4619F3 -> 7,4038461538E-003
7C74941627F3 -> 7,5961538465E-003
7E07E07E14F3 -> 7,6923076925E-003
User contributions licensed under CC BY-SA 3.0