/*
  FP8_E3M4 playground header file - functions

  Format: s.eee.mmmm (bias 3)
  Value = (-1)^s*(1+mmmm/16)*2^(eee-3)

  Special cases:
    eee==000 - unnormalized mantissa == (0+0.mmmm)
    0.000.0000 == 0
    s.111.0000 == +/-inf (signed infinity)
    s.111.1xxx == nan (not a number)

  Ranges of values :
    0.000.0001 - smallest unnormalized: (0+1/16)*2^(1-3) == 0.015625  [1/16/4]
    0.000.1111 - largest unnormalized:  (0+15/16)*2^(1-3) == 0.234375 [15/16/4]

    0.001.0000 - smallest normalized:   (1+0/16)*2^(1-3) == 0.25      [16/16/4]
    0.110.1111 - largest normalized:    (1+15/16)*2^(6-3) == 15.5     [31/16*8]

  ------------------------------------------------------------------------
  (C) 2026 Peeter Ellervee <peeter.ellervee@taltech.ee>
*/

#ifdef DEBUGGING
#include <stdio.h>
#endif /*DEBUGGING*/

#include "fp8_e3m4.h"


/*
 * Unpacking fp8_e3m4 into structure fp8_e3m4_str
 */
void fp8_e3m4_unpack(const fp8_e3m4 num, fp8_e3m4_str* data) {
  char sign=FP8_E3M4_SIGN(num);			/* Sign (8-bit int)     */
  char exp= FP8_E3M4_EXP(num);			/* Exponent (8-bit int) */
  char mant=FP8_E3M4_MANT(num);			/* Mantissa (8-bit int) */
  if (exp>0) {					/* Normalized mantissa? */
    mant|=FP8_E3M4_NORM;			/* Adjusting mantissa   */
    exp-=3;					/* Exponent (-bias)     */
  }
  else  exp=-2;					/* Exponent == -2       */
  data->sign=sign;  data->exp=exp;  data->mant=mant;
}


/*
 * Packing structure fp8_e3m4_str into fp8_e3m4
 * Checking & normalizing/unnormalizing...
 */
fp8_e3m4 fp8_e3m4_pack(const fp8_e3m4_str* data) {
  int sign_int, exp_int, mant_int;
  sign_int=data->sign;  exp_int=data->exp;  mant_int=data->mant;
  if (mant_int<0)	 			/* Absolute value & sign */
    { sign_int=0x80; mant_int=-mant_int; }
  else if (sign_int!=0) sign_int=0x80;
#ifdef DEBUGGING
  printf("**fp8_e3m4_pack*1* %d %d %d\n",sign_int,exp_int,mant_int);
#endif /*DEBUGGING*/
  /***** Adjusting mantissa: 1.0 <= mant < 2.0 *****/
  while (mant_int<16) {				/* Mantissa too small? */
    if (mant_int==0) { exp_int=-3; break; }	/* Zero? */
      mant_int=mant_int<<1; exp_int--;
  }
  while (mant_int>31)				/* Mantissa too large? */
    { mant_int=mant_int>>1; exp_int++; }
  /***** Adjusting exponent: -2 <= exp <= +3 *****/
  if (exp_int>3) { exp_int=4;  mant_int=0; }	/* +/- infinity */
  else if (exp_int<-6) { exp_int=-3; mant_int=0; }  /* Underflow? */
  else if (exp_int<-2) {			/* Unnormalized? */
    while (exp_int<-2)
      { mant_int=mant_int>>1; exp_int++; }	/* Unnormalizing... */
    exp_int=-3;
  }
#ifdef DEBUGGING
  printf("**fp8_e3m4_pack*2* %d %d %d\n",sign_int,exp_int,mant_int);
#endif /*DEBUGGING*/
  exp_int=((exp_int+3)<<4)&0x70;		/* Adjusting exponent (+bias) */
  return sign_int|exp_int|(mant_int&0x0f);	/* Packing... */
}


/*
 * Converting fp8_e3m4 to float
 */
float fp8_to_float(const fp8_e3m4 num) {
  fp8_e3m4_str val;
  fp8_e3m4_unpack(num,&val);
  float sign=(val.sign==1)?-1.0:1.0;		/* Sign (float)   */
  if (val.exp==4)				/* Special cases? */
    return val.mant&0x8?NAN:sign*HUGE_VALF;	/* NaN / infinity */
  float exp=val.exp-4;				/* Exponent (float) & /16.0 */
  return sign*(float)val.mant*pow(2.0,exp);	/* (-1)^s*(1+m/16)*2^(e-3)  */
}


/*
 * Converting float to fp8_e3m4
 */
fp8_e3m4 float_to_fp8(const float num) {
  if (num==HUGE_VALF) return FP8_INFP;		/* Plus infinity? */
  else if (num==-HUGE_VALF) return FP8_INFN;	/* Minus infinity? */
  else if (isnan(num)) return FP8_NAN;		/* NAN? */
  
  fp8_e3m4_str data;  data.sign=0;	/* Inital values... */
  float mant=frexp(num,&data.exp);	/* Extracting mantissa & exponent */
  mant*=2.0; data.exp--;		/* 1.0 <= mantissa < 2.0 */
  data.mant=(mant*32.0);		/* Fix point: 1.nnnnr (for rounding) */
  if (data.mant&1) {				/* Rounding? */
    if (data.mant<0) data.mant=(data.mant-1)>>1;
    else data.mant=(data.mant+1)>>1;
  }
  else data.mant=data.mant>>1;
#ifdef DEBUGGING
  printf("**float_to_fp8** %d %d %d %f\n",data.sign,data.exp,data.mant,mant);
#endif /*DEBUGGING*/
  return fp8_e3m4_pack(&data);
}