xref: /plan9/sys/src/cmd/page/rotate.c (revision a7529a1d594550a832531a572a2b23b6cce7f7a4)
17dd7cddfSDavid du Colombier /*
27dd7cddfSDavid du Colombier  * rotate an image 180° in O(log Dx + log Dy) /dev/draw writes,
37dd7cddfSDavid du Colombier  * using an extra buffer same size as the image.
47dd7cddfSDavid du Colombier  *
57dd7cddfSDavid du Colombier  * the basic concept is that you can invert an array by inverting
67dd7cddfSDavid du Colombier  * the top half, inverting the bottom half, and then swapping them.
77dd7cddfSDavid du Colombier  * the code does this slightly backwards to ensure O(log n) runtime.
87dd7cddfSDavid du Colombier  * (If you do it wrong, you can get O(log² n) runtime.)
97dd7cddfSDavid du Colombier  *
107dd7cddfSDavid du Colombier  * This is usually overkill, but it speeds up slow remote
117dd7cddfSDavid du Colombier  * connections quite a bit.
127dd7cddfSDavid du Colombier  */
137dd7cddfSDavid du Colombier 
147dd7cddfSDavid du Colombier #include <u.h>
157dd7cddfSDavid du Colombier #include <libc.h>
167dd7cddfSDavid du Colombier #include <bio.h>
177dd7cddfSDavid du Colombier #include <draw.h>
187dd7cddfSDavid du Colombier #include <event.h>
197dd7cddfSDavid du Colombier #include "page.h"
207dd7cddfSDavid du Colombier 
217dd7cddfSDavid du Colombier int ndraw = 0;
227dd7cddfSDavid du Colombier enum {
237dd7cddfSDavid du Colombier 	Xaxis = 0,
247dd7cddfSDavid du Colombier 	Yaxis = 1,
257dd7cddfSDavid du Colombier };
267dd7cddfSDavid du Colombier 
277dd7cddfSDavid du Colombier Image *mtmp;
287dd7cddfSDavid du Colombier 
297dd7cddfSDavid du Colombier void
writefile(char * name,Image * im,int gran)307dd7cddfSDavid du Colombier writefile(char *name, Image *im, int gran)
317dd7cddfSDavid du Colombier {
327dd7cddfSDavid du Colombier 	static int c = 100;
337dd7cddfSDavid du Colombier 	int fd;
347dd7cddfSDavid du Colombier 	char buf[200];
357dd7cddfSDavid du Colombier 
367dd7cddfSDavid du Colombier 	snprint(buf, sizeof buf, "%d%s%d", c++, name, gran);
377dd7cddfSDavid du Colombier 	fd = create(buf, OWRITE, 0666);
387dd7cddfSDavid du Colombier 	if(fd < 0)
397dd7cddfSDavid du Colombier 		return;
407dd7cddfSDavid du Colombier 	writeimage(fd, im, 0);
417dd7cddfSDavid du Colombier 	close(fd);
427dd7cddfSDavid du Colombier }
437dd7cddfSDavid du Colombier 
447dd7cddfSDavid du Colombier void
moveup(Image * im,Image * tmp,int a,int b,int c,int axis)457dd7cddfSDavid du Colombier moveup(Image *im, Image *tmp, int a, int b, int c, int axis)
467dd7cddfSDavid du Colombier {
477dd7cddfSDavid du Colombier 	Rectangle range;
487dd7cddfSDavid du Colombier 	Rectangle dr0, dr1;
497dd7cddfSDavid du Colombier 	Point p0, p1;
507dd7cddfSDavid du Colombier 
517dd7cddfSDavid du Colombier 	if(a == b || b == c)
527dd7cddfSDavid du Colombier 		return;
537dd7cddfSDavid du Colombier 
546b6b9ac8SDavid du Colombier 	drawop(tmp, tmp->r, im, nil, im->r.min, S);
557dd7cddfSDavid du Colombier 
567dd7cddfSDavid du Colombier 	switch(axis){
577dd7cddfSDavid du Colombier 	case Xaxis:
587dd7cddfSDavid du Colombier 		range = Rect(a, im->r.min.y,  c, im->r.max.y);
597dd7cddfSDavid du Colombier 		dr0 = range;
607dd7cddfSDavid du Colombier 		dr0.max.x = dr0.min.x+(c-b);
617dd7cddfSDavid du Colombier 		p0 = Pt(b, im->r.min.y);
627dd7cddfSDavid du Colombier 
637dd7cddfSDavid du Colombier 		dr1 = range;
647dd7cddfSDavid du Colombier 		dr1.min.x = dr1.max.x-(b-a);
657dd7cddfSDavid du Colombier 		p1 = Pt(a, im->r.min.y);
667dd7cddfSDavid du Colombier 		break;
677dd7cddfSDavid du Colombier 	case Yaxis:
687dd7cddfSDavid du Colombier 		range = Rect(im->r.min.x, a,  im->r.max.x, c);
697dd7cddfSDavid du Colombier 		dr0 = range;
707dd7cddfSDavid du Colombier 		dr0.max.y = dr0.min.y+(c-b);
717dd7cddfSDavid du Colombier 		p0 = Pt(im->r.min.x, b);
727dd7cddfSDavid du Colombier 
737dd7cddfSDavid du Colombier 		dr1 = range;
747dd7cddfSDavid du Colombier 		dr1.min.y = dr1.max.y-(b-a);
757dd7cddfSDavid du Colombier 		p1 = Pt(im->r.min.x, a);
767dd7cddfSDavid du Colombier 		break;
777dd7cddfSDavid du Colombier 	}
786b6b9ac8SDavid du Colombier 	drawop(im, dr0, tmp, nil, p0, S);
796b6b9ac8SDavid du Colombier 	drawop(im, dr1, tmp, nil, p1, S);
807dd7cddfSDavid du Colombier }
817dd7cddfSDavid du Colombier 
827dd7cddfSDavid du Colombier void
interlace(Image * im,Image * tmp,int axis,int n,Image * mask,int gran)837dd7cddfSDavid du Colombier interlace(Image *im, Image *tmp, int axis, int n, Image *mask, int gran)
847dd7cddfSDavid du Colombier {
857dd7cddfSDavid du Colombier 	Point p0, p1;
867dd7cddfSDavid du Colombier 	Rectangle r0, r1;
877dd7cddfSDavid du Colombier 
887dd7cddfSDavid du Colombier 	r0 = im->r;
897dd7cddfSDavid du Colombier 	r1 = im->r;
907dd7cddfSDavid du Colombier 	switch(axis) {
917dd7cddfSDavid du Colombier 	case Xaxis:
927dd7cddfSDavid du Colombier 		r0.max.x = n;
937dd7cddfSDavid du Colombier 		r1.min.x = n;
947dd7cddfSDavid du Colombier 		p0 = (Point){gran, 0};
957dd7cddfSDavid du Colombier 		p1 = (Point){-gran, 0};
967dd7cddfSDavid du Colombier 		break;
977dd7cddfSDavid du Colombier 	case Yaxis:
987dd7cddfSDavid du Colombier 		r0.max.y = n;
997dd7cddfSDavid du Colombier 		r1.min.y = n;
1007dd7cddfSDavid du Colombier 		p0 = (Point){0, gran};
1017dd7cddfSDavid du Colombier 		p1 = (Point){0, -gran};
1027dd7cddfSDavid du Colombier 		break;
1037dd7cddfSDavid du Colombier 	}
1047dd7cddfSDavid du Colombier 
1056b6b9ac8SDavid du Colombier 	drawop(tmp, im->r, im, display->opaque, im->r.min, S);
1066b6b9ac8SDavid du Colombier 	gendrawop(im, r0, tmp, p0, mask, mask->r.min, S);
1076b6b9ac8SDavid du Colombier 	gendrawop(im, r0, tmp, p1, mask, p1, S);
1087dd7cddfSDavid du Colombier }
1097dd7cddfSDavid du Colombier 
1107dd7cddfSDavid du Colombier /*
1117dd7cddfSDavid du Colombier  * Halve the grating period in the mask.
1127dd7cddfSDavid du Colombier  * The grating currently looks like
1137dd7cddfSDavid du Colombier  * ####____####____####____####____
1147dd7cddfSDavid du Colombier  * where #### is opacity.
1157dd7cddfSDavid du Colombier  *
1167dd7cddfSDavid du Colombier  * We want
1177dd7cddfSDavid du Colombier  * ##__##__##__##__##__##__##__##__
1187dd7cddfSDavid du Colombier  * which is achieved by shifting the mask
1197dd7cddfSDavid du Colombier  * and drawing on itself through itself.
1207dd7cddfSDavid du Colombier  * Draw doesn't actually allow this, so
1217dd7cddfSDavid du Colombier  * we have to copy it first.
1227dd7cddfSDavid du Colombier  *
1237dd7cddfSDavid du Colombier  *     ####____####____####____####____ (dst)
1247dd7cddfSDavid du Colombier  * +   ____####____####____####____#### (src)
1257dd7cddfSDavid du Colombier  * in  __####____####____####____####__ (mask)
1267dd7cddfSDavid du Colombier  * ===========================================
1277dd7cddfSDavid du Colombier  *     ##__##__##__##__##__##__##__##__
1287dd7cddfSDavid du Colombier  */
1297dd7cddfSDavid du Colombier int
nextmask(Image * mask,int axis,int maskdim)1307dd7cddfSDavid du Colombier nextmask(Image *mask, int axis, int maskdim)
1317dd7cddfSDavid du Colombier {
1327dd7cddfSDavid du Colombier 	Point δ;
1337dd7cddfSDavid du Colombier 
1347dd7cddfSDavid du Colombier 	δ = axis==Xaxis ? Pt(maskdim,0) : Pt(0,maskdim);
1356b6b9ac8SDavid du Colombier 	drawop(mtmp, mtmp->r, mask, nil, mask->r.min, S);
1366b6b9ac8SDavid du Colombier 	gendrawop(mask, mask->r, mtmp, δ, mtmp, divpt(δ,-2), S);
1377dd7cddfSDavid du Colombier //	writefile("mask", mask, maskdim/2);
1387dd7cddfSDavid du Colombier 	return maskdim/2;
1397dd7cddfSDavid du Colombier }
1407dd7cddfSDavid du Colombier 
1417dd7cddfSDavid du Colombier void
shuffle(Image * im,Image * tmp,int axis,int n,Image * mask,int gran,int lastnn)1427dd7cddfSDavid du Colombier shuffle(Image *im, Image *tmp, int axis, int n, Image *mask, int gran,
1437dd7cddfSDavid du Colombier 	int lastnn)
1447dd7cddfSDavid du Colombier {
1457dd7cddfSDavid du Colombier 	int nn, left;
1467dd7cddfSDavid du Colombier 
1477dd7cddfSDavid du Colombier 	if(gran == 0)
1487dd7cddfSDavid du Colombier 		return;
1497dd7cddfSDavid du Colombier 	left = n%(2*gran);
1507dd7cddfSDavid du Colombier 	nn = n - left;
1517dd7cddfSDavid du Colombier 
1527dd7cddfSDavid du Colombier 	interlace(im, tmp, axis, nn, mask, gran);
1537dd7cddfSDavid du Colombier //	writefile("interlace", im, gran);
1547dd7cddfSDavid du Colombier 
1557dd7cddfSDavid du Colombier 	gran = nextmask(mask, axis, gran);
1567dd7cddfSDavid du Colombier 	shuffle(im, tmp, axis, n, mask, gran, nn);
1577dd7cddfSDavid du Colombier //	writefile("shuffle", im, gran);
1587dd7cddfSDavid du Colombier 	moveup(im, tmp, lastnn, nn, n, axis);
1597dd7cddfSDavid du Colombier //	writefile("move", im, gran);
1607dd7cddfSDavid du Colombier }
1617dd7cddfSDavid du Colombier 
1627dd7cddfSDavid du Colombier void
rot180(Image * im)1637dd7cddfSDavid du Colombier rot180(Image *im)
1647dd7cddfSDavid du Colombier {
1655d459b5aSDavid du Colombier 	Image *tmp, *tmp0;
1667dd7cddfSDavid du Colombier 	Image *mask;
1677dd7cddfSDavid du Colombier 	Rectangle rmask;
1687dd7cddfSDavid du Colombier 	int gran;
1697dd7cddfSDavid du Colombier 
1705d459b5aSDavid du Colombier 	if(chantodepth(im->chan) < 8){
1715d459b5aSDavid du Colombier 		/* this speeds things up dramatically; draw is too slow on sub-byte pixel sizes */
1725d459b5aSDavid du Colombier 		tmp0 = xallocimage(display, im->r, CMAP8, 0, DNofill);
1736b6b9ac8SDavid du Colombier 		drawop(tmp0, tmp0->r, im, nil, im->r.min, S);
1745d459b5aSDavid du Colombier 	}else
1755d459b5aSDavid du Colombier 		tmp0 = im;
1767dd7cddfSDavid du Colombier 
1775d459b5aSDavid du Colombier 	tmp = xallocimage(display, tmp0->r, tmp0->chan, 0, DNofill);
1785d459b5aSDavid du Colombier 	if(tmp == nil){
1795d459b5aSDavid du Colombier 		if(tmp0 != im)
1805d459b5aSDavid du Colombier 			freeimage(tmp0);
1815d459b5aSDavid du Colombier 		return;
1825d459b5aSDavid du Colombier 	}
1837dd7cddfSDavid du Colombier 	for(gran=1; gran<Dx(im->r); gran *= 2)
1847dd7cddfSDavid du Colombier 		;
1857dd7cddfSDavid du Colombier 	gran /= 4;
1867dd7cddfSDavid du Colombier 
1877dd7cddfSDavid du Colombier 	rmask.min = ZP;
1887dd7cddfSDavid du Colombier 	rmask.max = (Point){2*gran, 100};
1897dd7cddfSDavid du Colombier 
1907dd7cddfSDavid du Colombier 	mask = xallocimage(display, rmask, GREY1, 1, DTransparent);
1917dd7cddfSDavid du Colombier 	mtmp = xallocimage(display, rmask, GREY1, 1, DTransparent);
1925316891fSDavid du Colombier 	if(mask == nil || mtmp == nil) {
1935316891fSDavid du Colombier 		fprint(2, "out of memory during rot180: %r\n");
1945316891fSDavid du Colombier 		wexits("memory");
1955316891fSDavid du Colombier 	}
1967dd7cddfSDavid du Colombier 	rmask.max.x = gran;
1976b6b9ac8SDavid du Colombier 	drawop(mask, rmask, display->opaque, nil, ZP, S);
1987dd7cddfSDavid du Colombier //	writefile("mask", mask, gran);
1997dd7cddfSDavid du Colombier 	shuffle(im, tmp, Xaxis, Dx(im->r), mask, gran, 0);
2007dd7cddfSDavid du Colombier 	freeimage(mask);
2017dd7cddfSDavid du Colombier 	freeimage(mtmp);
2027dd7cddfSDavid du Colombier 
2037dd7cddfSDavid du Colombier 	for(gran=1; gran<Dy(im->r); gran *= 2)
2047dd7cddfSDavid du Colombier 		;
2057dd7cddfSDavid du Colombier 	gran /= 4;
2067dd7cddfSDavid du Colombier 	rmask.max = (Point){100, 2*gran};
2077dd7cddfSDavid du Colombier 	mask = xallocimage(display, rmask, GREY1, 1, DTransparent);
2087dd7cddfSDavid du Colombier 	mtmp = xallocimage(display, rmask, GREY1, 1, DTransparent);
2095316891fSDavid du Colombier 	if(mask == nil || mtmp == nil) {
2105316891fSDavid du Colombier 		fprint(2, "out of memory during rot180: %r\n");
2115316891fSDavid du Colombier 		wexits("memory");
2125316891fSDavid du Colombier 	}
2137dd7cddfSDavid du Colombier 	rmask.max.y = gran;
2146b6b9ac8SDavid du Colombier 	drawop(mask, rmask, display->opaque, nil, ZP, S);
2157dd7cddfSDavid du Colombier 	shuffle(im, tmp, Yaxis, Dy(im->r), mask, gran, 0);
2167dd7cddfSDavid du Colombier 	freeimage(mask);
2177dd7cddfSDavid du Colombier 	freeimage(mtmp);
2187dd7cddfSDavid du Colombier 	freeimage(tmp);
2195d459b5aSDavid du Colombier 	if(tmp0 != im)
2205d459b5aSDavid du Colombier 		freeimage(tmp0);
2217dd7cddfSDavid du Colombier }
22208fd2d13SDavid du Colombier 
22308fd2d13SDavid du Colombier /* rotates an image 90 degrees clockwise */
22408fd2d13SDavid du Colombier Image *
rot90(Image * im)22508fd2d13SDavid du Colombier rot90(Image *im)
22608fd2d13SDavid du Colombier {
22708fd2d13SDavid du Colombier 	Image *tmp;
22808fd2d13SDavid du Colombier 	int i, j, dx, dy;
22908fd2d13SDavid du Colombier 
23008fd2d13SDavid du Colombier 	dx = Dx(im->r);
23108fd2d13SDavid du Colombier 	dy = Dy(im->r);
23208fd2d13SDavid du Colombier 	tmp = xallocimage(display, Rect(0, 0, dy, dx), im->chan, 0, DCyan);
2335316891fSDavid du Colombier 	if(tmp == nil) {
2345316891fSDavid du Colombier 		fprint(2, "out of memory during rot90: %r\n");
2355316891fSDavid du Colombier 		wexits("memory");
2365316891fSDavid du Colombier 	}
23708fd2d13SDavid du Colombier 
23808fd2d13SDavid du Colombier 	for(j = 0; j < dx; j++) {
23908fd2d13SDavid du Colombier 		for(i = 0; i < dy; i++) {
2405316891fSDavid du Colombier 			drawop(tmp, Rect(i, j, i+1, j+1), im, nil, Pt(j, dy-(i+1)), S);
24108fd2d13SDavid du Colombier 		}
24208fd2d13SDavid du Colombier 	}
24308fd2d13SDavid du Colombier 	freeimage(im);
24408fd2d13SDavid du Colombier 
24508fd2d13SDavid du Colombier 	return(tmp);
24608fd2d13SDavid du Colombier }
24708fd2d13SDavid du Colombier 
248*a7529a1dSDavid du Colombier /* rotates an image 270 degrees clockwise */
249*a7529a1dSDavid du Colombier Image *
rot270(Image * im)250*a7529a1dSDavid du Colombier rot270(Image *im)
251*a7529a1dSDavid du Colombier {
252*a7529a1dSDavid du Colombier 	Image *tmp;
253*a7529a1dSDavid du Colombier 	int i, j, dx, dy;
254*a7529a1dSDavid du Colombier 
255*a7529a1dSDavid du Colombier 	dx = Dx(im->r);
256*a7529a1dSDavid du Colombier 	dy = Dy(im->r);
257*a7529a1dSDavid du Colombier 	tmp = xallocimage(display, Rect(0, 0, dy, dx), im->chan, 0, DCyan);
258*a7529a1dSDavid du Colombier 	if(tmp == nil) {
259*a7529a1dSDavid du Colombier 		fprint(2, "out of memory during rot270: %r\n");
260*a7529a1dSDavid du Colombier 		wexits("memory");
261*a7529a1dSDavid du Colombier 	}
262*a7529a1dSDavid du Colombier 
263*a7529a1dSDavid du Colombier 	for(i = 0; i < dy; i++) {
264*a7529a1dSDavid du Colombier 		for(j = 0; j < dx; j++) {
265*a7529a1dSDavid du Colombier 			drawop(tmp, Rect(i, j, i+1, j+1), im, nil, Pt(dx-(j+1), i), S);
266*a7529a1dSDavid du Colombier 		}
267*a7529a1dSDavid du Colombier 	}
268*a7529a1dSDavid du Colombier 	freeimage(im);
269*a7529a1dSDavid du Colombier 
270*a7529a1dSDavid du Colombier 	return(tmp);
271*a7529a1dSDavid du Colombier }
272*a7529a1dSDavid du Colombier 
27308fd2d13SDavid du Colombier /* from resample.c -- resize from → to using interpolation */
27408fd2d13SDavid du Colombier 
27508fd2d13SDavid du Colombier 
27608fd2d13SDavid du Colombier #define K2 7	/* from -.7 to +.7 inclusive, meaning .2 into each adjacent pixel */
27708fd2d13SDavid du Colombier #define NK (2*K2+1)
27808fd2d13SDavid du Colombier double K[NK];
27908fd2d13SDavid du Colombier 
28008fd2d13SDavid du Colombier double
fac(int L)28108fd2d13SDavid du Colombier fac(int L)
28208fd2d13SDavid du Colombier {
28308fd2d13SDavid du Colombier 	int i, f;
28408fd2d13SDavid du Colombier 
28508fd2d13SDavid du Colombier 	f = 1;
28608fd2d13SDavid du Colombier 	for(i=L; i>1; --i)
28708fd2d13SDavid du Colombier 		f *= i;
28808fd2d13SDavid du Colombier 	return f;
28908fd2d13SDavid du Colombier }
29008fd2d13SDavid du Colombier 
29108fd2d13SDavid du Colombier /*
29208fd2d13SDavid du Colombier  * i0(x) is the modified Bessel function, Σ (x/2)^2L / (L!)²
29308fd2d13SDavid du Colombier  * There are faster ways to calculate this, but we precompute
29408fd2d13SDavid du Colombier  * into a table so let's keep it simple.
29508fd2d13SDavid du Colombier  */
29608fd2d13SDavid du Colombier double
i0(double x)29708fd2d13SDavid du Colombier i0(double x)
29808fd2d13SDavid du Colombier {
29908fd2d13SDavid du Colombier 	double v;
30008fd2d13SDavid du Colombier 	int L;
30108fd2d13SDavid du Colombier 
30208fd2d13SDavid du Colombier 	v = 1.0;
30308fd2d13SDavid du Colombier 	for(L=1; L<10; L++)
30408fd2d13SDavid du Colombier 		v += pow(x/2., 2*L)/pow(fac(L), 2);
30508fd2d13SDavid du Colombier 	return v;
30608fd2d13SDavid du Colombier }
30708fd2d13SDavid du Colombier 
30808fd2d13SDavid du Colombier double
kaiser(double x,doubleτ,doubleα)30908fd2d13SDavid du Colombier kaiser(double x, double τ, double α)
31008fd2d13SDavid du Colombier {
31108fd2d13SDavid du Colombier 	if(fabs(x) > τ)
31208fd2d13SDavid du Colombier 		return 0.;
31308fd2d13SDavid du Colombier 	return i0(α*sqrt(1-(x*x/(τ*τ))))/i0(α);
31408fd2d13SDavid du Colombier }
31508fd2d13SDavid du Colombier 
31608fd2d13SDavid du Colombier 
31708fd2d13SDavid du Colombier void
resamplex(uchar * in,int off,int d,int inx,uchar * out,int outx)31808fd2d13SDavid du Colombier resamplex(uchar *in, int off, int d, int inx, uchar *out, int outx)
31908fd2d13SDavid du Colombier {
32008fd2d13SDavid du Colombier 	int i, x, k;
32108fd2d13SDavid du Colombier 	double X, xx, v, rat;
32208fd2d13SDavid du Colombier 
32308fd2d13SDavid du Colombier 
32408fd2d13SDavid du Colombier 	rat = (double)inx/(double)outx;
32508fd2d13SDavid du Colombier 	for(x=0; x<outx; x++){
32608fd2d13SDavid du Colombier 		if(inx == outx){
32708fd2d13SDavid du Colombier 			/* don't resample if size unchanged */
32808fd2d13SDavid du Colombier 			out[off+x*d] = in[off+x*d];
32908fd2d13SDavid du Colombier 			continue;
33008fd2d13SDavid du Colombier 		}
33108fd2d13SDavid du Colombier 		v = 0.0;
33208fd2d13SDavid du Colombier 		X = x*rat;
33308fd2d13SDavid du Colombier 		for(k=-K2; k<=K2; k++){
33408fd2d13SDavid du Colombier 			xx = X + rat*k/10.;
33508fd2d13SDavid du Colombier 			i = xx;
33608fd2d13SDavid du Colombier 			if(i < 0)
33708fd2d13SDavid du Colombier 				i = 0;
33808fd2d13SDavid du Colombier 			if(i >= inx)
33908fd2d13SDavid du Colombier 				i = inx-1;
34008fd2d13SDavid du Colombier 			v += in[off+i*d] * K[K2+k];
34108fd2d13SDavid du Colombier 		}
34208fd2d13SDavid du Colombier 		out[off+x*d] = v;
34308fd2d13SDavid du Colombier 	}
34408fd2d13SDavid du Colombier }
34508fd2d13SDavid du Colombier 
34608fd2d13SDavid du Colombier void
resampley(uchar ** in,int off,int iny,uchar ** out,int outy)34708fd2d13SDavid du Colombier resampley(uchar **in, int off, int iny, uchar **out, int outy)
34808fd2d13SDavid du Colombier {
34908fd2d13SDavid du Colombier 	int y, i, k;
35008fd2d13SDavid du Colombier 	double Y, yy, v, rat;
35108fd2d13SDavid du Colombier 
35208fd2d13SDavid du Colombier 	rat = (double)iny/(double)outy;
35308fd2d13SDavid du Colombier 	for(y=0; y<outy; y++){
35408fd2d13SDavid du Colombier 		if(iny == outy){
35508fd2d13SDavid du Colombier 			/* don't resample if size unchanged */
35608fd2d13SDavid du Colombier 			out[y][off] = in[y][off];
35708fd2d13SDavid du Colombier 			continue;
35808fd2d13SDavid du Colombier 		}
35908fd2d13SDavid du Colombier 		v = 0.0;
36008fd2d13SDavid du Colombier 		Y = y*rat;
36108fd2d13SDavid du Colombier 		for(k=-K2; k<=K2; k++){
36208fd2d13SDavid du Colombier 			yy = Y + rat*k/10.;
36308fd2d13SDavid du Colombier 			i = yy;
36408fd2d13SDavid du Colombier 			if(i < 0)
36508fd2d13SDavid du Colombier 				i = 0;
36608fd2d13SDavid du Colombier 			if(i >= iny)
36708fd2d13SDavid du Colombier 				i = iny-1;
36808fd2d13SDavid du Colombier 			v += in[i][off] * K[K2+k];
36908fd2d13SDavid du Colombier 		}
37008fd2d13SDavid du Colombier 		out[y][off] = v;
37108fd2d13SDavid du Colombier 	}
37208fd2d13SDavid du Colombier 
37308fd2d13SDavid du Colombier }
37408fd2d13SDavid du Colombier 
37508fd2d13SDavid du Colombier Image*
resample(Image * from,Image * to)37608fd2d13SDavid du Colombier resample(Image *from, Image *to)
37708fd2d13SDavid du Colombier {
37808fd2d13SDavid du Colombier 	int i, j, bpl, nchan;
37908fd2d13SDavid du Colombier 	uchar **oscan, **nscan;
38008fd2d13SDavid du Colombier 	char tmp[20];
38108fd2d13SDavid du Colombier 	int xsize, ysize;
38208fd2d13SDavid du Colombier 	double v;
38308fd2d13SDavid du Colombier 	Image *t1, *t2;
38408fd2d13SDavid du Colombier 	ulong tchan;
38508fd2d13SDavid du Colombier 
38608fd2d13SDavid du Colombier 	for(i=-K2; i<=K2; i++){
38708fd2d13SDavid du Colombier 		K[K2+i] = kaiser(i/10., K2/10., 4.);
38808fd2d13SDavid du Colombier 	}
38908fd2d13SDavid du Colombier 
39008fd2d13SDavid du Colombier 	/* normalize */
39108fd2d13SDavid du Colombier 	v = 0.0;
39208fd2d13SDavid du Colombier 	for(i=0; i<NK; i++)
39308fd2d13SDavid du Colombier 		v += K[i];
39408fd2d13SDavid du Colombier 	for(i=0; i<NK; i++)
39508fd2d13SDavid du Colombier 		K[i] /= v;
39608fd2d13SDavid du Colombier 
39708fd2d13SDavid du Colombier 	switch(from->chan){
39808fd2d13SDavid du Colombier 	case GREY8:
39908fd2d13SDavid du Colombier 	case RGB24:
40008fd2d13SDavid du Colombier 	case RGBA32:
40108fd2d13SDavid du Colombier 	case ARGB32:
40208fd2d13SDavid du Colombier 	case XRGB32:
40308fd2d13SDavid du Colombier 		break;
40408fd2d13SDavid du Colombier 
40508fd2d13SDavid du Colombier 	case CMAP8:
40608fd2d13SDavid du Colombier 	case RGB15:
40708fd2d13SDavid du Colombier 	case RGB16:
40808fd2d13SDavid du Colombier 		tchan = RGB24;
40908fd2d13SDavid du Colombier 		goto Convert;
41008fd2d13SDavid du Colombier 
41108fd2d13SDavid du Colombier 	case GREY1:
41208fd2d13SDavid du Colombier 	case GREY2:
41308fd2d13SDavid du Colombier 	case GREY4:
41408fd2d13SDavid du Colombier 		tchan = GREY8;
41508fd2d13SDavid du Colombier 	Convert:
41608fd2d13SDavid du Colombier 		/* use library to convert to byte-per-chan form, then convert back */
4175316891fSDavid du Colombier 		t1 = xallocimage(display, Rect(0, 0, Dx(from->r), Dy(from->r)), tchan, 0, DNofill);
4185316891fSDavid du Colombier 		if(t1 == nil) {
4195316891fSDavid du Colombier 			fprint(2, "out of memory for temp image 1 in resample: %r\n");
4205316891fSDavid du Colombier 			wexits("memory");
4215316891fSDavid du Colombier 		}
4225316891fSDavid du Colombier 		drawop(t1, t1->r, from, nil, ZP, S);
4235316891fSDavid du Colombier 		t2 = xallocimage(display, to->r, tchan, 0, DNofill);
4245316891fSDavid du Colombier 		if(t2 == nil) {
4255316891fSDavid du Colombier 			fprint(2, "out of memory temp image 2 in resample: %r\n");
4265316891fSDavid du Colombier 			wexits("memory");
4275316891fSDavid du Colombier 		}
42808fd2d13SDavid du Colombier 		resample(t1, t2);
4295316891fSDavid du Colombier 		drawop(to, to->r, t2, nil, ZP, S);
43008fd2d13SDavid du Colombier 		freeimage(t1);
43108fd2d13SDavid du Colombier 		freeimage(t2);
43208fd2d13SDavid du Colombier 		return to;
43308fd2d13SDavid du Colombier 
43408fd2d13SDavid du Colombier 	default:
43508fd2d13SDavid du Colombier 		sysfatal("can't handle channel type %s", chantostr(tmp, from->chan));
43608fd2d13SDavid du Colombier 	}
43708fd2d13SDavid du Colombier 
43808fd2d13SDavid du Colombier 	xsize = Dx(to->r);
43908fd2d13SDavid du Colombier 	ysize = Dy(to->r);
44008fd2d13SDavid du Colombier 	oscan = malloc(Dy(from->r)*sizeof(uchar*));
44108fd2d13SDavid du Colombier 	nscan = malloc(max(ysize, Dy(from->r))*sizeof(uchar*));
44208fd2d13SDavid du Colombier 	if(oscan == nil || nscan == nil)
44308fd2d13SDavid du Colombier 		sysfatal("can't allocate: %r");
44408fd2d13SDavid du Colombier 
44508fd2d13SDavid du Colombier 	/* unload original image into scan lines */
44608fd2d13SDavid du Colombier 	bpl = bytesperline(from->r, from->depth);
44708fd2d13SDavid du Colombier 	for(i=0; i<Dy(from->r); i++){
44808fd2d13SDavid du Colombier 		oscan[i] = malloc(bpl);
44908fd2d13SDavid du Colombier 		if(oscan[i] == nil)
45008fd2d13SDavid du Colombier 			sysfatal("can't allocate: %r");
45108fd2d13SDavid du Colombier 		j = unloadimage(from, Rect(from->r.min.x, from->r.min.y+i, from->r.max.x, from->r.min.y+i+1), oscan[i], bpl);
45208fd2d13SDavid du Colombier 		if(j != bpl)
45308fd2d13SDavid du Colombier 			sysfatal("unloadimage");
45408fd2d13SDavid du Colombier 	}
45508fd2d13SDavid du Colombier 
45608fd2d13SDavid du Colombier 	/* allocate scan lines for destination. we do y first, so need at least Dy(from->r) lines */
45708fd2d13SDavid du Colombier 	bpl = bytesperline(Rect(0, 0, xsize, Dy(from->r)), from->depth);
45808fd2d13SDavid du Colombier 	for(i=0; i<max(ysize, Dy(from->r)); i++){
45908fd2d13SDavid du Colombier 		nscan[i] = malloc(bpl);
46008fd2d13SDavid du Colombier 		if(nscan[i] == nil)
46108fd2d13SDavid du Colombier 			sysfatal("can't allocate: %r");
46208fd2d13SDavid du Colombier 	}
46308fd2d13SDavid du Colombier 
46408fd2d13SDavid du Colombier 	/* resample in X */
46508fd2d13SDavid du Colombier 	nchan = from->depth/8;
46608fd2d13SDavid du Colombier 	for(i=0; i<Dy(from->r); i++){
46708fd2d13SDavid du Colombier 		for(j=0; j<nchan; j++){
46808fd2d13SDavid du Colombier 			if(j==0 && from->chan==XRGB32)
46908fd2d13SDavid du Colombier 				continue;
47008fd2d13SDavid du Colombier 			resamplex(oscan[i], j, nchan, Dx(from->r), nscan[i], xsize);
47108fd2d13SDavid du Colombier 		}
47208fd2d13SDavid du Colombier 		free(oscan[i]);
47308fd2d13SDavid du Colombier 		oscan[i] = nscan[i];
47408fd2d13SDavid du Colombier 		nscan[i] = malloc(bpl);
47508fd2d13SDavid du Colombier 		if(nscan[i] == nil)
47608fd2d13SDavid du Colombier 			sysfatal("can't allocate: %r");
47708fd2d13SDavid du Colombier 	}
47808fd2d13SDavid du Colombier 
47908fd2d13SDavid du Colombier 	/* resample in Y */
48008fd2d13SDavid du Colombier 	for(i=0; i<xsize; i++)
48108fd2d13SDavid du Colombier 		for(j=0; j<nchan; j++)
48208fd2d13SDavid du Colombier 			resampley(oscan, nchan*i+j, Dy(from->r), nscan, ysize);
48308fd2d13SDavid du Colombier 
48408fd2d13SDavid du Colombier 	/* pack data into destination */
48508fd2d13SDavid du Colombier 	bpl = bytesperline(to->r, from->depth);
48608fd2d13SDavid du Colombier 	for(i=0; i<ysize; i++){
48708fd2d13SDavid du Colombier 		j = loadimage(to, Rect(0, i, xsize, i+1), nscan[i], bpl);
48808fd2d13SDavid du Colombier 		if(j != bpl)
48908fd2d13SDavid du Colombier 			sysfatal("loadimage: %r");
49008fd2d13SDavid du Colombier 	}
49108fd2d13SDavid du Colombier 
49208fd2d13SDavid du Colombier 	for(i=0; i<Dy(from->r); i++){
49308fd2d13SDavid du Colombier 		free(oscan[i]);
49408fd2d13SDavid du Colombier 		free(nscan[i]);
49508fd2d13SDavid du Colombier 	}
49608fd2d13SDavid du Colombier 	free(oscan);
49708fd2d13SDavid du Colombier 	free(nscan);
49808fd2d13SDavid du Colombier 
49908fd2d13SDavid du Colombier 	return to;
50008fd2d13SDavid du Colombier }
501