plan9port

fork of plan9port with libvec, libstr and libsdb
Log | Files | Refs | README | LICENSE

rotate.c (10409B)


      1 /*
      2  * rotate an image 180° in O(log Dx + log Dy) /dev/draw writes,
      3  * using an extra buffer same size as the image.
      4  *
      5  * the basic concept is that you can invert an array by inverting
      6  * the top half, inverting the bottom half, and then swapping them.
      7  * the code does this slightly backwards to ensure O(log n) runtime.
      8  * (If you do it wrong, you can get O(log² n) runtime.)
      9  *
     10  * This is usually overkill, but it speeds up slow remote
     11  * connections quite a bit.
     12  */
     13 
     14 #include <u.h>
     15 #include <libc.h>
     16 #include <bio.h>
     17 #include <draw.h>
     18 #include <thread.h>
     19 #include <cursor.h>
     20 #include "page.h"
     21 
     22 int ndraw = 0;
     23 enum {
     24 	Xaxis = 0,
     25 	Yaxis = 1,
     26 };
     27 
     28 Image *mtmp;
     29 
     30 void
     31 writefile(char *name, Image *im, int gran)
     32 {
     33 	static int c = 100;
     34 	int fd;
     35 	char buf[200];
     36 
     37 	snprint(buf, sizeof buf, "%d%s%d", c++, name, gran);
     38 	fd = create(buf, OWRITE, 0666);
     39 	if(fd < 0)
     40 		return;
     41 	writeimage(fd, im, 0);
     42 	close(fd);
     43 }
     44 
     45 void
     46 moveup(Image *im, Image *tmp, int a, int b, int c, int axis)
     47 {
     48 	Rectangle range;
     49 	Rectangle dr0, dr1;
     50 	Point p0, p1;
     51 
     52 	if(a == b || b == c)
     53 		return;
     54 
     55 	drawop(tmp, tmp->r, im, nil, im->r.min, S);
     56 
     57 	switch(axis){
     58 	case Xaxis:
     59 		range = Rect(a, im->r.min.y,  c, im->r.max.y);
     60 		dr0 = range;
     61 		dr0.max.x = dr0.min.x+(c-b);
     62 		p0 = Pt(b, im->r.min.y);
     63 
     64 		dr1 = range;
     65 		dr1.min.x = dr1.max.x-(b-a);
     66 		p1 = Pt(a, im->r.min.y);
     67 		break;
     68 	case Yaxis:
     69 	default:
     70 		range = Rect(im->r.min.x, a,  im->r.max.x, c);
     71 		dr0 = range;
     72 		dr0.max.y = dr0.min.y+(c-b);
     73 		p0 = Pt(im->r.min.x, b);
     74 
     75 		dr1 = range;
     76 		dr1.min.y = dr1.max.y-(b-a);
     77 		p1 = Pt(im->r.min.x, a);
     78 		break;
     79 	}
     80 	drawop(im, dr0, tmp, nil, p0, S);
     81 	drawop(im, dr1, tmp, nil, p1, S);
     82 }
     83 
     84 void
     85 interlace(Image *im, Image *tmp, int axis, int n, Image *mask, int gran)
     86 {
     87 	Point p0, p1;
     88 	Rectangle r0;
     89 
     90 	r0 = im->r;
     91 	switch(axis) {
     92 	case Xaxis:
     93 		r0.max.x = n;
     94 		p0 = (Point){gran, 0};
     95 		p1 = (Point){-gran, 0};
     96 		break;
     97 	case Yaxis:
     98 	default:
     99 		r0.max.y = n;
    100 		p0 = (Point){0, gran};
    101 		p1 = (Point){0, -gran};
    102 		break;
    103 	}
    104 
    105 	drawop(tmp, im->r, im, display->opaque, im->r.min, S);
    106 	gendrawop(im, r0, tmp, p0, mask, mask->r.min, S);
    107 	gendrawop(im, r0, tmp, p1, mask, p1, S);
    108 }
    109 
    110 /*
    111  * Halve the grating period in the mask.
    112  * The grating currently looks like
    113  * ####____####____####____####____
    114  * where #### is opacity.
    115  *
    116  * We want
    117  * ##__##__##__##__##__##__##__##__
    118  * which is achieved by shifting the mask
    119  * and drawing on itself through itself.
    120  * Draw doesn't actually allow this, so
    121  * we have to copy it first.
    122  *
    123  *     ####____####____####____####____ (dst)
    124  * +   ____####____####____####____#### (src)
    125  * in  __####____####____####____####__ (mask)
    126  * ===========================================
    127  *     ##__##__##__##__##__##__##__##__
    128  */
    129 int
    130 nextmask(Image *mask, int axis, int maskdim)
    131 {
    132 	Point o;
    133 
    134 	o = axis==Xaxis ? Pt(maskdim,0) : Pt(0,maskdim);
    135 	drawop(mtmp, mtmp->r, mask, nil, mask->r.min, S);
    136 	gendrawop(mask, mask->r, mtmp, o, mtmp, divpt(o,-2), S);
    137 //	writefile("mask", mask, maskdim/2);
    138 	return maskdim/2;
    139 }
    140 
    141 void
    142 shuffle(Image *im, Image *tmp, int axis, int n, Image *mask, int gran,
    143 	int lastnn)
    144 {
    145 	int nn, left;
    146 
    147 	if(gran == 0)
    148 		return;
    149 	left = n%(2*gran);
    150 	nn = n - left;
    151 
    152 	interlace(im, tmp, axis, nn, mask, gran);
    153 //	writefile("interlace", im, gran);
    154 
    155 	gran = nextmask(mask, axis, gran);
    156 	shuffle(im, tmp, axis, n, mask, gran, nn);
    157 //	writefile("shuffle", im, gran);
    158 	moveup(im, tmp, lastnn, nn, n, axis);
    159 //	writefile("move", im, gran);
    160 }
    161 
    162 void
    163 rot180(Image *im)
    164 {
    165 	Image *tmp, *tmp0;
    166 	Image *mask;
    167 	Rectangle rmask;
    168 	int gran;
    169 
    170 	if(chantodepth(im->chan) < 8){
    171 		/* this speeds things up dramatically; draw is too slow on sub-byte pixel sizes */
    172 		tmp0 = xallocimage(display, im->r, CMAP8, 0, DNofill);
    173 		drawop(tmp0, tmp0->r, im, nil, im->r.min, S);
    174 	}else
    175 		tmp0 = im;
    176 
    177 	tmp = xallocimage(display, tmp0->r, tmp0->chan, 0, DNofill);
    178 	if(tmp == nil){
    179 		if(tmp0 != im)
    180 			freeimage(tmp0);
    181 		return;
    182 	}
    183 	for(gran=1; gran<Dx(im->r); gran *= 2)
    184 		;
    185 	gran /= 4;
    186 
    187 	rmask.min = ZP;
    188 	rmask.max = (Point){2*gran, 100};
    189 
    190 	mask = xallocimage(display, rmask, GREY1, 1, DTransparent);
    191 	mtmp = xallocimage(display, rmask, GREY1, 1, DTransparent);
    192 	if(mask == nil || mtmp == nil) {
    193 		fprint(2, "out of memory during rot180: %r\n");
    194 		wexits("memory");
    195 	}
    196 	rmask.max.x = gran;
    197 	drawop(mask, rmask, display->opaque, nil, ZP, S);
    198 //	writefile("mask", mask, gran);
    199 	shuffle(im, tmp, Xaxis, Dx(im->r), mask, gran, 0);
    200 	freeimage(mask);
    201 	freeimage(mtmp);
    202 
    203 	for(gran=1; gran<Dy(im->r); gran *= 2)
    204 		;
    205 	gran /= 4;
    206 	rmask.max = (Point){100, 2*gran};
    207 	mask = xallocimage(display, rmask, GREY1, 1, DTransparent);
    208 	mtmp = xallocimage(display, rmask, GREY1, 1, DTransparent);
    209 	if(mask == nil || mtmp == nil) {
    210 		fprint(2, "out of memory during rot180: %r\n");
    211 		wexits("memory");
    212 	}
    213 	rmask.max.y = gran;
    214 	drawop(mask, rmask, display->opaque, nil, ZP, S);
    215 	shuffle(im, tmp, Yaxis, Dy(im->r), mask, gran, 0);
    216 	freeimage(mask);
    217 	freeimage(mtmp);
    218 	freeimage(tmp);
    219 	if(tmp0 != im)
    220 		freeimage(tmp0);
    221 }
    222 
    223 /* rotates an image 90 degrees clockwise */
    224 Image *
    225 rot90(Image *im)
    226 {
    227 	Image *tmp;
    228 	int i, j, dx, dy;
    229 
    230 	dx = Dx(im->r);
    231 	dy = Dy(im->r);
    232 	tmp = xallocimage(display, Rect(0, 0, dy, dx), im->chan, 0, DCyan);
    233 	if(tmp == nil) {
    234 		fprint(2, "out of memory during rot90: %r\n");
    235 		wexits("memory");
    236 	}
    237 
    238 	for(j = 0; j < dx; j++) {
    239 		for(i = 0; i < dy; i++) {
    240 			drawop(tmp, Rect(i, j, i+1, j+1), im, nil, Pt(j, dy-(i+1)), S);
    241 		}
    242 	}
    243 	freeimage(im);
    244 
    245 	return(tmp);
    246 }
    247 
    248 /* rotates an image 270 degrees clockwise */
    249 Image *
    250 rot270(Image *im)
    251 {
    252 	Image *tmp;
    253 	int i, j, dx, dy;
    254 
    255 	dx = Dx(im->r);
    256 	dy = Dy(im->r);
    257 	tmp = xallocimage(display, Rect(0, 0, dy, dx), im->chan, 0, DCyan);
    258 	if(tmp == nil) {
    259 		fprint(2, "out of memory during rot270: %r\n");
    260 		wexits("memory");
    261 	}
    262 
    263 	for(i = 0; i < dy; i++) {
    264 		for(j = 0; j < dx; j++) {
    265 			drawop(tmp, Rect(i, j, i+1, j+1), im, nil, Pt(dx-(j+1), i), S);
    266 		}
    267 	}
    268 	freeimage(im);
    269 
    270 	return(tmp);
    271 }
    272 
    273 /* from resample.c -- resize from → to using interpolation */
    274 
    275 
    276 #define K2 7	/* from -.7 to +.7 inclusive, meaning .2 into each adjacent pixel */
    277 #define NK (2*K2+1)
    278 double K[NK];
    279 
    280 double
    281 fac(int L)
    282 {
    283 	int i, f;
    284 
    285 	f = 1;
    286 	for(i=L; i>1; --i)
    287 		f *= i;
    288 	return f;
    289 }
    290 
    291 /*
    292  * i0(x) is the modified Bessel function, Σ (x/2)^2L / (L!)²
    293  * There are faster ways to calculate this, but we precompute
    294  * into a table so let's keep it simple.
    295  */
    296 double
    297 i0(double x)
    298 {
    299 	double v;
    300 	int L;
    301 
    302 	v = 1.0;
    303 	for(L=1; L<10; L++)
    304 		v += pow(x/2., 2*L)/pow(fac(L), 2);
    305 	return v;
    306 }
    307 
    308 double
    309 kaiser(double x, double t, double a)
    310 {
    311 	if(fabs(x) > t)
    312 		return 0.;
    313 	return i0(a*sqrt(1-(x*x/(t*t))))/i0(a);
    314 }
    315 
    316 
    317 void
    318 resamplex(uchar *in, int off, int d, int inx, uchar *out, int outx)
    319 {
    320 	int i, x, k;
    321 	double X, xx, v, rat;
    322 
    323 
    324 	rat = (double)inx/(double)outx;
    325 	for(x=0; x<outx; x++){
    326 		if(inx == outx){
    327 			/* don't resample if size unchanged */
    328 			out[off+x*d] = in[off+x*d];
    329 			continue;
    330 		}
    331 		v = 0.0;
    332 		X = x*rat;
    333 		for(k=-K2; k<=K2; k++){
    334 			xx = X + rat*k/10.;
    335 			i = xx;
    336 			if(i < 0)
    337 				i = 0;
    338 			if(i >= inx)
    339 				i = inx-1;
    340 			v += in[off+i*d] * K[K2+k];
    341 		}
    342 		out[off+x*d] = v;
    343 	}
    344 }
    345 
    346 void
    347 resampley(uchar **in, int off, int iny, uchar **out, int outy)
    348 {
    349 	int y, i, k;
    350 	double Y, yy, v, rat;
    351 
    352 	rat = (double)iny/(double)outy;
    353 	for(y=0; y<outy; y++){
    354 		if(iny == outy){
    355 			/* don't resample if size unchanged */
    356 			out[y][off] = in[y][off];
    357 			continue;
    358 		}
    359 		v = 0.0;
    360 		Y = y*rat;
    361 		for(k=-K2; k<=K2; k++){
    362 			yy = Y + rat*k/10.;
    363 			i = yy;
    364 			if(i < 0)
    365 				i = 0;
    366 			if(i >= iny)
    367 				i = iny-1;
    368 			v += in[i][off] * K[K2+k];
    369 		}
    370 		out[y][off] = v;
    371 	}
    372 
    373 }
    374 
    375 Image*
    376 resample(Image *from, Image *to)
    377 {
    378 	int i, j, bpl, nchan;
    379 	uchar **oscan, **nscan;
    380 	char tmp[20];
    381 	int xsize, ysize;
    382 	double v;
    383 	Image *t1, *t2;
    384 	ulong tchan;
    385 
    386 	for(i=-K2; i<=K2; i++){
    387 		K[K2+i] = kaiser(i/10., K2/10., 4.);
    388 	}
    389 
    390 	/* normalize */
    391 	v = 0.0;
    392 	for(i=0; i<NK; i++)
    393 		v += K[i];
    394 	for(i=0; i<NK; i++)
    395 		K[i] /= v;
    396 
    397 	switch(from->chan){
    398 	case GREY8:
    399 	case RGB24:
    400 	case RGBA32:
    401 	case ARGB32:
    402 	case XRGB32:
    403 		break;
    404 
    405 	case CMAP8:
    406 	case RGB15:
    407 	case RGB16:
    408 		tchan = RGB24;
    409 		goto Convert;
    410 
    411 	case GREY1:
    412 	case GREY2:
    413 	case GREY4:
    414 		tchan = GREY8;
    415 	Convert:
    416 		/* use library to convert to byte-per-chan form, then convert back */
    417 		t1 = xallocimage(display, Rect(0, 0, Dx(from->r), Dy(from->r)), tchan, 0, DNofill);
    418 		if(t1 == nil) {
    419 			fprint(2, "out of memory for temp image 1 in resample: %r\n");
    420 			wexits("memory");
    421 		}
    422 		drawop(t1, t1->r, from, nil, ZP, S);
    423 		t2 = xallocimage(display, to->r, tchan, 0, DNofill);
    424 		if(t2 == nil) {
    425 			fprint(2, "out of memory temp image 2 in resample: %r\n");
    426 			wexits("memory");
    427 		}
    428 		resample(t1, t2);
    429 		drawop(to, to->r, t2, nil, ZP, S);
    430 		freeimage(t1);
    431 		freeimage(t2);
    432 		return to;
    433 
    434 	default:
    435 		sysfatal("can't handle channel type %s", chantostr(tmp, from->chan));
    436 	}
    437 
    438 	xsize = Dx(to->r);
    439 	ysize = Dy(to->r);
    440 	oscan = malloc(Dy(from->r)*sizeof(uchar*));
    441 	nscan = malloc(max(ysize, Dy(from->r))*sizeof(uchar*));
    442 	if(oscan == nil || nscan == nil)
    443 		sysfatal("can't allocate: %r");
    444 
    445 	/* unload original image into scan lines */
    446 	bpl = bytesperline(from->r, from->depth);
    447 	for(i=0; i<Dy(from->r); i++){
    448 		oscan[i] = malloc(bpl);
    449 		if(oscan[i] == nil)
    450 			sysfatal("can't allocate: %r");
    451 		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);
    452 		if(j != bpl)
    453 			sysfatal("unloadimage");
    454 	}
    455 
    456 	/* allocate scan lines for destination. we do y first, so need at least Dy(from->r) lines */
    457 	bpl = bytesperline(Rect(0, 0, xsize, Dy(from->r)), from->depth);
    458 	for(i=0; i<max(ysize, Dy(from->r)); i++){
    459 		nscan[i] = malloc(bpl);
    460 		if(nscan[i] == nil)
    461 			sysfatal("can't allocate: %r");
    462 	}
    463 
    464 	/* resample in X */
    465 	nchan = from->depth/8;
    466 	for(i=0; i<Dy(from->r); i++){
    467 		for(j=0; j<nchan; j++){
    468 			if(j==0 && from->chan==XRGB32)
    469 				continue;
    470 			resamplex(oscan[i], j, nchan, Dx(from->r), nscan[i], xsize);
    471 		}
    472 		free(oscan[i]);
    473 		oscan[i] = nscan[i];
    474 		nscan[i] = malloc(bpl);
    475 		if(nscan[i] == nil)
    476 			sysfatal("can't allocate: %r");
    477 	}
    478 
    479 	/* resample in Y */
    480 	for(i=0; i<xsize; i++)
    481 		for(j=0; j<nchan; j++)
    482 			resampley(oscan, nchan*i+j, Dy(from->r), nscan, ysize);
    483 
    484 	/* pack data into destination */
    485 	bpl = bytesperline(to->r, from->depth);
    486 	for(i=0; i<ysize; i++){
    487 		j = loadimage(to, Rect(0, i, xsize, i+1), nscan[i], bpl);
    488 		if(j != bpl)
    489 			sysfatal("loadimage: %r");
    490 	}
    491 
    492 	for(i=0; i<Dy(from->r); i++){
    493 		free(oscan[i]);
    494 		free(nscan[i]);
    495 	}
    496 	free(oscan);
    497 	free(nscan);
    498 
    499 	return to;
    500 }