sys/src/cmd/page/rotate.c

   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 <event.h>
  19 #include "page.h"
  20
  21 int ndraw = 0;
  22 enum {
  23         Xaxis = 0,
  24         Yaxis = 1,
  25 };
  26
  27 Image *mtmp;
  28
  29 void
  30 writefile(char *name, Image *im, int gran)
  31 {
  32         static int c = 100;
  33         int fd;
  34         char buf[200];
  35
  36         snprint(buf, sizeof buf, "%d%s%d", c++, name, gran);
  37         fd = create(buf, OWRITE, 0666);
  38         if(fd < 0)
  39                 return;
  40         writeimage(fd, im, 0);
  41         close(fd);
  42 }
  43
  44 void
  45 moveup(Image *im, Image *tmp, int a, int b, int c, int axis)
  46 {
  47         Rectangle range;
  48         Rectangle dr0, dr1;
  49         Point p0, p1;
  50
  51         if(a == b || b == c)
  52                 return;
  53
  54         drawop(tmp, tmp->r, im, nil, im->r.min, S);
  55
  56         switch(axis){
  57         case Xaxis:
  58                 range = Rect(a, im->r.min.y,  c, im->r.max.y);
  59                 dr0 = range;
  60                 dr0.max.x = dr0.min.x+(c-b);
  61                 p0 = Pt(b, im->r.min.y);
  62
  63                 dr1 = range;
  64                 dr1.min.x = dr1.max.x-(b-a);
  65                 p1 = Pt(a, im->r.min.y);
  66                 break;
  67         case Yaxis:
  68                 range = Rect(im->r.min.x, a,  im->r.max.x, c);
  69                 dr0 = range;
  70                 dr0.max.y = dr0.min.y+(c-b);
  71                 p0 = Pt(im->r.min.x, b);
  72
  73                 dr1 = range;
  74                 dr1.min.y = dr1.max.y-(b-a);
  75                 p1 = Pt(im->r.min.x, a);
  76                 break;
  77         }
  78         drawop(im, dr0, tmp, nil, p0, S);
  79         drawop(im, dr1, tmp, nil, p1, S);
  80 }
  81
  82 void
  83 interlace(Image *im, Image *tmp, int axis, int n, Image *mask, int gran)
  84 {
  85         Point p0, p1;
  86         Rectangle r0, r1;
  87
  88         r0 = im->r;
  89         r1 = im->r;
  90         switch(axis) {
  91         case Xaxis:
  92                 r0.max.x = n;
  93                 r1.min.x = n;
  94                 p0 = (Point){gran, 0};
  95                 p1 = (Point){-gran, 0};
  96                 break;
  97         case Yaxis:
  98                 r0.max.y = n;
  99                 r1.min.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 δ;
 133
 134         δ = 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, δ, mtmp, divpt(δ,-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 τ, double α)
 310 {
 311         if(fabs(x) > τ)
 312                 return 0.;
 313         return i0(α*sqrt(1-(x*x/(τ*τ))))/i0(α);
 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 }