MEDIAN.C
上传用户:alisonmail
上传日期:2013-02-28
资源大小:500k
文件大小:10k
- /*
** File: median.c Copyright (c) Truda Software
** Author: Anton Kruger 215 Marengo Rd, #2,
** Date: March 1992 Oxford, IA 52322-9383
** Revision: 1.0 March 1992
**
** Description: Contains an implementation of Heckbert's median-
** cut color quantization algorithm.
**
** Compilers: MSC 5.1, 6.0.
**
** Note: 1) Compile in large memory model.
** 2) Delete the "#define FAST_REMAP" statement below
** in order to deactivate fast remapping.
*/
#include <windows.h>
#define FAST_REMAP
#include <stdio.h>
#include <stddef.h> /* for NULL */
#include <stdlib.h> /* for "qsort" */
#include <float.h> /* for FLT_MAX, FLT_MIN */
#define MAXCOLORS 256 /* maximum # of output colors */
#define HSIZE 32768 /* size of image histogram */
//typedef unsigned char BYTE; /* range 0-255 */
//typedef unsigned short WORD; /* range 0-65,535 */
//typedef unsigned long DWORD; /* range 0-4,294,967,295 */
/* Macros for converting between (r,g,b)-colors and 15-bit */
/* colors follow. */
#define _RGB(r,g,b) (WORD)(((b)&~7)<<7)|(((g)&~7)<<2)|((r)>>3)
#define RED(x) (BYTE)(((x)&31)<<3)
#define GREEN(x) (BYTE)((((x)>>5)&255)<< 3)
#define BLUE(x) (BYTE)((((x)>>10)&255)<< 3)
typedef struct { /* structure for a cube in color space */
WORD lower; /* one corner's index in histogram */
WORD upper; /* another corner's index in histogram */
DWORD count; /* cube's histogram count */
int level; /* cube's level */
BYTE rmin,rmax;
BYTE gmin,gmax;
BYTE bmin,bmax;
} cube_t;
static cube_t list[MAXCOLORS]; /* list of cubes */
static int longdim; /* longest dimension of cube */
static WORD HistPtr[HSIZE]; /* points to colors in "Hist" */
void Shrink(cube_t * Cube);
void InvMap(WORD * Hist, BYTE ColMap[][3],WORD ncubes);
int compare(const void * a1, const void * a2);
WORD __declspec (dllexport) MedianCut(WORD Hist[],BYTE ColMap[][3], int maxcubes)
{
/*
** Accepts "Hist", a 32,768-element array that contains 15-bit
** color counts of the input image. Uses Heckbert's median-cut
** algorithm to divide the color space into "maxcubes" cubes,
** and returns the centroid (average value) of each cube in
** "ColMap". "Hist" is also updated so that it functions as an
** inverse color map. MedianCut returns the actual number of
** cubes, which may be less that "maxcubes".
*/
BYTE lr,lg,lb;
WORD i,median,color;
DWORD count;
int k,level,ncubes,splitpos;
void *base;
size_t num,width;
cube_t Cube,CubeA,CubeB;
/*
** Create the initial cube, which is the whole RGB-cube.
*/
ncubes = 0;
Cube.count = 0;
for (i=0,color=0;i<=HSIZE-1;i++){
if (Hist[i] != 0){
HistPtr[color++] = i;
Cube.count = Cube.count + Hist[i];
}
}
Cube.lower = 0; Cube.upper = color-1;
Cube.level = 0;
Shrink(&Cube);
list[ncubes++] = Cube;
/*
** The main loop follows. Search the list of cubes for the
** next cube to split, which is the lowest level cube. A
** special case is when a cube has only one color, so that it
** cannot be split.
*/
while (ncubes < maxcubes){
level = 255; splitpos = -1;
for (k=0;k<=ncubes-1;k++){
if (list[k].lower == list[k].upper)
; /* single color */
else if (list[k].level < level){
level = list[k].level;
splitpos = k;
}
}
if (splitpos == -1) /* no more cubes to split */
break;
/*
** Must split the cube "splitpos" in the list of cubes.
** Next find the longest dimension of this cube, and update
** the external variable "longdim", which is used by the
** sort routine so that it knows along which axis to sort.
*/
Cube = list[splitpos];
lr = Cube.rmax - Cube.rmin;
lg = Cube.gmax - Cube.gmin;
lb = Cube.bmax - Cube.bmin;
if (lr >= lg && lr >= lb) longdim = 0;
if (lg >= lr && lg >= lb) longdim = 1;
if (lb >= lr && lb >= lg) longdim = 2;
/*
** Sort along "longdim". This prepares for the next step,
** namely, finding the median. Use standard lib's "qsort".
*/
base = (void *)&HistPtr[Cube.lower];
num = (size_t)(Cube.upper - Cube.lower + 1);
width = (size_t)sizeof(HistPtr[0]);
qsort(base,num,width,compare);
/*
** Find median by scanning through the cube, and computing
** a running sum. When the running sum equals half the
** total for the cube, the median has been found.
*/
count = 0;
for (i=Cube.lower;i<=Cube.upper-1;i++){
if (count >= Cube.count/2) break;
color = HistPtr[i];
count = count + Hist[color];
}
median = i;
/*
** Now split "Cube" at the median. Then add the two new
** cubes to the list of cubes.
*/
CubeA = Cube; CubeA.upper = median-1;
CubeA.count = count;
CubeA.level = Cube.level + 1;
Shrink(&CubeA);
list[splitpos] = CubeA; /* add in old slot */
CubeB = Cube; CubeB.lower = median;
CubeB.count = Cube.count - count;
CubeB.level = Cube.level + 1;
Shrink(&CubeB);
list[ncubes++] = CubeB; /* add in new slot */
if ((ncubes % 10) == 0)
fprintf(stderr,"."); /* pacifier */
}
/*
** We have enough cubes, or we have split all we can. Now
** compute the color map, the inverse color map, and return
** the number of colors in the color map.
*/
InvMap(Hist, ColMap,ncubes);
return((WORD)ncubes);
}
void Shrink(cube_t * Cube)
{
/*
** Encloses "Cube" with a tight-fitting cube by updating the
** (rmin,gmin,bmin) and (rmax,gmax,bmax) members of "Cube".
*/
BYTE r,g,b;
WORD i,color;
Cube->rmin = 255; Cube->rmax = 0;
Cube->gmin = 255; Cube->gmax = 0;
Cube->bmin = 255; Cube->bmax = 0;
for (i=Cube->lower;i<=Cube->upper;i++){
color = HistPtr[i];
r = RED(color);
if (r > Cube->rmax) Cube->rmax = r;
if (r < Cube->rmin) Cube->rmin = r;
g = GREEN(color);
if (g > Cube->gmax) Cube->gmax = g;
if (g < Cube->gmin) Cube->gmin = g;
b = BLUE(color);
if (b > Cube->bmax) Cube->bmax = b;
if (b < Cube->bmin) Cube->bmin = b;
}
}
void InvMap(WORD * Hist, BYTE ColMap[][3],WORD ncubes)
{
/*
** For each cube in the list of cubes, computes the centroid
** (average value) of the colors enclosed by that cube, and
** then loads the centroids in the color map. Next loads the
** histogram with indices into the color map. A preprocessor
** directive "#define FAST_REMAP" controls whether the cube
** centroids become the output color for all the colors in a
** cube, or whether a "best remap" method is followed.
*/
BYTE r,g,b;
WORD i,k,color;
float rsum,gsum,bsum;
cube_t Cube;
for (k=0;k<=ncubes-1;k++){
Cube = list[k];
rsum = gsum = bsum = (float)0.0;
for (i=Cube.lower;i<=Cube.upper;i++){
color = HistPtr[i];
r = RED(color);
rsum += (float)r*(float)Hist[color];
g = GREEN(color);
gsum += (float)g*(float)Hist[color];
b = BLUE(color);
bsum += (float)b*(float)Hist[color];
}
/* Update the color map */
ColMap[k][0] = (BYTE)(rsum/(float)Cube.count);
ColMap[k][1] = (BYTE)(gsum/(float)Cube.count);
ColMap[k][2] = (BYTE)(bsum/(float)Cube.count);
}
#ifdef FAST_REMAP
/*
** Fast remap: for each color in each cube, load the corre-
** sponding slot in "Hist" with the centroid of the cube.
*/
for (k=0;k<=ncubes-1;k++){
Cube = list[k];
for (i=Cube.lower;i<=Cube.upper;i++){
color = HistPtr[i];
Hist[color] = k;
}
if ((k % 10) == 0) fprintf(stderr,"."); /* pacifier */
}
#else
/*
** Best remap: for each color in each cube, find the entry in
** "ColMap" that has the smallest Euclidian distance from the
** color. Record this information in "Hist".
*/
for (k=0;k<=ncubes-1;k++){
Cube = list[k];
for (i=Cube.lower;i<=Cube.upper;i++){
color = HistPtr[i];
r = RED(color); g = GREEN(color); b = BLUE(color);
/* Search for closest entry in "ColMap" */
dmin = (float)FLT_MAX;
for (j=0;j<=ncubes-1;j++){
dr = (float)ColMap[j][0] - (float)r;
dg = (float)ColMap[j][1] - (float)g;
db = (float)ColMap[j][2] - (float)b;
d = dr*dr + dg*dg + db*db;
if (d == (float)0.0){
index = j; break;
}
else if (d < dmin){
dmin = d; index = j;
}
}
Hist[color] = index;
}
if ((k % 10) == 0) fprintf(stderr,"."); /* pacifier */
}
#endif
return;
}
int compare(const void * a1, const void * a2)
{
/*
** Called by the sort routine in "MedianCut". Compares two
** colors based on the external variable "longdim".
*/
WORD color1,color2;
BYTE c1,c2;
color1 = (WORD)*(WORD *)a1;
color2 = (WORD)*(WORD *)a2;
switch (longdim){
case 0:
c1 = RED(color1), c2 = RED(color2);
break;
case 1:
c1 = GREEN(color1), c2 = GREEN(color2);
break;
case 2:
c1 = BLUE(color2), c2 = BLUE(color2);
break;
}
return ((int)(c1-c2));
}