#include "grib2.h"

void rdieee(g2int *rieee,g2float *a,g2int num)
//$$$  SUBPROGRAM DOCUMENTATION BLOCK
//                .      .    .                                       .
// SUBPROGRAM:    rdieee
//   PRGMMR: Gilbert         ORG: W/NP11    DATE: 2002-10-25
//
// ABSTRACT: This subroutine reads a list of real values in
//   32-bit IEEE floating point format.
//
// PROGRAM HISTORY LOG:
// 2002-10-25  Gilbert
//
// USAGE:    void rdieee(g2int *rieee,g2float *a,g2int num)
//   INPUT ARGUMENT LIST:
//     rieee    - g2int array of floating point values in 32-bit IEEE format.
//     num      - Number of floating point values to convert.
//
//   OUTPUT ARGUMENT LIST:
//     a        - float array of real values.  a must be allocated with at
//                least 4*num bytes of memory before calling this function.
//
// REMARKS: None
//
// ATTRIBUTES:
//   LANGUAGE: C
//   MACHINE:  IBM SP
//
//$$$
{

      g2int  j;
      g2int  isign,iexp,imant;

      g2float  sign,temp;
      static const double two23 = 1.1920928955078125e-07; // pow(2,-23)
      static const double two126 = 1.1754943508222875e-38; // pow(2,-126)


      g2intu msk1=0x80000000;        // 10000000000000000000000000000000 binary
      g2int msk2=0x7F800000;         // 01111111100000000000000000000000 binary
      g2int msk3=0x007FFFFF;         // 00000000011111111111111111111111 binary

      for (j=0;j<num;j++) {
//
//  Extract sign bit, exponent, and mantissa
//
        isign=(rieee[j]&msk1)>>31;
        iexp=(rieee[j]&msk2)>>23;
        imant=(rieee[j]&msk3);
        //printf("SAGieee= %ld %ld %ld\n",isign,iexp,imant);

        sign=1.0;
        if (isign == 1) sign=-1.0;

        if ( (iexp > 0) && (iexp < 255) ) {
          temp=(g2float)int_power(2.0,(iexp-127));
          a[j]=(float)(sign*temp*(1.0+(two23*(g2float)imant)));
        }
        else if ( iexp == 0 ) {
          if ( imant != 0 )
            a[j]=(float)(sign*two126*two23*(g2float)imant);
          else
            a[j]=(float)(sign*0.0);

        }
        else if ( iexp == 255 )
           a[j]=(float)(sign*(1E+37));


      }

}