#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <crblib/inc.h>

#define memcpyints(to,from,len)	memcpy(to,from,(len)*sizeof(int))

static void do_tdec_line(int * to,int *from,int len)
{
int x,*ptr,*low,*high,half,A,B;

	if ( len <= 2 ) errexit("special-cased to len > 2 , scumbag");
	if ( len&1 ) errexit("len shouldn't be odd");

	half = len>>1;

	high = to + half;
	ptr = from;
	for(x=0;x<(half-1);x++) {
		*high++ = ptr[1] - ((ptr[2] + ptr[0])>>1);
		ptr += 2;
	}
	*high = ptr[1] - ptr[0];

	low  = to;
	high = to + half;
	ptr = from;
	for(x=0;x<half;x++) {
		if ( x == 0 ) {
			A = ((high[0] )>>1) ;
		} else if ( x == 1 ||  x == (half-1)) {
			A = ((high[0] + high[-1])>>2) ;
		} else {
			A = ((19*(high[0] + high[-1]) - 3*(high[1] + high[-2]))>>6) ;
		}
		*low = *ptr + A;
		ptr += 2; high++; low++;
	}
}

static void un_tdec_line(int *to,int *from,int len)
{
int x,*ptr,*low,*high,half,A,B;

	half = len>>1;

	low  = from;
	ptr = to;
	high = from + half;
	for(x=0;x<half;x++) {
		if ( x == 0 ) {
			A = ((high[0] )>>1) ;
		} else if ( x == 1 ||  x == (half-1)) {
			A = ((high[0] + high[-1])>>2) ;
		} else {
			A = ((19*(high[0] + high[-1]) - 3*(high[1] + high[-2]))>>6) ;
		}
		*ptr = *low - A;
		ptr += 2; high++; low ++;
	}

	high = from + half;
	ptr = to;
	for(x=0;x<(half-1);x++) {
		ptr[1] = (*high++) + ((ptr[2] + ptr[0])>>1);
		ptr += 2;
	}
	ptr[1] = *high + ptr[0];
}

void cdf24_2D(int **rows, int width, int height, int levels,bool inverse)
{
int x, y, w, h, l;
int *buffer,*tempbuf,*temprow;

    if (width%(1 << (levels+1)) || height%(1 << (levels+1)))
		errexit("width and height must be divisible by 2^(levels+1)");
  
    /* Allocate a work array (for transposing columns) */
    
   	if ( (buffer = newarray(int,height+max(width,height)+height)) == NULL )
		errexit("malloc failed");
	temprow = buffer+height;
	tempbuf = buffer+height+height;

	if ( !inverse ) {

		for (l = 0; l < levels; l++) {
			w = width >> l;
			h = height >> l;
      
			/* Rows */
	
			do_tdec_line(temprow,rows[h-1],w);
			for (y = h-2; y >=0; y--) {
				do_tdec_line(rows[y+1],rows[y],w);
			}
	
			/* Columns */
	
			for (x = 0; x < w; x++) {
					for (y = 1; y < h; y++) buffer[y-1] = rows[y][x];
					buffer[h-1] = temprow[x];
				do_tdec_line(tempbuf,buffer,h);
					for (y = 0; y < h; y++) rows[y][x] = tempbuf[y];
			}
		}

    } else {

		for (l = levels-1; l >= 0; l--) { /** backwards in scale **/
			w = width >> l;
			h = height >> l;

			/* Columns */
	
			for (x = 0; x < w; x++) {
					for (y = 0; y < h; y++) buffer[y] = rows[y][x];
				un_tdec_line(tempbuf,buffer,h);
					for (y = 0; y < h-1; y++) rows[y+1][x] = tempbuf[y];
					temprow[x] = tempbuf[h-1];
			}

			/* Rows */
			for (y = 0; y < h-1; y++) {
				un_tdec_line(rows[y],rows[y+1],w);
			} 
			un_tdec_line(rows[h-1],temprow,w);
		}
	}

	free(buffer);
}