/* Copyright (C) 1994 Aladdin Enterprises.  All rights reserved.
  
  This file is part of Aladdin Ghostscript.
  
  Aladdin Ghostscript is distributed with NO WARRANTY OF ANY KIND.  No author
  or distributor accepts any responsibility for the consequences of using it,
  or for whether it serves any particular purpose or works at all, unless he
  or she says so in writing.  Refer to the Aladdin Ghostscript Free Public
  License (the "License") for full details.
  
  Every copy of Aladdin Ghostscript must include a copy of the License,
  normally in a plain ASCII text file named PUBLIC.  The License grants you
  the right to copy, modify and redistribute Aladdin Ghostscript, but only
  under certain conditions described in the License.  Among other things, the
  License requires that the copyright notice and this notice be preserved on
  all copies.
*/

/* shcgen.c */
/* Generate (bounded) Huffman code definitions from frequencies, */
/* and tables from definitions. */
#include "memory_.h"
#include "stdio_.h"
#include <stdlib.h>		/* for qsort */
#include "gdebug.h"
#include "gserror.h"
#include "gserrors.h"
#include "gsmemory.h"
#include "scommon.h"
#include "shc.h"
#include "shcgen.h"

/* ------ Frequency -> definition procedure ------ */

/* Define a node for the Huffman code tree. */
typedef struct count_node_s count_node;
struct count_node_s {
	long freq;		/* frequency of value */
	uint value;		/* data value being encoded */
	uint code_length;	/* length of Huffman code */
	count_node *next;	/* next node in freq-sorted list */
	count_node *left;	/* left child in tree (smaller code_length) */
	count_node *right;	/* right child in tree (greater code_length) */
};

#ifdef DEBUG
#  define debug_print_nodes(nodes, n, tag)\
     if ( gs_debug_c('W') ) print_nodes_proc(nodes, n, tag);
private void
print_nodes_proc(const count_node *nodes, int n, const char *tag)
{	int i;
	dprintf1("[w]---------------- %s ----------------\n", tag);
	for ( i = 0; i < n; ++i )
	  dprintf7("[w]node %d: f=%ld v=%d len=%d N=%d L=%d R=%d\n",
		   i, nodes[i].freq, nodes[i].value, nodes[i].code_length,
		   (nodes[i].next == 0 ? -1 : (int)(nodes[i].next - nodes)),
		   (nodes[i].left == 0 ? -1 : (int)(nodes[i].left - nodes)),
		   (nodes[i].right == 0 ? -1 : (int)(nodes[i].right - nodes)));
}
#else
#  define debug_print_nodes(nodes, n, tag) DO_NOTHING
#endif

/* Node comparison procedures for sorting. */
#define pn1 ((const count_node *)p1)
#define pn2 ((const count_node *)p2)
/* Sort by decreasing frequency. */
private int
compare_freqs(const void *p1, const void *p2)
{	long diff = pn2->freq - pn1->freq;
	return (diff < 0 ? -1 : diff > 0 ? 1 : 0);
}
/* Sort by increasing code length, and secondarily by decreasing frequency. */
private int
compare_code_lengths(const void *p1, const void *p2)
{	int diff = pn1->code_length - pn2->code_length;
	return (diff < 0 ? -1 : diff > 0 ? 1 : compare_freqs(p1, p2));
}
#undef pn1
#undef pn2

/* Adjust code lengths so that none of them exceeds max_length. */
/* We break this out just to help organize the code; it's only called */
/* from one place in hc_compute. */
private void
hc_limit_code_lengths(count_node *nodes, uint num_values, int max_length)
{	int needed;		/* # of max_length codes we need to free up */
	count_node *longest = nodes + num_values;

	{	/* Compute the number of additional max_length codes */
		/* we need to make available. */
		int length = longest[-1].code_length;
		int next_length;
		int avail = 0;
		while ( (next_length = longest[-1].code_length) > max_length )
		  {	avail >>= length - next_length;
			length = next_length;
			(--longest)->code_length = max_length;
			++avail;
		  }
		needed = (nodes + num_values - longest) -
		  (avail >> (length - max_length));
	}
	/* Skip over all max_length codes. */
	while ( longest[-1].code_length == max_length )
	  --longest;

	/*
	 * To free up a code of length N, suppose that the next shortest
	 * used code is of length M.  We take the lowest-frequency code
	 * of length M and change it to M+1; we then take the N-M-1
	 * highest-frequency codes of length N and change them to lengths
	 * M+1 ... N-1.  In the present situation, N = max_length.
	 */
	for ( ; needed > 0; --needed )
	  {	/* longest points to the first code of length max_length. */
		/* Since codes are sorted by increasing code length, */
		/* longest-1 is the desired code of length M. */
		int M = ++(longest[-1].code_length);
		if ( M == max_length )
		  --longest;
		else
		  while ( M < max_length )
		    longest++->code_length = M++;
	  }
}

/* Compute an optimal Huffman code from an input data set. */
/* The client must have set all the elements of *def. */
int
hc_compute(hc_definition *def, long *freqs, gs_memory_t *mem)
{	uint num_values = def->num_values;
	count_node *nodes =
	  (count_node *)gs_alloc_byte_array(mem, num_values * 2 - 1,
					    sizeof(count_node), "hc_compute");
	int i;
	count_node *lowest;
	count_node *comb;

	if ( nodes == 0 )
	  return_error(gs_error_VMerror);

	/* Create leaf nodes for the input data. */
	for ( i = 0; i < num_values; ++i )
	  nodes[i].freq = freqs[i], nodes[i].value = i;

	/* Create a list sorted by increasing frequency. */
	/* Also initialize the tree structure. */
	qsort(nodes, num_values, sizeof(count_node), compare_freqs);
	for ( i = 0; i < num_values; ++i )
	  nodes[i].next = &nodes[i - 1],
	  nodes[i].code_length = 0,
	  nodes[i].left = nodes[i].right = 0;
	nodes[0].next = 0;
	debug_print_nodes(nodes, num_values, "after sort");

	/* Construct the Huffman code tree. */
	for ( lowest = &nodes[num_values - 1], comb = &nodes[num_values]; ;
	      ++comb
	    )
	  {	count_node *pn1 = lowest;
		count_node *pn2 = pn1->next;
		long freq = pn1->freq + pn2->freq;

		lowest = pn2->next;
		comb->freq = freq;
		if ( pn1->code_length <= pn2->code_length )
		  comb->left = pn1, comb->right = pn2,
		  comb->code_length = pn2->code_length + 1;
		else
		  comb->left = pn2, comb->right = pn1,
		  comb->code_length = pn1->code_length + 1;
		if ( lowest == 0 )	/* no nodes left to combine */
		  break;
		/* Insert comb in the sorted list. */
		if ( freq < lowest->freq )
		  comb->next = lowest, lowest = comb;
		else
		  {	count_node *here = lowest;
			while ( here->next != 0 && freq >= here->freq )
			  here = here->next;
			comb->next = here->next;
			here->next = comb;
		  }
	  }

	/* comb (i.e., &nodes[num_values * 2 - 2] is the root of the tree. */
	/* Note that the left and right children of an interior node */
	/* were constructed before, and therefore have lower indices */
	/* in the nodes array than, the parent node.  Thus we can assign */
	/* the code lengths (node depths) in a single descending-order */
	/* sweep. */
	comb->code_length = 0;
	for ( ++comb; comb > nodes + num_values; )
	  {	--comb;
		comb->left->code_length = comb->right->code_length =
		  comb->code_length + 1;
	  }
	debug_print_nodes(nodes, num_values * 2 - 1, "after combine");

	/* Sort the leaves again by code length. */
	qsort(nodes, num_values, sizeof(count_node), compare_code_lengths);
	debug_print_nodes(nodes, num_values, "after re-sort");

	/* Limit the code length to def->num_counts. */
	hc_limit_code_lengths(nodes, num_values, def->num_counts);
	debug_print_nodes(nodes, num_values, "after limit");

	memset(def->counts, 0, sizeof(*def->counts) * (def->num_counts + 1));
	for ( i = 0; i < num_values; ++i )
	  {	def->values[i] = nodes[i].value;
		def->counts[nodes[i].code_length]++;
	  }

	/* All done, release working storage. */
	gs_free_object(mem, nodes, "hc_compute");
	return 0;
}

/* ------ Definition -> table procedures ------ */

/* Generate the encoding table from the definition. */
/* The size of the encode array is def->num_values. */
void
hc_make_encoding(hce_code *encode, const hc_definition *def)
{	uint next = 0;
	const ushort *pvalue = def->values;
	uint i, k;
	for ( i = 1; i <= def->num_counts; i++ )
	  {	for ( k = 0; k < def->counts[i]; k++, pvalue++, next++ )
		  {	hce_code *pce = encode + *pvalue;
			pce->code = next;
			pce->code_length = i;
		  }
		next <<= 1;
	  }
}

/* We decode in two steps, first indexing into a table with */
/* a fixed number of bits from the source, and then indexing into */
/* an auxiliary table if necessary.  (See shc.h for details.) */

/* Calculate the size of the decoding table. */
uint
hc_sizeof_decoding(const hc_definition *def, int initial_bits)
{	uint size = 1 << initial_bits;
	uint carry = 0, mask = ~1;
	uint i;
	for ( i = initial_bits + 1; i <= def->num_counts;
	      i++, carry <<= 1, mask <<= 1
	    )
	  {	carry += def->counts[i];
		size += carry & mask;
		carry &= ~mask;
	  }
	return size;
}

/* Generate the decoding tables. */
void
hc_make_decoding(hcd_code *decode, const hc_definition *def,
  int initial_bits)
{	/* Make entries for single-dispatch codes. */
	{	hcd_code *pcd = decode;
		const ushort *pvalue = def->values;
		uint i, k, d;
		for ( i = 0; i <= initial_bits; i++ )
		  {	for ( k = 0; k < def->counts[i]; k++, pvalue++ )
			  {	for ( d = 1 << (initial_bits - i); d > 0;
				      d--, pcd++
				    )
				  pcd->value = *pvalue,
				  pcd->code_length = i;
			  }
		  }
	}
	/* Make entries for two-dispatch codes. */
	/* By working backward, we can do this more easily */
	/* in a single pass. */
	{	uint dsize = hc_sizeof_decoding(def, initial_bits);
		hcd_code *pcd = decode + (1 << initial_bits);
		hcd_code *pcd2 = decode + dsize;
		const ushort *pvalue = def->values + def->num_values;
		uint entries_left = 0, slots_left = 0, mult_shift = 0;
		uint i = def->num_counts + 1, j;
		for ( ; ; )
		  {	if ( slots_left == 0 )
			  {	if ( entries_left != 0 )
				  {	slots_left = 1 << (i - initial_bits);
					mult_shift = 0;
					continue;
				  }
				if ( --i <= initial_bits )
					break;
				entries_left = def->counts[i];
				continue;
			  }
			if ( entries_left == 0 )
			  {	entries_left = def->counts[--i];
				mult_shift++;
				continue;
			  }
			--entries_left, --pvalue;
			for ( j = 1 << mult_shift; j > 0; j-- )
			  {	--pcd2;
				pcd2->value = *pvalue;
				pcd2->code_length = i - initial_bits;
			  }
			if ( (slots_left -= 1 << mult_shift) == 0 )
			  {	--pcd;
				pcd->value = pcd2 - decode;
				pcd->code_length = i + mult_shift;
			  }
		  }
	}
}
