readjpg.c (33031B)
1 #include <u.h> 2 #include <libc.h> 3 #include <bio.h> 4 #include <draw.h> 5 #include "imagefile.h" 6 7 enum { 8 /* Constants, all preceded by byte 0xFF */ 9 SOF =0xC0, /* Start of Frame */ 10 SOF2=0xC2, /* Start of Frame; progressive Huffman */ 11 JPG =0xC8, /* Reserved for JPEG extensions */ 12 DHT =0xC4, /* Define Huffman Tables */ 13 DAC =0xCC, /* Arithmetic coding conditioning */ 14 RST =0xD0, /* Restart interval termination */ 15 RST7 =0xD7, /* Restart interval termination (highest value) */ 16 SOI =0xD8, /* Start of Image */ 17 EOI =0xD9, /* End of Image */ 18 SOS =0xDA, /* Start of Scan */ 19 DQT =0xDB, /* Define quantization tables */ 20 DNL =0xDC, /* Define number of lines */ 21 DRI =0xDD, /* Define restart interval */ 22 DHP =0xDE, /* Define hierarchical progression */ 23 EXP =0xDF, /* Expand reference components */ 24 APPn =0xE0, /* Reserved for application segments */ 25 JPGn =0xF0, /* Reserved for JPEG extensions */ 26 COM =0xFE, /* Comment */ 27 28 CLAMPOFF = 300, 29 NCLAMP = CLAMPOFF+700 30 }; 31 32 typedef struct Framecomp Framecomp; 33 typedef struct Header Header; 34 typedef struct Huffman Huffman; 35 36 struct Framecomp /* Frame component specifier from SOF marker */ 37 { 38 int C; 39 int H; 40 int V; 41 int Tq; 42 }; 43 44 struct Huffman 45 { 46 int *size; /* malloc'ed */ 47 int *code; /* malloc'ed */ 48 int *val; /* malloc'ed */ 49 int mincode[17]; 50 int maxcode[17]; 51 int valptr[17]; 52 /* fast lookup */ 53 int value[256]; 54 int shift[256]; 55 }; 56 57 58 struct Header 59 { 60 Biobuf *fd; 61 char err[256]; 62 jmp_buf errlab; 63 /* variables in i/o routines */ 64 int sr; /* shift register, right aligned */ 65 int cnt; /* # bits in right part of sr */ 66 uchar *buf; 67 int nbuf; 68 int peek; 69 70 int Nf; 71 72 Framecomp comp[3]; 73 uchar mode; 74 int X; 75 int Y; 76 int qt[4][64]; /* quantization tables */ 77 Huffman dcht[4]; 78 Huffman acht[4]; 79 int **data[3]; 80 int ndata[3]; 81 82 uchar *sf; /* start of frame; do better later */ 83 uchar *ss; /* start of scan; do better later */ 84 int ri; /* restart interval */ 85 86 /* progressive scan */ 87 Rawimage *image; 88 Rawimage **array; 89 int *dccoeff[3]; 90 int **accoeff[3]; /* only need 8 bits plus quantization */ 91 int naccoeff[3]; 92 int nblock[3]; 93 int nacross; 94 int ndown; 95 int Hmax; 96 int Vmax; 97 }; 98 99 static uchar clamp[NCLAMP]; 100 101 static Rawimage *readslave(Header*, int); 102 static int readsegment(Header*, int*); 103 static void quanttables(Header*, uchar*, int); 104 static void huffmantables(Header*, uchar*, int); 105 static void soiheader(Header*); 106 static int nextbyte(Header*, int); 107 static int int2(uchar*, int); 108 static void nibbles(int, int*, int*); 109 static int receive(Header*, int); 110 static int receiveEOB(Header*, int); 111 static int receivebit(Header*); 112 static void restart(Header*, int); 113 static int decode(Header*, Huffman*); 114 static Rawimage* baselinescan(Header*, int); 115 static void progressivescan(Header*, int); 116 static Rawimage* progressiveIDCT(Header*, int); 117 static void idct(int*); 118 static void colormap1(Header*, int, Rawimage*, int*, int, int); 119 static void colormapall1(Header*, int, Rawimage*, int*, int*, int*, int, int); 120 static void colormap(Header*, int, Rawimage*, int**, int**, int**, int, int, int, int, int*, int*); 121 static void jpgerror(Header*, char*, ...); 122 123 static char readerr[] = "ReadJPG: read error: %r"; 124 static char memerr[] = "ReadJPG: malloc failed: %r"; 125 126 static int zig[64] = { 127 0, 1, 8, 16, 9, 2, 3, 10, 17, /* 0-7 */ 128 24, 32, 25, 18, 11, 4, 5, /* 8-15 */ 129 12, 19, 26, 33, 40, 48, 41, 34, /* 16-23 */ 130 27, 20, 13, 6, 7, 14, 21, 28, /* 24-31 */ 131 35, 42, 49, 56, 57, 50, 43, 36, /* 32-39 */ 132 29, 22, 15, 23, 30, 37, 44, 51, /* 40-47 */ 133 58, 59, 52, 45, 38, 31, 39, 46, /* 48-55 */ 134 53, 60, 61, 54, 47, 55, 62, 63 /* 56-63 */ 135 }; 136 137 static 138 void 139 jpginit(void) 140 { 141 int k; 142 static int inited; 143 144 if(inited) 145 return; 146 inited = 1; 147 for(k=0; k<CLAMPOFF; k++) 148 clamp[k] = 0; 149 for(; k<CLAMPOFF+256; k++) 150 clamp[k] = k-CLAMPOFF; 151 for(; k<NCLAMP; k++) 152 clamp[k] = 255; 153 } 154 155 static 156 void* 157 jpgmalloc(Header *h, int n, int clear) 158 { 159 void *p; 160 161 p = malloc(n); 162 if(p == nil) 163 jpgerror(h, memerr); 164 if(clear) 165 memset(p, 0, n); 166 return p; 167 } 168 169 static 170 void 171 clear(void *pp) 172 { 173 void **p = (void**)pp; 174 175 if(*p){ 176 free(*p); 177 *p = nil; 178 } 179 } 180 181 static 182 void 183 jpgfreeall(Header *h, int freeimage) 184 { 185 int i, j; 186 187 clear(&h->buf); 188 if(h->dccoeff[0]) 189 for(i=0; i<3; i++) 190 clear(&h->dccoeff[i]); 191 if(h->accoeff[0]) 192 for(i=0; i<3; i++){ 193 if(h->accoeff[i]) 194 for(j=0; j<h->naccoeff[i]; j++) 195 clear(&h->accoeff[i][j]); 196 clear(&h->accoeff[i]); 197 } 198 for(i=0; i<4; i++){ 199 clear(&h->dcht[i].size); 200 clear(&h->acht[i].size); 201 clear(&h->dcht[i].code); 202 clear(&h->acht[i].code); 203 clear(&h->dcht[i].val); 204 clear(&h->acht[i].val); 205 } 206 if(h->data[0]) 207 for(i=0; i<3; i++){ 208 if(h->data[i]) 209 for(j=0; j<h->ndata[i]; j++) 210 clear(&h->data[i][j]); 211 clear(&h->data[i]); 212 } 213 if(freeimage && h->image!=nil){ 214 clear(&h->array); 215 clear(&h->image->cmap); 216 for(i=0; i<3; i++) 217 clear(&h->image->chans[i]); 218 clear(&h->image); 219 } 220 } 221 222 static 223 void 224 jpgerror(Header *h, char *fmt, ...) 225 { 226 va_list arg; 227 228 va_start(arg, fmt); 229 vseprint(h->err, h->err+sizeof h->err, fmt, arg); 230 va_end(arg); 231 232 werrstr(h->err); 233 jpgfreeall(h, 1); 234 longjmp(h->errlab, 1); 235 } 236 237 Rawimage** 238 Breadjpg(Biobuf *b, int colorspace) 239 { 240 Rawimage *r, **array; 241 Header *h; 242 char buf[ERRMAX]; 243 244 buf[0] = '\0'; 245 if(colorspace!=CYCbCr && colorspace!=CRGB){ 246 errstr(buf, sizeof buf); /* throw it away */ 247 werrstr("ReadJPG: unknown color space"); 248 return nil; 249 } 250 jpginit(); 251 h = malloc(sizeof(Header)); 252 array = malloc(sizeof(Header)); 253 if(h==nil || array==nil){ 254 free(h); 255 free(array); 256 return nil; 257 } 258 h->array = array; 259 memset(h, 0, sizeof(Header)); 260 h->fd = b; 261 errstr(buf, sizeof buf); /* throw it away */ 262 if(setjmp(h->errlab)) 263 r = nil; 264 else 265 r = readslave(h, colorspace); 266 jpgfreeall(h, 0); 267 free(h); 268 array[0] = r; 269 array[1] = nil; 270 return array; 271 } 272 273 Rawimage** 274 readjpg(int fd, int colorspace) 275 { 276 Rawimage** a; 277 Biobuf b; 278 279 if(Binit(&b, fd, OREAD) < 0) 280 return nil; 281 a = Breadjpg(&b, colorspace); 282 Bterm(&b); 283 return a; 284 } 285 286 static 287 Rawimage* 288 readslave(Header *header, int colorspace) 289 { 290 Rawimage *image; 291 int nseg, i, H, V, m, n; 292 uchar *b; 293 294 soiheader(header); 295 nseg = 0; 296 image = nil; 297 298 header->buf = jpgmalloc(header, 4096, 0); 299 header->nbuf = 4096; 300 while(header->err[0] == '\0'){ 301 nseg++; 302 n = readsegment(header, &m); 303 b = header->buf; 304 switch(m){ 305 case -1: 306 return image; 307 308 case APPn+0: 309 if(nseg==1 && strncmp((char*)b, "JFIF", 4)==0) /* JFIF header; check version */ 310 if(b[5]>1 || b[6]>2) 311 sprint(header->err, "ReadJPG: can't handle JFIF version %d.%2d", b[5], b[6]); 312 break; 313 314 case APPn+1: case APPn+2: case APPn+3: case APPn+4: case APPn+5: 315 case APPn+6: case APPn+7: case APPn+8: case APPn+9: case APPn+10: 316 case APPn+11: case APPn+12: case APPn+13: case APPn+14: case APPn+15: 317 break; 318 319 case DQT: 320 quanttables(header, b, n); 321 break; 322 323 case SOF: 324 case SOF2: 325 header->Y = int2(b, 1); 326 header->X = int2(b, 3); 327 header->Nf =b[5]; 328 for(i=0; i<header->Nf; i++){ 329 header->comp[i].C = b[6+3*i+0]; 330 nibbles(b[6+3*i+1], &H, &V); 331 if(H<=0 || V<=0) 332 jpgerror(header, "non-positive sampling factor (Hsamp or Vsamp)"); 333 header->comp[i].H = H; 334 header->comp[i].V = V; 335 header->comp[i].Tq = b[6+3*i+2]; 336 } 337 header->mode = m; 338 header->sf = b; 339 break; 340 341 case SOS: 342 header->ss = b; 343 switch(header->mode){ 344 case SOF: 345 image = baselinescan(header, colorspace); 346 break; 347 case SOF2: 348 progressivescan(header, colorspace); 349 break; 350 default: 351 sprint(header->err, "unrecognized or unspecified encoding %d", header->mode); 352 break; 353 } 354 break; 355 356 case DHT: 357 huffmantables(header, b, n); 358 break; 359 360 case DRI: 361 header->ri = int2(b, 0); 362 break; 363 364 case COM: 365 break; 366 367 case EOI: 368 if(header->mode == SOF2) 369 image = progressiveIDCT(header, colorspace); 370 return image; 371 372 default: 373 sprint(header->err, "ReadJPG: unknown marker %.2x", m); 374 break; 375 } 376 } 377 return image; 378 } 379 380 /* readsegment is called after reading scan, which can have */ 381 /* read ahead a byte. so we must check peek here */ 382 static 383 int 384 readbyte(Header *h) 385 { 386 uchar x; 387 388 if(h->peek >= 0){ 389 x = h->peek; 390 h->peek = -1; 391 }else if(Bread(h->fd, &x, 1) != 1) 392 jpgerror(h, readerr); 393 return x; 394 } 395 396 static 397 int 398 marker(Header *h) 399 { 400 int c; 401 402 while((c=readbyte(h)) == 0) 403 fprint(2, "ReadJPG: skipping zero byte at offset %lld\n", Boffset(h->fd)); 404 if(c != 0xFF) 405 jpgerror(h, "ReadJPG: expecting marker; found 0x%x at offset %lld\n", c, Boffset(h->fd)); 406 while(c == 0xFF) 407 c = readbyte(h); 408 return c; 409 } 410 411 static 412 int 413 int2(uchar *buf, int n) 414 { 415 return (buf[n]<<8) + buf[n+1]; 416 } 417 418 static 419 void 420 nibbles(int b, int *p0, int *p1) 421 { 422 *p0 = (b>>4) & 0xF; 423 *p1 = b & 0xF; 424 } 425 426 static 427 void 428 soiheader(Header *h) 429 { 430 h->peek = -1; 431 if(marker(h) != SOI) 432 jpgerror(h, "ReadJPG: unrecognized marker in header"); 433 h->err[0] = '\0'; 434 h->mode = 0; 435 h->ri = 0; 436 } 437 438 static 439 int 440 readsegment(Header *h, int *markerp) 441 { 442 int m, n; 443 uchar tmp[2]; 444 445 m = marker(h); 446 switch(m){ 447 case EOI: 448 *markerp = m; 449 return 0; 450 case 0: 451 jpgerror(h, "ReadJPG: expecting marker; saw %.2x at offset %lld", m, Boffset(h->fd)); 452 } 453 if(Bread(h->fd, tmp, 2) != 2) 454 Readerr: 455 jpgerror(h, readerr); 456 n = int2(tmp, 0); 457 if(n < 2) 458 goto Readerr; 459 n -= 2; 460 if(n > h->nbuf){ 461 free(h->buf); 462 h->buf = jpgmalloc(h, n+1, 0); /* +1 for sentinel */ 463 h->nbuf = n; 464 } 465 if(Bread(h->fd, h->buf, n) != n) 466 goto Readerr; 467 *markerp = m; 468 return n; 469 } 470 471 static 472 int 473 huffmantable(Header *h, uchar *b) 474 { 475 Huffman *t; 476 int Tc, th, n, nsize, i, j, k, v, cnt, code, si, sr, m; 477 int *maxcode; 478 479 nibbles(b[0], &Tc, &th); 480 if(Tc > 1) 481 jpgerror(h, "ReadJPG: unknown Huffman table class %d", Tc); 482 if(th>3 || (h->mode==SOF && th>1)) 483 jpgerror(h, "ReadJPG: unknown Huffman table index %d", th); 484 if(Tc == 0) 485 t = &h->dcht[th]; 486 else 487 t = &h->acht[th]; 488 489 /* flow chart C-2 */ 490 nsize = 0; 491 for(i=0; i<16; i++) 492 nsize += b[1+i]; 493 t->size = jpgmalloc(h, (nsize+1)*sizeof(int), 1); 494 k = 0; 495 for(i=1; i<=16; i++){ 496 n = b[i]; 497 for(j=0; j<n; j++) 498 t->size[k++] = i; 499 } 500 t->size[k] = 0; 501 502 /* initialize HUFFVAL */ 503 t->val = jpgmalloc(h, nsize*sizeof(int), 1); 504 for(i=0; i<nsize; i++) 505 t->val[i] = b[17+i]; 506 507 /* flow chart C-3 */ 508 t->code = jpgmalloc(h, (nsize+1)*sizeof(int), 1); 509 k = 0; 510 code = 0; 511 si = t->size[0]; 512 for(;;){ 513 do 514 t->code[k++] = code++; 515 while(t->size[k] == si); 516 if(t->size[k] == 0) 517 break; 518 do{ 519 code <<= 1; 520 si++; 521 }while(t->size[k] != si); 522 } 523 524 /* flow chart F-25 */ 525 i = 0; 526 j = 0; 527 for(;;){ 528 for(;;){ 529 i++; 530 if(i > 16) 531 goto outF25; 532 if(b[i] != 0) 533 break; 534 t->maxcode[i] = -1; 535 } 536 t->valptr[i] = j; 537 t->mincode[i] = t->code[j]; 538 j += b[i]-1; 539 t->maxcode[i] = t->code[j]; 540 j++; 541 } 542 outF25: 543 544 /* create byte-indexed fast path tables */ 545 maxcode = t->maxcode; 546 /* stupid startup algorithm: just run machine for each byte value */ 547 for(v=0; v<256; ){ 548 cnt = 7; 549 m = 1<<7; 550 code = 0; 551 sr = v; 552 i = 1; 553 for(;;i++){ 554 if(sr & m) 555 code |= 1; 556 if(code <= maxcode[i]) 557 break; 558 code <<= 1; 559 m >>= 1; 560 if(m == 0){ 561 t->shift[v] = 0; 562 t->value[v] = -1; 563 goto continueBytes; 564 } 565 cnt--; 566 } 567 t->shift[v] = 8-cnt; 568 t->value[v] = t->val[t->valptr[i]+(code-t->mincode[i])]; 569 570 continueBytes: 571 v++; 572 } 573 574 return nsize; 575 } 576 577 static 578 void 579 huffmantables(Header *h, uchar *b, int n) 580 { 581 int l, mt; 582 583 for(l=0; l<n; l+=17+mt) 584 mt = huffmantable(h, &b[l]); 585 } 586 587 static 588 int 589 quanttable(Header *h, uchar *b) 590 { 591 int i, pq, tq, *q; 592 593 nibbles(b[0], &pq, &tq); 594 if(pq > 1) 595 jpgerror(h, "ReadJPG: unknown quantization table class %d", pq); 596 if(tq > 3) 597 jpgerror(h, "ReadJPG: unknown quantization table index %d", tq); 598 q = h->qt[tq]; 599 for(i=0; i<64; i++){ 600 if(pq == 0) 601 q[i] = b[1+i]; 602 else 603 q[i] = int2(b, 1+2*i); 604 } 605 return 64*(1+pq); 606 } 607 608 static 609 void 610 quanttables(Header *h, uchar *b, int n) 611 { 612 int l, m; 613 614 for(l=0; l<n; l+=1+m) 615 m = quanttable(h, &b[l]); 616 } 617 618 static 619 Rawimage* 620 baselinescan(Header *h, int colorspace) 621 { 622 int Ns, z, k, m, Hmax, Vmax, comp; 623 int allHV1, nblock, ri, mcu, nacross, nmcu; 624 Huffman *dcht, *acht; 625 int block, t, diff, *qt; 626 uchar *ss; 627 Rawimage *image; 628 int Td[3], Ta[3], H[3], V[3], DC[3]; 629 int ***data, *zz; 630 631 ss = h->ss; 632 Ns = ss[0]; 633 if((Ns!=3 && Ns!=1) || Ns!=h->Nf) 634 jpgerror(h, "ReadJPG: can't handle scan not 3 components"); 635 636 image = jpgmalloc(h, sizeof(Rawimage), 1); 637 h->image = image; 638 image->r = Rect(0, 0, h->X, h->Y); 639 image->cmap = nil; 640 image->cmaplen = 0; 641 image->chanlen = h->X*h->Y; 642 image->fields = 0; 643 image->gifflags = 0; 644 image->gifdelay = 0; 645 image->giftrindex = 0; 646 if(Ns == 3) 647 image->chandesc = colorspace; 648 else 649 image->chandesc = CY; 650 image->nchans = h->Nf; 651 for(k=0; k<h->Nf; k++) 652 image->chans[k] = jpgmalloc(h, h->X*h->Y, 0); 653 654 /* compute maximum H and V */ 655 Hmax = 0; 656 Vmax = 0; 657 for(comp=0; comp<Ns; comp++){ 658 if(h->comp[comp].H > Hmax) 659 Hmax = h->comp[comp].H; 660 if(h->comp[comp].V > Vmax) 661 Vmax = h->comp[comp].V; 662 } 663 664 /* initialize data structures */ 665 allHV1 = 1; 666 data = h->data; 667 for(comp=0; comp<Ns; comp++){ 668 /* JPEG requires scan components to be in same order as in frame, */ 669 /* so if both have 3 we know scan is Y Cb Cr and there's no need to */ 670 /* reorder */ 671 nibbles(ss[2+2*comp], &Td[comp], &Ta[comp]); 672 H[comp] = h->comp[comp].H; 673 V[comp] = h->comp[comp].V; 674 nblock = H[comp]*V[comp]; 675 if(nblock != 1) 676 allHV1 = 0; 677 data[comp] = jpgmalloc(h, nblock*sizeof(int*), 0); 678 h->ndata[comp] = nblock; 679 DC[comp] = 0; 680 for(m=0; m<nblock; m++) 681 data[comp][m] = jpgmalloc(h, 8*8*sizeof(int), 0); 682 } 683 684 ri = h->ri; 685 686 h->cnt = 0; 687 h->sr = 0; 688 h->peek = -1; 689 nacross = ((h->X+(8*Hmax-1))/(8*Hmax)); 690 nmcu = ((h->Y+(8*Vmax-1))/(8*Vmax))*nacross; 691 for(mcu=0; mcu<nmcu; ){ 692 for(comp=0; comp<Ns; comp++){ 693 dcht = &h->dcht[Td[comp]]; 694 acht = &h->acht[Ta[comp]]; 695 qt = h->qt[h->comp[comp].Tq]; 696 697 for(block=0; block<H[comp]*V[comp]; block++){ 698 /* F-22 */ 699 t = decode(h, dcht); 700 diff = receive(h, t); 701 DC[comp] += diff; 702 703 /* F-23 */ 704 zz = data[comp][block]; 705 memset(zz, 0, 8*8*sizeof(int)); 706 zz[0] = qt[0]*DC[comp]; 707 k = 1; 708 709 for(;;){ 710 t = decode(h, acht); 711 if((t&0x0F) == 0){ 712 if((t&0xF0) != 0xF0) 713 break; 714 k += 16; 715 }else{ 716 k += t>>4; 717 z = receive(h, t&0xF); 718 zz[zig[k]] = z*qt[k]; 719 if(k == 63) 720 break; 721 k++; 722 } 723 } 724 725 idct(zz); 726 } 727 } 728 729 /* rotate colors to RGB and assign to bytes */ 730 if(Ns == 1) /* very easy */ 731 colormap1(h, colorspace, image, data[0][0], mcu, nacross); 732 else if(allHV1) /* fairly easy */ 733 colormapall1(h, colorspace, image, data[0][0], data[1][0], data[2][0], mcu, nacross); 734 else /* miserable general case */ 735 colormap(h, colorspace, image, data[0], data[1], data[2], mcu, nacross, Hmax, Vmax, H, V); 736 /* process restart marker, if present */ 737 mcu++; 738 if(ri>0 && mcu<nmcu && mcu%ri==0){ 739 restart(h, mcu); 740 for(comp=0; comp<Ns; comp++) 741 DC[comp] = 0; 742 } 743 } 744 return image; 745 } 746 747 static 748 void 749 restart(Header *h, int mcu) 750 { 751 int rest, rst, nskip; 752 753 rest = mcu/h->ri-1; 754 nskip = 0; 755 do{ 756 do{ 757 rst = nextbyte(h, 1); 758 nskip++; 759 }while(rst>=0 && rst!=0xFF); 760 if(rst == 0xFF){ 761 rst = nextbyte(h, 1); 762 nskip++; 763 } 764 }while(rst>=0 && (rst&~7)!=RST); 765 if(nskip != 2) 766 sprint(h->err, "ReadJPG: skipped %d bytes at restart %d\n", nskip-2, rest); 767 if(rst < 0) 768 jpgerror(h, readerr); 769 if((rst&7) != (rest&7)) 770 jpgerror(h, "ReadJPG: expected RST%d got %d", rest&7, rst&7); 771 h->cnt = 0; 772 h->sr = 0; 773 } 774 775 static 776 Rawimage* 777 progressiveIDCT(Header *h, int colorspace) 778 { 779 int k, m, comp, block, Nf, bn; 780 int allHV1, nblock, mcu, nmcu; 781 int H[3], V[3], blockno[3]; 782 int *dccoeff, **accoeff; 783 int ***data, *zz; 784 785 Nf = h->Nf; 786 allHV1 = 1; 787 data = h->data; 788 789 for(comp=0; comp<Nf; comp++){ 790 H[comp] = h->comp[comp].H; 791 V[comp] = h->comp[comp].V; 792 nblock = h->nblock[comp]; 793 if(nblock != 1) 794 allHV1 = 0; 795 h->ndata[comp] = nblock; 796 data[comp] = jpgmalloc(h, nblock*sizeof(int*), 0); 797 for(m=0; m<nblock; m++) 798 data[comp][m] = jpgmalloc(h, 8*8*sizeof(int), 0); 799 } 800 801 memset(blockno, 0, sizeof blockno); 802 nmcu = h->nacross*h->ndown; 803 for(mcu=0; mcu<nmcu; mcu++){ 804 for(comp=0; comp<Nf; comp++){ 805 dccoeff = h->dccoeff[comp]; 806 accoeff = h->accoeff[comp]; 807 bn = blockno[comp]; 808 for(block=0; block<h->nblock[comp]; block++){ 809 zz = data[comp][block]; 810 memset(zz, 0, 8*8*sizeof(int)); 811 zz[0] = dccoeff[bn]; 812 813 for(k=1; k<64; k++) 814 zz[zig[k]] = accoeff[bn][k]; 815 816 idct(zz); 817 bn++; 818 } 819 blockno[comp] = bn; 820 } 821 822 /* rotate colors to RGB and assign to bytes */ 823 if(Nf == 1) /* very easy */ 824 colormap1(h, colorspace, h->image, data[0][0], mcu, h->nacross); 825 else if(allHV1) /* fairly easy */ 826 colormapall1(h, colorspace, h->image, data[0][0], data[1][0], data[2][0], mcu, h->nacross); 827 else /* miserable general case */ 828 colormap(h, colorspace, h->image, data[0], data[1], data[2], mcu, h->nacross, h->Hmax, h->Vmax, H, V); 829 } 830 831 return h->image; 832 } 833 834 static 835 void 836 progressiveinit(Header *h, int colorspace) 837 { 838 int Nf, Ns, j, k, nmcu, comp; 839 uchar *ss; 840 Rawimage *image; 841 842 ss = h->ss; 843 Ns = ss[0]; 844 Nf = h->Nf; 845 if((Ns!=3 && Ns!=1) || Ns!=Nf) 846 jpgerror(h, "ReadJPG: image must have 1 or 3 components"); 847 848 image = jpgmalloc(h, sizeof(Rawimage), 1); 849 h->image = image; 850 image->r = Rect(0, 0, h->X, h->Y); 851 image->cmap = nil; 852 image->cmaplen = 0; 853 image->chanlen = h->X*h->Y; 854 image->fields = 0; 855 image->gifflags = 0; 856 image->gifdelay = 0; 857 image->giftrindex = 0; 858 if(Nf == 3) 859 image->chandesc = colorspace; 860 else 861 image->chandesc = CY; 862 image->nchans = h->Nf; 863 for(k=0; k<Nf; k++){ 864 image->chans[k] = jpgmalloc(h, h->X*h->Y, 0); 865 h->nblock[k] = h->comp[k].H*h->comp[k].V; 866 } 867 868 /* compute maximum H and V */ 869 h->Hmax = 0; 870 h->Vmax = 0; 871 for(comp=0; comp<Nf; comp++){ 872 if(h->comp[comp].H > h->Hmax) 873 h->Hmax = h->comp[comp].H; 874 if(h->comp[comp].V > h->Vmax) 875 h->Vmax = h->comp[comp].V; 876 } 877 h->nacross = ((h->X+(8*h->Hmax-1))/(8*h->Hmax)); 878 h->ndown = ((h->Y+(8*h->Vmax-1))/(8*h->Vmax)); 879 nmcu = h->nacross*h->ndown; 880 881 for(k=0; k<Nf; k++){ 882 h->dccoeff[k] = jpgmalloc(h, h->nblock[k]*nmcu * sizeof(int), 1); 883 h->accoeff[k] = jpgmalloc(h, h->nblock[k]*nmcu * sizeof(int*), 1); 884 h->naccoeff[k] = h->nblock[k]*nmcu; 885 for(j=0; j<h->nblock[k]*nmcu; j++) 886 h->accoeff[k][j] = jpgmalloc(h, 64*sizeof(int), 1); 887 } 888 889 } 890 891 static 892 void 893 progressivedc(Header *h, int comp, int Ah, int Al) 894 { 895 int Ns, z, ri, mcu, nmcu; 896 int block, t, diff, qt, *dc, bn; 897 Huffman *dcht; 898 uchar *ss; 899 int Td[3], DC[3], blockno[3]; 900 901 ss= h->ss; 902 Ns = ss[0]; 903 if(Ns!=h->Nf) 904 jpgerror(h, "ReadJPG: can't handle progressive with Nf!=Ns in DC scan"); 905 906 /* initialize data structures */ 907 h->cnt = 0; 908 h->sr = 0; 909 h->peek = -1; 910 for(comp=0; comp<Ns; comp++){ 911 /* 912 * JPEG requires scan components to be in same order as in frame, 913 * so if both have 3 we know scan is Y Cb Cr and there's no need to 914 * reorder 915 */ 916 nibbles(ss[2+2*comp], &Td[comp], &z); /* z is ignored */ 917 DC[comp] = 0; 918 } 919 920 ri = h->ri; 921 922 nmcu = h->nacross*h->ndown; 923 memset(blockno, 0, sizeof blockno); 924 for(mcu=0; mcu<nmcu; ){ 925 for(comp=0; comp<Ns; comp++){ 926 dcht = &h->dcht[Td[comp]]; 927 qt = h->qt[h->comp[comp].Tq][0]; 928 dc = h->dccoeff[comp]; 929 bn = blockno[comp]; 930 931 for(block=0; block<h->nblock[comp]; block++){ 932 if(Ah == 0){ 933 t = decode(h, dcht); 934 diff = receive(h, t); 935 DC[comp] += diff; 936 dc[bn] = qt*DC[comp]<<Al; 937 }else 938 dc[bn] |= qt*receivebit(h)<<Al; 939 bn++; 940 } 941 blockno[comp] = bn; 942 } 943 944 /* process restart marker, if present */ 945 mcu++; 946 if(ri>0 && mcu<nmcu && mcu%ri==0){ 947 restart(h, mcu); 948 for(comp=0; comp<Ns; comp++) 949 DC[comp] = 0; 950 } 951 } 952 } 953 954 static 955 void 956 progressiveac(Header *h, int comp, int Al) 957 { 958 int Ns, Ss, Se, z, k, eobrun, x, y, nver, tmcu, blockno, *acc, rs; 959 int ri, mcu, nacross, ndown, nmcu, nhor; 960 Huffman *acht; 961 int *qt, rrrr, ssss, q; 962 uchar *ss; 963 int Ta, H, V; 964 965 ss = h->ss; 966 Ns = ss[0]; 967 if(Ns != 1) 968 jpgerror(h, "ReadJPG: illegal Ns>1 in progressive AC scan"); 969 Ss = ss[1+2]; 970 Se = ss[2+2]; 971 H = h->comp[comp].H; 972 V = h->comp[comp].V; 973 974 nacross = h->nacross*H; 975 ndown = h->ndown*V; 976 q = 8*h->Hmax/H; 977 nhor = (h->X+q-1)/q; 978 q = 8*h->Vmax/V; 979 nver = (h->Y+q-1)/q; 980 981 /* initialize data structures */ 982 h->cnt = 0; 983 h->sr = 0; 984 h->peek = -1; 985 nibbles(ss[1+1], &z, &Ta); /* z is thrown away */ 986 987 ri = h->ri; 988 989 eobrun = 0; 990 acht = &h->acht[Ta]; 991 qt = h->qt[h->comp[comp].Tq]; 992 nmcu = nacross*ndown; 993 mcu = 0; 994 for(y=0; y<nver; y++){ 995 for(x=0; x<nhor; x++){ 996 /* Figure G-3 */ 997 if(eobrun > 0){ 998 --eobrun; 999 continue; 1000 } 1001 1002 /* arrange blockno to be in same sequence as original scan calculation. */ 1003 tmcu = x/H + (nacross/H)*(y/V); 1004 blockno = tmcu*H*V + H*(y%V) + x%H; 1005 acc = h->accoeff[comp][blockno]; 1006 k = Ss; 1007 for(;;){ 1008 rs = decode(h, acht); 1009 /* XXX remove rrrr ssss as in baselinescan */ 1010 nibbles(rs, &rrrr, &ssss); 1011 if(ssss == 0){ 1012 if(rrrr < 15){ 1013 eobrun = 0; 1014 if(rrrr > 0) 1015 eobrun = receiveEOB(h, rrrr)-1; 1016 break; 1017 } 1018 k += 16; 1019 }else{ 1020 k += rrrr; 1021 z = receive(h, ssss); 1022 acc[k] = z*qt[k]<<Al; 1023 if(k == Se) 1024 break; 1025 k++; 1026 } 1027 } 1028 } 1029 1030 /* process restart marker, if present */ 1031 mcu++; 1032 if(ri>0 && mcu<nmcu && mcu%ri==0){ 1033 restart(h, mcu); 1034 eobrun = 0; 1035 } 1036 } 1037 } 1038 1039 static 1040 void 1041 increment(Header *h, int acc[], int k, int Pt) 1042 { 1043 if(acc[k] == 0) 1044 return; 1045 if(receivebit(h) != 0) 1046 if(acc[k] < 0) 1047 acc[k] -= Pt; 1048 else 1049 acc[k] += Pt; 1050 } 1051 1052 static 1053 void 1054 progressiveacinc(Header *h, int comp, int Al) 1055 { 1056 int Ns, i, z, k, Ss, Se, Ta, **ac, H, V; 1057 int ri, mcu, nacross, ndown, nhor, nver, eobrun, nzeros, pending, x, y, tmcu, blockno, q, nmcu; 1058 Huffman *acht; 1059 int *qt, rrrr, ssss, *acc, rs; 1060 uchar *ss; 1061 1062 ss = h->ss; 1063 Ns = ss[0]; 1064 if(Ns != 1) 1065 jpgerror(h, "ReadJPG: illegal Ns>1 in progressive AC scan"); 1066 Ss = ss[1+2]; 1067 Se = ss[2+2]; 1068 H = h->comp[comp].H; 1069 V = h->comp[comp].V; 1070 1071 nacross = h->nacross*H; 1072 ndown = h->ndown*V; 1073 q = 8*h->Hmax/H; 1074 nhor = (h->X+q-1)/q; 1075 q = 8*h->Vmax/V; 1076 nver = (h->Y+q-1)/q; 1077 1078 /* initialize data structures */ 1079 h->cnt = 0; 1080 h->sr = 0; 1081 h->peek = -1; 1082 nibbles(ss[1+1], &z, &Ta); /* z is thrown away */ 1083 ri = h->ri; 1084 1085 eobrun = 0; 1086 ac = h->accoeff[comp]; 1087 acht = &h->acht[Ta]; 1088 qt = h->qt[h->comp[comp].Tq]; 1089 nmcu = nacross*ndown; 1090 mcu = 0; 1091 pending = 0; 1092 nzeros = -1; 1093 for(y=0; y<nver; y++){ 1094 for(x=0; x<nhor; x++){ 1095 /* Figure G-7 */ 1096 1097 /* arrange blockno to be in same sequence as original scan calculation. */ 1098 tmcu = x/H + (nacross/H)*(y/V); 1099 blockno = tmcu*H*V + H*(y%V) + x%H; 1100 acc = ac[blockno]; 1101 if(eobrun > 0){ 1102 if(nzeros > 0) 1103 jpgerror(h, "ReadJPG: zeros pending at block start"); 1104 for(k=Ss; k<=Se; k++) 1105 increment(h, acc, k, qt[k]<<Al); 1106 --eobrun; 1107 continue; 1108 } 1109 1110 for(k=Ss; k<=Se; ){ 1111 if(nzeros >= 0){ 1112 if(acc[k] != 0) 1113 increment(h, acc, k, qt[k]<<Al); 1114 else if(nzeros-- == 0) 1115 acc[k] = pending; 1116 k++; 1117 continue; 1118 } 1119 rs = decode(h, acht); 1120 nibbles(rs, &rrrr, &ssss); 1121 if(ssss == 0){ 1122 if(rrrr < 15){ 1123 eobrun = 0; 1124 if(rrrr > 0) 1125 eobrun = receiveEOB(h, rrrr)-1; 1126 while(k <= Se){ 1127 increment(h, acc, k, qt[k]<<Al); 1128 k++; 1129 } 1130 break; 1131 } 1132 for(i=0; i<16; k++){ 1133 increment(h, acc, k, qt[k]<<Al); 1134 if(acc[k] == 0) 1135 i++; 1136 } 1137 continue; 1138 }else if(ssss != 1) 1139 jpgerror(h, "ReadJPG: ssss!=1 in progressive increment"); 1140 nzeros = rrrr; 1141 pending = receivebit(h); 1142 if(pending == 0) 1143 pending = -1; 1144 pending *= qt[k]<<Al; 1145 } 1146 } 1147 1148 /* process restart marker, if present */ 1149 mcu++; 1150 if(ri>0 && mcu<nmcu && mcu%ri==0){ 1151 restart(h, mcu); 1152 eobrun = 0; 1153 nzeros = -1; 1154 } 1155 } 1156 } 1157 1158 static 1159 void 1160 progressivescan(Header *h, int colorspace) 1161 { 1162 uchar *ss; 1163 int Ns, Ss, Ah, Al, c, comp, i; 1164 1165 if(h->dccoeff[0] == nil) 1166 progressiveinit(h, colorspace); 1167 1168 ss = h->ss; 1169 Ns = ss[0]; 1170 Ss = ss[1+2*Ns]; 1171 nibbles(ss[3+2*Ns], &Ah, &Al); 1172 c = ss[1]; 1173 comp = -1; 1174 for(i=0; i<h->Nf; i++) 1175 if(h->comp[i].C == c) 1176 comp = i; 1177 if(comp == -1) 1178 jpgerror(h, "ReadJPG: bad component index in scan header"); 1179 1180 if(Ss == 0){ 1181 progressivedc(h, comp, Ah, Al); 1182 return; 1183 } 1184 if(Ah == 0){ 1185 progressiveac(h, comp, Al); 1186 return; 1187 } 1188 progressiveacinc(h, comp, Al); 1189 } 1190 1191 enum { 1192 c1 = 2871, /* 1.402 * 2048 */ 1193 c2 = 705, /* 0.34414 * 2048 */ 1194 c3 = 1463, /* 0.71414 * 2048 */ 1195 c4 = 3629 /* 1.772 * 2048 */ 1196 }; 1197 1198 static 1199 void 1200 colormap1(Header *h, int colorspace, Rawimage *image, int data[8*8], int mcu, int nacross) 1201 { 1202 uchar *pic; 1203 int x, y, dx, dy, minx, miny; 1204 int r, k, pici; 1205 1206 USED(colorspace); 1207 pic = image->chans[0]; 1208 minx = 8*(mcu%nacross); 1209 dx = 8; 1210 if(minx+dx > h->X) 1211 dx = h->X-minx; 1212 miny = 8*(mcu/nacross); 1213 dy = 8; 1214 if(miny+dy > h->Y) 1215 dy = h->Y-miny; 1216 pici = miny*h->X+minx; 1217 k = 0; 1218 for(y=0; y<dy; y++){ 1219 for(x=0; x<dx; x++){ 1220 r = clamp[(data[k+x]+128)+CLAMPOFF]; 1221 pic[pici+x] = r; 1222 } 1223 pici += h->X; 1224 k += 8; 1225 } 1226 } 1227 1228 static 1229 void 1230 colormapall1(Header *h, int colorspace, Rawimage *image, int data0[8*8], int data1[8*8], int data2[8*8], int mcu, int nacross) 1231 { 1232 uchar *rpic, *gpic, *bpic, *rp, *gp, *bp; 1233 int *p0, *p1, *p2; 1234 int x, y, dx, dy, minx, miny; 1235 int r, g, b, k, pici; 1236 int Y, Cr, Cb; 1237 1238 rpic = image->chans[0]; 1239 gpic = image->chans[1]; 1240 bpic = image->chans[2]; 1241 minx = 8*(mcu%nacross); 1242 dx = 8; 1243 if(minx+dx > h->X) 1244 dx = h->X-minx; 1245 miny = 8*(mcu/nacross); 1246 dy = 8; 1247 if(miny+dy > h->Y) 1248 dy = h->Y-miny; 1249 pici = miny*h->X+minx; 1250 k = 0; 1251 for(y=0; y<dy; y++){ 1252 p0 = data0+k; 1253 p1 = data1+k; 1254 p2 = data2+k; 1255 rp = rpic+pici; 1256 gp = gpic+pici; 1257 bp = bpic+pici; 1258 if(colorspace == CYCbCr) 1259 for(x=0; x<dx; x++){ 1260 *rp++ = clamp[*p0++ + 128 + CLAMPOFF]; 1261 *gp++ = clamp[*p1++ + 128 + CLAMPOFF]; 1262 *bp++ = clamp[*p2++ + 128 + CLAMPOFF]; 1263 } 1264 else 1265 for(x=0; x<dx; x++){ 1266 Y = (*p0++ + 128) << 11; 1267 Cb = *p1++; 1268 Cr = *p2++; 1269 r = Y+c1*Cr; 1270 g = Y-c2*Cb-c3*Cr; 1271 b = Y+c4*Cb; 1272 *rp++ = clamp[(r>>11)+CLAMPOFF]; 1273 *gp++ = clamp[(g>>11)+CLAMPOFF]; 1274 *bp++ = clamp[(b>>11)+CLAMPOFF]; 1275 } 1276 pici += h->X; 1277 k += 8; 1278 } 1279 } 1280 1281 static 1282 void 1283 colormap(Header *h, int colorspace, Rawimage *image, int *data0[8*8], int *data1[8*8], int *data2[8*8], int mcu, int nacross, int Hmax, int Vmax, int *H, int *V) 1284 { 1285 uchar *rpic, *gpic, *bpic; 1286 int x, y, dx, dy, minx, miny; 1287 int r, g, b, pici, H0, H1, H2; 1288 int t, b0, b1, b2, y0, y1, y2, x0, x1, x2; 1289 int Y, Cr, Cb; 1290 1291 rpic = image->chans[0]; 1292 gpic = image->chans[1]; 1293 bpic = image->chans[2]; 1294 minx = 8*Hmax*(mcu%nacross); 1295 dx = 8*Hmax; 1296 if(minx+dx > h->X) 1297 dx = h->X-minx; 1298 miny = 8*Vmax*(mcu/nacross); 1299 dy = 8*Vmax; 1300 if(miny+dy > h->Y) 1301 dy = h->Y-miny; 1302 pici = miny*h->X+minx; 1303 H0 = H[0]; 1304 H1 = H[1]; 1305 H2 = H[2]; 1306 for(y=0; y<dy; y++){ 1307 t = y*V[0]; 1308 b0 = H0*(t/(8*Vmax)); 1309 y0 = 8*((t/Vmax)&7); 1310 t = y*V[1]; 1311 b1 = H1*(t/(8*Vmax)); 1312 y1 = 8*((t/Vmax)&7); 1313 t = y*V[2]; 1314 b2 = H2*(t/(8*Vmax)); 1315 y2 = 8*((t/Vmax)&7); 1316 x0 = 0; 1317 x1 = 0; 1318 x2 = 0; 1319 for(x=0; x<dx; x++){ 1320 if(colorspace == CYCbCr){ 1321 rpic[pici+x] = clamp[data0[b0][y0+x0++*H0/Hmax] + 128 + CLAMPOFF]; 1322 gpic[pici+x] = clamp[data1[b1][y1+x1++*H1/Hmax] + 128 + CLAMPOFF]; 1323 bpic[pici+x] = clamp[data2[b2][y2+x2++*H2/Hmax] + 128 + CLAMPOFF]; 1324 }else{ 1325 Y = (data0[b0][y0+x0++*H0/Hmax]+128)<<11; 1326 Cb = data1[b1][y1+x1++*H1/Hmax]; 1327 Cr = data2[b2][y2+x2++*H2/Hmax]; 1328 r = Y+c1*Cr; 1329 g = Y-c2*Cb-c3*Cr; 1330 b = Y+c4*Cb; 1331 rpic[pici+x] = clamp[(r>>11)+CLAMPOFF]; 1332 gpic[pici+x] = clamp[(g>>11)+CLAMPOFF]; 1333 bpic[pici+x] = clamp[(b>>11)+CLAMPOFF]; 1334 } 1335 if(x0*H0/Hmax >= 8){ 1336 x0 = 0; 1337 b0++; 1338 } 1339 if(x1*H1/Hmax >= 8){ 1340 x1 = 0; 1341 b1++; 1342 } 1343 if(x2*H2/Hmax >= 8){ 1344 x2 = 0; 1345 b2++; 1346 } 1347 } 1348 pici += h->X; 1349 } 1350 } 1351 1352 /* 1353 * decode next 8-bit value from entropy-coded input. chart F-26 1354 */ 1355 static 1356 int 1357 decode(Header *h, Huffman *t) 1358 { 1359 int code, v, cnt, m, sr, i; 1360 int *maxcode; 1361 static int badcode; 1362 1363 maxcode = t->maxcode; 1364 if(h->cnt < 8) 1365 nextbyte(h, 0); 1366 /* fast lookup */ 1367 code = (h->sr>>(h->cnt-8))&0xFF; 1368 v = t->value[code]; 1369 if(v >= 0){ 1370 h->cnt -= t->shift[code]; 1371 return v; 1372 } 1373 1374 h->cnt -= 8; 1375 if(h->cnt == 0) 1376 nextbyte(h, 0); 1377 h->cnt--; 1378 cnt = h->cnt; 1379 m = 1<<cnt; 1380 sr = h->sr; 1381 code <<= 1; 1382 i = 9; 1383 for(;;i++){ 1384 if(sr & m) 1385 code |= 1; 1386 if(code <= maxcode[i]) 1387 break; 1388 code <<= 1; 1389 m >>= 1; 1390 if(m == 0){ 1391 sr = nextbyte(h, 0); 1392 m = 0x80; 1393 cnt = 8; 1394 } 1395 cnt--; 1396 } 1397 if(i >= 17){ 1398 if(badcode == 0) 1399 fprint(2, "badly encoded %dx%d JPEG file; ignoring bad value\n", h->X, h->Y); 1400 badcode = 1; 1401 i = 0; 1402 } 1403 h->cnt = cnt; 1404 return t->val[t->valptr[i]+(code-t->mincode[i])]; 1405 } 1406 1407 /* 1408 * load next byte of input 1409 */ 1410 static 1411 int 1412 nextbyte(Header *h, int marker) 1413 { 1414 int b, b2; 1415 1416 if(h->peek >= 0){ 1417 b = h->peek; 1418 h->peek = -1; 1419 }else{ 1420 b = Bgetc(h->fd); 1421 if(b == Beof) 1422 jpgerror(h, "truncated file"); 1423 b &= 0xFF; 1424 } 1425 1426 if(b == 0xFF){ 1427 if(marker) 1428 return b; 1429 b2 = Bgetc(h->fd); 1430 if(b2 != 0){ 1431 if(b2 == Beof) 1432 jpgerror(h, "truncated file"); 1433 b2 &= 0xFF; 1434 if(b2 == DNL) 1435 jpgerror(h, "ReadJPG: DNL marker unimplemented"); 1436 /* decoder is reading into marker; satisfy it and restore state */ 1437 Bungetc(h->fd); 1438 h->peek = b; 1439 } 1440 } 1441 h->cnt += 8; 1442 h->sr = (h->sr<<8) | b; 1443 return b; 1444 } 1445 1446 /* 1447 * return next s bits of input, MSB first, and level shift it 1448 */ 1449 static 1450 int 1451 receive(Header *h, int s) 1452 { 1453 int v, m; 1454 1455 while(h->cnt < s) 1456 nextbyte(h, 0); 1457 h->cnt -= s; 1458 v = h->sr >> h->cnt; 1459 m = (1<<s); 1460 v &= m-1; 1461 /* level shift */ 1462 if(v < (m>>1)) 1463 v += ~(m-1)+1; 1464 return v; 1465 } 1466 1467 /* 1468 * return next s bits of input, decode as EOB 1469 */ 1470 static 1471 int 1472 receiveEOB(Header *h, int s) 1473 { 1474 int v, m; 1475 1476 while(h->cnt < s) 1477 nextbyte(h, 0); 1478 h->cnt -= s; 1479 v = h->sr >> h->cnt; 1480 m = (1<<s); 1481 v &= m-1; 1482 /* level shift */ 1483 v += m; 1484 return v; 1485 } 1486 1487 /* 1488 * return next bit of input 1489 */ 1490 static 1491 int 1492 receivebit(Header *h) 1493 { 1494 if(h->cnt < 1) 1495 nextbyte(h, 0); 1496 h->cnt--; 1497 return (h->sr >> h->cnt) & 1; 1498 } 1499 1500 /* 1501 * Scaled integer implementation. 1502 * inverse two dimensional DCT, Chen-Wang algorithm 1503 * (IEEE ASSP-32, pp. 803-816, Aug. 1984) 1504 * 32-bit integer arithmetic (8 bit coefficients) 1505 * 11 mults, 29 adds per DCT 1506 * 1507 * coefficients extended to 12 bit for IEEE1180-1990 compliance 1508 */ 1509 1510 enum { 1511 W1 = 2841, /* 2048*sqrt(2)*cos(1*pi/16)*/ 1512 W2 = 2676, /* 2048*sqrt(2)*cos(2*pi/16)*/ 1513 W3 = 2408, /* 2048*sqrt(2)*cos(3*pi/16)*/ 1514 W5 = 1609, /* 2048*sqrt(2)*cos(5*pi/16)*/ 1515 W6 = 1108, /* 2048*sqrt(2)*cos(6*pi/16)*/ 1516 W7 = 565, /* 2048*sqrt(2)*cos(7*pi/16)*/ 1517 1518 W1pW7 = 3406, /* W1+W7*/ 1519 W1mW7 = 2276, /* W1-W7*/ 1520 W3pW5 = 4017, /* W3+W5*/ 1521 W3mW5 = 799, /* W3-W5*/ 1522 W2pW6 = 3784, /* W2+W6*/ 1523 W2mW6 = 1567, /* W2-W6*/ 1524 1525 R2 = 181 /* 256/sqrt(2)*/ 1526 }; 1527 1528 static 1529 void 1530 idct(int b[8*8]) 1531 { 1532 int x, y, eighty, v; 1533 int x0, x1, x2, x3, x4, x5, x6, x7, x8; 1534 int *p; 1535 1536 /* transform horizontally*/ 1537 for(y=0; y<8; y++){ 1538 eighty = y<<3; 1539 /* if all non-DC components are zero, just propagate the DC term*/ 1540 p = b+eighty; 1541 if(p[1]==0) 1542 if(p[2]==0 && p[3]==0) 1543 if(p[4]==0 && p[5]==0) 1544 if(p[6]==0 && p[7]==0){ 1545 v = p[0]<<3; 1546 p[0] = v; 1547 p[1] = v; 1548 p[2] = v; 1549 p[3] = v; 1550 p[4] = v; 1551 p[5] = v; 1552 p[6] = v; 1553 p[7] = v; 1554 continue; 1555 } 1556 /* prescale*/ 1557 x0 = (p[0]<<11)+128; 1558 x1 = p[4]<<11; 1559 x2 = p[6]; 1560 x3 = p[2]; 1561 x4 = p[1]; 1562 x5 = p[7]; 1563 x6 = p[5]; 1564 x7 = p[3]; 1565 /* first stage*/ 1566 x8 = W7*(x4+x5); 1567 x4 = x8 + W1mW7*x4; 1568 x5 = x8 - W1pW7*x5; 1569 x8 = W3*(x6+x7); 1570 x6 = x8 - W3mW5*x6; 1571 x7 = x8 - W3pW5*x7; 1572 /* second stage*/ 1573 x8 = x0 + x1; 1574 x0 -= x1; 1575 x1 = W6*(x3+x2); 1576 x2 = x1 - W2pW6*x2; 1577 x3 = x1 + W2mW6*x3; 1578 x1 = x4 + x6; 1579 x4 -= x6; 1580 x6 = x5 + x7; 1581 x5 -= x7; 1582 /* third stage*/ 1583 x7 = x8 + x3; 1584 x8 -= x3; 1585 x3 = x0 + x2; 1586 x0 -= x2; 1587 x2 = (R2*(x4+x5)+128)>>8; 1588 x4 = (R2*(x4-x5)+128)>>8; 1589 /* fourth stage*/ 1590 p[0] = (x7+x1)>>8; 1591 p[1] = (x3+x2)>>8; 1592 p[2] = (x0+x4)>>8; 1593 p[3] = (x8+x6)>>8; 1594 p[4] = (x8-x6)>>8; 1595 p[5] = (x0-x4)>>8; 1596 p[6] = (x3-x2)>>8; 1597 p[7] = (x7-x1)>>8; 1598 } 1599 /* transform vertically*/ 1600 for(x=0; x<8; x++){ 1601 /* if all non-DC components are zero, just propagate the DC term*/ 1602 p = b+x; 1603 if(p[8*1]==0) 1604 if(p[8*2]==0 && p[8*3]==0) 1605 if(p[8*4]==0 && p[8*5]==0) 1606 if(p[8*6]==0 && p[8*7]==0){ 1607 v = (p[8*0]+32)>>6; 1608 p[8*0] = v; 1609 p[8*1] = v; 1610 p[8*2] = v; 1611 p[8*3] = v; 1612 p[8*4] = v; 1613 p[8*5] = v; 1614 p[8*6] = v; 1615 p[8*7] = v; 1616 continue; 1617 } 1618 /* prescale*/ 1619 x0 = (p[8*0]<<8)+8192; 1620 x1 = p[8*4]<<8; 1621 x2 = p[8*6]; 1622 x3 = p[8*2]; 1623 x4 = p[8*1]; 1624 x5 = p[8*7]; 1625 x6 = p[8*5]; 1626 x7 = p[8*3]; 1627 /* first stage*/ 1628 x8 = W7*(x4+x5) + 4; 1629 x4 = (x8+W1mW7*x4)>>3; 1630 x5 = (x8-W1pW7*x5)>>3; 1631 x8 = W3*(x6+x7) + 4; 1632 x6 = (x8-W3mW5*x6)>>3; 1633 x7 = (x8-W3pW5*x7)>>3; 1634 /* second stage*/ 1635 x8 = x0 + x1; 1636 x0 -= x1; 1637 x1 = W6*(x3+x2) + 4; 1638 x2 = (x1-W2pW6*x2)>>3; 1639 x3 = (x1+W2mW6*x3)>>3; 1640 x1 = x4 + x6; 1641 x4 -= x6; 1642 x6 = x5 + x7; 1643 x5 -= x7; 1644 /* third stage*/ 1645 x7 = x8 + x3; 1646 x8 -= x3; 1647 x3 = x0 + x2; 1648 x0 -= x2; 1649 x2 = (R2*(x4+x5)+128)>>8; 1650 x4 = (R2*(x4-x5)+128)>>8; 1651 /* fourth stage*/ 1652 p[8*0] = (x7+x1)>>14; 1653 p[8*1] = (x3+x2)>>14; 1654 p[8*2] = (x0+x4)>>14; 1655 p[8*3] = (x8+x6)>>14; 1656 p[8*4] = (x8-x6)>>14; 1657 p[8*5] = (x0-x4)>>14; 1658 p[8*6] = (x3-x2)>>14; 1659 p[8*7] = (x7-x1)>>14; 1660 } 1661 }