Loading...
Searching...
No Matches
simpack.c
Go to the documentation of this file.
#include "grib2_int.h"
#include <math.h>
#include <stdlib.h>
void
simpack(float *fld, g2int ndpts, g2int *idrstmpl,
unsigned char *cpack, g2int *lcpack)
{
static g2int zero = 0;
g2int *ifld;
g2int j, nbits, imin, imax, maxdif, nbittot, left;
float bscale, dscale, rmax, rmin, temp;
double maxnum;
static float alog2 = ALOG2; /* ln(2.0) */
LOG((3, "simpack ndpts %ld", ndpts));
bscale = int_power(2.0, -idrstmpl[1]);
dscale = int_power(10.0, idrstmpl[2]);
if (idrstmpl[3] <= 0 || idrstmpl[3] > 31)
nbits = 0;
else
nbits = idrstmpl[3];
/* Find max and min values in the data. */
rmax = fld[0];
rmin = fld[0];
for (j = 1; j < ndpts; j++)
{
if (fld[j] > rmax)
rmax = fld[j];
if (fld[j] < rmin)
rmin = fld[j];
}
ifld = calloc(ndpts, sizeof(g2int));
/* If max and min values are not equal, pack up field. If they are
* equal, we have a constant field, and the reference value (rmin)
* is the value for each point in the field and set nbits to 0. */
if (rmin != rmax)
{
/* Determine which algorithm to use based on user-supplied
* binary scale factor and number of bits. */
if (nbits == 0 && idrstmpl[1] == 0)
{
/* No binary scaling and calculate minumum number of bits
* in which the data will fit. */
imin = (g2int)rint(rmin * dscale);
imax = (g2int)rint(rmax * dscale);
maxdif = imax - imin;
temp = log((double)(maxdif + 1)) / alog2;
nbits = (g2int)ceil(temp);
rmin = (float)imin;
/* scale data */
for (j = 0; j < ndpts; j++)
ifld[j] = (g2int)rint(fld[j] * dscale) - imin;
}
else if (nbits != 0 && idrstmpl[1] == 0)
{
/* Use minimum number of bits specified by user and adjust
* binary scaling factor to accomodate data. */
rmin = rmin * dscale;
rmax = rmax * dscale;
maxnum = int_power(2.0, nbits) - 1;
temp = log(maxnum / (rmax - rmin)) / alog2;
idrstmpl[1] = (g2int)ceil(-1.0 * temp);
bscale = int_power(2.0, -idrstmpl[1]);
/* scale data */
for (j = 0; j < ndpts; j++)
ifld[j] = (g2int)rint(((fld[j] * dscale) - rmin) * bscale);
}
else if (nbits == 0 && idrstmpl[1] != 0)
{
/* Use binary scaling factor and calculate minumum number
* of bits in which the data will fit. */
rmin = rmin * dscale;
rmax = rmax * dscale;
maxdif = (g2int)rint((rmax - rmin) * bscale);
temp = log((double)(maxdif + 1)) / alog2;
nbits = (g2int)ceil(temp);
/* scale data */
for (j = 0; j < ndpts; j++)
ifld[j] = (g2int)rint(((fld[j] * dscale) - rmin) * bscale);
}
else if (nbits != 0 && idrstmpl[1] != 0)
{
/* Use binary scaling factor and use minumum number of
* bits specified by user. Dangerous - may loose
* information if binary scale factor and nbits not set
* properly by user. */
rmin = rmin * dscale;
/* scale data */
for (j = 0; j < ndpts; j++)
ifld[j] = (g2int)rint(((fld[j] * dscale) - rmin) * bscale);
}
/* Pack data, Pad last octet with Zeros, if necessary, and
* calculate the length of the packed data in bytes. */
sbits(cpack, ifld, 0, nbits, 0, ndpts);
nbittot = nbits * ndpts;
left = 8 - (nbittot % 8);
if (left != 8)
{
sbit(cpack, &zero, nbittot, left); /* Pad with zeros to fill Octet. */
nbittot = nbittot + left;
}
*lcpack = nbittot / 8;
}
else
{
nbits = 0;
*lcpack = 0;
}
/* Fill in ref value and number of bits in Template 5.0. */
mkieee(&rmin, idrstmpl, 1); /* ensure reference value is IEEE format. */
idrstmpl[3] = nbits;
idrstmpl[4] = 0; /* original data were reals. */
free(ifld);
}
sbits
void sbits(unsigned char *out, g2int *in, g2int iskip, g2int nbits, g2int nskip, g2int n)
Store arbitrary size values into a packed bit string, taking the low order bits from each value in th...
Definition gbits.c:178
sbit
void sbit(unsigned char *out, g2int *in, g2int iskip, g2int nbits)
Store arbitrary size values into a packed bit string, taking the low order bits from each value in th...
Definition gbits.c:38
g2int
int64_t g2int
Long integer type.
Definition grib2.h:32
grib2_int.h
Header file with internal function prototypes NCEPLIBS-g2c library.
mkieee
void mkieee(float *a, g2int *rieee, g2int num)
Store a list of real values in 32-bit IEEE floating point format.
Definition mkieee.c:22
int_power
double int_power(double x, g2int y)
Function similar to C pow() power function.
Definition int_power.c:18
LOG
#define LOG(e)
Ignore logging to stdout.
Definition grib2_int.h:428
ALOG2
#define ALOG2
ln(2.0)
Definition grib2_int.h:30
simpack
void simpack(float *fld, g2int ndpts, g2int *idrstmpl, unsigned char *cpack, g2int *lcpack)
Packs a data field using the simple packing algorithm.
Definition simpack.c:36