This is pyrlib.c in view mode; [Download] [Up]
/*
* pyrlib.c - Library routines for pyramids
*
* Author: Rod Bogart
* Computer Science Dept.
* University of Utah
* Date: Thu Mar 12 1987
* Copyright (c) 1987 Rod Bogart
*
*/
#ifndef lint
static char rcs_ident[] = "$Id: pyrlib.c,v 3.0 90/08/03 15:22:11 spencer Exp $";
#endif
#include <stdio.h>
#include "rle.h"
#include "pyramid.h"
#ifdef USE_STDLIB_H
#include <stdlib.h>
#else
#ifdef VOID_STAR
extern void *malloc();
#else
extern char *malloc();
#endif
extern void free();
#endif /* USE_STDLIB_H */
void extrap_level(), extrap_int_level();
void alloc_pyramid(), gauss_level(), band_level(), expand_level();
void dump_level_line(), band_to_int_level(), add_int_level();
int
rle_to_pyramids(infile, gausspyr, bandpyr,
in_hdr, levellimit, mask_mult_table)
FILE * infile;
pyramid * gausspyr;
pyramid * bandpyr;
rle_hdr * in_hdr;
int levellimit;
float * mask_mult_table;
{
rle_pixel *rastptr, *rasterbase;
rle_pixel *firstbase, *outbase;
int xlen, ylen, nchan;
int xsize, ysize, levels, xlinewidth;
int i, chan, x, y, minlen;
rle_pixel *firstrastptr, *outrastptr;
rle_pixel *firstpxl, *outpxl, *signpxl;
rle_pixel **rows;
int rle_err;
*in_hdr = rle_dflt_hdr;
in_hdr->rle_file = infile;
if ( (rle_err = rle_get_setup( in_hdr )) < 0 )
{
return rle_err;
}
if ( in_hdr->alpha )
RLE_SET_BIT( *in_hdr, RLE_ALPHA );
in_hdr->xmax -= in_hdr->xmin;
in_hdr->xmin = 0;
xlen = in_hdr->xmax + 1;
ylen = in_hdr->ymax - in_hdr->ymin + 1;
xlinewidth = xlen + MASKSIZE - 1;
/* assumes that alpha is 1 or 0 */
nchan = in_hdr->ncolors + in_hdr->alpha;
if (nchan > 8)
{
fprintf(stderr,"pyrlib: %d is too many input channels, 8 max\n",nchan);
exit(-2);
}
/*************************************************************
* Figure out how many levels
* Note that the smallest level will be between 8 and 4 pixels
* in its minimum direction.
*************************************************************/
minlen = (xlen < ylen) ? xlen : ylen;
levels = levellimit;
if (!levellimit)
for(levels = 1; minlen > 8; minlen = minlen >> 1, levels++);
else if ((minlen >> (levels - 1)) <= 4)
{
fprintf(stderr,
"pyrlib: %d is too many levels to make a pyramid\n",levels);
exit(-2);
}
if (levels < 2)
{
fprintf(stderr,
"pyrlib: %d is too few levels to make a pyramid\n",levels);
exit(-2);
}
gausspyr->levels = levels;
gausspyr->nchan = nchan;
gausspyr->xlen = (int *) malloc( levels * sizeof( int ) );
gausspyr->ylen = (int *) malloc( levels * sizeof( int ) );
gausspyr->corners = (rle_pixel **) malloc( levels * sizeof( rle_pixel *) );
if ( (!gausspyr->corners) or (!gausspyr->xlen) or (!gausspyr->ylen) )
{
fprintf(stderr,
"pyrlib: could not allocate space for gaussian pyramid\n");
exit(-2);
}
gausspyr->xlen[0] = xlen;
gausspyr->ylen[0] = ylen;
/* fprintf(stderr,"Allocating gaussian pyramid structure\n"); */
alloc_pyramid(gausspyr);
if ( bandpyr ) /* zero means do gauss only */
{
bandpyr->levels = levels;
bandpyr->nchan = nchan + 1; /* plus for sign chan */
bandpyr->xlen = (int *) malloc( levels * sizeof( int ) );
bandpyr->ylen = (int *) malloc( levels * sizeof( int ) );
bandpyr->corners = (rle_pixel **) malloc(levels * sizeof(rle_pixel *));
if ( (!bandpyr->corners) or (!bandpyr->xlen) or (!bandpyr->ylen) )
{
fprintf(stderr,
"pyrlib: could not allocate space for band pyramid\n");
exit(-2);
}
bandpyr->xlen[0] = xlen;
bandpyr->ylen[0] = ylen;
/* fprintf(stderr,"Allocating band pyramid structure\n"); */
alloc_pyramid(bandpyr);
}
rows = (rle_pixel **) malloc( nchan * sizeof( rle_pixel * ) );
if (! rows)
{
fprintf(stderr, "pyramid: malloc failed\n");
exit(-2);
}
rasterbase = gausspyr->corners[0];
rastptr = &(rasterbase[MASKBELOW * xlinewidth * nchan + MASKBELOW]);
/****************************************************
* Read in all of the pixels
****************************************************/
/* fprintf(stderr,"Reading input image\n"); */
for (i = in_hdr->ymin; i <= in_hdr->ymax; i++)
{
for (chan=0; chan < nchan; chan++)
{
rows[chan] = rastptr;
/* Increment pointer by xlinewidth */
rastptr = &(rastptr[xlinewidth]);
}
/* assumes 1 or 0 for alpha */
rle_getrow( in_hdr, &rows[in_hdr->alpha] );
}
/****************************************************
* Build each gauss level
****************************************************/
for(i=1; i < levels; i++)
{
/* fprintf(stderr,"Building gauss level %d\n", i); */
/* extrapolate the level we have */
extrap_level(i-1, gausspyr);
/* make gauss image in this level */
gauss_level(i, gausspyr, mask_mult_table);
}
/* extrapolate the highest level */
extrap_level(levels-1, gausspyr);
/****************************************************
* Build each band level from the gauss levels
****************************************************/
if (bandpyr)
{
for(i=levels - 2; i >= 0; i--)
{
/* fprintf(stderr,"Building band level %d\n", i); */
/* make band image from gauss images */
band_level(i, bandpyr, gausspyr, mask_mult_table);
}
/****************************************************
* The highest band is just the highest gauss (remnant)
****************************************************/
firstbase = gausspyr->corners[levels - 1];
outbase = bandpyr->corners[levels - 1];
xsize = xlen >> (levels - 1);
ysize = ylen >> (levels - 1);
xlinewidth = xsize + MASKSIZE - 1;
for (y = 0; y < ysize; y++)
{
firstrastptr = &(firstbase[MASKBELOW + (y+MASKBELOW) *
xlinewidth * nchan]);
outrastptr = &(outbase[MASKBELOW + (y+MASKBELOW) *
xlinewidth * (nchan + 1)]);
signpxl = &(outrastptr[nchan * xlinewidth]);
for(x=0; x < xsize; x++)
{
/* zero out sign bits for each scanline */
*signpxl = 0;
signpxl++;
}
for (chan = 0; chan < nchan; chan++)
{
firstpxl = &(firstrastptr[chan * xlinewidth]);
outpxl = &(outrastptr[chan * xlinewidth]);
for(x=0; x < xsize; x++)
{
*outpxl = *firstpxl;
firstpxl++;
outpxl++;
}
}
}
}
return RLE_SUCCESS;
}
/* alloc_pyramid -
* Allocates a pyramid structure with "channels" channels per level
* and level zero of size "xlen" by "ylen". The pointers to each
* corner are returned in the array "corners", so it is assumed that
* space has been set aside for the corners to be returned to.
*/
void
alloc_pyramid(pyr)
pyramid * pyr;
{
int i, xsize = pyr->xlen[0], ysize = pyr->ylen[0];
for(i=0; i < pyr->levels; i++)
{
pyr->corners[i] = (rle_pixel *) malloc( (xsize + MASKSIZE - 1) *
(ysize + MASKSIZE - 1) * pyr->nchan );
if ( ! pyr->corners[i] )
{
fprintf(stderr, "pyramid: could not allocate pyramid\n");
exit(-2);
}
pyr->xlen[i] = xsize;
pyr->ylen[i] = ysize;
xsize = xsize >> 1;
ysize = ysize >> 1;
}
}
/* extrap_level -
* "corners" is assumed to be an array of pointers to allocated
* levels of a gaussian pyramid (band pyramids have sign bits). The
* level indicated by "level" (0 to levels - 1) is extrapolated. This
* currently means that boundary pixels are reflected thru the edge
* pixel to find the outer boundary pixels. In actuality, there
* should be some nifty mask thingy that flags known and unknown
* pixels. Then one could extrapolate any unknown pixel that has
* known neighbors.
*/
void
extrap_level(level,pyr)
int level;
pyramid * pyr;
{
int chan, x, y, i, newvalue, xsize, ysize;
int xlinewidth, nchan;
rle_pixel *rastptr, *rasterbase;
rle_pixel *bot_edge_pxl, *bot_data_pxl, *bot_extrap_pxl;
rle_pixel *top_edge_pxl, *top_data_pxl, *top_extrap_pxl;
rle_pixel *left_edge_pxl, *left_data_pxl, *left_extrap_pxl;
rle_pixel *right_edge_pxl, *right_data_pxl, *right_extrap_pxl;
rasterbase = pyr->corners[level];
xsize = pyr->xlen[level];
ysize = pyr->ylen[level];
nchan = pyr->nchan;
xlinewidth = xsize + MASKSIZE - 1;
rastptr = &(rasterbase[MASKBELOW + MASKBELOW * xlinewidth * nchan]);
for (chan = 0; chan < nchan; chan++)
{
bot_edge_pxl = &(rastptr[chan * xlinewidth]);
top_edge_pxl = &(rastptr[(ysize-1) * nchan * xlinewidth
+ chan * xlinewidth]);
for(x=0; x < xsize; x++)
{
for (i=1; i <= MASKBELOW; i++)
{
bot_data_pxl = bot_edge_pxl + i * xlinewidth * nchan;
bot_extrap_pxl = bot_edge_pxl - i * xlinewidth * nchan;
newvalue = 2 * (*bot_edge_pxl) - (*bot_data_pxl);
*bot_extrap_pxl = (newvalue < 0) ? 0 :
((newvalue > 255) ? 255 : newvalue);
}
for (i=1; i <= MASKABOVE; i++)
{
top_data_pxl = top_edge_pxl - i * xlinewidth * nchan;
top_extrap_pxl = top_edge_pxl + i * xlinewidth * nchan;
newvalue = 2 * (*top_edge_pxl) - (*top_data_pxl);
*top_extrap_pxl = (newvalue < 0) ? 0 :
((newvalue > 255) ? 255 : newvalue);
}
bot_edge_pxl++;
top_edge_pxl++;
}
}
left_edge_pxl = &(rastptr[(-MASKBELOW) * nchan * xlinewidth]);
right_edge_pxl = &(rastptr[(-MASKBELOW) * nchan * xlinewidth
+ xlinewidth - MASKSIZE]);
for (chan = 0; chan < nchan; chan++)
{
for(y=0; y < ysize + MASKSIZE - 1; y++)
{
for (i=1; i <= MASKBELOW; i++)
{
left_data_pxl = left_edge_pxl + i;
left_extrap_pxl = left_edge_pxl - i;
newvalue = 2 * (*left_edge_pxl) - (*left_data_pxl);
*left_extrap_pxl = (newvalue < 0) ? 0 :
((newvalue > 255) ? 255 : newvalue);
}
for (i=1; i <= MASKABOVE; i++)
{
right_data_pxl = right_edge_pxl - i;
right_extrap_pxl = right_edge_pxl + i;
newvalue = 2 * (*right_edge_pxl) - (*right_data_pxl);
*right_extrap_pxl = (newvalue < 0) ? 0 :
((newvalue > 255) ? 255 : newvalue);
}
left_edge_pxl += xlinewidth;
right_edge_pxl += xlinewidth;
}
}
}
/* gauss_mask -
* Create the gaussian 5x5 mask. It gets allocated here, a pointer
* is returned.
*/
float *
gauss_mask(siz)
int siz;
{
int i,j;
float w[5];
float *mptr;
w[0] = 0.05;
w[1] = 0.25;
w[2] = 0.4;
w[3] = 0.25;
w[4] = 0.05;
mptr = (float *) malloc( sizeof(float) * siz * siz );
if (siz != 5)
{
fprintf(stderr,"mask size not 5\n");
exit(-1);
}
for(i=0; i< 5; i++)
{
for(j=0;j<5;j++)
{
mptr[j+i*siz] = w[j]*w[i];
}
}
return mptr;
}
/* gauss_level -
* Build gaussian level from level - 1. It is assumed that the
* previous level is full of data, and that it has been extrapolated
* at its boundaries. Then every other pixel is used for the center
* of the gaussian mask convolution. The result is half the size of
* level - 1.
*/
void
gauss_level(level,pyr,mask_mult_table)
int level;
pyramid * pyr;
float *mask_mult_table;
{
int xsize, ysize, x, y, xlinewidth, lilxlinewidth;
int chan,i,j,nchan;
rle_pixel *inbase, *outbase, *inrastptr, *outrastptr;
rle_pixel *centerpxl, *outpxl, *calcpxl;
float total;
float *maskmult;
/* this assumes that the level-1 image is extrapolated */
inbase = pyr->corners[level-1];
outbase = pyr->corners[level];
xsize = pyr->xlen[level - 1];
ysize = pyr->ylen[level - 1];
nchan = pyr->nchan;
xlinewidth = xsize + MASKSIZE - 1;
lilxlinewidth = (xsize >> 1) + MASKSIZE - 1;
for (y = 0; y < ysize; y+=2)
{
inrastptr = &(inbase[MASKBELOW + (y+MASKBELOW) *
xlinewidth * nchan]);
outrastptr = &(outbase[MASKBELOW + ((y/2)+MASKBELOW) *
lilxlinewidth * nchan]);
for (chan = 0; chan < nchan; chan++)
{
centerpxl = &(inrastptr[chan * xlinewidth]);
outpxl = &(outrastptr[chan * lilxlinewidth]);
for(x=0; x < xsize; x+=2)
{
total = 0.0;
for (i=0; i < 5; i++)
{
calcpxl = centerpxl + (i - MASKBELOW)
* xlinewidth * nchan - MASKBELOW;
maskmult = &(mask_mult_table[(i*5) << 8]);
for (j=0; j < 5; j++)
{
/* look up multiply in maskmult table */
total += maskmult[ calcpxl[j] ];
maskmult += 256;
}
}
/* bound result to 0 to 255 */
*outpxl = (total < 0.0) ? 0 :
((total > 255.0) ? 255 : (rle_pixel) total );
centerpxl++; centerpxl++;
outpxl++;
}
}
}
}
/* band_level -
* Given a complete gaussian pyramid, this creates the band pass
* levels, which are ((Gn) - (Gn+1)). The "level + 1" gaussian level is
* expanded, put in temp location, and then is subtracted pixel by
* pixel from the "level" gaussian level. Since this can cause
* negative values, the band level has an additional channel which is
* the sign bits (1 = negative, 0 = positive) for up to eight
* channels. There is a limit of eight channels for an input image.
* Hence, a band pyramid has nchan + 1 channels in it.
*/
void
band_level(level, bandpyr, gausspyr, mask_mult_table)
int level;
pyramid * bandpyr;
pyramid * gausspyr;
float *mask_mult_table;
{
int xsize, ysize, x, y, xlinewidth;
int chan,subtractval,nchan;
rle_pixel *firstbase, *secondbase, *outbase;
rle_pixel *firstrastptr, *secondrastptr, *outrastptr;
rle_pixel *firstpxl, *secondpxl, *outpxl, *signpxl;
/* this assumes that the gauss images are extrapolated */
/* the result pyramid will NOT be extrapolated */
/* second level will be expanded and subtracted from first
* to produce out */
firstbase = gausspyr->corners[level];
expand_level(level+1, gausspyr, &secondbase, mask_mult_table);
outbase = bandpyr->corners[level];
xsize = gausspyr->xlen[level];
ysize = gausspyr->ylen[level];
nchan = gausspyr->nchan;
xlinewidth = xsize + MASKSIZE - 1;
for (y = 0; y < ysize; y++)
{
firstrastptr = &(firstbase[MASKBELOW + (y+MASKBELOW) *
xlinewidth * nchan]);
secondrastptr = &(secondbase[MASKBELOW + (y+MASKBELOW) *
xlinewidth * nchan]);
outrastptr = &(outbase[MASKBELOW + (y+MASKBELOW) *
xlinewidth * (nchan + 1)]);
signpxl = &(outrastptr[nchan * xlinewidth]);
for(x=0; x < xsize; x++)
{
/* zero out sign bits for this scanline */
*signpxl = 0;
signpxl++;
}
for (chan = 0; chan < nchan; chan++)
{
firstpxl = &(firstrastptr[chan * xlinewidth]);
secondpxl = &(secondrastptr[chan * xlinewidth]);
outpxl = &(outrastptr[chan * xlinewidth]);
signpxl = &(outrastptr[nchan * xlinewidth]);
for(x=0; x < xsize; x++)
{
subtractval = (int) (*firstpxl) - (int) (*secondpxl);
if (subtractval >= 0)
*outpxl = (rle_pixel) subtractval;
else
{
*outpxl = (rle_pixel) -subtractval;
*signpxl = (1 << chan) | (*signpxl);
}
firstpxl++;
secondpxl++;
outpxl++;
signpxl++;
}
}
}
free( *secondbase );
}
/* expand_level -
* This is used for expanding a level of a gaussian pyramid. It
* allocates the storage for the result, and the pointer is returned
* in corner. The mask centering process is way weird. Actually, it
* might even, be slightly wrong, but don't tell anyone, K?
*/
void
expand_level(level, gausspyr, corner, mask_mult_table)
int level;
pyramid * gausspyr;
rle_pixel ** corner;
float *mask_mult_table;
{
int xsize, ysize, x, y, xlinewidth, lilxlinewidth;
int chan,i,j,nchan;
rle_pixel *inbase, *outbase;
rle_pixel *inrastptr, *outrastptr;
rle_pixel *centerpxl, *outpxl, *calcpxl, expandedval;
float total;
float *maskmult;
/* this assumes that the gauss images are extrapolated */
xsize = gausspyr->xlen[level - 1];
ysize = gausspyr->ylen[level - 1];
nchan = gausspyr->nchan;
/* level is the gauss level to be expanded */
outbase = (rle_pixel *) malloc( (xsize + MASKSIZE - 1) *
(ysize + MASKSIZE - 1) * nchan );
if (!outbase)
{
fprintf(stderr, "pyrhalf: could not allocate expand level\n");
exit(-2);
}
inbase = gausspyr->corners[level];
*corner = outbase;
xlinewidth = xsize + MASKSIZE - 1;
lilxlinewidth = (xsize >> 1) + MASKSIZE - 1;
for (y = 0; y < ysize; y++)
{
inrastptr = &(inbase[MASKBELOW + ((y/2)+MASKBELOW) *
lilxlinewidth * nchan]);
outrastptr = &(outbase[MASKBELOW + (y+MASKBELOW) *
xlinewidth * nchan]);
for (chan = 0; chan < nchan; chan++)
{
centerpxl = &(inrastptr[chan * lilxlinewidth]);
outpxl = &(outrastptr[chan * xlinewidth]);
for(x=0; x < xsize; x++)
{
total = 0.0;
for (i=0; i < 5; i++)
{
if ((y + i - MASKBELOW) % 2 == 0)
{
if ((y % 2) == 0)
calcpxl = centerpxl + ((i - MASKBELOW) / 2)
* lilxlinewidth * nchan - ( MASKBELOW / 2);
else
calcpxl = centerpxl + ((i + 1 - MASKBELOW) / 2)
* lilxlinewidth * nchan - ( MASKBELOW / 2);
maskmult = &(mask_mult_table[(i*5) << 8]);
for (j=0; j < 5; j++)
{
if ((x + j - MASKBELOW) % 2 == 0)
{
/* look up multiply in maskmult table */
if ((x % 2) == 0)
total += maskmult[ calcpxl[j / 2] ];
else
total += maskmult[ calcpxl[(j + 1) / 2] ];
}
maskmult += 256;
}
}
}
/* bound result to 0 to 255 */
total *= 4.0;
expandedval = (total < 0.0) ? 0 :
((total > 255.0) ? 255 : (rle_pixel) total );
*outpxl = (rle_pixel) expandedval;
if ((x % 2) == 1) centerpxl++;
outpxl++;
}
}
}
}
/* dump_pyramid -
* Given a file pointer (an open file, hmmm...), this dumps a
* pyramid out on multiple channels. If the input pyramid is only one
* level, an attempt is made to do the right thing with alpha, and
* just dump a 'normal' RLE file. If it has more channels ( > 1 ),
* then the levels of the pyramid end up in higher channels, in groups
* of "channels" channels.
*/
void
dump_pyramid(outfile,levels,corners,xlen,ylen,channels,in_hdr)
FILE * outfile;
int levels;
rle_pixel **corners;
int xlen,ylen,channels;
rle_hdr in_hdr;
{
rle_hdr out_hdr;
int chan, level, x, y, activelevels, adjust_alpha = 0;
rle_pixel **scanout;
rle_pixel **outrows;
/* fprintf(stderr,"Dumping %d level pyramid with %d channels per level\n",
levels, channels);
*/
if ((levels == 1) && (in_hdr.alpha))
adjust_alpha = 1; /* for dumping a normal image */
out_hdr = in_hdr;
out_hdr.alpha = adjust_alpha;
out_hdr.ncolors = channels * levels - adjust_alpha;
out_hdr.rle_file = outfile;
for ( chan=0; chan < channels*levels - adjust_alpha; chan++)
{
RLE_SET_BIT( out_hdr, chan );
}
if (rle_row_alloc( &out_hdr, &scanout ))
fprintf(stderr, "pyrhalf: Ran out of malloc space\n");
outrows = (rle_pixel **) malloc( channels * levels * sizeof( rle_pixel * ) );
if (! outrows)
{
fprintf(stderr, "pyrhalf: malloc failed\n");
exit( -2 );
}
if (adjust_alpha) scanout--;
for ( chan=0; chan < channels*levels; chan++)
{
outrows[chan] = scanout[chan];
/* fill the scanline with black */
for (x = 0; x < xlen; x++)
scanout[chan][x] = 0;
}
rle_put_setup( &out_hdr );
activelevels = levels;
for (y = 0; y < ylen; y++)
{
/* for each active level, dump the data into the right part
* of the scanline.
*/
for (level = 0; level < activelevels; level++)
{
dump_level_line(scanout,y,level,corners,xlen,channels);
}
rle_putrow( &(outrows[adjust_alpha]), xlen, &out_hdr );
/* check if top level becomes inactive, dump black
* to the scanline for next time.
*/
if (y == (ylen >> (activelevels -1)))
{
for (chan = 0; chan < channels; chan++)
for (x = 0; x < (xlen >> (activelevels -1)); x++)
scanout[chan + channels * (activelevels-1)][x] = 0;
activelevels--;
}
}
rle_puteof( &out_hdr );
}
void
dump_level_line(scanout,y,level,corners,xlen,nchan)
rle_pixel ** scanout;
int y, level;
rle_pixel ** corners;
int xlen,nchan;
{
int xsize, x, chan, xlinewidth;
rle_pixel * outpxl, * outrastptr, *outbase, *scanoutptr;
xsize = xlen >> level;
xlinewidth = xsize + MASKSIZE - 1;
outbase = corners[level];
outrastptr = &(outbase[MASKBELOW + (y + MASKBELOW) * xlinewidth * nchan]);
/* for each channel in the level
* copy the pixel into scanout
*/
for ( chan=0; chan < nchan; chan++)
{
outpxl = &(outrastptr[chan * xlinewidth]);
scanoutptr = scanout[chan + nchan * level];
for (x=0; x < xsize; x++)
{
*scanoutptr = *outpxl;
scanoutptr++;
outpxl++;
}
}
}
void
rebuild_image(imgcorner, bandpyr, mask_mult_table)
rle_pixel ** imgcorner;
pyramid * bandpyr;
float *mask_mult_table;
{
int i, xsize, ysize, level, nchan, xlen, ylen;
int x, y, chan, xlinewidth,levels;
int ** intcorners;
int *intbase, *intrastptr, *intpxl;
rle_pixel *outbase, *outrastptr, *outpxl;
/* band pyramid is probably not extrapolated
* allocate int_pyramid without sign channel
* copy band pyramid to int pyramid (do correct thing with sign)
* for each level (highest to lowest+1)
* extrapolate level
* expand level and add result directly to next level
* free temp expanded level
* allocate rle_pixel image of level 1 size
* copy from int level 1 to image
* free int pyramid
*/
/*************************************************************
* Allocate int pyramid
*************************************************************/
levels = bandpyr->levels;
nchan = bandpyr->nchan - 1;
xlen = bandpyr->xlen[0]; ylen = bandpyr->ylen[0];
intcorners = (int **) malloc( levels * sizeof( int *) );
if ( ! intcorners )
{
fprintf(stderr, "pyrhalf: could not allocate pyramid\n");
exit(-2);
}
xsize = xlen; ysize = ylen;
for(i=0; i < levels; i++)
{
intcorners[i] = (int *) malloc( (xsize + MASKSIZE - 1)
* (ysize + MASKSIZE - 1) * nchan * sizeof(int) );
if ( ! intcorners[i] )
{
fprintf(stderr, "pyrhalf: could not allocate pyramid\n");
exit(-2);
}
xsize = xsize >> 1;
ysize = ysize >> 1;
}
/*************************************************************
* Copy band pyramid to int (with sign correct)
*************************************************************/
for(level = levels - 1; level >= 0; level--)
{
band_to_int_level(level,bandpyr->corners,intcorners,
xlen,ylen,nchan);
}
/*************************************************************
* For each level in the int pyramid (highest to lowest+1)
*************************************************************/
for(level = levels - 1; level >= 1; level--)
{
/* extrapolate level */
extrap_int_level(level,intcorners,xlen,ylen,nchan);
/* expand level and add result directly to next level */
add_int_level(level,intcorners,xlen,ylen,
nchan,mask_mult_table);
}
/*************************************************************
* Allocate single image buffer
*************************************************************/
*imgcorner = (rle_pixel *) malloc( (xlen + MASKSIZE - 1) *
(ylen + MASKSIZE - 1) * nchan );
if ( ! (*imgcorner) )
{
fprintf(stderr, "pyrhalf: could not allocate pyramid\n");
exit(-2);
}
/*************************************************************
* Copy lowest level of int pyramid to image buffer (clamp 0 - 255)
*************************************************************/
intbase = intcorners[0];
outbase = *imgcorner;
xsize = xlen;
ysize = ylen;
xlinewidth = xsize + MASKSIZE - 1;
for (y = 0; y < ysize; y++)
{
intrastptr = &(intbase[MASKBELOW + (y+MASKBELOW) *
xlinewidth * nchan]);
outrastptr = &(outbase[MASKBELOW + (y+MASKBELOW) *
xlinewidth * nchan]);
for (chan = 0; chan < nchan; chan++)
{
intpxl = &(intrastptr[chan * xlinewidth]);
outpxl = &(outrastptr[chan * xlinewidth]);
for(x=0; x < xsize; x++)
{
*outpxl = (rle_pixel) ((*intpxl) > 255) ? 255 :
(((*intpxl) < 0) ? 0 : (*intpxl));
intpxl++;
outpxl++;
}
}
}
/*************************************************************
* Free levels of int pyramid
*************************************************************/
for(i=0; i < levels; i++)
free(intcorners[i]);
}
void
band_to_int_level(level,bandcorners,intcorners, xlen,ylen,nchan)
int level;
rle_pixel ** bandcorners;
int ** intcorners;
int xlen,ylen,nchan;
{
int xsize, ysize, x, y, chan, xlinewidth, sign;
rle_pixel *firstbase, *firstrastptr, *firstpxl, *signpxl;
int *outbase, *outrastptr, *outpxl;
/* copy from band level to int level keeping track of the sign */
firstbase = bandcorners[level];
outbase = intcorners[level];
xsize = xlen >> level;
ysize = ylen >> level;
xlinewidth = xsize + MASKSIZE - 1;
for (y = 0; y < ysize; y++)
{
firstrastptr = &(firstbase[MASKBELOW + (y+MASKBELOW) *
xlinewidth * (nchan + 1)]);
outrastptr = &(outbase[MASKBELOW + (y+MASKBELOW) *
xlinewidth * nchan]);
for (chan = 0; chan < nchan; chan++)
{
firstpxl = &(firstrastptr[chan * xlinewidth]);
outpxl = &(outrastptr[chan * xlinewidth]);
signpxl = &(firstrastptr[nchan * xlinewidth]);
for(x=0; x < xsize; x++)
{
sign = ((1 << chan) & (*signpxl)) ? -1 : 1;
*outpxl = (int) (*firstpxl) * sign;
firstpxl++;
outpxl++;
signpxl++;
}
}
}
}
void
extrap_int_level(level,corners,xlen,ylen,nchan)
int level;
int **corners;
int xlen,ylen,nchan;
{
int chan, x, y, i, newvalue, xsize, ysize;
int xlinewidth;
int *rastptr, *rasterbase;
int *bot_edge_pxl, *bot_data_pxl, *bot_extrap_pxl;
int *top_edge_pxl, *top_data_pxl, *top_extrap_pxl;
int *left_edge_pxl, *left_data_pxl, *left_extrap_pxl;
int *right_edge_pxl, *right_data_pxl, *right_extrap_pxl;
rasterbase = corners[level];
xsize = xlen >> level;
ysize = ylen >> level;
xlinewidth = xsize + MASKSIZE - 1;
rastptr = &(rasterbase[MASKBELOW + MASKBELOW * xlinewidth * nchan]);
for (chan = 0; chan < nchan; chan++)
{
bot_edge_pxl = &(rastptr[chan * xlinewidth]);
top_edge_pxl = &(rastptr[(ysize-1) * nchan * xlinewidth
+ chan * xlinewidth]);
for(x=0; x < xsize; x++)
{
for (i=1; i <= MASKBELOW; i++)
{
bot_data_pxl = bot_edge_pxl + i * xlinewidth * nchan;
bot_extrap_pxl = bot_edge_pxl - i * xlinewidth * nchan;
newvalue = 2 * (*bot_edge_pxl) - (*bot_data_pxl);
*bot_extrap_pxl = newvalue;
}
for (i=1; i <= MASKABOVE; i++)
{
top_data_pxl = top_edge_pxl - i * xlinewidth * nchan;
top_extrap_pxl = top_edge_pxl + i * xlinewidth * nchan;
newvalue = 2 * (*top_edge_pxl) - (*top_data_pxl);
*top_extrap_pxl = newvalue;
}
bot_edge_pxl++;
top_edge_pxl++;
}
}
left_edge_pxl = &(rastptr[(-MASKBELOW) * nchan * xlinewidth]);
right_edge_pxl = &(rastptr[(-MASKBELOW) * nchan * xlinewidth
+ xlinewidth - MASKSIZE]);
for (chan = 0; chan < nchan; chan++)
{
for(y=0; y < ysize + MASKSIZE - 1; y++)
{
for (i=1; i <= MASKBELOW; i++)
{
left_data_pxl = left_edge_pxl + i;
left_extrap_pxl = left_edge_pxl - i;
newvalue = 2 * (*left_edge_pxl) - (*left_data_pxl);
*left_extrap_pxl = newvalue;
}
for (i=1; i <= MASKABOVE; i++)
{
right_data_pxl = right_edge_pxl - i;
right_extrap_pxl = right_edge_pxl + i;
newvalue = 2 * (*right_edge_pxl) - (*right_data_pxl);
*right_extrap_pxl = newvalue;
}
left_edge_pxl += xlinewidth;
right_edge_pxl += xlinewidth;
}
}
}
void
add_int_level(level,intcorners,xlen,ylen,nchan,mask_mult_table)
int level;
int **intcorners;
int xlen,ylen,nchan;
float *mask_mult_table;
{
int xsize, ysize, x, y, xlinewidth, lilxlinewidth;
int chan,i,j;
int *inbase, *outbase;
int *inrastptr, *outrastptr;
int *centerpxl, *outpxl, *calcpxl, expandedval;
float total;
float *maskmult;
/* this assumes that the gauss images are extrapolated */
xsize = xlen >> (level - 1);
ysize = ylen >> (level - 1);
/* level is the int level to be expanded */
inbase = intcorners[level];
outbase = intcorners[level-1];
xlinewidth = xsize + MASKSIZE - 1;
lilxlinewidth = (xsize >> 1) + MASKSIZE - 1;
for (y = 0; y < ysize; y++)
{
inrastptr = &(inbase[MASKBELOW + ((y/2)+MASKBELOW) *
lilxlinewidth * nchan]);
outrastptr = &(outbase[MASKBELOW + (y+MASKBELOW) *
xlinewidth * nchan]);
for (chan = 0; chan < nchan; chan++)
{
centerpxl = &(inrastptr[chan * lilxlinewidth]);
outpxl = &(outrastptr[chan * xlinewidth]);
for(x=0; x < xsize; x++)
{
total = 0.0;
for (i=0; i < 5; i++)
{
if ((y + i - MASKBELOW) % 2 == 0)
{
if ((y % 2) == 0)
calcpxl = centerpxl + ((i - MASKBELOW) / 2)
* lilxlinewidth * nchan - ( MASKBELOW / 2);
else
calcpxl = centerpxl + ((i + 1 - MASKBELOW) / 2)
* lilxlinewidth * nchan - ( MASKBELOW / 2);
maskmult = &(mask_mult_table[(i*5) << 8]);
for (j=0; j < 5; j++)
{
if ((x + j - MASKBELOW) % 2 == 0)
{
/* since we gots ints, we can't use lookup */
/* so just multiply by float value (yucko) */
if ((x % 2) == 0)
{
total += maskmult[ 1 ] * calcpxl[j / 2];
}
else
{
total += maskmult[ 1 ] *
calcpxl[(j + 1) / 2];
}
}
maskmult += 256;
}
}
}
total *= 4.0;
expandedval = (int) total;
*outpxl += expandedval;
if ((x % 2) == 1) centerpxl++;
outpxl++;
}
}
}
}
These are the contents of the former NiCE NeXT User Group NeXTSTEP/OpenStep software archive, currently hosted by Netfuture.ch.