home *** CD-ROM | disk | FTP | other *** search
/ Microsoft Programmer's Library 1.3 / Microsoft-Programers-Library-v1.3.iso / sampcode / os2sdk / os2sdk10 / apps / mandel / mandel.c < prev    next >
Encoding:
C/C++ Source or Header  |  1988-08-11  |  10.0 KB  |  403 lines

  1. /***    Mandel.c - compute and display Mandelbrot set
  2.  *    R. A. Garmoe  86/12/10
  3.  */
  4.  
  5.  
  6.  
  7. /**    Mandel compute the elements of the Mandelbrot set
  8.  *
  9.  *        z = z**2 + c
  10.  *
  11.  *    where z and c are in the complex plane and z = 0 + 0i for the first
  12.  *    iteration.  The elements of the set along with the description of the
  13.  *    computation are written to the .cnt file.
  14.  */
  15.  
  16.  
  17.  
  18.  
  19. /**    Call
  20.  *
  21.  *    mandel -v {-u a+bi} {-l c+di} {-r num} {-c num} {-i itr} {-a f} {-f name}
  22.  *        {-d display}
  23.  *
  24.  *    where
  25.  *
  26.  *        a    compute each point with an aspect ratio of f where
  27.  *            f is :
  28.  *                f.ffff    floating point number describing the
  29.  *                    ratio of width vs height of each pixel
  30.  *
  31.  *        c    number of points across the real coordinate.  This
  32.  *            value cannot exceed 2000 and is defaulted to 640.
  33.  *
  34.  *        d    set display to device.    The default display is EGA.
  35.  *                CGA    Color Graphics Adapter
  36.  *                    columns = 350
  37.  *                    rows    = 200
  38.  *                    aspect    = 2.4
  39.  *
  40.  *                EGA    Extended Graphics Adapter
  41.  *                    columns = 640
  42.  *                    rows    = 350
  43.  *                    aspect    = 2.4
  44.  *
  45.  *                PGA    Professional Graphics Adapter
  46.  *                    columns = 640
  47.  *                    rows    = 480
  48.  *                    aspect    = 1.0
  49.  *
  50.  *        f    name for the .cnt file
  51.  *            If name is not specified, "mandel" is used.
  52.  *
  53.  *        i    iteration limit.  This value cannot exceed 1000 and
  54.  *            the default is 256
  55.  *
  56.  *        l    lower right corner as c+di
  57.  *
  58.  *
  59.  *        r    number of points down the imaginary coordinate.  This
  60.  *            value is defaulted to 350.
  61.  *
  62.  *        u    upper left corner as a+bi
  63.  *
  64.  *        v    display information about each scan line as processed
  65.  *
  66.  *    The .cnt file has the format
  67.  *        int    number of points along the real axis
  68.  *        int    number of points along the imaginary axis
  69.  *        int    maximum iteration point for each point
  70.  *        double    real coordinate of upper left
  71.  *        double    imaginary coordinate of upper left
  72.  *        double    real coordinate of lower right
  73.  *        double    imaginary coordinate lower rightft
  74.  *        long    (loop + 1) counters for histogram values
  75.  */
  76.  
  77.  
  78.  
  79.  
  80. #include <stdio.h>
  81. #include <doscalls.h>
  82.  
  83. int mandinter ();
  84.  
  85. #define FALSE    0
  86. #define TRUE    ~FALSE
  87.  
  88. #define MAXREAL 2000        /* maximum number of points in line */
  89. #define MAXLOOP 1000        /* maximum number of iterations */
  90.  
  91. #define CGA_DEV 0        /* Color Graphics Adapter */
  92. #define CGA_COL 350        /* number of display columns for CGA */
  93. #define CGA_ROW 200        /* number of display rows for CGA */
  94. #define CGA_ASP 12/5        /* aspect ratio (width/height) for CGA */
  95.  
  96. #define EGA_DEV 1        /* Extended Graphics Adapter */
  97. #define EGA_COL 640        /* number of display columns for EGA */
  98. #define EGA_ROW 350        /* number of display rows for EGA */
  99. #define EGA_ASP 1.33        /* aspect ratio (width/height) for EGA */
  100.  
  101. #define PGA_DEV 2        /* Professional Graphics Adapter */
  102. #define PGA_COL 640        /* number of display columns for PGA */
  103. #define PGA_ROW 480        /* number of display rows for PGA */
  104. #define PGA_ASP 1.0        /* aspect ratio (width/height) for PGA */
  105.  
  106.  
  107. struct cmplx {
  108.     double     realp;     /* real part of number */
  109.     double     imagp;     /* imaginary part of number */
  110. };
  111.  
  112.  
  113.  
  114. char    fp287;            /* 287 processor status from DOSDEVCONFIG */
  115. char     pmand[60] = "mandel.cnt"; /* file name for loop counts */
  116. FILE    *fmand;
  117.  
  118. double    rinc;            /* increment in real coordinate */
  119. double    iinc;            /* increment in imaginary coordinate */
  120.  
  121. int    device = EGA_DEV;    /* default device type */
  122. int    arg_device = -1;    /* device type from arguments */
  123. int    mcount[MAXREAL];    /* array holding interation counters */
  124. int    rlcount[MAXREAL + 2];    /* run length encoded interation counters */
  125. int    nreal = EGA_COL;    /* number of points on real coordinate */
  126. int    arg_nreal = 0;        /* number real coordinate from arguments */
  127. int    nimag = EGA_ROW;    /* number of points on imaginary coordinate */
  128. int    arg_nimag = 0;        /* number imaginary coordinate from arguments */
  129. int    loop = 256;        /* maximum loop count for each point */
  130. int    verbose = FALSE;    /* display scan line information if true */
  131. double    aspect = EGA_ASP;    /* default screen aspect ratio */
  132. double    arg_aspect = 0.0;    /* aspect from arguments */
  133.  
  134. long    hist[MAXLOOP + 1] = {0}; /* histogram counters */
  135.  
  136. struct    cmplx c;        /* value of mandelbrot seed */
  137. struct    cmplx ul;        /* coordinates of upper left corner */
  138. struct    cmplx lr;        /* coordinates of lower right corner */
  139.  
  140.  
  141. main (argc, argv)
  142. int    argc;
  143. char    **argv;
  144. {
  145.     int    nr;
  146.     int    ni;
  147.     int    i;
  148.     int    temp;
  149.     int    *rl;
  150.     long    position;
  151.  
  152.     if (DOSDEVCONFIG ((char far *)&fp287, 3, 0) != 0) {
  153.         printf ("Unable to get fp processor presence flag\n");
  154.         exit (1);
  155.     }
  156.     if (fp287 != 1) {
  157.         printf ("There is not a 287 fp processor present\n");
  158.         exit (1);
  159.     }
  160.     nextarg (argc, argv);
  161.     if ((fmand = fopen (pmand, "wb")) == NULL) {
  162.         printf ("Unable to open count file %s\n", pmand);
  163.         exit (3);
  164.     }
  165.     if (ul.realp == lr.realp && ul.imagp == lr.imagp) {
  166.         printf ("%15.13lf+%15.13lfi, %15.13lf+%15.13lfi\n",ul.realp,ul.imagp,lr.realp,lr.imagp);
  167.         printf ("Upper left real = lower right real\n");
  168.         exit (3);
  169.     }
  170.  
  171.     /* reset nreal, nimag, aspect from command line device type */
  172.  
  173.     switch (arg_device) {
  174.         case CGA_DEV:
  175.             nreal = CGA_COL;
  176.             nimag = CGA_ROW;
  177.             aspect = CGA_ASP;
  178.             break;
  179.  
  180.         case EGA_DEV:
  181.             nreal = EGA_COL;
  182.             nimag = EGA_ROW;
  183.             aspect = EGA_ASP;
  184.             break;
  185.  
  186.         case PGA_DEV:
  187.             nreal = PGA_COL;
  188.             nimag = PGA_ROW;
  189.             aspect = PGA_ASP;
  190.             break;
  191.     }
  192.  
  193.     /* reset individual parameters */
  194.  
  195.     if (arg_nreal != 0)
  196.         nreal = arg_nreal;
  197.     if (arg_nimag != 0)
  198.         nimag = arg_nimag;
  199.     if (arg_aspect != 0)
  200.         aspect = arg_aspect;
  201.  
  202.     /* compute the step increments for the specifed aspect */
  203.  
  204.     rinc = (lr.realp - ul.realp) / (double)nreal;
  205.     iinc = rinc * aspect;
  206.  
  207.     lr.imagp = ul.imagp - iinc * nimag;
  208.  
  209.     fwrite ((char *)&nreal, sizeof (int), 1, fmand);
  210.     fwrite ((char *)&nimag, sizeof (int), 1, fmand);
  211.     fwrite ((char *)&loop, sizeof (int), 1, fmand);
  212.     fwrite ((char *)&ul, sizeof (ul), 1, fmand);
  213.     fwrite ((char *)&lr, sizeof (lr), 1, fmand);
  214.     fwrite ((char *)&rinc, sizeof (rinc), 1, fmand);
  215.     fwrite ((char *)&iinc, sizeof (iinc), 1, fmand);
  216.     fwrite ((char *)&aspect, sizeof (aspect), 1, fmand);
  217.     position = ftell (fmand);
  218.     fwrite ((char *)hist, sizeof (long), loop + 1, fmand);
  219.  
  220.     c.imagp = ul.imagp;
  221.     for (ni = 0; ni < nimag; ni++) {
  222.         c.realp = ul.realp;
  223.         if (verbose)
  224.             printf ("%4d, %15.13lf + %15.13lfi", ni, c.realp, c.imagp);
  225.         manditer (c.realp, c.imagp, loop, nreal, mcount, rinc);
  226.         for (nr = 0; nr < nreal; nr++) {
  227.             i = mcount[nr];
  228.             hist[i]++;
  229.         }
  230.  
  231.         /* run length encode counts for this scan line */
  232.  
  233.         temp = *mcount;
  234.         nr = 1;
  235.         rl = rlcount + 1;
  236.         i = -1;
  237.         for (; nr < nreal; nr++) {
  238.             i++;
  239.             if (mcount[nr] != temp) {
  240.                 /* process change in value */
  241.                 if (i != 0)
  242.                     /* output count of duplicate values */
  243.                     *rl++ = -(i + 1);
  244.                 /* output value */
  245.                 *rl++ = temp;
  246.                 /* reset compression values */
  247.                 temp = mcount[nr];
  248.                 i = -1;
  249.             }
  250.         }
  251.         if (i != -1) {
  252.             if (i != 0)
  253.                 /* output count of duplicate values */
  254.                 *rl++ = -(i + 2);
  255.             /* output value */
  256.             *rl++ = temp;
  257.         }
  258.         *rlcount = rl - rlcount -1;
  259.         if (verbose)
  260.             printf ("%4d\n", *rlcount);
  261.         if (fwrite ((char *)rlcount, sizeof (int), *rlcount + 1, fmand) !=
  262.             *rlcount + 1)  {
  263.             printf ("error writing count file %s\n", pmand);
  264.             exit (4);
  265.         }
  266.         c.imagp -= iinc;
  267.     }
  268.     fseek (fmand, position, 0);
  269.     if (fwrite ((char *)hist, sizeof (long), loop + 1, fmand) != loop + 1) {
  270.         printf ("error writing histogram to file %s\n", pmand);
  271.         exit (4);
  272.     }
  273. }
  274.  
  275.  
  276.  
  277. /**    nextarg - process arguments
  278.  *
  279.  */
  280.  
  281.  
  282. nextarg (argc, argv)
  283. int    argc;
  284. char    **argv;
  285. {
  286.  
  287.     argv++;
  288.     while ((argc-- > 0) && (**argv == '-')) {
  289.         switch (*(*argv+1)) {
  290.             case 'a':
  291.                 /* specify aspect ratio */
  292.                 argv++;
  293.                 if ((argc < 1) ||
  294.                     (sscanf (*argv++, "%lf", &arg_aspect) != 1)) {
  295.                     printf ("Error specifing aspect\n");
  296.                     exit (1);
  297.                 }
  298.                 argc--;
  299.                 break;
  300.  
  301.             case 'c':
  302.                 /* specify number of columns and limit to maximum value */
  303.                 argv++;
  304.                 if ((argc < 1) ||
  305.                     (sscanf (*argv++, " %d", &arg_nreal) != 1)) {
  306.                     printf ("Error specifying number of real points\n");
  307.                     exit (1);
  308.                 }
  309.                 argc--;
  310.                 if (arg_nreal > MAXREAL) {
  311.                     arg_nreal = MAXREAL;
  312.                     printf ("Number of columns limited to %d\n", MAXREAL);
  313.                 }
  314.                 break;
  315.  
  316.             case 'd':
  317.                 /* specify default device paramters */
  318.                 argv++;
  319.                 if (strcmp (*argv, "CGA") == 0)
  320.                     arg_device = CGA_DEV;
  321.                 else if (strcmp (*argv, "EGA") == 0)
  322.                     arg_device = EGA_DEV;
  323.                 else if (strcmp (*argv, "PGA") == 0)
  324.                     arg_device = PGA_DEV;
  325.                 else {
  326.                     printf ("Error specifing device\n");
  327.                     exit (1);
  328.                 }
  329.                 argv++;
  330.                 argc--;
  331.                 break;
  332.  
  333.             case 'f':
  334.                 /* specify output file name */
  335.                 argv++;
  336.                 /* set file name */
  337.                 strcpy (pmand, *argv);
  338.                 strcat (pmand, ".cnt");
  339.                 break;
  340.  
  341.             case 'i':
  342.                 /* specify iteration count and limit to maximum */
  343.                 argv++;
  344.                 if ((argc < 1) ||
  345.                     (sscanf (*argv++, " %d", &loop) != 1)) {
  346.                     printf ("Error specifying iteration limit\n");
  347.                     exit (1);
  348.                 }
  349.                 argc--;
  350.                 if (loop > MAXLOOP) {
  351.                     loop = MAXREAL;
  352.                     printf ("Iteration count limited to %d\n", MAXLOOP);
  353.                 }
  354.                 break;
  355.  
  356.             case 'l':
  357.                 /* specify lower right corner as c+di */
  358.                 argv++;
  359.                 if ((argc < 1) || (sscanf (*argv++, "%lf%lfi",
  360.                     &lr.realp, &lr.imagp) != 2)) {
  361.                     printf ("Error specifying lower right\n");
  362.                     exit (1);
  363.                 }
  364.                 argc--;
  365.                 break;
  366.  
  367.  
  368.             case 'r':
  369.                 /* specify number of rows */
  370.                 argv++;
  371.                 if ((argc < 1) ||
  372.                     (sscanf (*argv++, " %d", &arg_nimag) != 1)) {
  373.                     printf ("Error specifying number of imaginary points\n");
  374.                     exit (1);
  375.                 }
  376.                 argc--;
  377.                 break;
  378.  
  379.             case 'u':
  380.                 /* specify upper right corner as a+bi */
  381.                 argv++;
  382.                 if ((argc < 1) || (sscanf (*argv++, "%lf%lfi",
  383.                     &ul.realp, &ul.imagp) != 2)) {
  384.                     printf ("Error specifying upper left\n");
  385.                     exit (1);
  386.                 }
  387.                 argc--;
  388.                 break;
  389.  
  390.             case 'v':
  391.                 /* specify verbose output */
  392.                 argv++;
  393.                 verbose = TRUE;
  394.                 break;
  395.  
  396.             default:
  397.                 printf ("Unknown argument %s\n", *argv);
  398.                 exit (1);
  399.  
  400.         }
  401.     }
  402. }
  403.