/*
C=============================================================================C
C Comprehensive Ocean-Atmosphere Data Set (COADS):                  C Library C
C Filename:level: gsbytes.c:v 2.1 86/11/03 10:28:23               18 May 1998 C
C Function: Bit-string packing/unpacking            Author: S.Woodruff et al. C
C=============================================================================C
C Software documentation for the (invariant) user-interface routines {gbytes,
C gbyte,sbytes,sbyte}:
C
C A C language implementation of the Fortran subroutines designed by Jenne and
C Joseph (1974).  {gbytes} unpacks n b bit bytes from p into u, starting by
C skipping q bits in p, then skipping s bits between bytes.  {gbyte} unpacks
C one such byte.  {sbytes} packs n b bit bytes from u into p, starting by
C skipping q bits in p, then skipping s bits between bytes. {sbyte} packs one
C such byte.  SWORD must be adjusted if necessary for the computer word size
C (by default 32).  Following are the constraints on q, b, s, and n:
C      q >= 0     number of bits to be skipped preceeding first byte in p
C  0 < b < sword  byte size
C      s >= 0     number of bits to be skipped between bytes
C      n >= 0     number of bytes to be packed/unpacked
C
C Reference: Jenne, R.L., and D.H. Joseph, 1974: Techniques for the processing,
C storage, and exchange of data.  NCAR Tech. Note IA-93, National Center for
C Atmospheric Research, Boulder, Colo., 46 pp.
C
C Note: Owing to authorship of an original code version, we distribute without
C copyright; a copyright notice has been attached by UCAR.  Implementations of
C these routines in Fortran and for different computer systems, plus additional
C documentation, is available on the Web at ftp://ncardata.ucar.edu/libraries/
C or via anonymous ftp at ncardata.ucar.edu, cd libraries/gbytes.
C
C Machine dependencies: Word size (SWORD) is by default 32 bits.
C For more information: See <soft_info> (electronic document).
C-----------------------------------------------------------------------3456789
 */
#define MASK	( (unsigned) 0xffffffff)
#define	SWORD	32		/* Word size in bits  */
/*
C=============================================================================C
C WARNING:  Code beyond this point should not require any modification.       C
C=============================================================================C
 */
#define G1BYTE(p,q,b) ((unsigned) p << q >> (SWORD - b))
/* right justify b bits, located q bits from left of p.  */

gbytes_(p, u, q, b, s, n)
	long           *p, *u, *q, *b, *s, *n;
{
	register long   i = 0, jp = 0;
	long            jq = *q;
	if (*n > 0) {
		while (1) {
			gbyte_(p + jp, u + i, &jq, b);
			if (++i == *n)
				break;
			jq += *b + *s;
			jp += jq / SWORD;
			jq %= SWORD;
		}
	}
}


sbytes_(p, u, q, b, s, n)
	long           *p, *u, *q, *b, *s, *n;
{
	register long   i = 0, jp = 0;
	long            jq = *q;
	if (*n > 0) {
		while (1) {
			sbyte_(p + jp, u + i, &jq, b);
			if (++i == *n)
				break;
			jq += *b + *s;
			jp += jq / SWORD;
			jq %= SWORD;
		}
	}
}


gbyte_(p, u, q, b)
	long           *p, *u, *q, *b;
{
	register long   j, jq = *q, jb = *b, lb, qb;

	if (jq < SWORD)
		j = 0;
	else {
		j = jq / SWORD;	/* number of words offset */
		jq %= SWORD;	/* odd bits of offset     */
	}

	if ((jq + jb) > SWORD) {
		qb = SWORD - jq;
		jb -= qb;
		lb = ((~(MASK << qb)) & (*(p + j))) << jb;
		jq = 0;
		j++;		/* increment to next word */
		*u = lb + (G1BYTE(*(p + j), jq, jb));
	} else
		*u = (G1BYTE(*(p + j), jq, jb));
}


sbyte_(p, u, q, b)
	long           *p, *u, *q, *b;
{
	register long  *t, jq = *q, jb = *b, qb;

	if (jq < SWORD)
		t = p;
	else {
		t = p + jq / SWORD;	/* number of words offset */
		jq %= SWORD;	/* odd bit offset         */
	}

	if ((jq + jb) > SWORD) {
		qb = SWORD - jq;
		jq = SWORD - jb;
		jb -= qb;
		*t = ((*t >> qb) << qb) + (G1BYTE(*u, jq, qb));
		jq = 0;
		t++;		/* point to next word */
	}
	/*
	 * *t = (*t & ~(MASK << (SWORD - jb) >> jq)) + ((~(MASK << jb) & *u)
	 * << SWORD-(jb+jq));
	 */
	*t = (*t & ~(MASK << (SWORD - jb) >> jq)) + (((unsigned) *u) <<
						     (SWORD - jb) >> jq);
}