Path: utzoo!utgpu!news-server.csri.toronto.edu!rutgers!uwm.edu!ogicse!plains!bakke
From: bakke@plains.NoDak.edu (Jeffrey P. Bakke)
Newsgroups: comp.graphics
Subject: Image Scaling (Replies)
Message-ID: <6364@plains.NoDak.edu>
Date: 19 Oct 90 16:18:32 GMT
Organization: North Dakota State University, Fargo
Lines: 319

I had a request to post all the replies and information that I received
regarding Image Scaling (shrinking an image).  Here it is.  Thanks for
all of the replies and help.


From ted@aps1.spa.umn.edu Wed Oct 10 13:33:54 1990

Take a look at the pbmplus package, available via ftp (see Freq. asked question
posting for how to get  it).  It includes code for scaling color images
(ppmscale).  The code should be easy to adapt to your purpose.


From sloan@cs.washington.edu Wed Oct 10 22:19:56 1990

Subsampling an image (what you asked about) is really quite easy - once you
build the right tools.  Supersampling (blowing up the image) is a bit harder
- but you didn't ask about that.

The first try is to simply subsample.  Given a high resolution image with
dimensions NxN, that you want to resample at MxM, with M < N, all you have
to do is:

for(x=0;x<M;x++)
 for(y=0;y<M;y++)
  Output[x][y] = Input[x*N/M][y*N/M];

Streaming the input image in ad the output image out is left as an exercise.

The problem with this is, of course, aliasing.  High frequency energy in the
original, high resolution image will "alias" as low frequency energy in the
low resolution image.  The classic solution is to low pass filter the high
resolution image first, BEFORE you take the samples.  

One way to do this is to take the Fourier transform of the high resolution
image,  mask out the high frequency terms, and take the inverse Fourier
transform.  Now, the sampling program above will work just fine.

The problem now is that the Fourier Transform, and it's inverse, are
computationally intensive.  You may not have the horsepower to take the
transform and the inverse fast enough.

No problem!  It turns out that you can APPROXIMATE the low pass filter by
"convolving the image with a kernel" (i.e., take weighted averages).
Example kernels are:

    0.125 0.125 0.125
    0.125 0.125 0.125        BOX filter ("average over 3x3 neighborhood")
    0.125 0.125 0.125

or, perhaps something like:
    0.050 0.100 0.050
    0.100 0.400 0.100        TENT filter ("weighted average")         
    0.050 0.100 0.050

The bigger the neighborhood, the better you can approximate the low-pass
filter that you want.  Discussion of this topic, and design of these
kernels, can be found in any text on "image processing".
But, the second kernel above will do just fine for your application, as long
as you only want to go from NxN -> N/2xN/2.  If you want to make the image
much smaller than half the linear size of the original, you need to average
over larger neighborhoods (or, filter several times).

Now, my preference is to build a program which performs the 3x3 convolution,
and a separate subsampling program.  so, to "properly" downsize an image, you
do:
  cat NxNImage | Convolve | SubSample >MxMImage

But, you can also write one program which convolves and samples
simultaneously.  It looks something like:
for(x=0;x<M;x++)
 for(y=0;y<M;y++)
  Output[x][y] = AverageOverNeighborhood(Input,x*N/M,y*N/M]);

Again, you can spend a lot of programmer time making this efficient.

Finally, two advanced points:

*for isotropic kernels, you can get the exact same effect by filtering
first with a 3x1 kernel (say, 0.25 0.50 0.25) and then a 1x3 kernel,
     0.25
say, 0.50  (these numbers are just illustrations).  This can make your
     0.50
program marginally more efficient for 3x3 kernels (only 6 multiplications
instead of 9) and MUCH MORE efficient for larger kernels (only 2N
multiplications instead of N*N)

* if you write one monolithic "filter and sample" program, you can custome
design your filter "on the fly" to do just the right amount of averaging

I suggest that you start by writing one program to do a 3x3 convolution with
an arbitrary 3x3 kernel, and another to subsampl by exactly half.
Experiment with different kernels (start with the two I gave above) until
you get a feel for the effects involved.  Compare the filtered versions
against each other, and against the simple subsampling with no filtering.

Good luck.  If you still have questions - just ask.

-Ken Sloan


From cristy%ESSUN3%ESVAX%dupont.com@RELAY.CS.NET Wed Oct 10 18:00:40 1990

Get contrib/ImageMagick.tar.Z on expo.lcs.mit.edu.  It has two algorithms
for image scaling; an antialias digital filter and straight pixel
replication.

From wolberg@cs.columbia.edu Wed Oct 10 15:59:15 1990

Hi,
There's a new book available on this subject called
"Digital Image Warping." It goes over all issues of
geometric transformations on digital images, i.e.,
scaling, rotation, warping, etc. It includes C code.

Title: Digital Image Warping
Author: George Wolberg
Publisher: IEEE Computer Society Press (1990)
IEEE order #: 1944
Price: $45 for IEEE members, $60 for nonmembers
1-800-CS-BOOKS


From hawley%adobe.uucp%decwrl.uucp@ogicse.UUCP Thu Oct 11 16:37:26 1990

The task that you want to solve can be generalized to the task of mapping
an rectangular bitmap to an arbitrary quadrilateral.  This usually gets
solved by ray tracers for pattern mapping.  I suggest you look through
siggraph papers for this more general topic.

The typical technique for bitmap remapping is to loop over the destination
pixels (NOT the source pixels) and apply a reverse transformation.

Scaling turns out to be fairly straightforward as it is just 2 linear
transformations, but the tricky part (which I can't help you with) is to
do a degree of antialiasing so that mappings that fall on the cracks will
look good.


From dal%syntel.uucp%tcnet.uucp@jhereg.osa.com Wed Oct 17 10:14:21 1990

/*
 * rescale -- rescale an image to the target size (brute force method)
 *
 * Copyright (C) 1990 by Dale Schumacher.
 *
 * Permission to use, copy, modify, and distribute this software and its
 * documentation for any purpose and without fee is hereby granted, provided
 * that the above copyright notice appear in all copies and that both that
 * copyright notice and this permission notice appear in supporting
 * documentation.  This software is provided "as is" without express or
 * implied warranty.
 *
 */

static char	_Program[] = "rescale";
static char	_Version[] = "1.1";
static char	_Copyright[] = "Copyright 1990 by Dale Schumacher";

#include <stdio.h>
#include "pxm.h"

#define	DEBUG(x)	/* x */
#define	TRACE(x)	/* x */

long	sx = 0;		/* source x size */
long	sy = 0;		/* source y size */
long	tx = 0;		/* target x size */
long	ty = 0;		/* target y size */
int	verbose = 0;	/* verbose flag */

banner()
{
	printf("%s v%s -- %s\n", _Program, _Version, _Copyright);
}

usage()
{
	fprintf(stderr,
		"usage: %s [-v] [-x newsize] [-y newsize] inpxm outpxm\n",
		_Program);
	exit(1);
}

progress(a, b)
long a, b;
{
	printf("%3ld%%\r", ((a * 100) / b));
}

main(argc, argv)
int argc;
char *argv[];
{
	register PX_DEF *ipx, *opx;
	FILE *f;
	register int c;
	register char *p;
	extern int optind;
	extern char *optarg;
	int w, h;

	while((c = getopt(argc, argv, "x:y:vV")) != EOF) {
		switch(c) {
		case 'x':
			tx = atoi(optarg);
			break;
		case 'y':
			ty = atoi(optarg);
			break;
		case 'v':
			verbose = !verbose;
			break;
		case 'V':
			banner();
			exit(0);
		case '?':
		default:
			usage();
		}
	}
	if((argc - optind) < 2) {
		usage();
	} else {
		setbuf(stdout, NULL);
TRACE(fprintf(stderr, "rescale: opening pxms\n"));
		ipx = px_ropen(argv[optind++]);
		sx = (long) ipx->px_width;
		sy = (long) ipx->px_height;
		if(tx == 0) {
			tx = sx;
		}
		if(ty == 0) {
			ty = sy;
		}
		if(verbose) {
			printf("inpxm (%ldx%ld) --> outpxm (%ldx%ld)\n",
				sx, sy, tx, ty);
		}
		opx = px_wopen(argv[optind++], ipx->px_type,
			(int) tx, (int) ty, ipx->px_depth, ipx->px_map);
TRACE(fprintf(stderr, "rescale: scaling pxms\n"));
		rescale(ipx, opx);
TRACE(fprintf(stderr, "rescale: closing pxms\n"));
		px_close(ipx);
		opx->px_map = NULL;	/* prevent double free() */
		px_close(opx);
	}
	exit(0);
}

#define	XSCALE(x)	(int)((((long)(x)) * sx) / tx)
#define	YSCALE(y)	(int)((((long)(y)) * sy) / ty)

rescale(ipx, opx)
PX_DEF *ipx, *opx;
{
	int y0, yy, y;
	register int x, xx, n, i;
	register PX_BYTE *ibuf, *obuf;

	ibuf = px_rowalloc(ipx->px_width, ipx->px_psize);
	obuf = px_rowalloc(opx->px_width, opx->px_psize);
	n = opx->px_psize;
	y0 = -1;
	for(y = 0; y < ty; ++y) {
		if(verbose) {
			progress((long) y, ty);
		}
		yy = YSCALE(y);
		while(yy > y0) {
TRACE(fprintf(stderr, "rescale: y=%d yy=%d y0=%d\n", y, yy, y0));
			px_rrow(ipx, ibuf);
			++y0;
		}
		if(n == 1) {
			for(x = 0; x < tx; ++x) {
				xx = XSCALE(x);
				obuf[x] = ibuf[xx];
			}
		} else {
			for(x = 0; x < tx; ++x) {
				xx = XSCALE(x);
				for(i = 0; i < n; ++i) {
					obuf[(x * n) + i] = ibuf[(xx * n) + i];
				}
			}
		}
		px_wrow(opx, obuf);
	}
	free(ibuf);
	free(obuf);
}


From the NET:
-------------

Watson, Andrew B.
"Ideal shrinking and expansion of discrete sequences"
NASA Technical Memorandum #88202

The best way to get ahold of this is simply to write the author
and ask for a copy:  A.B. Watson -> beau@gauss.arc.nasa.gov
*OR*

Andrew B. Watson
NASA Ames Research Center
Mail Stop 262-2
Moffett Field CA 94035-1000





Thanks again for all replies and information.
-- 
Jeffrey P. Bakke               Internet: bakke@plains.NoDak.edu 
                      UUCP    : ...!uunet!plains!bakke
           BITNET  : bakke@plains.bitnet