home *** CD-ROM | disk | FTP | other *** search
/ Collection of Hack-Phreak Scene Programs / cleanhpvac.zip / cleanhpvac / FRCOMP.ZIP / fracpack.c next >
C/C++ Source or Header  |  1992-11-12  |  13KB  |  327 lines
  1. Here is some C source for fractal image compression.
  2. Notes:
  3. 1)  You will need DJGPP (a C compiler for DOS ) to handle the arrays.
  4. 2)  It is available for anonymous ftp at
  5. ftp.math.niu.edu:/pub/msdos/djgpp.stuff/v109   (i think, look around)
  6. 3)  A 386+387 or higher is REALLY required.
  7. 4)  This should compile on other systems... graphics stuff will need to 
  8. change though...
  9. 5)  Images are 256x256x256 shades of gray.  Totals 65536 bytes. 1 byte/pixel
  10. 6)  I'm not a C programmer, so don't complain about the code or lack of
  11. documentation.  
  12. 7)  Get the paper /pub/inls-ucsd/fractal-2.0.tar from lyapunov.ucsd.edu
  13. It is a very good description of the process by yuval fisher.
  14. 8)  Decompression code will be sent in another message.
  15. 9)  Enjoy... there is also some later work of mine in /pub/young-fractal
  16. at the same anonymous ftp site.
  17. 10) It will be updated soon... (2+ weeks)
  18.  
  19.  
  20. ----------------- Cut Here -------------------
  21.  
  22. /*****************************************************************************
  23.                          FRACPACK - a program for
  24.                          FRACTAL IMAGE COMPRESSION
  25.                           image.krd ==> image.ifs
  26. *****************************************************************************/
  27.  
  28. #include <stdio.h>
  29. #include <stdlib.h>
  30. #include <graphics.h>
  31. #include <math.h>
  32. int main(int argc, char **argv)
  33. {  
  34.     unsigned char IMAGE [256][256], Range[8][8][8], Domain[8][8];
  35.     int REDGE[32][32], DEDGE[241][241];
  36.     FILE *in, *out;
  37.     char *inf, *outf, rmsstr[13];
  38.     int xsize=256, ysize=256, x, y, dx, dy, rx, ry, bestdx, bestdy, inrx, inry;
  39.     int class1, class2, class3, class4, i, besti, nflips=8, rsize, rsize2,
  40. plotdx, plotdy, usegraph;
  41.     long int suma,  sumb,  sumab,  suma2,  sumb2, sumx, sumx2;
  42.     long int DSUMS[241][241][2];
  43.     float   fsuma, fsumb, fsumab, fsuma2, fsumb2, fmagica;
  44.     float bests, besto, scale, offset, rms, bestrms, ifs_tolerance;
  45.     struct trans_t {
  46.       unsigned char dx;
  47.       unsigned char dy;
  48.       signed char scale;
  49.       short int offset : 13;
  50.       unsigned short int flip : 3;
  51.       } trans[1];
  52.  
  53.    
  54.     if (argc < 4)
  55.         {
  56.         printf("usage: fracpack rms infile.ext outfile.ext [graphics flag]");
  57.         return 1;
  58.         }
  59.     ifs_tolerance = atof(argv[1]);
  60.     inf = argv[2]; outf = argv[3];
  61.     if (argc == 5) usegraph = -1;
  62.     else usegraph = 0;
  63.     if ((in = fopen(inf, "rb")) == NULL)
  64.         {
  65.         fprintf(stderr, "Cannot open input file.\n");
  66.         return 1;
  67.         }
  68.     if ((out = fopen(outf, "wb")) == NULL)
  69.         {
  70.         fprintf(stderr, "Cannot open output file.\n");
  71.         return 1;
  72.         }
  73.     fclose(out);
  74.     printf("Compressing %s to %s with an rms of
  75. %4.2f\n",inf,outf,ifs_tolerance);
  76.  
  77.     ifs_tolerance=ifs_tolerance*ifs_tolerance;
  78.     if (usegraph)
  79.         {
  80.         GrSetMode(GR_default_graphics);
  81.         for (x = 0; x < xsize; x++)
  82.         GrSetColor(x,x,x,x);
  83.         }
  84.  
  85. /* Get .KRD bitmap 8-bit-grey-scale pixels, put in IMAGE array.  */
  86.     for (y = 0; y < ysize; y++)
  87.     for (x = 0; x < xsize; x++)
  88.         {
  89.         IMAGE[y][x] = fgetc(in);
  90.         if (usegraph)
  91.             {
  92.             GrPlot(x,y,IMAGE[y][x]);
  93.             GrPlot(x+300,y,IMAGE[y][x]);
  94.             }
  95.         }
  96.    fclose(in);
  97.  
  98. /*************************************************************************/
  99. /*                            Classify Range's                           */
  100. /*************************************************************************/
  101.     if (!usegraph) printf("Classifying RangeS\n");
  102.     for (ry = 0; ry < ysize; ry+=8)
  103.     for (rx = 0; rx < xsize; rx+=8)
  104.         {
  105.         for (y = 0; y < 8; y++)
  106.         for (x = 0; x < 8; x++)
  107.             Range[0][y][x] = IMAGE [ ry  +y ] [ rx  +x ];
  108.  
  109.         class1 = (Range[0][0][0]==Range[0][0][1]) +
  110. (Range[0][0][1]==Range[0][0][2])
  111.                + (Range[0][0][2]==Range[0][0][3]) +
  112. (Range[0][0][3]==Range[0][0][4])
  113.                + (Range[0][0][4]==Range[0][0][5]) +
  114. (Range[0][0][5]==Range[0][0][6])
  115.                + (Range[0][0][6]==Range[0][0][7]);
  116.  
  117.         class2 = (Range[0][7][0]==Range[0][7][1]) +
  118. (Range[0][7][1]==Range[0][7][2])
  119.                + (Range[0][7][2]==Range[0][7][3]) +
  120. (Range[0][7][3]==Range[0][7][4])
  121.                + (Range[0][7][4]==Range[0][7][5]) +
  122. (Range[0][7][5]==Range[0][7][6])
  123.                + (Range[0][7][6]==Range[0][7][7]);
  124.  
  125.         class3 = (Range[0][0][0]==Range[0][1][0]) +
  126. (Range[0][1][0]==Range[0][2][0])
  127.                + (Range[0][2][0]==Range[0][3][0]) +
  128. (Range[0][3][0]==Range[0][4][0])
  129.                + (Range[0][4][0]==Range[0][5][0]) +
  130. (Range[0][5][0]==Range[0][6][0])
  131.                + (Range[0][6][0]==Range[0][7][0]);
  132.  
  133.         class4 = (Range[0][0][7]==Range[0][1][7]) +
  134. (Range[0][1][7]==Range[0][2][7])
  135.                + (Range[0][2][7]==Range[0][3][7]) +
  136. (Range[0][3][7]==Range[0][4][7])
  137.                + (Range[0][4][7]==Range[0][5][7]) +
  138. (Range[0][5][7]==Range[0][6][7])
  139.                + (Range[0][6][7]==Range[0][7][7]);
  140.         class1 = class1*1000 + class2*100 + class3*10 + class4;
  141.         REDGE[ry>>3][rx>>3] = class1;
  142.         }
  143.  
  144.  
  145. /*************************************************************************/
  146. /*               Classify Domain's and compute DSUMS array               */
  147. /*************************************************************************/
  148.         if (!usegraph) printf("Classifying DomainS\n");
  149.         for (dy = 0; dy < ysize - 15; dy++)
  150.         for (dx = 0; dx < xsize - 15; dx++)
  151.             {
  152.             DSUMS[dy][dx][0] = DSUMS[dy][dx][1] = 0;
  153.             if (usegraph) GrPlot(dx+308,dy+8,384);
  154.             for (y = 0; y < 16; y+=2)
  155.             for (x = 0; x < 16; x+=2)
  156.                 {
  157.                 Domain[y>>1][x>>1] = (IMAGE[dy+y  ][dx+x  ]
  158.                                     + IMAGE[dy+y  ][dx+x+1]
  159.                                     + IMAGE[dy+y+1][dx+x  ]
  160.                                     + IMAGE[dy+y+1][dx+x+1]) >> 2;
  161.  
  162.                 DSUMS[dy][dx][0] += Domain[y>>1][x>>1];
  163.                 DSUMS[dy][dx][1] += Domain[y>>1][x>>1]*Domain[y>>1][x>>1];
  164.                 }
  165.  
  166.         class1 = (Domain[0][0]==Domain[0][1]) + (Domain[0][1]==Domain[0][2])
  167.                + (Domain[0][2]==Domain[0][3]) + (Domain[0][3]==Domain[0][4])
  168.                + (Domain[0][4]==Domain[0][5]) + (Domain[0][5]==Domain[0][6])
  169.                + (Domain[0][6]==Domain[0][7]);
  170.  
  171.         class2 = (Domain[7][0]==Domain[7][1]) + (Domain[7][1]==Domain[7][2])
  172.                + (Domain[7][2]==Domain[7][3]) + (Domain[7][3]==Domain[7][4])
  173.                + (Domain[7][4]==Domain[7][5]) + (Domain[7][5]==Domain[7][6])
  174.                + (Domain[7][6]==Domain[7][7]);
  175.  
  176.         class3 = (Domain[0][0]==Domain[1][0]) + (Domain[1][0]==Domain[2][0])
  177.                + (Domain[2][0]==Domain[3][0]) + (Domain[3][0]==Domain[4][0])
  178.                + (Domain[4][0]==Domain[5][0]) + (Domain[5][0]==Domain[6][0])
  179.                + (Domain[6][0]==Domain[7][0]);
  180.  
  181.         class4 = (Domain[0][7]==Domain[1][7]) + (Domain[1][7]==Domain[2][7])
  182.                + (Domain[2][7]==Domain[3][7]) + (Domain[3][7]==Domain[4][7])
  183.                + (Domain[4][7]==Domain[5][7]) + (Domain[5][7]==Domain[6][7])
  184.                + (Domain[6][7]==Domain[7][7]);
  185.  
  186.         class1 = class1*1000 + class2*100 + class3*10 + class4;
  187.         DEDGE[dy][dx] = class1;
  188.         if (usegraph) GrPlot(dx+308,dy+8,384);
  189.         }
  190.  
  191. /* Run through all non-overlapping 8x8 "R" blocks in the image  */
  192.     for (ry = 0; ry < ysize; ry+=8)
  193.     for (rx = 0; rx < xsize; rx+=8)
  194.         {
  195.         if (usegraph)
  196.             {
  197.             GrLine(rx  ,ry  ,rx  ,ry+7,384);
  198.             GrLine(rx  ,ry  ,rx+7,ry  ,384);
  199.             GrLine(rx+7,ry  ,rx+7,ry+7,384);
  200.             GrLine(rx  ,ry+7,rx+7,ry+7,384);
  201.             }
  202.         else printf("Entering Range Loop\n");
  203.         sumb = sumb2 = 0;
  204. /*************************************************************************/
  205. /*                   Grab Range   (and its "flips").                     */
  206. /*************************************************************************/
  207.         for (y = 0; y < 8; y++)
  208.         for (x = 0; x < 8; x++)
  209.             {
  210.             Range[0][y][x] = IMAGE [ ry  +y ] [ rx  +x ];
  211.             Range[1][y][x] = IMAGE [ ry+7-x ] [ rx  +y ];
  212.             Range[2][y][x] = IMAGE [ ry+7-y ] [ rx+7-x ];
  213.             Range[3][y][x] = IMAGE [ ry  +x ] [ rx+7-y ];
  214.  
  215.             Range[4][y][x] = IMAGE [ ry+7-y ] [ rx  +x ];
  216.             Range[5][y][x] = IMAGE [ ry+7-x ] [ rx+7-y ];
  217.             Range[6][y][x] = IMAGE [ ry  +y ] [ rx+7-x ];
  218.             Range[7][y][x] = IMAGE [ ry  +x ] [ rx  +y ];
  219.             sumb+=Range[0][y][x]; sumb2+=Range[0][y][x]*Range[0][y][x];
  220.             }
  221.  
  222. /*************************************************************************/
  223. /*                  Grab Domain, compare with Range(s)                   */
  224. /*************************************************************************/
  225.         bestrms = 10000000000.;
  226.         for (dy = 0; dy < ysize - 15; dy++)
  227.         for (dx = 0; dx < xsize - 15; dx++)
  228.             if ((DEDGE[dy][dx] == REDGE[ry>>3][rx>>3])||
  229.                 (((dx - 16 < rx)&&(rx < dx + 16))&&
  230.                 ((dy - 16 < ry)&&(ry < dy + 16))))
  231.                 {
  232.                 if (usegraph)
  233.                     {
  234.                     plotdy = dy; plotdx = dx;
  235.                     GrPlot(dx+308,dy+8,384);
  236.                     }
  237.                 for (y = 0; y < 16; y+=2)
  238.                 for (x = 0; x < 16; x+=2)
  239.                     Domain[y>>1][x>>1] = (IMAGE[dy+y  ][dx+x  ]
  240.                                         + IMAGE[dy+y  ][dx+x+1]
  241.                                         + IMAGE[dy+y+1][dx+x  ]
  242.                                         + IMAGE[dy+y+1][dx+x+1]) >> 2;
  243.  
  244.                 suma = DSUMS[dy][dx][0]; suma2 = DSUMS[dy][dx][1];
  245.                 fsuma=suma; fsuma2=suma2;
  246.                 fsumb=sumb; fsumb2=sumb2;
  247.                 fmagica = (float) (suma2 - suma*suma/64.);
  248.  
  249.                 for (i = 0; i < nflips; i++)
  250.                     {
  251.                     sumab = 0;
  252.                     for (y = 0; y < 8; y++)
  253.                     for (x = 0; x < 8; x++)
  254.                         sumab += Domain[y][x] * Range[i][y][x];
  255.                     fsumab = sumab;
  256.                     if (fmagica != 0.)
  257.                         scale = (fsumab - fsuma*fsumb/64.)/fmagica;
  258.                     else scale = 0;
  259.                     if (scale*scale < 1.44)
  260.                         {
  261.                         scale = (signed char) 127. * scale / 1.2;
  262.                         scale = 1.2 * scale / 127.;
  263.                         offset = (short int) (fsumb - scale*fsuma)/64.;
  264.                         rms = (fsumb2 + scale*(scale*fsuma2 - 2*fsumab +
  265. 2*offset*fsuma)
  266.                             + offset*(offset*64. - 2.*fsumb)) / 64.;
  267.                         if (rms > 10000) (i = nflips);
  268.                         if (rms < bestrms)
  269.                             {
  270.                             if (!usegraph) printf("%3i%4i%4i%4i%4i%2i
  271. %12.8f\n",REDGE[ry][rx],dx, dy, rx, ry, i, rms);
  272.                             besti = i;  bestrms = rms; bestdx = dx; bestdy =
  273. dy;
  274.                             bests = scale;  besto = offset;
  275.                             }
  276.                         if (rms < ifs_tolerance) goto gotbest;
  277.                         }
  278.                     }
  279.                 } /* End of Domain LOOP */
  280.  
  281. gotbest:
  282.         if (usegraph)
  283.             {
  284.             for (dy = 0; dy <= plotdy; dy++)
  285.             for (dx = 0; (dx <= plotdx) || ((dy < plotdy)&&(dx < xsize - 15));
  286. dx++)
  287.                 if ((DEDGE[dy][dx] == REDGE[ry>>3][rx>>3])||
  288.                 (((dx - 16 < rx)&&(rx < dx + 16))&&
  289.                 ((dy - 16 < ry)&&(ry < dy + 16))))
  290.                     GrPlot(dx+308,dy+8,384);
  291.             sprintf(rmsstr,"%8.4f",sqrt((double)bestrms));
  292.             GrTextXY(560,20,rmsstr,255,0);
  293.             }
  294.         if ((out = fopen(outf, "ab")) == NULL)
  295.             {
  296.             fprintf(stderr, "Cannot open output file.\n");
  297.             return 1;
  298.             }
  299.  
  300.         if (!usegraph)
  301.             {
  302.             printf("best rms=%12.8f,dx=%5i,dy=%5i,flip=%i for
  303.             rx=%5i,ry=%5i,s=%12.5f\n",  bestrms,
  304.             bestdx, bestdy, besti, rx, ry,bests);
  305.             }
  306.         trans[0].dx = bestdx;
  307.         trans[0].dy = bestdy;
  308.         trans[0].flip = besti;
  309.         trans[0].scale = 127. * bests / 1.2;
  310.         trans[0].offset = besto;
  311.         fwrite(trans, sizeof(struct trans_t), 1, out);  /* write out compressed
  312. file */
  313.         fclose(out);
  314.         if (usegraph)
  315.             {
  316.             GrLine(rx  ,ry  ,rx  ,ry+7,384);
  317.             GrLine(rx  ,ry  ,rx+7,ry  ,384);
  318.             GrLine(rx+7,ry  ,rx+7,ry+7,384);
  319.             GrLine(rx  ,ry+7,rx+7,ry+7,384);
  320.             }
  321.       }
  322.    if (usegraph) GrSetMode(GR_default_text);
  323. /* All done. Whew... */
  324.    return 0;
  325. }
  326.  
  327.