fractaldecode.c

分形图像压缩源代码

/*   FractalDecode.c is a modified version of dec.c by Yuval Fisher  */

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#include <png.h>

#include "FractalSqueeze.h"

static IMAGE_TYPE **image,*imptr,**image1; 
                             /* The input image data and a dummy         */ 

static double max_scale = 1.0;      /* maximum allowable grey level scale factor */

static int    s_bits = 5,           /* number of bits used to store scale factor */
       o_bits = 7,           /* number of bits used to store offset       */
       hsize = -1,           /* The horizontal size of the input image    */
       vsize = -1,           /* The vertical size                         */
       scaledhsize,          /* hsize*scalefactor                         */
       scaledvsize,          /* vsize*scalefactor                         */
       size,                 /* largest square image that fits in image   */
       min_part = 3,         /* min and max _part determine a range of    */
       max_part = 4,         /* Range sizes from hsize>>min to hsize>>max */
       dom_step = 4,         /* Density of domains relative to size       */
       dom_step_type = 0,    /* Flag for dom_step a multiplier or divisor */
       dom_type = 0,         /* Method of generating domain pool 0,1,2..  */
       post_process = 1,     /* Flag for postprocessing.                  */
       num_iterations= 10,   /* Number of decoding iterations used.       */
       *no_h_domains,        /* Number of horizontal domains.             */
       *domain_hstep,        /* Domain density step size.                 */
       *domain_vstep,        /* Domain density step size.                 */
       *bits_needed,         /* Number of bits to encode domain position. */
       zero_ialpha,          /* the const ialpha when alpha = 0           */
       output_partition=0,   /* A flag for outputing the partition        */
       max_exponent;         /* The max power of 2 side of square image   */
                             /* that fits in our input image.             */

struct transformation_node {
       int rx,ry,            /* The range position and size in a trans.   */
              xsize, ysize,  
              rrx,rry,
              dx,dy;         /* The domain position.                      */
       int sym_op;           /* The symmetry operation used in the trans. */
       int depth;            /* The depth in the quadtree partition.      */
       double scale, offset; /* scalling and offset values.               */
       struct transformation_node *next;       /* The next trans. in list */
} transformations, *trans; 

static FILE *input;

//Prototypes:
void scale_transformations(double scalefactor);
void partition_image_read(int atx, int aty, int hsize, int vsize);
void apply_transformations(void);
void smooth_image(void);
void read_transformations(int atx, int aty, int xsize, int ysize, int depth);

/* ************************************************************ */
/* unpack value using size bits read from fin.                  */
/* ************************************************************ */
long unpack(int size, FILE* fin)
{
    int i;
    int value = 0;
    static int ptr = 1; /* how many bits are packed in sum so far */
    static int sum; 


    /* size == -2 means we initialize things */
    if (size == -2) {
        sum = fgetc(fin);
        sum <<= 1;
        return((long)0);
    }

    /* size == -1 means we want to peek at the next bit without */
    /* advancing the pointer */
    if (size == -1)
        return((long)((sum&256)>>8));

     for (i=0; i<size; ++i, ++ptr,  sum <<= 1) {
        if (sum & 256) value |= 1<<i;

         if (ptr == 8) {
            sum = getc(fin);
            ptr=0;
        }
     }
    return((long)value);
}

void fractal_decode(png_infop image_data, char* inputfilename)
{
    /* Defaults are set initially */
    double         scalefactor = 1.0;           /* Scale factor for output */
    int            i, x_exponent, y_exponent;
    int            domain_size, no_domains;
 
    if ((input = fopen(inputfilename, "rb")) == NULL)
        ERROROUT("Can't open input file.");
        
    unpack(-2,input); /* initialize the unpacking routine */

	/* read the header data from the input file. This should probably */
    /* be put into one read which reads a structure with the info     */
    min_part = (int)unpack(4,input);
    max_part = (int)unpack(4,input);
    dom_step = (int)unpack(4,input);
    dom_step_type = (int)unpack(1,input);
    dom_type = (int)unpack(2,input);
    hsize = (int)unpack(12,input);
    vsize = (int)unpack(12,input);

    /* we now compute size */
    x_exponent = (int)floor(log((double)hsize)/log(2.0));
    y_exponent = (int)floor(log((double)vsize)/log(2.0));
   
    /* exponent is min of x_ and y_ exponent */
    max_exponent = (x_exponent > y_exponent ? y_exponent : x_exponent);
    /* size is the size of the largest square that fits in the image */
    /* It is used to compute the domain and range sizes.             */
    size =  1<<max_exponent; 

    /* This is the quantized value of zero scaling */
    zero_ialpha = (int)(0.5 + (max_scale)/(2.0*max_scale)*(1<<s_bits));

    /* allocate memory for the output image. Allocating one chunck saves  */
    /* work and time later.                                              */
	scaledhsize = (int)(scalefactor*hsize);
	scaledvsize = (int)(scalefactor*vsize);
    matrix_allocate(image, scaledhsize,scaledvsize, IMAGE_TYPE);
    matrix_allocate(image1, scaledhsize, scaledvsize, IMAGE_TYPE);

    /* since max_ and min_ part are variable, these must be allocated */
    i = max_part - min_part + 1;
    bits_needed = (int *)malloc(sizeof(int)*i);
    no_h_domains = (int *)malloc(sizeof(int)*i);
    domain_hstep = (int *)malloc(sizeof(int)*i);
    domain_vstep = (int *)malloc(sizeof(int)*i);
    
    /* compute bits needed to read each domain type */
    for (i=0; i <= max_part-min_part; ++i) {
       /* first compute how many domains there are horizontally */
       domain_size = size >> (min_part+i-1);
       if (dom_type == 2) 
             domain_hstep[i] = dom_step; 
       else if (dom_type == 1)
             if (dom_step_type ==1)
                domain_hstep[i] = (size >> (max_part - i-1))*dom_step;
             else 
                domain_hstep[i] = (size >> (max_part - i-1))/dom_step;
       else 
             if (dom_step_type ==1)
                domain_hstep[i] = domain_size*dom_step;
             else 
                domain_hstep[i] = domain_size/dom_step;

       no_h_domains[i] = 1+(hsize-domain_size)/domain_hstep[i];
       /* bits_needed[i][0] = ceil(log(no_domains)/log(2));  */

       /* now compute how many domains there are vertically */
       if (dom_type == 2) 
             domain_vstep[i] = dom_step; 
       else if (dom_type == 1)
             if (dom_step_type ==1)
             domain_vstep[i] = (size >> (max_part - i-1))*dom_step;
             else 
             domain_vstep[i] = (size >> (max_part - i-1))/dom_step;
       else 
             if (dom_step_type ==1)
                domain_vstep[i] = domain_size*dom_step;
             else 
                domain_vstep[i] = domain_size/dom_step;

       no_domains = 1+(vsize-domain_size)/domain_vstep[i];
       bits_needed[i] = (int)(ceil(log((double)no_domains*(double)no_h_domains[i])/log(2.0))); 
     }

    /* Read in the transformation data */
    trans = &transformations;
    printf("\nReading transformations.....");
    fflush(stdout);
    partition_image_read(0, 0, hsize,vsize );
    fclose(input);
    printf("Done.");
    fflush(stdout);
   
	if (scalefactor != 1.0) {
		printf("\nScaling image to %d x %d.", scaledhsize,scaledvsize);
		scale_transformations(scalefactor); 
	}
 
    for (i=0; i<num_iterations; ++i) {
      apply_transformations();
	}

    if (post_process) 	
      smooth_image();

	// copy resultant image to the main buffer
	{
		unsigned int i,j;
		for (i = 0; i < image_data->width; ++i) {
			for (j = 0; j < image_data->height; ++j) {
				image_data->row_pointers[j][i] = image[i][j];
			}
		}
	}


    printf("\n%d pixels written to output file.\n", i);
}

/* ************************************************************ */
/* Read in the transformation data from *input.                 */
/* This is a recursive routine whose recursion tree follows the */
/* recursion done by the encoding program.                      */
/* ************************************************************ */
void read_transformations(int atx, int aty, int xsize, int ysize, int depth)
{
    /* Having all these locals in a recursive procedure is hard on the  */
    /* stack.. but it is more readable.                                 */
	int ialpha,                   /* Intgerized scalling factor         */
        ibeta;                    /* Intgerized offset                  */

    long domain_ref; 

    double alpha, beta;