27 files changed, 2596 insertions, 1216 deletions

diff --git a/generator/Makefile b/generator/Makefile
index d54a6cc5..761181bd 100644
--- a/generator/Makefile
+++ b/generator/Makefile
@@ -18,7 +18,7 @@ SUBDIRS = pamtris
 
 PORTBINARIES = pamcrater pamgauss pamgradient \
 	       pamseq pamshadedrelief pamstereogram \
-	       pbmpage pbmmake pbmtext pbmupc \
+	       pbmpage pbmmake pbmnoise pbmtext pbmupc \
 	       pgmkernel pgmmake pgmnoise pgmramp \
 	       ppmcie ppmcolors ppmforge ppmmake ppmpat ppmrough ppmwheel \
 
diff --git a/generator/pamcrater.c b/generator/pamcrater.c
index 43c27dbc..8f0e422d 100644
--- a/generator/pamcrater.c
+++ b/generator/pamcrater.c
@@ -48,6 +48,7 @@
 
 #include "pm_c_util.h"
 #include "mallocvar.h"
+#include "rand.h"
 #include "shhopt.h"
 #include "nstring.h"
 #include "pam.h"
@@ -162,18 +163,18 @@ parseCommandLine(int argc, const char ** const argv,
 
 
 
-static double const arand       = 32767.0;  /* Random number parameters */
 static double const CdepthPower = 1.5;      /* Crater depth power factor */
 static double const DepthBias2  = 0.5;      /* Square of depth bias */
 
 
 
 static double const
-cast(double const high) {
+cast(double             const high,
+     struct pm_randSt * const randStP) {
 /*----------------------------------------------------------------------------
    A random number in the range [0, 'high'].
 -----------------------------------------------------------------------------*/
-    return high * ((rand() & 0x7FFF) / arand);
+  return high * ((double) pm_rand(randStP) / randStP->max);
 }
 
 
@@ -220,7 +221,6 @@ terrainModP(struct pam * const pamP,
 
 
 
-
 static sample
 terrainMod(struct pam * const pamP,
            tuple **     const terrain,
@@ -252,11 +252,12 @@ setElev(struct pam * const pamP,
 
 
 static void
-smallCrater(struct pam * const pamP,
-            tuple **     const terrain,
-            int          const cx,
-            int          const cy,
-            double       const radius) {
+smallCrater(struct pam *       const pamP,
+            tuple **           const terrain,
+            int                const cx,
+            int                const cy,
+            double             const radius,
+            struct pm_randSt * const randStP) {
 /*----------------------------------------------------------------------------
    Generate a crater with a special method for tiny craters.
 
@@ -283,10 +284,10 @@ smallCrater(struct pam * const pamP,
             /* The mean elevation of the Moore neighborhood (9 pixels
                centered on the crater location).
             */
-        
+
         /* Perturb the mean elevation by a small random factor. */
 
-        int const x = radius >= 1 ? ((rand() >> 8) & 0x3) - 1 : 0;
+        int const x = radius >= 1 ? ((pm_rand(randStP) >> 8) & 0x3) - 1 : 0;
 
         assert(axelev > 0);
 
@@ -374,7 +375,7 @@ normalCrater(struct pam * const pamP,
                 av = (axelev + cz) * (1 - roll) +
                     (terrainMod(pamP, terrain, x, y) + cz) * roll;
                 av = MAX(1000, MIN(64000, av));
-                
+
                 setElev(pamP, terrain, x, y, av);
             }
         }
@@ -388,19 +389,20 @@ normalCrater(struct pam * const pamP,
 
 
 static void
-plopCrater(struct pam * const pamP,
-           tuple **     const terrain,
-           int          const cx,
-           int          const cy,
-           double       const radius,
-           bool         const verbose) {
+plopCrater(struct pam *       const pamP,
+           tuple **           const terrain,
+           int                const cx,
+           int                const cy,
+           double             const radius,
+           bool               const verbose,
+           struct pm_randSt * const randStP) {
 
     if (verbose && pm_have_float_format())
         pm_message("Plopping crater at (%4d, %4d) with radius %g",
                    cx, cy, radius);
 
     if (radius < 3)
-        smallCrater (pamP, terrain, cx, cy, radius);
+        smallCrater (pamP, terrain, cx, cy, radius, randStP);
     else
         normalCrater(pamP, terrain, cx, cy, radius);
 }
@@ -448,6 +450,7 @@ genCraters(struct CmdlineInfo const cmdline) {
 -----------------------------------------------------------------------------*/
     tuple ** terrain;    /* elevation array */
     struct pam pam;
+    struct pm_randSt randSt;
 
     /* Allocate the elevation array and initialize it to mean surface
        elevation.
@@ -455,18 +458,22 @@ genCraters(struct CmdlineInfo const cmdline) {
 
     initCanvas(cmdline.width, cmdline.height, &pam, &terrain);
 
+    pm_randinit(&randSt);
+    pm_srand2(&randSt, cmdline.randomseedSpec, cmdline.randomseed);
+
     if (cmdline.test)
         plopCrater(&pam, terrain,
                    pam.width/2 + cmdline.offset,
                    pam.height/2 + cmdline.offset,
-                   (double) cmdline.radius, cmdline.verbose);
+                   (double) cmdline.radius, cmdline.verbose,
+                   &randSt);
     else {
         unsigned int const ncraters = cmdline.number; /* num of craters */
         unsigned int l;
 
         for (l = 0; l < ncraters; ++l) {
-            int const cx = cast((double) pam.width  - 1);
-            int const cy = cast((double) pam.height - 1);
+            int const cx = cast((double) pam.width  - 1, &randSt);
+            int const cy = cast((double) pam.height - 1, &randSt);
 
             /* Thanks, Rudy, for this equation that maps the uniformly
                distributed numbers from cast() into an area-law distribution
@@ -475,9 +482,11 @@ genCraters(struct CmdlineInfo const cmdline) {
                Produces values within the interval:
                0.56419 <= radius <= 56.419
             */
-            double const radius = sqrt(1 / (M_PI * (1 - cast(0.9999))));
+            double const radius =
+                sqrt(1 / (M_PI * (1 - cast(0.9999, &randSt))));
 
-            plopCrater(&pam, terrain, cx, cy, radius, cmdline.verbose);
+            plopCrater(&pam, terrain, cx, cy, radius,
+                       cmdline.verbose, &randSt);
 
             if (((l + 1) % 100000) == 0)
                 pm_message("%u craters generated of %u (%u%% done)",
@@ -485,6 +494,8 @@ genCraters(struct CmdlineInfo const cmdline) {
         }
     }
 
+    pm_randterm(&randSt);
+
     pnm_writepam(&pam, terrain);
 
     pnm_freepamarray(terrain, &pam);
@@ -503,8 +514,6 @@ main(int argc, const char ** argv) {
 
     parseCommandLine(argc, argv, &cmdline);
 
-    srand(cmdline.randomseedSpec ? cmdline.randomseed : pm_randseed());
-
     genCraters(cmdline);
 
     return 0;
diff --git a/generator/pamgauss.c b/generator/pamgauss.c
index 9656708b..5553e757 100644
--- a/generator/pamgauss.c
+++ b/generator/pamgauss.c
@@ -31,11 +31,11 @@ static void
 parseCommandLine(int argc, const char ** argv,
                  struct CmdlineInfo * const cmdlineP) {
 /*----------------------------------------------------------------------------
-  Convert program invocation arguments (argc,argv) into a format the 
+  Convert program invocation arguments (argc,argv) into a format the
   program can use easily, struct cmdlineInfo.  Validate arguments along
   the way and exit program with message if invalid.
 
-  Note that some string information we return as *cmdlineP is in the storage 
+  Note that some string information we return as *cmdlineP is in the storage
   argv[] points to.
 -----------------------------------------------------------------------------*/
     optEntry * option_def;
@@ -49,11 +49,11 @@ parseCommandLine(int argc, const char ** argv,
     MALLOCARRAY_NOFAIL(option_def, 100);
 
     option_def_index = 0;   /* incremented by OPTENTRY */
-    OPTENT3(0,   "tupletype",  OPT_STRING, &cmdlineP->tupletype, 
+    OPTENT3(0,   "tupletype",  OPT_STRING, &cmdlineP->tupletype,
             &tupletypeSpec,      0);
-    OPTENT3(0,   "maxval",     OPT_UINT,   &cmdlineP->maxval, 
+    OPTENT3(0,   "maxval",     OPT_UINT,   &cmdlineP->maxval,
             &maxvalSpec,         0);
-    OPTENT3(0,   "sigma",      OPT_FLOAT,  &cmdlineP->sigma, 
+    OPTENT3(0,   "sigma",      OPT_FLOAT,  &cmdlineP->sigma,
             &sigmaSpec,          0);
     OPTENT3(0,   "maximize",   OPT_FLAG,   NULL,
             &cmdlineP->maximize, 0);
@@ -75,12 +75,12 @@ parseCommandLine(int argc, const char ** argv,
             pm_error("The tuple type you specified is too long.  "
                      "Maximum %u characters.",
                      (unsigned)sizeof(pam.tuple_type)-1);
-    }        
+    }
 
     if (!sigmaSpec)
         pm_error("You must specify the -sigma option.");
     else if (cmdlineP->sigma <= 0.0)
-        pm_error("-sigma must be positive.  You specified %f", 
+        pm_error("-sigma must be positive.  You specified %f",
                  cmdlineP->sigma);
 
     if (!maxvalSpec)
@@ -91,7 +91,7 @@ parseCommandLine(int argc, const char ** argv,
                      "Maximum is %u", cmdlineP->maxval, PNM_OVERALLMAXVAL);
         if (cmdlineP->maxval < 1)
             pm_error("-maxval must be at least 1");
-    }    
+    }
 
     if (oversampleSpec) {
         if (cmdlineP->oversample < 1)
@@ -106,8 +106,8 @@ parseCommandLine(int argc, const char ** argv,
         pm_error("Only two arguments allowed: width and height.  "
                  "You specified %d", argc-1);
     else {
-        cmdlineP->width = atoi(argv[1]);
-        cmdlineP->height = atoi(argv[2]);
+        cmdlineP->width  = pm_parse_width(argv[1]);
+        cmdlineP->height = pm_parse_height(argv[2]);
         if (cmdlineP->width <= 0)
             pm_error("width argument must be a positive number.  You "
                      "specified '%s'", argv[1]);
@@ -164,7 +164,7 @@ pixelValue(unsigned int const width,
   The gaussian function has standard deviation 'sigma' and amplitude 1.
 -----------------------------------------------------------------------------*/
     double const offset = 1.0 / (subpixDivision * 2);
-    double const y0     = (double)row + offset; 
+    double const y0     = (double)row + offset;
     double const x0     = (double)col + offset;
 
     double const subpixSize = 1.0 / subpixDivision;
@@ -236,7 +236,7 @@ maximumKernelValue(double **    const kernel,
        we know the maximum value is at the center.
     */
     return kernel[height/2][width/2];
-}        
+}
 
 
 
@@ -250,13 +250,13 @@ totalKernelValue(double **    const kernel,
 
     for (row = 0, total = 0.0; row < height; ++row) {
         unsigned int col;
-        
+
         for (col = 0; col < width; ++col)
             total += kernel[row][col];
     }
 
     return total;
-}        
+}
 
 
 
@@ -311,7 +311,7 @@ writePam(double **    const kernel,
         unsigned int col;
         for (col = 0; col < pam.width; ++col) {
             tuplerown[col][0] = kernel[row][col] / normalizer;
-        
+
             assert(tuplerown[col][0] <= 1.0);
         }
         pnm_writepamrown(&pam, tuplerown);
@@ -319,7 +319,7 @@ writePam(double **    const kernel,
 
     pnm_freepamrown(tuplerown);
 }
-    
+
 
 
 int
@@ -327,7 +327,7 @@ main(int argc, const char **argv) {
     struct CmdlineInfo cmdline;
     double ** kernel;
     double    normalizer;
-    
+
     pm_proginit(&argc, argv);
 
     parseCommandLine(argc, argv, &cmdline);
@@ -344,6 +344,9 @@ main(int argc, const char **argv) {
              normalizer, stdout);
 
     pm_freearray2((void **)kernel);
-    
+
     return 0;
 }
+
+
+
diff --git a/generator/pamgradient.c b/generator/pamgradient.c
index 526efdae..bbbaad30 100644
--- a/generator/pamgradient.c
+++ b/generator/pamgradient.c
@@ -53,8 +53,8 @@ parseCommandLine(int argc, const char **argv,
     else {
         if (cmdlineP->maxval > PAM_OVERALL_MAXVAL)
             pm_error("The value you specified for -maxval (%u) is too big.  "
-                     "Max allowed is %u", cmdlineP->maxval,
-                     PAM_OVERALL_MAXVAL);
+                     "Max allowed is %lu", cmdlineP->maxval,
+                     (unsigned long int) PAM_OVERALL_MAXVAL);
 
         if (cmdlineP->maxval < 1)
             pm_error("You cannot specify 0 for -maxval");
@@ -68,8 +68,8 @@ parseCommandLine(int argc, const char **argv,
         cmdlineP->colorTopRight    = pnm_parsecolor(argv[2], cmdlineP->maxval);
         cmdlineP->colorBottomLeft  = pnm_parsecolor(argv[3], cmdlineP->maxval);
         cmdlineP->colorBottomRight = pnm_parsecolor(argv[4], cmdlineP->maxval);
-        cmdlineP->cols = atoi(argv[5]);
-        cmdlineP->rows = atoi(argv[6]);
+        cmdlineP->cols = pm_parse_width(argv[5]);
+        cmdlineP->rows = pm_parse_height(argv[6]);
         if (cmdlineP->cols <= 0)
             pm_error("width argument must be a positive number.  You "
                      "specified '%s'", argv[5]);
@@ -211,3 +211,6 @@ main(int argc, const char *argv[]) {
 
     return 0;
 }
+
+
+
diff --git a/generator/pamseq.c b/generator/pamseq.c
index 1af5252a..71905a67 100644
--- a/generator/pamseq.c
+++ b/generator/pamseq.c
@@ -1,3 +1,4 @@
+#include <assert.h>
 #include <string.h>
 #include <unistd.h>
 #include <stdlib.h>
@@ -6,49 +7,173 @@
 #include "pm_c_util.h"
 #include "pam.h"
 #include "shhopt.h"
+#include "mallocvar.h"
+#include "nstring.h"
 
 
 
-struct cmdlineInfo {
+struct CmdlineInfo {
     /* All the information the user supplied in the command line,
        in a form easy for the program to use.
     */
     unsigned int depth;
     sample maxval;
     const char * tupletype;
+    sample * min;   /* array of size 'depth' */
+    sample * max;   /* array of size 'depth' */
+    sample * step;  /* array of size 'depth' */
 };
 
 
 
 static void
-parseCommandLine(int argc, char ** argv,
-                 struct cmdlineInfo * const cmdlineP) {
+destroyCmdline(struct CmdlineInfo * const cmdlineP)  {
+
+    if (cmdlineP->min)
+        free(cmdlineP->min);
+    if (cmdlineP->max)
+        free(cmdlineP->max);
+    if (cmdlineP->step)
+        free(cmdlineP->step);
+}
+
+
+
+static unsigned int
+entryCt(char ** const stringList) {
+
+    unsigned int i;
+
+    for (i = 0; stringList[i]; ++i) {}
+
+    return i;
+}
+
+
+
+static void
+parseOptList(bool         const isSpec,
+             char **      const stringList,
+             unsigned int const depth,
+             sample       const maxval,
+             const char * const optNm,
+             sample **    const sampleListP) {
+
+    if (!isSpec)
+        *sampleListP = NULL;
+    else {
+        unsigned int i;
+        sample * sampleList;
+        const char * memberError;
+
+        if (entryCt(stringList) != depth) {
+            pm_error("Wrong number of values for -%s: %u.  Need %u",
+                     optNm, entryCt(stringList), depth);
+        }
+
+        MALLOCARRAY(sampleList, depth);
+
+        for (i = 0, memberError = NULL; i < depth && !memberError; ++i) {
+            char * endPtr;
+            long const n = strtol(stringList[i], &endPtr, 10);
+
+            if (strlen(stringList[i]) == 0)
+                pm_asprintf(&memberError, "is null string");
+            else if (*endPtr != '\0')
+                pm_asprintf(&memberError,
+                            "contains non-numeric character '%c'",
+                            *endPtr);
+            else if (n < 0)
+                pm_asprintf(&memberError, "is negative");
+            else if (n > maxval)
+                pm_asprintf(&memberError, "is greater than maxval %lu",
+                            maxval);
+            else
+                sampleList[i] = n;
+        }
+        if (memberError) {
+            free(sampleList);
+            pm_errormsg("Value in -%s %s", optNm, memberError);
+            pm_longjmp();
+        }
+        *sampleListP = sampleList;
+    }
+}
+
+
+
+static void
+validateMinIsAtMostMax(sample *     const min,
+                       sample *     const max,
+                       unsigned int const depth) {
+
+    unsigned int plane;
+
+    for (plane = 0; plane < depth; ++plane) {
+        if (min[plane] > max[plane])
+            pm_error("-min for plane %u (%lu) is greater than -max (%lu)",
+                     plane, min[plane], max[plane]);
+    }
+}
+
+
+
+static void
+validateStepIsPositive(sample *     const step,
+                       unsigned int const depth) {
+
+    unsigned int plane;
+
+    for (plane = 0; plane < depth; ++plane) {
+        if (step[plane] <= 0)
+            pm_error("-step for plane %u (%lu) is not positive",
+                     plane, step[plane]);
+    }
+}
+
+
+
+static void
+parseCommandLine(int argc, const char ** argv,
+                 struct CmdlineInfo * const cmdlineP) {
 /*----------------------------------------------------------------------------
-  Convert program invocation arguments (argc,argv) into a format the 
+  Convert program invocation arguments (argc,argv) into a format the
   program can use easily, struct cmdlineInfo.  Validate arguments along
   the way and exit program with message if invalid.
 
-  Note that some string information we return as *cmdlineP is in the storage 
+  Note that some string information we return as *cmdlineP is in the storage
   argv[] points to.
 -----------------------------------------------------------------------------*/
-    optEntry *option_def = malloc(100*sizeof(optEntry));
-        /* Instructions to OptParseOptions2 on how to parse our options.
-         */
+    optEntry *option_def;
     optStruct3 opt;
-
+        /* Instructions to pm_optParseOptions3 on how to parse our options. */
+
+    unsigned int maxSpec;
+    char ** max;
+    unsigned int minSpec;
+    char ** min;
+    unsigned int stepSpec;
+    char ** step;
     unsigned int tupletypeSpec;
     unsigned int option_def_index;
 
+    MALLOCARRAY_NOFAIL(option_def, 100);
+
     option_def_index = 0;   /* incremented by OPTENT3 */
-    OPTENT3(0,   "tupletype",  OPT_STRING, &cmdlineP->tupletype, 
+    OPTENT3(0,   "tupletype",  OPT_STRING, &cmdlineP->tupletype,
             &tupletypeSpec,     0);
-
+    OPTENT3(0,   "min",         OPT_STRINGLIST, &min,
+            &minSpec,           0);
+    OPTENT3(0,   "max",         OPT_STRINGLIST, &max,
+            &maxSpec,           0);
+    OPTENT3(0,   "step",        OPT_STRINGLIST, &step,
+            &stepSpec,          0);
 
     opt.opt_table = option_def;
     opt.short_allowed = FALSE;  /* We have no short (old-fashioned) options */
     opt.allowNegNum = FALSE;  /* We have no parms that are negative numbers */
 
-    pm_optParseOptions3(&argc, argv, opt, sizeof(opt), 0);
+    pm_optParseOptions3(&argc, (char **)argv, opt, sizeof(opt), 0);
         /* Uses and sets argc, argv, and some of *cmdlineP and others. */
 
     if (!tupletypeSpec)
@@ -57,9 +182,9 @@ parseCommandLine(int argc, char ** argv,
         struct pam pam;
         if (strlen(cmdlineP->tupletype)+1 > sizeof(pam.tuple_type))
             pm_error("The tuple type you specified is too long.  "
-                     "Maximum %u characters.", 
+                     "Maximum %u characters.",
                      (unsigned)sizeof(pam.tuple_type)-1);
-    }        
+    }
 
     if (argc-1 < 2)
         pm_error("Need two arguments: depth and maxval.");
@@ -67,83 +192,175 @@ parseCommandLine(int argc, char ** argv,
         pm_error("Only two argumeents allowed: depth and maxval.  "
                  "You specified %d", argc-1);
     else {
-        cmdlineP->depth = atoi(argv[1]);
-        cmdlineP->maxval = atoi(argv[2]);
-        if (cmdlineP->depth <= 0)
+        const char * error;
+        unsigned int depth, maxval;
+
+        pm_string_to_uint(argv[1], &depth, &error);
+        if (error) {
+            pm_error("'%s' is invalid as an image depth.  %s", argv[1], error);
+            pm_strfree(error);
+        }
+        else if (depth <= 0)
             pm_error("depth argument must be a positive number.  You "
                      "specified '%s'", argv[1]);
-        if (cmdlineP->maxval <= 0)
-            pm_error("maxval argument must be a positive number.  You "
-                     "specified '%s'", argv[2]);
-        if (cmdlineP->maxval > PNM_OVERALLMAXVAL)
+        else
+            cmdlineP->depth = depth;
+
+        maxval = pm_parse_maxval(argv[2]);
+
+        if (maxval > PAM_OVERALL_MAXVAL)
             pm_error("The maxval you specified (%u) is too big.  "
-                     "Maximum is %u", (unsigned int) cmdlineP->maxval, 
-                     PNM_OVERALLMAXVAL);
-        if (pm_maxvaltobits(cmdlineP->maxval) + 
+                     "Maximum is %lu", maxval, PAM_OVERALL_MAXVAL);
+        else
+            cmdlineP->maxval = maxval;
+
+
+        if (pm_maxvaltobits(cmdlineP->maxval) +
             pm_maxvaltobits(cmdlineP->depth-1) > sizeof(unsigned int)*8)
             pm_error("The maxval (%u) and depth (%u) you specified result "
                      "in a larger number of tuples than this program can "
-                     "handle (roughly %u)", 
+                     "handle (roughly %u)",
                      (unsigned int) cmdlineP->maxval, cmdlineP->depth,
                      (unsigned int) -1);
     }
+    parseOptList(minSpec, min,  cmdlineP->depth, cmdlineP->maxval, "min",
+                 &cmdlineP->min);
+    parseOptList(maxSpec, max,  cmdlineP->depth, cmdlineP->maxval, "max",
+                 &cmdlineP->max);
+    parseOptList(stepSpec, step, cmdlineP->depth, cmdlineP->maxval, "step",
+                 &cmdlineP->step);
+
+    if (cmdlineP->min && cmdlineP->max)
+        validateMinIsAtMostMax(cmdlineP->min, cmdlineP->max, cmdlineP->depth);
+
+    if (cmdlineP->step)
+        validateStepIsPositive(cmdlineP->step, cmdlineP->depth);
+
+    if (minSpec)
+        free(min);
+    if (maxSpec)
+        free(max);
+    if (stepSpec)
+        free(step);
 }
 
 
 
-static unsigned int
-powint(unsigned int base, unsigned int exponent) {
+static void
+computeMinMaxStep(unsigned int   const depth,
+                  sample         const maxval,
+                  const sample * const min,
+                  const sample * const max,
+                  const sample * const step,
+                  sample **      const minP,
+                  sample **      const maxP,
+                  sample **      const stepP) {
+
+    unsigned int plane;
+
+    MALLOCARRAY(*minP,  depth);
+    MALLOCARRAY(*maxP,  depth);
+    MALLOCARRAY(*stepP, depth);
+
+    for (plane = 0; plane < depth; ++plane) {
+        (*minP)[plane]  = min  ? min[plane]  : 0;
+        (*maxP)[plane]  = max  ? max[plane]  : maxval;
+        (*stepP)[plane] = step ? step[plane] : 1;
+    }
+}
+
+
+
+static int
+imageWidth(unsigned int   const depth,
+           const sample * const min,
+           const sample * const max,
+           const sample * const step) {
 /*----------------------------------------------------------------------------
-   This is standard pow(), but for integers and optimized for small
-   exponents.
+   The width of the output image (i.e. the number of pixels in the image),
+   given that the minimum and maximum sample values in Plane P are min[P] and
+   max[P] and the samples step by step[P].
+
+   E.g. in the example case of min 0, max 4, and step 1 everywhere, with
+   depth 2,  We return 5*5 = 25.
 -----------------------------------------------------------------------------*/
-    unsigned int result;
-    unsigned int i;
+    unsigned int product;
+    unsigned int plane;
 
-    result = 1;  /* initial value */
-    for (i = 0; i < exponent; ++i) 
-        result *= base;
+    for (plane = 0, product=1; plane < depth; ++plane) {
+        assert(max[plane] >= min[plane]);
 
-    return(result);
+        unsigned int const valueCtThisPlane =
+            ROUNDUP((max[plane] - min[plane] + 1), step[plane])/step[plane];
+
+        if (INT_MAX / valueCtThisPlane < product)
+            pm_error("Uncomputably large number of pixels (greater than %u)",
+                     INT_MAX);
+
+        product *= valueCtThisPlane;
+    }
+    assert(product < INT_MAX);
+
+    return product;
 }
 
 
+
 static void
-permuteHigherPlanes(struct pam const pam, int const nextplane, 
-                    tuple * const tuplerow, int * const colP, 
-                    tuple const lowerPlanes) {
+permuteHigherPlanes(unsigned int   const depth,
+                    const sample * const min,
+                    const sample * const max,
+                    const sample * const step,
+                    unsigned int   const nextPlane,
+                    tuple *        const tuplerow,
+                    int *          const colP,
+                    tuple          const lowerPlanes) {
 /*----------------------------------------------------------------------------
    Create all the possible permutations of tuples whose lower-numbered planes
    contain the values from 'lowerPlanes'.  I.e. vary the higher-numbered
-   planes between zero and maxval.
+   planes according to min[], max[], and step[].
 
    Write them sequentially into *tuplerow, starting at *colP.  Adjust
    *colP to next the column after the ones we write.
 
-   lower-numbered means with plane numbers less than 'nextplane'.
+   lower-numbered means with plane numbers less than 'nextPlane'.
 
    We modify 'lowerPlanes' in the higher planes to undefined values.
 -----------------------------------------------------------------------------*/
-    if (nextplane == pam.depth - 1) {
+    if (nextPlane == depth - 1) {
+        /* lowerPlanes[] contains values for all the planes except the
+           highest, so we just vary the highest plane and combine that
+           with lowerPlanes[] and output that to tuplerow[].
+        */
         sample value;
-        for (value = 0; value <= pam.maxval; ++value) {
+
+        for (value = min[nextPlane];
+             value <= max[nextPlane];
+             value += step[nextPlane]) {
+
             unsigned int plane;
-            for (plane = 0; plane < nextplane; ++plane)
+
+            for (plane = 0; plane < nextPlane; ++plane)
                 tuplerow[*colP][plane] = lowerPlanes[plane];
-            tuplerow[*colP][nextplane] = value;
+
+            tuplerow[*colP][nextPlane] = value;
+
             ++(*colP);
         }
     } else {
         sample value;
 
-        for (value = 0; value <= pam.maxval; ++value) {
+        for (value = min[nextPlane];
+             value <= max[nextPlane];
+             value += step[nextPlane]) {
             /* We do something sleazy here and use Caller's lowerPlanes[]
                variable as a local variable, modifying it in the higher
                plane positions.  That's just for speed.
             */
-            lowerPlanes[nextplane] = value;
+            lowerPlanes[nextPlane] = value;
 
-            permuteHigherPlanes(pam, nextplane+1, tuplerow, colP, lowerPlanes);
+            permuteHigherPlanes(depth, min, max, step,
+                                nextPlane+1, tuplerow, colP, lowerPlanes);
         }
     }
 }
@@ -151,52 +368,65 @@ permuteHigherPlanes(struct pam const pam, int const nextplane,
 
 
 int
-main(int argc, char **argv) {
+main(int argc, const char **argv) {
 
-    struct cmdlineInfo cmdline;
+    struct CmdlineInfo cmdline;
     struct pam pam;
     int col;
     tuple lowerPlanes;
         /* This is working storage passed to permuteHigherPlanes(),
            which we call.  Note that because we always pass zero as the
-           "planes" argument to permuteHigherPlanes(), none of the 
-           "lower planes" value is defined as an input to 
+           "planes" argument to permuteHigherPlanes(), none of the
+           "lower planes" value is defined as an input to
            permuteHigherPlanes().
         */
     tuple * tuplerow;
-    
-    pnm_init(&argc, argv);
-   
+    sample * min;   /* malloc'ed array */
+    sample * max;   /* malloc'ed array */
+    sample * step;  /* malloc'ed array */
+
+    pm_proginit(&argc, argv);
+
     parseCommandLine(argc, argv, &cmdline);
 
+    computeMinMaxStep(cmdline.depth, cmdline.maxval,
+                      cmdline.min, cmdline.max, cmdline.step,
+                      &min, &max, &step);
+
     pam.size = sizeof(pam);
     pam.len = PAM_STRUCT_SIZE(tuple_type);
     pam.file = stdout;
     pam.format = PAM_FORMAT;
     pam.plainformat = 0;
-    pam.width = powint(cmdline.maxval+1, cmdline.depth);
+    pam.width = imageWidth(cmdline.depth, min, max, step);
     pam.height = 1;
     pam.depth = cmdline.depth;
     pam.maxval = cmdline.maxval;
     strcpy(pam.tuple_type, cmdline.tupletype);
 
     pnm_writepaminit(&pam);
-   
+
     tuplerow = pnm_allocpamrow(&pam);
 
     lowerPlanes = pnm_allocpamtuple(&pam);
 
     col = 0;
 
-    permuteHigherPlanes(pam, 0, tuplerow, &col, lowerPlanes);
+    permuteHigherPlanes(pam.depth, min, max, step,
+                        0, tuplerow, &col, lowerPlanes);
 
     if (col != pam.width)
         pm_error("INTERNAL ERROR: Wrote %d columns; should have written %d.",
                  col, pam.width);
 
     pnm_writepamrow(&pam, tuplerow);
-    
+
     pnm_freepamrow(tuplerow);
 
+    destroyCmdline(&cmdline);
+
     return 0;
 }
+
+
+
diff --git a/generator/pamshadedrelief.c b/generator/pamshadedrelief.c
index 35d1e3e0..80b34a0d 100644
--- a/generator/pamshadedrelief.c
+++ b/generator/pamshadedrelief.c
@@ -8,7 +8,7 @@
 
   The input array is a one-channel PAM image.  The sample values are
   elevations of terrain.
-  
+
   This is derived from John Walker's 'pgmcrater' which not only does this
   shading, but first generates a terrain map of fractal craters on which to
   run it.
@@ -95,7 +95,7 @@ parseCommandLine(int argc, const char ** const argv,
     if (cmdlineP->gamma <= 0.0)
         pm_error("gamma correction must be greater than 0");
 
-    if (argc-1 == 0) 
+    if (argc-1 == 0)
         cmdlineP->inputFileName = "-";
     else if (argc-1 != 1)
         pm_error("Program takes zero or one argument (filename).  You "
@@ -201,7 +201,7 @@ writeShadedRelief(struct pam * const terrainPamP,
         }
         {
             /* Wrap around to determine shade of pixel on right edge */
-            int const slope = 
+            int const slope =
                 terrain[row][0][0] - terrain[row][outpam.width-1][0];
             outrow[outpam.width - 1][0] =
                 brightnessOfSlope(slope, slopeGrayMap);
@@ -222,8 +222,8 @@ readTerrain(FILE *       const ifP,
 
     *tuplesP = pnm_readpam(ifP, pamP, PAM_STRUCT_SIZE(tuple_type));
 }
-            
-            
+
+
 
 int
 main(int argc, const char ** argv) {
@@ -243,7 +243,7 @@ main(int argc, const char ** argv) {
     readTerrain(ifP, &terrainPam, &terrain);
 
     writeShadedRelief(&terrainPam, terrain, cmdline.gamma, stdout);
-    
+
     return 0;
 }
 
diff --git a/generator/pamstereogram.c b/generator/pamstereogram.c
index 9b861b51..4bc1ea92 100644
--- a/generator/pamstereogram.c
+++ b/generator/pamstereogram.c
@@ -16,7 +16,7 @@
  *
  * ----------------------------------------------------------------------
  *
- * Copyright (C) 2006-2020 Scott Pakin <scott+pbm@pakin.org>
+ * Copyright (C) 2006-2021 Scott Pakin <scott+pbm@pakin.org>
  *
  * All rights reserved.
  *
@@ -58,11 +58,12 @@
 #include "pm_c_util.h"
 #include "mallocvar.h"
 #include "nstring.h"
+#include "rand.h"
 #include "shhopt.h"
 #include "pam.h"
 
 
-enum outputType {OUTPUT_BW, OUTPUT_GRAYSCALE, OUTPUT_COLOR};
+enum OutputType {OUTPUT_BW, OUTPUT_GRAYSCALE, OUTPUT_COLOR};
 
 /* ---------------------------------------------------------------------- */
 
@@ -85,13 +86,14 @@ struct cmdlineInfo {
     unsigned int magnifypat;     /* -magnifypat option */
     int xshift;                  /* -xshift option */
     int yshift;                  /* -yshift option */
+    int yfillshift;              /* -yfillshift option */
     const char * patfile;        /* -patfile option.  Null if none */
     const char * texfile;        /* -texfile option.  Null if none */
     const char * bgcolor;        /* -bgcolor option */
     unsigned int smoothing;      /* -smoothing option */
     unsigned int randomseed;     /* -randomseed option */
     unsigned int randomseedSpec; /* -randomseed option count */
-    enum outputType outputType;  /* Type of output file */
+    enum OutputType outputType;  /* Type of output file */
     unsigned int xbegin;         /* -xbegin option */
     unsigned int xbeginSpec;     /* -xbegin option count */
     unsigned int tileable;       /* -tileable option */
@@ -153,8 +155,8 @@ parseCommandLine(int                  argc,
     unsigned int option_def_index;
 
     unsigned int patfileSpec, texfileSpec, dpiSpec, eyesepSpec, depthSpec,
-        guidesizeSpec, magnifypatSpec, xshiftSpec, yshiftSpec,
-      bgcolorSpec, smoothingSpec, planesSpec;
+        guidesizeSpec, magnifypatSpec, xshiftSpec, yshiftSpec, yfillshiftSpec,
+        bgcolorSpec, smoothingSpec, planesSpec;
 
     unsigned int blackandwhite, grayscale, color;
     const char ** planes;
@@ -194,6 +196,8 @@ parseCommandLine(int                  argc,
             &xshiftSpec,              0);
     OPTENT3(0, "yshift",          OPT_INT,    &cmdlineP->yshift,
             &yshiftSpec,              0);
+    OPTENT3(0, "yfillshift",      OPT_INT,    &cmdlineP->yfillshift,
+            &yfillshiftSpec,          0);
     OPTENT3(0, "patfile",         OPT_STRING, &cmdlineP->patfile,
             &patfileSpec,             0);
     OPTENT3(0, "texfile",         OPT_STRING, &cmdlineP->texfile,
@@ -282,9 +286,10 @@ parseCommandLine(int                  argc,
 
     if (!xshiftSpec)
         cmdlineP->xshift = 0;
-
     if (!yshiftSpec)
         cmdlineP->yshift = 0;
+    if (!yfillshiftSpec)
+        cmdlineP->yfillshift = 0;
 
     if (xshiftSpec && !cmdlineP->patfile)
         pm_error("-xshift is valid only with -patfile");
@@ -402,11 +407,12 @@ typedef struct outGenerator {
 
 
 
-struct randomState {
+struct RandomState {
     /* The state of a randomColor generator. */
     unsigned int magnifypat;
     tuple *      currentRow;
     unsigned int prevy;
+    struct pm_randSt * randStP;
 };
 
 
@@ -421,7 +427,7 @@ randomColor(outGenerator * const outGenP,
 /*----------------------------------------------------------------------------
    Return a random RGB value.
 -----------------------------------------------------------------------------*/
-    struct randomState * const stateP = outGenP->stateP;
+    struct RandomState * const stateP = outGenP->stateP;
 
     /* Every time we start a new row, we select a new sequence of random
        colors.
@@ -436,7 +442,7 @@ randomColor(outGenerator * const outGenP,
             unsigned int plane;
 
             for (plane = 0; plane < outGenP->pam.depth; ++plane) {
-                unsigned int const randval = rand();
+                unsigned int const randval = pm_rand(stateP->randStP);
                 thisTuple[plane] = randval % modulus;
             }
         }
@@ -456,9 +462,13 @@ static outGenStateTerm termRandomColor;
 static void
 termRandomColor(outGenerator * const outGenP) {
 
-    struct randomState * const stateP = outGenP->stateP;
+    struct RandomState * const stateP = outGenP->stateP;
 
     pnm_freepamrow(stateP->currentRow);
+
+    pm_randterm(stateP->randStP);
+
+    free(stateP->randStP);
 }
 
 
@@ -468,7 +478,7 @@ initRandomColor(outGenerator *     const outGenP,
                 const struct pam * const inPamP,
                 struct cmdlineInfo const cmdline) {
 
-    struct randomState * stateP;
+    struct RandomState * stateP;
 
     outGenP->pam.format      = PAM_FORMAT;
     outGenP->pam.plainformat = 0;
@@ -499,6 +509,10 @@ initRandomColor(outGenerator *     const outGenP,
     stateP->magnifypat = cmdline.magnifypat;
     stateP->prevy      = (unsigned int)(-cmdline.magnifypat);
 
+    MALLOCVAR_NOFAIL(stateP->randStP);
+    pm_randinit(stateP->randStP);
+    pm_srand2(stateP->randStP, cmdline.randomseedSpec, cmdline.randomseed);
+
     outGenP->stateP         = stateP;
     outGenP->getTuple       = &randomColor;
     outGenP->terminateState = &termRandomColor;
@@ -1178,6 +1192,8 @@ static void
 makeImageRow(outGenerator *       const outGenP,
              unsigned int         const row,
              unsigned int         const xbegin,
+             unsigned int         const farWidth,
+             int                  const yfillshift,
              const unsigned int * const sameL,
              const unsigned int * const sameR,
              const tuple *        const outRow) {
@@ -1200,42 +1216,59 @@ makeImageRow(outGenerator *       const outGenP,
 
   sameL[N] > N is not allowed.
 -----------------------------------------------------------------------------*/
-    int col;
-    int lastLinked;
+    unsigned int const width  = outGenP->pam.width;
+    unsigned int const height = outGenP->pam.height;
 
-    for (col = (int)xbegin, lastLinked = INT_MIN;
-         col < outGenP->pam.width;
-         ++col) {
+    unsigned int col;
+    bool colHasBeenLinked;
+    unsigned int lastLinkedCol;
+        /* Last column to have been linked; meaningless if 'colHasBeenLinked'
+           is false.
+        */
+
+    for (col = xbegin, colHasBeenLinked = false; col < width; ++col) {
 
         tuple newtuple;
 
-        if (sameL[col] == col || sameL[col] < (int)xbegin) {
-            if (lastLinked == col - 1)
-                newtuple = outRow[col - 1];
-            else
-                newtuple = outGenP->getTuple(outGenP, col, row);
+        if (sameL[col] == col || sameL[col] < xbegin) {
+            if (colHasBeenLinked && lastLinkedCol == col - 1)
+                newtuple = outRow[lastLinkedCol];
+            else {
+                if (col < xbegin + farWidth)
+                    newtuple = outGenP->getTuple(outGenP, col, row);
+                else
+                    newtuple = outGenP->getTuple(
+                        outGenP, col, (row + height - yfillshift) % height);
+            }
         } else {
-          newtuple = outRow[sameL[col]];
-          lastLinked = col;
-              /* Keep track of the last pixel to be constrained. */
+            newtuple = outRow[sameL[col]];
+            colHasBeenLinked = true;
+            lastLinkedCol = col;
+                /* Keep track of the last pixel to be constrained. */
         }
         pnm_assigntuple(&outGenP->pam, outRow[col], newtuple);
     }
 
-    for (col = (int)xbegin - 1, lastLinked = INT_MIN; col >= 0; --col) {
+    for (col = xbegin, colHasBeenLinked = false; col > 0; --col) {
         tuple newtuple;
 
-        if (sameR[col] == col) {
-            if (lastLinked == col + 1)
-                newtuple = outRow[col + 1];
-            else
-                newtuple = outGenP->getTuple(outGenP, col, row);
+        if (sameR[col-1] == col-1) {
+            if (colHasBeenLinked && lastLinkedCol == col)
+                newtuple = outRow[lastLinkedCol];
+            else {
+                if (col > xbegin - farWidth)
+                    newtuple = outGenP->getTuple(outGenP, col-1, row);
+                else
+                    newtuple = outGenP->getTuple(
+                        outGenP, col-1, (row + height - yfillshift) % height);
+            }
         } else {
-            newtuple = outRow[sameR[col]];
-            lastLinked = col;
+            newtuple = outRow[sameR[col-1]];
+            colHasBeenLinked = true;
+            lastLinkedCol = col - 1;
                 /* Keep track of the last pixel to be constrained. */
         }
-        pnm_assigntuple(&outGenP->pam, outRow[col], newtuple);
+        pnm_assigntuple(&outGenP->pam, outRow[col-1], newtuple);
     }
 }
 
@@ -1261,6 +1294,8 @@ makeOneImageRow(unsigned int         const row,
                 const unsigned int * const sameR,
                 unsigned int *       const colNumBuffer,
                 unsigned int         const xbegin,
+                unsigned int         const farWidth,
+                int                  const yfillshift,
                 tuple *              const rowBuffer,
                 tuple *              const outRow) {
 
@@ -1272,7 +1307,7 @@ makeOneImageRow(unsigned int         const row,
                             rowBuffer, outRow);
         else
             makeImageRow(outputGeneratorP, row,
-                         xbegin, sameL, sameR, outRow);
+                         xbegin, farWidth, yfillshift, sameL, sameR, outRow);
     }
 }
 
@@ -1284,6 +1319,8 @@ constructRowTileable(const struct pam *   const inPamP,
                      bool                 const makeMask,
                      const unsigned int * const sameL,
                      const unsigned int * const sameR,
+                     unsigned int         const farWidth,
+                     int                  const yfillshift,
                      unsigned int         const row,
                      tuple *              const outRow,
                      tuple *              const outRowBuf1,
@@ -1300,10 +1337,10 @@ constructRowTileable(const struct pam *   const inPamP,
        xbegin=maximum case.
     */
     makeOneImageRow(row, outputGeneratorP, makeMask, sameL, sameR, colNumBuf2,
-                    0, outRowBuf2, outRow);
+                    farWidth, yfillshift, 0, outRowBuf2, outRow);
 
     makeOneImageRow(row, outputGeneratorP, makeMask, sameL, sameR, colNumBuf2,
-                    inPamP->width - 1, outRowBuf2, outRowMax);
+                    farWidth, yfillshift, inPamP->width - 1, outRowBuf2, outRowMax);
 
     for (col = 0; col < inPamP->width; ++col) {
         unsigned int plane;
@@ -1362,6 +1399,8 @@ doRow(const struct pam * const inPamP,
       unsigned int       const magnifypat,
       unsigned int       const smoothing,
       unsigned int       const xbegin,
+      unsigned int       const farWidth,
+      int                const yfillshift,
       unsigned int       const row,
       tuple *            const inRow,
       tuple *            const outRowBuf0,
@@ -1397,13 +1436,16 @@ doRow(const struct pam * const inPamP,
     /* Construct a single row. */
     if (tileable) {
         constructRowTileable(inPamP, outputGeneratorP, makeMask,
-                             sameL, sameR, row, outRow,
+                             sameL, sameR,
+                             farWidth, yfillshift,
+                             row, outRow,
                              outRowBuf1, outRowBuf2, colNumBuf2);
 
     } else {
         makeOneImageRow(row, outputGeneratorP, makeMask,
                         sameL, sameR, colNumBuf2,
-                        xbegin, outRowBuf1, outRow);
+                        xbegin, farWidth, yfillshift,
+                        outRowBuf1, outRow);
     }
 
     /* Write the resulting row. */
@@ -1423,7 +1465,8 @@ makeImageRows(const struct pam * const inPamP,
               bool               const tileable,
               unsigned int       const magnifypat,
               unsigned int       const smoothing,
-              unsigned int       const xbegin) {
+              unsigned int       const xbegin,
+              int                const yfillshift) {
 
     tuple * inRow;    /* Buffer for use in reading from the height-map image */
     tuple * outRowBuf0; /* Buffer for use in generating output rows */
@@ -1433,6 +1476,8 @@ makeImageRows(const struct pam * const inPamP,
     unsigned int * colNumBuf1;
     unsigned int * colNumBuf2;
     unsigned int row;      /* Current row in the input and output files */
+    unsigned int pixelEyesep;
+    unsigned int farWidth;
 
     inRow = pnm_allocpamrow(inPamP);
     outRowBuf0 = pnm_allocpamrow(&outputGeneratorP->pam);
@@ -1451,10 +1496,13 @@ makeImageRows(const struct pam * const inPamP,
         pm_error("Unable to allocate space for %u column buffer",
                  inPamP->width);
 
+    pixelEyesep = ROUNDU(eyesep * dpi);
+    farWidth = pixelEyesep/(magnifypat * 2);
+
     for (row = 0; row < inPamP->height; ++row) {
         doRow(inPamP, outputGeneratorP,  depthOfField, eyesep, dpi,
               crossEyed, makeMask, tileable, magnifypat, smoothing, xbegin,
-              row,
+              farWidth, yfillshift, row,
               inRow, outRowBuf0, outRowBuf1, outRowBuf2,
               colNumBuf0, colNumBuf1, colNumBuf2);
     }
@@ -1510,7 +1558,8 @@ produceStereogram(FILE *             const ifP,
     makeImageRows(&inPam, outputGeneratorP,
                   cmdline.depth, cmdline.eyesep, cmdline.dpi,
                   cmdline.crosseyed, cmdline.makemask, cmdline.tileable,
-                  cmdline.magnifypat, cmdline.smoothing, xbegin);
+                  cmdline.magnifypat, cmdline.smoothing, xbegin,
+                  cmdline.yfillshift);
 
     if (cmdline.guidebottom)
         drawguides(cmdline.guidesize, &outputGeneratorP->pam,
@@ -1567,8 +1616,6 @@ main(int argc, const char *argv[]) {
     if (cmdline.verbose)
         reportParameters(cmdline);
 
-    srand(cmdline.randomseedSpec ? cmdline.randomseed : pm_randseed());
-
     ifP = pm_openr(cmdline.inputFilespec);
 
     /* Produce a stereogram. */
@@ -1580,4 +1627,3 @@ main(int argc, const char *argv[]) {
 }
 
 
-
diff --git a/generator/pamstereogram.test b/generator/pamstereogram.test
deleted file mode 100644
index 80f70ee0..00000000
--- a/generator/pamstereogram.test
+++ /dev/null
@@ -1,70 +0,0 @@
-# Make some input files
-pamdepth -quiet 255 ../testgrid.pbm >/tmp/testgrid.pgm
-
-
-# Random pattern
-
-echo Test 01.  Should print 610673698 293:
-./pamstereogram -randomseed=1 ../testgrid.pbm | cksum 
-echo Test 02.  Should print 610673698 293:
-./pamstereogram -randomseed=1 -blackandwhite ../testgrid.pbm | cksum 
-echo Test 03.  Should print 3439084201 170:
-pamseq -tupletype=GRAYSCALE 1 100 | ./pamstereogram -randomseed=1 | cksum 
-echo Test 04.  Should print 2484923390 1070:
-pamgauss 100 10 -maxval=10000 -sigma 20 | pamfunc -multiplier=500 | \
-  ./pamstereogram -randomseed=1 -dpi=10 | cksum
-
-# Makemask
-
-echo Test 10.  Should print 1266273778 293:
-./pamstereogram -randomseed=1 -makemask ../testgrid.pbm | cksum 
-
-echo Test 11.  Should print 3034751595 1070:
-pamgauss 100 10 -maxval=10000 -sigma 20 | pamfunc -multiplier=500 | \
-  ./pamstereogram -randomseed=1 -dpi=10 -makemask | cksum
-
-# Grayscale
-
-echo Test 20.  Should print 2468969328 289:
-./pamstereogram -randomseed=1 -grayscale ../testgrid.pbm | cksum 
-echo Test 21.  Should print 1946982115 4068:
-pamseq 1 100 | pnmtile 200 20 | \
-  ./pamstereogram -randomseed=1 -dpi=10 -grayscale | \
-  cksum
-echo Test 22.  Should print 2078013430 4068:
-pamseq 1 100 | pnmtile 200 20 | \
-  ./pamstereogram -randomseed=1 -dpi=10 -grayscale -maxval 255 | \
-  cksum
-
-# Color
-
-echo Test 30.  Should print 1319392622 731:
-./pamstereogram -randomseed=1 -color ../testgrid.pbm | cksum 
-echo Test 31.  Should print 389886159 12062:
-pamseq 1 100 | pnmtile 200 20 | \
-  ./pamstereogram -randomseed=1 -dpi=10 -color | \
-  cksum
-
-# Pattern file
-
-echo Test 40.  Should print 1834916830 660:
-pamgradient black gray50 white gray50 100 50 | \
-  ./pamstereogram -patfile ../testgrid.pbm -eyesep=.1 -crosseyed | cksum
-
-echo Test 41.  Should print 4016818756 5014:
-pamgradient black gray50 white gray50 100 50 | \
-  ./pamstereogram -patfile /tmp/testgrid.pgm -eyesep=.1 -crosseyed | cksum
-
-# drawguides
-
-echo Test 51.  Should print 2365956562 11071:
-pamgradient black gray50 white gray50 100 50 | \
-  ./pamstereogram -randomseed=1 -dpi 10 -guidesize=20 | cksum
-
-echo Test 51.  Should print 3502025270 1441:
-pamgradient black gray50 white gray50 100 50 | \
-  ./pamstereogram -patfile=../testgrid.pbm -dpi 10 -guidesize=20 | cksum
-
-
-# Clean up files we created
-rm /tmp/testgrid.pgm
diff --git a/generator/pamtris/framebuffer.c b/generator/pamtris/framebuffer.c
index 5665e9a4..4318b1fc 100644
--- a/generator/pamtris/framebuffer.c
+++ b/generator/pamtris/framebuffer.c
@@ -320,7 +320,7 @@ draw_span(uint32_t           const base,
         fbi->z.buffer[j] = (fbi->z.buffer[j] & ~d_mask) | (d & d_mask);
 
         for (l = 0; l < z; l++) {
-	    varying const newval = multiply_varyings(attribs[l], inverse_w);
+            varying const newval = multiply_varyings(attribs[l], inverse_w);
 
             fbi->img.buffer[k + l] =
                 (fbi->img.buffer[k + l] & ~d_mask) |
diff --git a/generator/pamtris/pamtris.c b/generator/pamtris/pamtris.c
index e0becf7a..477b657e 100644
--- a/generator/pamtris/pamtris.c
+++ b/generator/pamtris/pamtris.c
@@ -54,7 +54,7 @@ parse_command_line(int *         const argc_ptr,
     if (!width_spec || !height_spec || (!attribs_spec && !(rgb_spec || grayscale_spec))) {
         pm_errormsg(
             "you must at least specify -width, -height and "
-	    "either -num_attribs, -rgb or -grayscale.");
+            "either -num_attribs, -rgb or -grayscale.");
 
         return 0;
     }
diff --git a/generator/pamtris/triangle.c b/generator/pamtris/triangle.c
index 5143f9ee..264fc221 100644
--- a/generator/pamtris/triangle.c
+++ b/generator/pamtris/triangle.c
@@ -190,8 +190,8 @@ draw_triangle(Xy                 const xy_input,
 
     for (i = 0; i < 3; i++) {
         int32_to_varying_array(attribs_input._[i], attribs[i], n);
-	attribs[i][z] = compute_varying_z(attribs_input._[i][z]);
-	attribs[i][w] = inverse_varying(attribs[i][w]);
+        attribs[i][z] = compute_varying_z(attribs_input._[i][z]);
+        attribs[i][w] = inverse_varying(attribs[i][w]);
         multiply_varying_array_by_varying(attribs[i], attribs[i][w], z);
     }
 
diff --git a/generator/pamtris/utils.c b/generator/pamtris/utils.c
index a6b6e7d4..53f1ea9e 100644
--- a/generator/pamtris/utils.c
+++ b/generator/pamtris/utils.c
@@ -61,7 +61,7 @@ multiply_varying_array_by_varying(varying * const vars,
 
     for (i = 0; i < elements; i++) {
         vars[i].v *= multiplier.v;
-	vars[i].s  = 0.0;
+        vars[i].s  = 0.0;
     }
 }
 
@@ -77,7 +77,7 @@ divide_varying_array_by_varying(varying * const vars,
 
     for (i = 0; i < elements; i++) {
         vars[i].v *= inverse_divisor;
-	vars[i].s  = 0.0;
+        vars[i].s  = 0.0;
     }
 }
 
diff --git a/generator/pbmnoise.c b/generator/pbmnoise.c
new file mode 100644
index 00000000..541ce92e
--- /dev/null
+++ b/generator/pbmnoise.c
@@ -0,0 +1,485 @@
+/* pbmnoise.c - create a random bitmap of a specified size
+                with a specified ratio of white/black pixels
+
+   Written by Akira F Urushibata and contributed to the public domain
+   December 2021
+*/
+
+#include <math.h>
+#include <string.h>
+#include <assert.h>
+
+#include "pm_c_util.h"
+#include "shhopt.h"
+#include "mallocvar.h"
+#include "rand.h"
+#include "nstring.h"
+
+#include "pbm.h"
+
+
+
+static void
+parseFraction(const char *   const fraction,
+              unsigned int * const numeratorP,
+              unsigned int * const precisionP) {
+
+    unsigned int numerator, denominator;
+
+    sscanf(fraction, "%u/%u", &numerator, &denominator);
+
+    if (denominator > 65536)
+        pm_error("Denominator (%u) too large.", denominator);
+    else if (denominator == 0 || (denominator & (denominator - 1)) != 0)
+        pm_error("Denominator must be a power of two.  You specified %u.",
+                 denominator);
+    else if (numerator > denominator)
+        pm_error("Invalid fraction (%s).  Denominator must be larger than "
+                 "numerator.", fraction);
+    else {
+        /* Reduce fraction to lowest terms */
+        unsigned int numerator2, denominator2;
+            /* The fraction reduced to lowest terms */
+        unsigned int precision2;
+
+        if (numerator == 0) { /* all white image */
+            numerator2   = 0;
+            denominator2 = 1;
+            precision2   = 0;
+        } else if (numerator == denominator) { /* all black */
+            numerator2   = 1;
+            denominator2 = 1;
+            precision2   = 0;
+        } else {
+            numerator2   = numerator;   /* initial value */
+            denominator2 = denominator; /* initial value */
+
+            while ((numerator2 & 0x01) != 0x01) {
+                denominator2 /= 2;
+                numerator2   /= 2;
+            }
+            precision2 = 1;
+        }
+
+      if (denominator != denominator2)
+          pm_message("Ratio %u/%u = %u/%u",
+                     numerator, denominator, numerator2, denominator2);
+
+      *precisionP = (precision2 == 0) ? 0 : pm_maxvaltobits(denominator2 - 1);
+          /* pm_maxvaltobits(N):  Max of N is 65535 */
+
+      *numeratorP = numerator2;
+    }
+}
+
+
+
+static void
+setRatio(const char *   const ratioArg,
+         unsigned int * const numeratorP,
+         unsigned int * const precisionP) {
+/*----------------------------------------------------------------------------
+    Convert string "ratioArg" to ratio: numerator / (2 ^ precision) The input
+    string must be in fraction "n/d" form and the denominator must be a power
+    of 2.
+
+    Ratio is the probability of one binary digit being "1".  The ratio of "1"
+    (=PBM black) pixels in the entire output image will be close to this
+    value.
+
+    Most invalid strings are rejected here.
+----------------------------------------------------------------------------*/
+    if (strspn(ratioArg, "0123456789/") == strlen(ratioArg) &&
+             ratioArg[0] != '/' &&
+             ratioArg[strlen(ratioArg) - 1] != '/' &&
+             strchr(ratioArg, '/') != NULL &&
+             strchr(ratioArg, '/') == strrchr(ratioArg, '/'))
+        parseFraction(ratioArg, numeratorP, precisionP);
+    else
+        pm_error("Invalid ratio: '%s'", ratioArg);
+}
+
+
+
+struct CmdlineInfo {
+    /* All the information the user supplied in the command line,
+       in a form easy for the program to use.
+    */
+    unsigned int width;
+    unsigned int height;
+    unsigned int numerator;
+    unsigned int precision;
+    unsigned int randomseed;
+    unsigned int randomseedSpec;
+    unsigned int bswap;   /* boolean */
+    unsigned int pack;    /* boolean */
+};
+
+
+
+static void
+parseCommandLine(int argc, const char ** argv,
+                 struct CmdlineInfo * const cmdlineP) {
+/*----------------------------------------------------------------------------
+   Note that the file spec array we return is stored in the storage that
+   was passed to us as the argv array.
+-----------------------------------------------------------------------------*/
+    optEntry *option_def;
+        /* Instructions to pm_optParseOptions3 on how to parse our options.
+         */
+    optStruct3 opt;
+
+    unsigned int option_def_index;
+    const char * ratioArg;
+    unsigned int ratioSpec;
+    const char * endianArg;
+    unsigned int endianSpec;
+
+    MALLOCARRAY_NOFAIL(option_def, 100);
+
+    option_def_index = 0;   /* incremented by OPTENTRY */
+    OPTENT3(0, "ratio",         OPT_STRING, &ratioArg,
+            &ratioSpec,                     0);
+    OPTENT3(0, "randomseed",    OPT_UINT,   &cmdlineP->randomseed,
+            &cmdlineP->randomseedSpec,      0);
+    OPTENT3(0, "endian",        OPT_STRING, &endianArg,
+            &endianSpec,                    0);
+    OPTENT3(0, "pack",          OPT_FLAG,   NULL,
+            &cmdlineP->pack,                0);
+
+    opt.opt_table = option_def;
+    opt.short_allowed = FALSE;  /* We have no short (old-fashioned) options */
+    opt.allowNegNum = FALSE;  /* We may have parms that are negative numbers */
+
+    pm_optParseOptions3(&argc, (char **)argv, opt, sizeof(opt), 0);
+    /* Uses and sets argc, argv, and some of *cmdlineP and others. */
+    free(option_def);
+
+    if (ratioSpec)
+        setRatio(ratioArg, &cmdlineP->numerator, &cmdlineP->precision);
+    else {
+        /* Set ratio to default: 1/2 */
+        cmdlineP->numerator = 1;
+        cmdlineP->precision = 1;
+    }
+
+    if (!endianSpec)
+        cmdlineP->bswap = false;
+    else {
+        if      (streq(endianArg, "native"))
+            cmdlineP->bswap = false;
+        else if (streq(endianArg, "swap"))
+            cmdlineP->bswap = true;
+        else if (streq(endianArg, "big"))
+            cmdlineP->bswap = (BYTE_ORDER == LITTLE_ENDIAN);
+        else if (streq(endianArg, "little"))
+            cmdlineP->bswap = (BYTE_ORDER != LITTLE_ENDIAN);
+        else
+            pm_error("Invalid value '%s' for -endian argument.", endianArg);
+    }
+
+    if (argc-1 != 2)
+        pm_error("Wrong number of arguments (%d).  There are two "
+                 "non-option arguments: width and height in pixels",
+                 argc-1);
+    else {
+        cmdlineP->width  = pm_parse_width (argv[1]);
+        cmdlineP->height = pm_parse_height(argv[2]);
+    }
+}
+
+
+
+static void
+writeSingleColorRaster(unsigned int const cols,
+                       unsigned int const rows,
+                       bit          const color,
+                       FILE *       const ofP) {
+/*-----------------------------------------------------------------------------
+  Generate a single-color raster of color 'color', dimensions
+  'cols' by 'rows', to output file *ofP.
+-----------------------------------------------------------------------------*/
+    unsigned int const lastColChar = (cols - 1) / 8;
+
+    unsigned char * bitrow0;
+    unsigned int i;
+
+    bitrow0 = pbm_allocrow_packed(cols + 32);
+
+    for (i = 0; i <= lastColChar; ++i)
+        bitrow0[i] = color * 0xff;
+
+    if (color != 0)
+        pbm_cleanrowend_packed(bitrow0, cols);
+
+    /* row end trimming, not necessary with white */
+
+    {
+        unsigned int row;
+        for (row = 0; row < rows; ++row)
+            pbm_writepbmrow_packed(ofP, bitrow0, cols, 0);
+    }
+    pbm_freerow(bitrow0);
+}
+
+
+
+static uint32_t
+randombits(unsigned int       const precision,
+           unsigned int       const numerator,
+           struct pm_randSt * const randStP) {
+/*----------------------------------------------------------------------------
+  Generate 32 random bits so that for each bit the probability of "1"
+  being generated is numerator / (2 ^ precision).
+
+  How this works:
+
+  Ratios such as 1/8, 7/8, 1/16, 15/16, 1/32, 31/32 are straightforward.  How
+  do you get intermediate values such as 3/8, 5/8, 3/16, 5/16, 7/16?
+
+  Imagine a set of 10 bits which are 90% 1, 10% 0 and a random number source
+  which produces 1 and 0 in even proportions.
+
+  Conduct "and" and "or" on these bits:
+
+           0011111111 (90%)       0011111111 (90%)
+      and) 0101010101 (50%)   or) 0101010101 (50%)
+      ---------------------   --------------------
+           0001010101 (45%)       0111111111 (95%)
+
+  It can be seen from this example that an "and" operation gives a new ratio
+  which is halfway between the old ratio (90% in this example) and 0%, while
+  "or" gives one at the middle of the old ratio and 100% The corresponding
+  binary operations for fixed-point fractions are: "right-shift by one and
+  insert a 0 behind the fixed point" and "right-shift by one and insert a 1
+  behind the fixed point" respecatbly.
+
+  115/128 = 89.84% (near 90%)  In binary fix-point: 0.1110011
+  0.01110011 = 115/256 = 44.92%
+  0.11110011 = 243/256 = 94.92%
+
+  So to achieve the desired ratio, start at the LSB (right end) of
+  'numerator'.  Initialize the output bits to zero.  Conduct an "and" for each
+  0 and an "or" for each 1 with a freshly drawn random number until the fixed
+  point is reached.
+
+  An "and" operation of a random number and zero always yields zero.  To avoid
+  waste, we reduce terms to eliminate the trailing zeroes in 'numerator' and
+  indicate the fixed point with 'precision'.  Each time the program is
+  executed the location of the fixed point is set anew, but it stays constant
+  until the program exits.
+----------------------------------------------------------------------------*/
+    unsigned int i;
+    uint32_t mask;
+    uint32_t retval;
+
+    for (i = 0, mask = 0x01, retval=0x00000000;
+         i < precision;
+         ++i, mask <<= 1) {
+
+        uint32_t const randval = pm_rand32(randStP);
+
+        if ((numerator & mask) != 0 )
+            retval |= randval;
+        else
+            retval &= randval;
+    }
+    return retval;
+}
+
+
+
+static uint32_t
+swapWord(uint32_t const word) {
+/*----------------------------------------------------------------------------
+  Swap four bytes.
+  This swap method works regardless of native system endianness.
+----------------------------------------------------------------------------*/
+    uint32_t const retval =
+        ((word >> 24) & 0xff) |
+        ((word >>  8) & 0xff00) |
+        ((word <<  8) & 0xff0000) |
+        ((word << 24) & 0xff000000)
+        ;
+
+    return retval;
+}
+
+
+
+static void
+swapBitrow(unsigned char * const bitrow,
+           unsigned int    const words,
+           bool            const bswap) {
+/*----------------------------------------------------------------------------
+  Modify bits in 'bitrow', swapping as indicated.
+----------------------------------------------------------------------------*/
+    uint32_t * const bitrowByWord = (uint32_t *) bitrow;
+
+    unsigned int wordCnt;
+
+    for (wordCnt=0; wordCnt < words; ++wordCnt) {
+        uint32_t const inWord = bitrowByWord[wordCnt];
+
+        bitrowByWord[wordCnt] = bswap ? swapWord(inWord) : inWord;
+    }
+}
+
+
+
+static void
+pbmnoise(FILE *             const ofP,
+         unsigned int       const cols,
+         unsigned int       const rows,
+         unsigned int       const numerator,
+         unsigned int       const precision,
+         bool               const bswap,
+         struct pm_randSt * const randStP) {
+/*----------------------------------------------------------------------------
+  Default method of constructing rows.
+
+  Generate pixels in units of 32 bits.
+
+  If cols is not a multiple of 32, discard pixels beyond row end.
+-----------------------------------------------------------------------------*/
+    unsigned int const words = (cols + 31) / 32;
+
+    unsigned char * bitrow;
+    unsigned int row;
+    unsigned int wordCnt;
+
+    bitrow = pbm_allocrow_packed(cols + 32);
+
+    for (row = 0; row < rows; ++row) {
+        uint32_t * const bitrowByWord = (uint32_t *) bitrow;
+
+        for (wordCnt = 0; wordCnt < words; ++wordCnt)
+            bitrowByWord[wordCnt] = randombits(precision, numerator, randStP);
+
+        if (bswap)
+            swapBitrow(bitrow, words, bswap);
+
+        pbm_cleanrowend_packed(bitrow, cols);
+        pbm_writepbmrow_packed(ofP, bitrow, cols, 0);
+    }
+    pbm_freerow(bitrow);
+}
+
+
+
+static void
+pbmnoise_packed(FILE *             const ofP,
+                unsigned int       const cols,
+                unsigned int       const rows,
+                unsigned int       const numerator,
+                unsigned int       const precision,
+                bool               const bswap,
+                struct pm_randSt * const randStP) {
+/*----------------------------------------------------------------------------
+  Alternate method of constructing rows.
+  Like the default pbmnoise(), generate pixels in units of 32 bits
+  but carry over unused pixel data at row end to the next row.
+-----------------------------------------------------------------------------*/
+    unsigned char * bitrow0;
+    uint32_t * bitrowByWord;
+    unsigned int offset;
+    unsigned int row;
+    uint32_t wordSave;    /* Pixels carried over to next row */
+
+    bitrow0 = pbm_allocrow_packed(cols + 63);
+    bitrowByWord = (uint32_t *) bitrow0;
+
+    for (row = 0, offset = 0; row < rows; ++row) {
+        if (offset == 0) {
+            unsigned int const words = (cols + 31 ) / 32;
+
+            unsigned int wordCnt;
+
+            for (wordCnt = 0; wordCnt< words; ++wordCnt) {
+                bitrowByWord[wordCnt] =
+                    randombits(precision, numerator, randStP);
+            }
+
+            if (bswap)
+                swapBitrow(bitrow0, words, bswap);
+
+            wordSave = bitrowByWord[words - 1];
+
+            pbm_cleanrowend_packed(bitrow0, cols);
+            pbm_writepbmrow_packed(ofP, bitrow0, cols, 0);
+            offset = cols % 32;
+        } else {
+            unsigned int const wordsToFetch = (cols - (32 - offset) + 31) / 32;
+            unsigned int const lastWord = wordsToFetch;
+
+            unsigned int wordCnt;
+
+            bitrowByWord[0] = wordSave;
+
+            for (wordCnt = 0; wordCnt < wordsToFetch; ++wordCnt) {
+                bitrowByWord[wordCnt + 1] =
+                    randombits(precision, numerator, randStP);
+            }
+
+            if (bswap)
+                swapBitrow((unsigned char *) & bitrowByWord[1],
+                           wordsToFetch, bswap);
+
+            wordSave = bitrowByWord [lastWord];
+
+            pbm_writepbmrow_bitoffset(ofP, bitrow0, cols, 0, offset);
+            offset = (offset + cols) % 32;
+        }
+    }
+    pbm_freerow(bitrow0);
+}
+
+
+int
+main(int argc, const char *argv[]) {
+
+    struct CmdlineInfo cmdline;
+
+    pm_proginit(&argc, argv);
+
+    parseCommandLine(argc, argv, &cmdline);
+
+    pbm_writepbminit(stdout, cmdline.width, cmdline.height, 0);
+
+    if (cmdline.precision == 0) {
+        bit color;
+
+        if (cmdline.numerator == 0)
+            color = PBM_WHITE;
+        else {
+            assert (cmdline.numerator == 1);
+            color = PBM_BLACK;
+        }
+        writeSingleColorRaster(cmdline.width, cmdline.height, color, stdout);
+    } else if (cmdline.width % 32 == 0 || !cmdline.pack) {
+        struct pm_randSt randSt;
+
+        pm_randinit(&randSt);
+        pm_srand2(&randSt, cmdline.randomseedSpec, cmdline.randomseed);
+
+        pbmnoise(stdout, cmdline.width, cmdline.height,
+                 cmdline.numerator, cmdline.precision,
+                 cmdline.bswap, &randSt);
+
+        pm_randterm(&randSt);
+    } else {
+        struct pm_randSt randSt;
+
+        pm_randinit(&randSt);
+        pm_srand2(&randSt, cmdline.randomseedSpec, cmdline.randomseed);
+
+        pbmnoise_packed(stdout, cmdline.width, cmdline.height,
+                        cmdline.numerator, cmdline.precision,
+                        cmdline.bswap, &randSt);
+
+        pm_randterm(&randSt);
+    }
+    return 0;
+}
+
+
diff --git a/generator/pbmtext.c b/generator/pbmtext.c
index 7ea64857..6d4ab8c5 100644
--- a/generator/pbmtext.c
+++ b/generator/pbmtext.c
@@ -1,4 +1,4 @@
-/* pbmtext.c - render text into a bitmap
+/* pbmtext.c - render text into a PBM
 **
 ** Copyright (C) 1991 by Jef Poskanzer.
 **
@@ -81,8 +81,13 @@ textFmCmdLine(int argc, const char ** argv) {
     text[0] = '\0';
 
     for (i = 1, totaltextsize = 0; i < argc; ++i) {
-        if (i > 1) {
+        if (i > 1)
             strcat(text, " ");
+
+        if (strlen(argv[i]) > MAXLINECHARS) { /* avoid arithmetic overflow */
+            pm_error("Command line argument %u is %u characters.  "
+                     "Cannot process longer than %u",
+                     i, (unsigned) strlen(argv[i]), (unsigned) MAXLINECHARS);
         }
         totaltextsize += strlen(argv[i]) + (i > 1 ? 1 : 0);
         if (totaltextsize > MAXLINECHARS)
@@ -1105,7 +1110,7 @@ getText(PM_WCHAR             const cmdlineText[],
 
         unsigned int const lineBufTerm = LINEBUFSIZE - 1;
 
-        unsigned int maxlines;
+        unsigned int textArraySz;
             /* Maximum number of lines for which we currently have space in
                the text array
             */
@@ -1122,12 +1127,12 @@ getText(PM_WCHAR             const cmdlineText[],
         buf[lineBufTerm] = L'\1';  /* Initialize to non-zero value */
                                    /* to detect input overrun */
 
-        maxlines = 50;  /* initial value */
-        MALLOCARRAY(textArray, maxlines);
+        textArraySz = 50;  /* initial value */
+        MALLOCARRAY(textArray, textArraySz);
 
         if (!textArray)
             pm_error("Unable to allocate memory for a buffer for up to %u "
-                     "lines of text", maxlines);
+                     "lines of text", textArraySz);
 
         for (lineCount = 0, eof = false; !eof; ) {
             const char * error;
@@ -1143,9 +1148,13 @@ getText(PM_WCHAR             const cmdlineText[],
                             "is longer than %u characters. "
                             "Cannot process",
                             lineCount, (unsigned int) MAXLINECHARS);
-                    if (lineCount >= maxlines) {
-                        maxlines *= 2;
-                        REALLOCARRAY(textArray, maxlines);
+                    if (lineCount >= textArraySz) {
+                        if (textArraySz > UINT_MAX/2)
+                            pm_error("Too many lines of input for "
+                                     "computation (more than %u)",
+                                     textArraySz);
+                        textArraySz *= 2;
+                        REALLOCARRAY(textArray, textArraySz);
                         if (textArray == NULL)
                             pm_error("out of memory");
                     }
@@ -1282,6 +1291,8 @@ computeImageWidth(struct Text          const formattedText,
                   unsigned int *       const colsP,
                   unsigned int *       const maxleftbP) {
 
+    assert (pbm_maxfontwidth() < (INT_MAX - 10) / LINEBUFSIZE);
+
     if (intercharacterSpace < 0 && fontP->maxwidth < -intercharacterSpace)
         pm_error("negative -space value %.2f exceeds font width",
                  intercharacterSpace);
@@ -1562,4 +1573,3 @@ main(int argc, const char *argv[]) {
 }
 
 
-
diff --git a/generator/pbmtextps.c b/generator/pbmtextps.c
index 0c3656e7..2a1173de 100644
--- a/generator/pbmtextps.c
+++ b/generator/pbmtextps.c
@@ -27,6 +27,7 @@
 
 #define _XOPEN_SOURCE 500
   /* Make sure popen() is in stdio.h, strdup() is in string.h */
+
 #include <unistd.h>
 #include <stdio.h>
 #include <stdlib.h>
@@ -41,6 +42,8 @@
 #include "pm_system.h"
 #include "pbm.h"
 
+enum InputFmt {INPUT_LITERAL, INPUT_ASCIIHEX, INPUT_ASCII85};
+
 static void
 validateFontName(const char * const name) {
 /*-----------------------------------------------------------------------------
@@ -52,6 +55,9 @@ validateFontName(const char * const name) {
 -----------------------------------------------------------------------------*/
     unsigned int idx;
 
+    if (name[0] == '\0')
+        pm_error("Font name is empty string");
+
     for (idx = 0; name[idx] != '\0'; ++idx) {
         char const c = name[idx];
 
@@ -71,8 +77,8 @@ validateFontName(const char * const name) {
 
 
 static void
-asciiHexEncode(char *          const inbuff,
-               char *          const outbuff) {
+asciiHexEncode(const char * const inbuff,
+               char *       const outbuff) {
 /*-----------------------------------------------------------------------------
   Convert the input text string to ASCII-Hex encoding.
 
@@ -83,69 +89,357 @@ asciiHexEncode(char *          const inbuff,
 
     unsigned int idx;
 
+    outbuff[0] = '<';
+
     for (idx = 0; inbuff[idx] != '\0'; ++idx) {
         unsigned int const item = (unsigned char) inbuff[idx];
 
-        outbuff[idx*2]   = hexits[item >> 4];
-        outbuff[idx*2+1] = hexits[item & 0xF];
+        outbuff[1 + idx * 2 + 0] = hexits[item >> 4];
+        outbuff[1 + idx * 2 + 1] = hexits[item & 0xF];
     }
 
-    outbuff[idx * 2] = '\0';
+    if (idx == 0)
+        pm_message("Empty input string");
+
+    outbuff[1 + idx * 2] = '>';
+    outbuff[1 + idx * 2 + 1] = '\0';
 }
 
 
 
 static void
-buildTextFromArgs(int           const argc,
-                  const char ** const argv,
-                  const char ** const asciiHexTextP) {
+failForInvalidChar(char         const c,
+                   const char * const type) {
+
+    if (c >= 0x32 || c < 0x7F)
+        pm_error("Invalid character '%c' in '%s' input string",   c, type);
+    else
+        pm_error("Invalid character 0x%02x in '%s' input string", c, type);
+}
+
+
+
+static void
+formatNonemptyAsciiHex(const char * const inbuff,
+                       char       * const outbuff,
+                       unsigned int const inLen) {
+
+    unsigned int inIdx, outIdx;
+        /* Cursors in input and output buffers, respectively */
+    unsigned int validCharCt;
+        /* Number of valid hex characters we've processed so far in input */
+    unsigned int startIdx, endIdx;
+        /* Limits of input buffer cursor ('inIdx') */
+
+     if (inbuff[0] == '<') {
+         startIdx = 1;
+         endIdx = inLen - 2;
+     } else {
+         startIdx = 0;
+         endIdx = inLen - 1;
+     }
+
+     validCharCt = 0;  /* No valid characters seen yet */
+
+     outIdx = 0;  /* Start at beginning of output buffer */
+
+     outbuff[outIdx++] = '<';
+
+     for (inIdx = startIdx; inIdx <= endIdx; ++inIdx) {
+         switch (inbuff[inIdx]) {
+         case '<':
+         case '>':
+             pm_error("Misplaced character '%c' in Ascii Hex input string",
+                      inbuff[inIdx]);
+           break;
+         case '\f':
+         case '\n':
+         case '\r':
+         case ' ':
+         case '\t':
+             /* ignore whitespace chars */
+             break;
+         default:
+             if ((inbuff [inIdx] >='0' && inbuff [inIdx] <='9') ||
+                 (inbuff [inIdx] >='A' && inbuff [inIdx] <='F') ||
+                 (inbuff [inIdx] >='a' && inbuff [inIdx] <='f')) {
+
+                 outbuff[outIdx++] = inbuff [inIdx];
+                 ++validCharCt;
+             } else
+                 failForInvalidChar(inbuff[inIdx], "Ascii Hex");
+             break;
+         }
+     }
+
+     if (validCharCt == 0)
+         pm_message("Empty Ascii Hex input string");
+     else if (validCharCt % 2 != 0)
+         pm_error("Number of characters in Ascii Hex input string "
+                  "is not even");
+
+     outbuff[outIdx++] = '>';
+     outbuff[outIdx++] = '\0';
+}
+
+
+
+static void
+formatAsciiHexString(const char * const inbuff,
+                     char       * const outbuff,
+                     unsigned int const inLen) {
 /*----------------------------------------------------------------------------
-   Build the array of text to be included in the Postscript program to
-   be rendered, from the arguments of this program.
-
-   We encode it in ASCII-Hex notation as opposed to using the plain text from
-   the command line because 1) the command line might have Postscript control
-   characters in it; and 2) the command line might have text in 8-bit or
-   multibyte code, but a Postscript program is supposed to be entirely
-   printable ASCII characters.
------------------------------------------------------------------------------*/
-    char * text;
+  Format the ASCII Hex input 'inbuff' as a Postscript ASCII Hex string,
+  e.g. "<313233>".  Input can be just the ASCII Hex (e.g. "313233") or already
+  formatted (e.g. "<313233>").  Input may also contain white space, which we
+  ignore -- our output never contains white space.  Though in Postscript, an
+  ASCII NUL character counts as white space, we consider it the end of the
+  input.
+
+  We consider white space outside of the <> delimiters to be an error.
+
+  Abort with error message if there is anything other than valid hex digits in
+  the ASCII hex string proper.  This is necessary to prevent code injection.
+----------------------------------------------------------------------------*/
+    if (inLen == 0 ||
+        (inLen == 2 && inbuff[0] == '<' && inbuff[inLen-1] == '>' )) {
+        pm_message("Empty Ascii Hex input string");
+        strncpy(outbuff, "<>", 3);
+    } else {
+        if (inbuff[0] == '<' && inbuff[inLen-1] != '>' )
+            pm_error("Ascii Hex input string starts with '<' "
+                     "but does not end with '>'");
+        else if (inbuff[0] != '<' && inbuff[inLen-1] == '>' )
+            pm_error("Ascii Hex input string ends with '>' "
+                     "but does not start with '<'");
+
+        formatNonemptyAsciiHex(inbuff, outbuff, inLen);
+    }
+}
+
+
+
+static void
+formatNonemptyAscii85(const char * const inbuff,
+                      char       * const outbuff,
+                      unsigned int const inLen) {
+
+    unsigned int inIdx, outIdx;
+        /* Cursors in input and output buffers, respectively */
+    unsigned int seqPos;
+        /* Position in 5-character Ascii-85 sequence where 'inIdx' points */
+    unsigned int startIdx, endIdx;
+        /* Limits of input buffer cursor ('inIdx') */
+
+    if (inLen > 4 && inbuff[0] == '<' && inbuff[1] == '~' &&
+        inbuff[inLen-2] == '~' && inbuff[inLen-1] == '>') {
+        startIdx = 2;
+        endIdx = inLen - 3;
+    } else {
+        startIdx = 0;
+        endIdx = inLen - 1;
+    }
+
+    seqPos = 0;  /* No 5-character Ascii-85 sequence encountered yet */
+    outIdx = 0;  /* Start filling output buffer from beginning */
+
+    outbuff[outIdx++] = '<';
+    outbuff[outIdx++] = '~';
+
+    for (inIdx = startIdx; inIdx <= endIdx; ++inIdx) {
+      switch (inbuff[inIdx]) {
+      case '~':
+          pm_error("Misplaced character '~' in Ascii 85 input string");
+          break;
+      case '\f':
+      case '\n':
+      case '\r':
+      case ' ':
+      case '\t':
+          break;
+      case 'z':
+          /* z extension */
+          if (seqPos > 0)
+              pm_error("Special 'z' character appears in the middle of a "
+                       "5-character Ascii-85 sequence, which is invalid");
+          else
+              outbuff[outIdx++] = inbuff[inIdx];
+          break;
+        default:
+          /* valid Ascii 85 char */
+          if (inbuff [inIdx] >='!' && inbuff [inIdx] <='u') {
+              outbuff[outIdx++] = inbuff[inIdx];
+              seqPos = (seqPos + 1) % 5;
+          } else
+              failForInvalidChar(inbuff[inIdx], "Ascii 85");
+          break;
+      }
+    }
+
+    if (outIdx == 2) {
+        pm_message("Empty Ascii 85 input string");
+    }
+
+    outbuff[outIdx++]   = '~';
+    outbuff[outIdx++] = '>';
+    outbuff[outIdx++] = '\0';
+}
+
+
+
+static void
+formatAscii85String(const char * const inbuff,
+                    char       * const outbuff,
+                    unsigned int const inLen) {
+/*----------------------------------------------------------------------------
+  Format the Ascii-85 input 'inbuff' as a Postscript Ascii-85 string,
+  e.g. "<~313233~>".  Input can be just the Ascii-85 (e.g. "313233") or
+  already formatted (e.g. "<~313233~>").  Input may also contain white space,
+  which we ignore -- our output never contains white space.  Though in
+  Postscript, an ASCII NUL character counts as white space, we consider it the
+  end of the input.
+
+  We consider white space outside of the <~~> delimiters to be an error.
+
+  Abort with error message if we encounter anything other than valid Ascii-85
+  encoding characters in the string proper.  Note that the Adobe variant
+  does not support the "y" extention.
+----------------------------------------------------------------------------*/
+    if (inLen == 0 || (inLen == 4 && strncmp (inbuff, "<~~>", 4) == 0)) {
+        pm_message("Empty Ascii 85 input string");
+        strncpy(outbuff,"<~~>", 5);
+    } else {
+        if (inLen >= 2) {
+            if (inbuff[0] == '<' && inbuff[1] == '~' &&
+                (inLen < 4 || inbuff[inLen-2] != '~'
+                 || inbuff[inLen-1] != '>' ))
+                pm_error("Ascii 85 input string starts with '<~' "
+                         "but does not end with '~>'");
+            else if (inbuff[inLen-2] == '~' && inbuff[inLen-1] == '>' &&
+                     (inLen < 4 || inbuff[0] != '<' || inbuff[1] != '~'))
+                pm_error("Ascii 85 input string ends with '~>' "
+                         "but does not start with '<~'");
+        }
+        formatNonemptyAscii85(inbuff, outbuff, inLen);
+    }
+}
+
+
+
+static void
+combineArgs(int            const argc,
+            const char **  const argv,
+            const char **  const textP,
+            unsigned int * const textSizeP) {
+
     unsigned int totalTextSize;
+    char * text;
+    size_t * argSize;
     unsigned int i;
+    size_t idx;
+
+    MALLOCARRAY_NOFAIL(argSize, argc);
+        /* argv[0] not accessed; argSize[0] not used */
+
+    for (i = 1, totalTextSize = 0; i < argc; ++i) {
+        argSize[i] = strlen(argv[i]);
+        totalTextSize += argSize[i];
+    }
+
+    totalTextSize = totalTextSize + (argc - 2);  /* adjust for spaces */
+
+    MALLOCARRAY(text, totalTextSize + 1); /* add one for \0 at end */
+    if (text == NULL)
+        pm_error("out of memory allocating buffer for "
+                 "%u characters of text", totalTextSize);
+
+    strncpy(text, argv[1], argSize[1]);
+    for (i = 2, idx = argSize[1]; i < argc; ++i) {
+        text[idx++] = ' ';
+        strncpy(&text[idx], argv[i], argSize[i]);
+        idx += argSize[i];
+    }
+
+    assert(idx == totalTextSize);
+
+    text[idx++] = '\0';
+
+    assert(strlen(text) == totalTextSize);
+
+    *textP     = text;
+    *textSizeP = totalTextSize;
+
+    free(argSize);
+}
+
+
 
-    text = strdup("");
-    totalTextSize = 1;
+static void
+buildTextFromArgs(int           const argc,
+                  const char ** const argv,
+                  const char ** const inputTextP,
+                  enum InputFmt const inputFmt) {
+/*----------------------------------------------------------------------------
+   Build the string of text to be included in the Postscript program to be
+   rendered, from the arguments of this program.
+
+   We encode it in either ASCII-Hex or ASCII-85 as opposed to using the plain
+   text from the command line because 1) the command line might have
+   Postscript control characters in it; and 2) the command line might have
+   text in 8-bit or multibyte code, but a Postscript program is supposed to be
+   entirely printable ASCII characters.
+
+   Official Postscript specifications have a limit on the length of a program
+   line, which means there is a limit on the length of text string such as we
+   generate.  But we don't worry about that because Ghostscript actually
+   accepts very long program lines.
+-----------------------------------------------------------------------------*/
+    const char * text;
+        /* All the arguments ('argv') concatenated */
+    unsigned int textSize;
+        /* Length of 'text' */
 
     if (argc-1 < 1)
         pm_error("No text");
 
-    for (i = 1; i < argc; ++i) {
-        if (i > 1) {
-            totalTextSize += 1;
-            text = realloc(text, totalTextSize);
-            if (text == NULL)
-                pm_error("out of memory");
-            strcat(text, " ");
-        }
-        totalTextSize += strlen(argv[i]);
-        text = realloc(text, totalTextSize);
-        if (text == NULL)
-            pm_error("out of memory");
-        strcat(text, argv[i]);
-    }
+    combineArgs(argc, argv, &text, &textSize);
 
-    {
+    switch (inputFmt) {
+    case INPUT_LITERAL: {
         char * asciiHexText;
 
-        MALLOCARRAY(asciiHexText, totalTextSize * 2);
-
+        MALLOCARRAY(asciiHexText, textSize * 2 + 3);
         if (!asciiHexText)
             pm_error("Unable to allocate memory for hex encoding of %u "
-                     "characters of text", totalTextSize);
+                     "characters of text", textSize);
 
         asciiHexEncode(text, asciiHexText);
-        *asciiHexTextP = asciiHexText;
+        *inputTextP = asciiHexText;
+    } break;
+    case INPUT_ASCIIHEX: {
+        char * asciiHexText;
+
+        MALLOCARRAY(asciiHexText, textSize + 3);
+        if (!asciiHexText)
+            pm_error("Unable to allocate memory for hex encoding of %u "
+                     "characters of text", textSize);
+
+        formatAsciiHexString(text, asciiHexText, textSize);
+        *inputTextP = asciiHexText;
+    } break;
+    case INPUT_ASCII85: {
+        char * ascii85Text;
+
+        MALLOCARRAY(ascii85Text, textSize + 5);
+        if (!ascii85Text)
+            pm_error("Unable to allocate memory for hex encoding of %u "
+                     "characters of text", textSize);
+
+        formatAscii85String(text, ascii85Text, textSize);
+        *inputTextP = ascii85Text;
+    } break;
     }
+
     pm_strfree(text);
 }
 
@@ -155,20 +449,23 @@ struct CmdlineInfo {
     /* All the information the user supplied in the command line,
        in a form easy for the program to use.
     */
-    unsigned int res;
-    float        fontsize;
-    const char * font;
-    float        stroke;
-    float        ascent;
-    float        descent;
-    float        leftmargin;
-    float        rightmargin;
-    float        topmargin;
-    float        bottommargin;
-    unsigned int pad;
-    unsigned int verbose;
-    unsigned int dump;
-    const char * text;
+    unsigned int  res;
+    float         fontsize;
+    const char *  font;
+    float         stroke;
+    float         ascent;
+    float         descent;
+    float         leftmargin;
+    float         rightmargin;
+    float         topmargin;
+    float         bottommargin;
+    unsigned int  pad;
+    unsigned int  verbose;
+    unsigned int  dump;
+    const char *  text;
+        /* Text to render, in Postscript format, either Ascii-hex
+           (e.g. <313233>) or Ascii-85 (e.g. <~aBc-~>)
+        */
 };
 
 
@@ -189,6 +486,7 @@ parseCommandLine(int argc, const char ** argv,
     unsigned int cropSpec, ascentSpec, descentSpec;
     unsigned int leftmarginSpec, rightmarginSpec;
     unsigned int topmarginSpec, bottommarginSpec;
+    unsigned int asciihexSpec, ascii85Spec;
 
     MALLOCARRAY_NOFAIL(option_def, 100);
 
@@ -217,6 +515,10 @@ parseCommandLine(int argc, const char ** argv,
             NULL,                    &cropSpec,                 0);
     OPTENT3(0, "pad",           OPT_FLAG,
             NULL,                    &cmdlineP->pad,            0);
+    OPTENT3(0, "asciihex",      OPT_FLAG,
+            NULL,                    &asciihexSpec,             0);
+    OPTENT3(0, "ascii85",       OPT_FLAG,
+            NULL,                    &ascii85Spec,              0);
     OPTENT3(0, "verbose",       OPT_FLAG,
             NULL,                    &cmdlineP->verbose,        0);
     OPTENT3(0, "dump-ps",       OPT_FLAG,
@@ -233,8 +535,6 @@ parseCommandLine(int argc, const char ** argv,
     cmdlineP->leftmargin  = 0;
     cmdlineP->topmargin   = 0;
     cmdlineP->bottommargin = 0;
-    cropSpec       = FALSE;
-    cmdlineP->pad  = FALSE;
 
     opt.opt_table = option_def;
     opt.short_allowed = FALSE;  /* We have no short (old-fashioned) options */
@@ -244,6 +544,8 @@ parseCommandLine(int argc, const char ** argv,
 
     validateFontName(cmdlineP->font);
 
+    if (cmdlineP->res <= 0)
+        pm_error("-resolution must be positive");
     if (cmdlineP->fontsize <= 0)
         pm_error("-fontsize must be positive");
     if (cmdlineP->ascent < 0)
@@ -259,7 +561,7 @@ parseCommandLine(int argc, const char ** argv,
     if (cmdlineP->bottommargin <0)
         pm_error("-bottommargin must not be negative");
 
-    if (cropSpec == TRUE) {
+    if (cropSpec) {
         if (ascentSpec || descentSpec ||
             leftmarginSpec || rightmarginSpec ||
             topmarginSpec || bottommarginSpec ||
@@ -275,8 +577,24 @@ parseCommandLine(int argc, const char ** argv,
              cmdlineP->leftmargin = cmdlineP->fontsize / 2;
     }
 
-    buildTextFromArgs(argc, argv, &cmdlineP->text);
+    {
+        enum InputFmt inputFmt;
+
+        if (asciihexSpec) {
+            if (ascii85Spec)
+                pm_error("You cannot specify both -asciihex and -ascii85");
+            else
+                inputFmt = INPUT_ASCIIHEX;
+        } else if (ascii85Spec)
+            inputFmt = INPUT_ASCII85;
+        else
+            inputFmt = INPUT_LITERAL;
+
+        if (argc-1 < 1)
+            pm_error("No text");
 
+        buildTextFromArgs(argc, argv, &cmdlineP->text, inputFmt);
+    }
     free(option_def);
 }
 
@@ -325,7 +643,7 @@ postscriptProgram(struct CmdlineInfo const cmdline) {
         "/FindFont {/%s findfont} def\n"
         "/fontsize %f def\n"
         "/pensize %f def\n"
-        "/textstring <%s> def\n"
+        "/textstring %s def\n"
         "/ascent %f def\n"
         "/descent %f def\n"
         "/leftmargin %f def\n"
@@ -386,6 +704,13 @@ postscriptProgram(struct CmdlineInfo const cmdline) {
     const char * retval;
     const char * psVariable;
 
+    /* According to Adobe, maximum line length in a Postscript program is
+       255 characters excluding the newline.  The following may generated a
+       line longer than that if 'cmdline.text' or 'cmdline.font' is long.
+       However, Ghostscript accepts much longer lines, so we don't worry
+       about that.
+    */
+
     pm_asprintf(&psVariable, psTemplate, cmdline.font,
                 cmdline.fontsize, cmdline.stroke, cmdline.text,
                 cmdline.ascent, cmdline.descent,
@@ -413,11 +738,8 @@ gsArgList(const char *       const outputFilename,
     const char ** retval;
     unsigned int argCt;  /* Number of arguments in 'retval' so far */
 
-    if (cmdline.res <= 0)
-         pm_error("Resolution (dpi) must be positive.");
-
-    if (cmdline.fontsize <= 0)
-         pm_error("Font size must be positive.");
+    assert(cmdline.res > 0);
+    assert(cmdline.fontsize > 0);
 
     MALLOCARRAY_NOFAIL(retval, maxArgCt+2);
 
diff --git a/generator/pbmupc.c b/generator/pbmupc.c
index 6ef75654..28ccffb7 100644
--- a/generator/pbmupc.c
+++ b/generator/pbmupc.c
@@ -26,146 +26,29 @@
 #define SHORT_HEIGHT ( 8 * LINES_WIDTH )
 #define TALL_HEIGHT ( SHORT_HEIGHT + DIGIT_HEIGHT / 2 )
 
-static int alldig ARGS(( char* cp ));
-static void putdigit ARGS(( int d, bit** bits, int row0, int col0 ));
-static int addlines ARGS(( int d, bit** bits, int row0, int col0, int height, bit color ));
-static int rect ARGS(( bit** bits, int row0, int col0, int height, int width, bit color ));
 
-int
-main( argc, argv )
-    int argc;
-    char* argv[];
-    {
-    register bit** bits;
-    int argn, style, rows, cols, row, digrow, col, digcolofs;
-    char* typecode;
-    char* manufcode;
-    char* prodcode;
-    int sum, p, lc0, lc1, lc2, lc3, lc4, rc0, rc1, rc2, rc3, rc4;
-    const char* const usage = "[-s1|-s2] <type> <manufac> <product>";
 
+static int
+alldig(const char * const cp) {
 
-    pbm_init( &argc, argv );
-
-    argn = 1;
-    style = 1;
-
-    /* Check for flags. */
-    while ( argn < argc && argv[argn][0] == '-' && argv[argn][1] != '\0' )
-	{
-	if ( pm_keymatch( argv[argn], "-s1", 3 ) )
-	    style = 1;
-	else if ( pm_keymatch( argv[argn], "-s2", 3 ) )
-	    style = 2;
-	else
-	    pm_usage( usage );
-	argn++;
-	}
-
-    if ( argn + 3 < argc )
-	pm_usage( usage );
-    typecode = argv[argn];
-    manufcode = argv[argn + 1];
-    prodcode = argv[argn + 2];
-    argn += 3;
-
-    if ( argn != argc )
-	pm_usage( usage );
+    unsigned int i;
 
-    if ( strlen( typecode ) != 1 || ( ! alldig( typecode ) ) ||
-	 strlen( manufcode ) != 5 || ( ! alldig ( manufcode ) ) ||
-	 strlen( prodcode ) != 5 || ( ! alldig ( prodcode ) ) )
-	pm_error(
-	    "type code must be one digit, and\n    manufacturer and product codes must be five digits" );
-    p = typecode[0] - '0';
-    lc0 = manufcode[0] - '0';
-    lc1 = manufcode[1] - '0';
-    lc2 = manufcode[2] - '0';
-    lc3 = manufcode[3] - '0';
-    lc4 = manufcode[4] - '0';
-    rc0 = prodcode[0] - '0';
-    rc1 = prodcode[1] - '0';
-    rc2 = prodcode[2] - '0';
-    rc3 = prodcode[3] - '0';
-    rc4 = prodcode[4] - '0';
-    sum = ( 10 - ( ( ( p + lc1 + lc3 + rc0 + rc2 + rc4 ) * 3 + lc0 + lc2 + lc4 + rc1 + rc3 ) % 10 ) ) % 10;
+    for (i = 0; cp[i] != '\0'; ++i)
+        if (cp[i] < '0' || cp[i] > '9')
+            return 0;
 
-    rows = 2 * MARGIN + SHORT_HEIGHT + DIGIT_HEIGHT;
-    cols = 2 * MARGIN + 12 * LINES_WIDTH + 11 * LINE1_WIDTH;
-    bits = pbm_allocarray( cols, rows );
+    return 1;
+}
 
-    (void) rect( bits, 0, 0, rows, cols, PBM_WHITE );
-    
-    row = MARGIN;
-    digrow = row + SHORT_HEIGHT;
-    col = MARGIN;
-    digcolofs = ( LINES_WIDTH - DIGIT_WIDTH ) / 2;
-
-    if ( style == 1 )
-	putdigit( p, bits, digrow, col - DIGIT_WIDTH - LINE1_WIDTH );
-    else if ( style == 2 )
-	putdigit(
-	    p, bits, row + SHORT_HEIGHT / 2, col - DIGIT_WIDTH - LINE1_WIDTH );
-    col = rect( bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_BLACK );
-    col = rect( bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_WHITE );
-    col = rect( bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_BLACK );
-    col = addlines( p, bits, row, col, TALL_HEIGHT, PBM_WHITE );
-    putdigit( lc0, bits, digrow, col + digcolofs );
-    col = addlines( lc0, bits, row, col, SHORT_HEIGHT, PBM_WHITE );
-    putdigit( lc1, bits, digrow, col + digcolofs );
-    col = addlines( lc1, bits, row, col, SHORT_HEIGHT, PBM_WHITE );
-    putdigit( lc2, bits, digrow, col + digcolofs );
-    col = addlines( lc2, bits, row, col, SHORT_HEIGHT, PBM_WHITE );
-    putdigit( lc3, bits, digrow, col + digcolofs );
-    col = addlines( lc3, bits, row, col, SHORT_HEIGHT, PBM_WHITE );
-    putdigit( lc4, bits, digrow, col + digcolofs );
-    col = addlines( lc4, bits, row, col, SHORT_HEIGHT, PBM_WHITE );
-    col = rect( bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_WHITE );
-    col = rect( bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_BLACK );
-    col = rect( bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_WHITE );
-    col = rect( bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_BLACK );
-    col = rect( bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_WHITE );
-    putdigit( rc0, bits, digrow, col + digcolofs );
-    col = addlines( rc0, bits, row, col, SHORT_HEIGHT, PBM_BLACK );
-    putdigit( rc1, bits, digrow, col + digcolofs );
-    col = addlines( rc1, bits, row, col, SHORT_HEIGHT, PBM_BLACK );
-    putdigit( rc2, bits, digrow, col + digcolofs );
-    col = addlines( rc2, bits, row, col, SHORT_HEIGHT, PBM_BLACK );
-    putdigit( rc3, bits, digrow, col + digcolofs );
-    col = addlines( rc3, bits, row, col, SHORT_HEIGHT, PBM_BLACK );
-    putdigit( rc4, bits, digrow, col + digcolofs );
-    col = addlines( rc4, bits, row, col, SHORT_HEIGHT, PBM_BLACK );
-    col = addlines( sum, bits, row, col, TALL_HEIGHT, PBM_BLACK );
-    col = rect( bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_BLACK );
-    col = rect( bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_WHITE );
-    col = rect( bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_BLACK );
-    if ( style == 1 )
-	putdigit( sum, bits, digrow, col + LINE1_WIDTH );
-
-    pbm_writepbm( stdout, bits, cols, rows, 0 );
-    pm_close( stdout );
-
-    exit( 0 );
-    }
 
-static int
-alldig( cp )
-    char* cp;
-    {
-    for ( ; *cp != '\0'; cp++ )
-	if ( *cp < '0' || *cp > '9' )
-	    return 0;
-    return 1;
-    }
 
 static void
-putdigit( d, bits, row0, col0 )
-    int d;
-    bit** bits;
-    int row0, col0;
-    {
-    int row, col;
-    static bit digits[10][DIGIT_HEIGHT][DIGIT_WIDTH] = {
+putDigit(int          const d,
+         bit **       const bits,
+         unsigned int const row0,
+         unsigned int const col0) {
+
+    static bit const digits[10][DIGIT_HEIGHT][DIGIT_WIDTH] = {
         /* 0 */
         {
             {0,0,0,0,0,0,0,0,0,0,0,0,0,0},
@@ -399,7 +282,7 @@ putdigit( d, bits, row0, col0 )
             {0,0,0,1,1,1,1,1,1,1,1,0,0,0},
             {0,0,0,0,0,0,0,0,0,0,0,0,0,0},
             {0,0,0,0,0,0,0,0,0,0,0,0,0,0}
-        }, 
+        },
         /* 9 */
         {
             {0,0,0,0,0,0,0,0,0,0,0,0,0,0},
@@ -425,103 +308,255 @@ putdigit( d, bits, row0, col0 )
             {0,0,0,0,0,1,1,0,0,0,0,0,0,0},
             {0,0,0,0,0,0,0,0,0,0,0,0,0,0},
             {0,0,0,0,0,0,0,0,0,0,0,0,0,0}
-        } 
+        }
     };
 
-    for ( row = 0; row < DIGIT_HEIGHT; row++ )
-	for ( col = 0; col < DIGIT_WIDTH; col++ )
-	    bits[row0 + row][col0 + col] = digits[d][row][col];
-    }
+    unsigned int row;
 
-static int
-addlines( int d, bit** bits, int row0, int col0, int height, bit color )
-    {
-    switch ( d )
-	{
-	case 0:
-	col0 = rect( bits, row0, col0, height, LINE3_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE2_WIDTH, 1 - color );
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, 1 - color );
-	break;
-
-	case 1:
-	col0 = rect( bits, row0, col0, height, LINE2_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE2_WIDTH, 1 - color );
-	col0 = rect( bits, row0, col0, height, LINE2_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, 1 - color );
-	break;
-
-	case 2:
-	col0 = rect( bits, row0, col0, height, LINE2_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, 1 - color );
-	col0 = rect( bits, row0, col0, height, LINE2_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE2_WIDTH, 1 - color );
-	break;
-
-	case 3:
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE4_WIDTH, 1 - color );
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, 1 - color );
-	break;
-
-	case 4:
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, 1 - color );
-	col0 = rect( bits, row0, col0, height, LINE3_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE2_WIDTH, 1 - color );
-	break;
-
-	case 5:
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE2_WIDTH, 1 - color );
-	col0 = rect( bits, row0, col0, height, LINE3_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, 1 - color );
-	break;
-
-	case 6:
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, 1 - color );
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE4_WIDTH, 1 - color );
-	break;
-
-	case 7:
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE3_WIDTH, 1 - color );
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE2_WIDTH, 1 - color );
-	break;
-
-	case 8:
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE2_WIDTH, 1 - color );
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE3_WIDTH, 1 - color );
-	break;
-
-	case 9:
-	col0 = rect( bits, row0, col0, height, LINE3_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, 1 - color );
-	col0 = rect( bits, row0, col0, height, LINE1_WIDTH, color );
-	col0 = rect( bits, row0, col0, height, LINE2_WIDTH, 1 - color );
-	break;
-
-	default:
-	pm_error( "can't happen" );
-	}
+    for (row = 0; row < DIGIT_HEIGHT; ++row) {
+        unsigned int col;
 
-    return col0;
+        for (col = 0; col < DIGIT_WIDTH; ++col)
+            bits[row0 + row][col0 + col] = digits[d][row][col];
     }
+}
 
-static int
-rect( bit** bits, int row0, int col0, int height, int width, bit color )
-    {
-    int row, col;
 
-    for ( row = row0; row < row0 + height; row++ )
-	for ( col = col0; col < col0 + width; col++ )
-	    bits[row][col] = color;
+
+static unsigned int
+rect(bit ** const bits,
+     unsigned int const row0,
+     unsigned int const col0,
+     unsigned int const height,
+     unsigned int const width,
+     bit          const color) {
+
+    unsigned int row;
+
+    for (row = row0; row < row0 + height; ++row) {
+        unsigned int col;
+
+        for (col = col0; col < col0 + width; ++col)
+            bits[row][col] = color;
+    }
     return col0 + width;
+}
+
+
+
+static unsigned int
+addLines(int          const d,
+         bit **       const bits,
+         unsigned int const row0,
+         unsigned int const startCol,
+         unsigned int const height,
+         bit          const color) {
+
+    unsigned int col0;
+
+    col0 = startCol;  /* initial value */
+
+    switch (d) {
+    case 0:
+        col0 = rect(bits, row0, col0, height, LINE3_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE2_WIDTH, 1 - color);
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, 1 - color);
+        break;
+
+    case 1:
+        col0 = rect(bits, row0, col0, height, LINE2_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE2_WIDTH, 1 - color);
+        col0 = rect(bits, row0, col0, height, LINE2_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, 1 - color);
+        break;
+
+    case 2:
+        col0 = rect(bits, row0, col0, height, LINE2_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, 1 - color);
+        col0 = rect(bits, row0, col0, height, LINE2_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE2_WIDTH, 1 - color);
+        break;
+
+    case 3:
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE4_WIDTH, 1 - color);
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, 1 - color);
+        break;
+
+    case 4:
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, 1 - color);
+        col0 = rect(bits, row0, col0, height, LINE3_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE2_WIDTH, 1 - color);
+        break;
+
+    case 5:
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE2_WIDTH, 1 - color);
+        col0 = rect(bits, row0, col0, height, LINE3_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, 1 - color);
+        break;
+
+    case 6:
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, 1 - color);
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE4_WIDTH, 1 - color);
+        break;
+
+    case 7:
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE3_WIDTH, 1 - color);
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE2_WIDTH, 1 - color);
+        break;
+
+    case 8:
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE2_WIDTH, 1 - color);
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE3_WIDTH, 1 - color);
+        break;
+
+    case 9:
+        col0 = rect(bits, row0, col0, height, LINE3_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, 1 - color);
+        col0 = rect(bits, row0, col0, height, LINE1_WIDTH, color);
+        col0 = rect(bits, row0, col0, height, LINE2_WIDTH, 1 - color);
+        break;
+
+    default:
+        pm_error("INTERNAL ERROR: invalid digit passed to 'addlines'");
     }
+
+    return col0;
+}
+
+
+
+int
+main(int argc, const char ** argv) {
+
+    const char* const usage = "[-s1|-s2] <type> <manufac> <product>";
+
+    bit ** bits;
+    int argn, style, rows, cols, row, digrow, col, digcolofs;
+    const char * typecode;
+    const char * manufcode;
+    const char * prodcode;
+    int sum, p, lc0, lc1, lc2, lc3, lc4, rc0, rc1, rc2, rc3, rc4;
+
+    pm_proginit(&argc, argv);
+
+    argn = 1;
+    style = 1;
+
+    /* Check for flags. */
+    while ( argn < argc && argv[argn][0] == '-' && argv[argn][1] != '\0' )
+        {
+        if ( pm_keymatch( argv[argn], "-s1", 3 ) )
+            style = 1;
+        else if ( pm_keymatch( argv[argn], "-s2", 3 ) )
+            style = 2;
+        else
+            pm_usage( usage );
+        argn++;
+        }
+
+    if ( argn + 3 < argc )
+        pm_usage( usage );
+    typecode = argv[argn];
+    manufcode = argv[argn + 1];
+    prodcode = argv[argn + 2];
+    argn += 3;
+
+    if ( argn != argc )
+        pm_usage( usage );
+
+    if ( strlen( typecode ) != 1 || ( ! alldig( typecode ) ) ||
+         strlen( manufcode ) != 5 || ( ! alldig ( manufcode ) ) ||
+         strlen( prodcode ) != 5 || ( ! alldig ( prodcode ) ) )
+        pm_error(
+            "type code must be one digit, and\n    manufacturer and product codes must be five digits" );
+    p = typecode[0] - '0';
+    lc0 = manufcode[0] - '0';
+    lc1 = manufcode[1] - '0';
+    lc2 = manufcode[2] - '0';
+    lc3 = manufcode[3] - '0';
+    lc4 = manufcode[4] - '0';
+    rc0 = prodcode[0] - '0';
+    rc1 = prodcode[1] - '0';
+    rc2 = prodcode[2] - '0';
+    rc3 = prodcode[3] - '0';
+    rc4 = prodcode[4] - '0';
+    sum = (10 -
+           (((p + lc1 + lc3 + rc0 + rc2 + rc4 ) * 3 +
+             lc0 + lc2 + lc4 + rc1 + rc3)
+            % 10)
+        )
+        % 10;
+
+    rows = 2 * MARGIN + SHORT_HEIGHT + DIGIT_HEIGHT;
+    cols = 2 * MARGIN + 12 * LINES_WIDTH + 11 * LINE1_WIDTH;
+    bits = pbm_allocarray(cols, rows);
+
+    rect(bits, 0, 0, rows, cols, PBM_WHITE);
+
+    row = MARGIN;
+    digrow = row + SHORT_HEIGHT;
+    col = MARGIN;
+    digcolofs = (LINES_WIDTH - DIGIT_WIDTH) / 2;
+
+    if (style == 1)
+        putDigit(p, bits, digrow, col - DIGIT_WIDTH - LINE1_WIDTH);
+    else if (style == 2)
+        putDigit(
+            p, bits, row + SHORT_HEIGHT / 2, col - DIGIT_WIDTH - LINE1_WIDTH);
+
+    col = rect(bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_BLACK);
+    col = rect(bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_WHITE);
+    col = rect(bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_BLACK);
+    col = addLines(p, bits, row, col, TALL_HEIGHT, PBM_WHITE);
+    putDigit(lc0, bits, digrow, col + digcolofs);
+    col = addLines(lc0, bits, row, col, SHORT_HEIGHT, PBM_WHITE);
+    putDigit(lc1, bits, digrow, col + digcolofs);
+    col = addLines(lc1, bits, row, col, SHORT_HEIGHT, PBM_WHITE);
+    putDigit(lc2, bits, digrow, col + digcolofs);
+    col = addLines(lc2, bits, row, col, SHORT_HEIGHT, PBM_WHITE);
+    putDigit(lc3, bits, digrow, col + digcolofs);
+    col = addLines(lc3, bits, row, col, SHORT_HEIGHT, PBM_WHITE);
+    putDigit(lc4, bits, digrow, col + digcolofs);
+    col = addLines(lc4, bits, row, col, SHORT_HEIGHT, PBM_WHITE);
+    col = rect(bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_WHITE);
+    col = rect(bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_BLACK);
+    col = rect(bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_WHITE);
+    col = rect(bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_BLACK);
+    col = rect(bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_WHITE);
+    putDigit(rc0, bits, digrow, col + digcolofs);
+    col = addLines(rc0, bits, row, col, SHORT_HEIGHT, PBM_BLACK);
+    putDigit(rc1, bits, digrow, col + digcolofs);
+    col = addLines(rc1, bits, row, col, SHORT_HEIGHT, PBM_BLACK);
+    putDigit(rc2, bits, digrow, col + digcolofs);
+    col = addLines(rc2, bits, row, col, SHORT_HEIGHT, PBM_BLACK);
+    putDigit(rc3, bits, digrow, col + digcolofs);
+    col = addLines(rc3, bits, row, col, SHORT_HEIGHT, PBM_BLACK);
+    putDigit(rc4, bits, digrow, col + digcolofs);
+    col = addLines(rc4, bits, row, col, SHORT_HEIGHT, PBM_BLACK);
+    col = addLines(sum, bits, row, col, TALL_HEIGHT, PBM_BLACK);
+    col = rect(bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_BLACK);
+    col = rect(bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_WHITE);
+    col = rect(bits, row, col, TALL_HEIGHT, LINE1_WIDTH, PBM_BLACK);
+    if (style == 1)
+        putDigit(sum, bits, digrow, col + LINE1_WIDTH);
+
+    pbm_writepbm(stdout, bits, cols, rows, 0);
+
+    pm_close(stdout );
+
+    exit(0);
+}
+
+
+
diff --git a/generator/pgmcrater b/generator/pgmcrater
index c66f5576..6d81df39 100755
--- a/generator/pgmcrater
+++ b/generator/pgmcrater
@@ -36,6 +36,16 @@ exec perl -w -x -S -- "$0" "$@"
 use strict;
 
 use Getopt::Long;
+use IO::Handle;
+
+
+
+sub pm_message($) {
+    STDERR->print("pgmcrater: $_[0]\n");
+}
+
+
+
 
 sub doVersionHack($) {
     my ($argvR) = @_;
@@ -69,7 +79,7 @@ my $validOptions = GetOptions(
     'randomseed=i' => \my $randomseedOpt);
 
 if (!$validOptions) {
-    print STDERR "Invalid syntax\n";
+    pm_message("Invalid syntax");
     exit(100);
 }
 
diff --git a/generator/pgmkernel.c b/generator/pgmkernel.c
index 37072c38..8d99a76d 100644
--- a/generator/pgmkernel.c
+++ b/generator/pgmkernel.c
@@ -20,7 +20,7 @@
 #include "shhopt.h"
 #include "mallocvar.h"
 #include "pgm.h"
-
+#include "nstring.h"
 
 
 struct CmdlineInfo {
@@ -39,11 +39,11 @@ static void
 parseCommandLine(int argc, const char ** argv,
                  struct CmdlineInfo * const cmdlineP) {
 /*----------------------------------------------------------------------------
-  Convert program invocation arguments (argc,argv) into a format the 
+  Convert program invocation arguments (argc,argv) into a format the
   program can use easily, struct cmdlineInfo.  Validate arguments along
   the way and exit program with message if invalid.
 
-  Note that some string information we return as *cmdlineP is in the storage 
+  Note that some string information we return as *cmdlineP is in the storage
   argv[] points to.
 -----------------------------------------------------------------------------*/
     optEntry *option_def;
@@ -57,9 +57,9 @@ parseCommandLine(int argc, const char ** argv,
     MALLOCARRAY_NOFAIL(option_def, 100);
 
     option_def_index = 0;   /* incremented by OPTENTRY */
-    OPTENT3(0,   "weight",  OPT_FLOAT, &cmdlineP->weight, 
+    OPTENT3(0,   "weight",  OPT_FLOAT, &cmdlineP->weight,
             &weightSpec,     0);
-    OPTENT3(0,   "maxval",  OPT_UINT, &cmdlineP->maxval, 
+    OPTENT3(0,   "maxval",  OPT_UINT, &cmdlineP->maxval,
             &maxvalSpec,     0);
 
     opt.opt_table = option_def;
@@ -89,20 +89,30 @@ parseCommandLine(int argc, const char ** argv,
     if (argc-1 < 1)
         pm_error("Need at least one argument: size of (square) kernel");
     else if (argc-1 == 1) {
-        if (atoi(argv[1]) <= 0)
+        unsigned int dimension;
+        const char * error;
+
+        pm_string_to_uint(argv[1], &dimension, &error);
+        if (error) {
+            pm_error("'%s' is invalid as an image width/height.  %s", argv[1], error);
+            pm_strfree(error);
+        }
+        if (dimension <= 0)
             pm_error("Dimension must be a positive number.  "
                      "You specified '%s'", argv[1]);
-        cmdlineP->cols = atoi(argv[1]);
-        cmdlineP->rows = atoi(argv[1]);
+        cmdlineP->cols = cmdlineP->rows = dimension;
+
     } else if (argc-1 == 2) {
-        if (atoi(argv[1]) <= 0)
+        unsigned int const width  = pm_parse_width(argv[1]);
+        unsigned int const height = pm_parse_height(argv[2]);
+        if (width <= 0)
             pm_error("Width must be a positive number.  "
                      "You specified '%s'", argv[1]);
-        if (atoi(argv[2]) <= 0)
+        if (height <= 0)
             pm_error("Height must be a positive number.  "
                      "You specified '%s'", argv[2]);
-        cmdlineP->cols = atoi(argv[1]);
-        cmdlineP->rows = atoi(argv[2]);
+        cmdlineP->cols = width;
+        cmdlineP->rows = height;
     } else
         pm_error("At most two arguments allowed.  "
                  "You specified %u", argc-1);
@@ -243,3 +253,6 @@ main(int argc, const char * argv[]) {
 
     return 0;
 }
+
+
+
diff --git a/generator/pgmmake.c b/generator/pgmmake.c
index ae706639..6fd8d426 100644
--- a/generator/pgmmake.c
+++ b/generator/pgmmake.c
@@ -52,11 +52,11 @@ static void
 parseCommandLine(int argc, const char ** argv,
                  struct CmdlineInfo * const cmdlineP) {
 /*----------------------------------------------------------------------------
-  Convert program invocation arguments (argc,argv) into a format the 
+  Convert program invocation arguments (argc,argv) into a format the
   program can use easily, struct cmdlineInfo.  Validate arguments along
   the way and exit program with message if invalid.
 
-  Note that some string information we return as *cmdlineP is in the storage 
+  Note that some string information we return as *cmdlineP is in the storage
   argv[] points to.
 -----------------------------------------------------------------------------*/
     optEntry * option_def;
diff --git a/generator/pgmnoise.c b/generator/pgmnoise.c
index e3e4eca6..ea35956b 100644
--- a/generator/pgmnoise.c
+++ b/generator/pgmnoise.c
@@ -7,12 +7,17 @@
 #include "pm_c_util.h"
 #include "mallocvar.h"
 #include "nstring.h"
+#include "rand.h"
 #include "shhopt.h"
 #include "pgm.h"
 
+/* constants */
+static unsigned long int const ceil31bits = 0x7fffffffUL;
+static unsigned long int const ceil32bits = 0xffffffffUL;
 
 
-struct cmdlineInfo {
+
+struct CmdlineInfo {
     /* All the information the user supplied in the command line,
        in a form easy for the program to use.
     */
@@ -21,13 +26,15 @@ struct cmdlineInfo {
     unsigned int maxval;
     unsigned int randomseed;
     unsigned int randomseedSpec;
+    unsigned int verbose;
 };
 
 
 
 static void
-parseCommandLine(int argc, const char ** const argv,
-                 struct cmdlineInfo * const cmdlineP) {
+parseCommandLine(int argc,
+                 const char **        const argv,
+                 struct CmdlineInfo * const cmdlineP) {
 /*----------------------------------------------------------------------------
    Note that the file spec array we return is stored in the storage that
    was passed to us as the argv array.
@@ -46,6 +53,8 @@ parseCommandLine(int argc, const char ** const argv,
             &cmdlineP->randomseedSpec,      0);
     OPTENT3(0,   "maxval",       OPT_UINT,    &cmdlineP->maxval,
             &maxvalSpec,                    0);
+    OPTENT3(0,   "verbose",      OPT_FLAG,    NULL,
+            &cmdlineP->verbose,             0);
 
     opt.opt_table = option_def;
     opt.short_allowed = FALSE;  /* We have no short (old-fashioned) options */
@@ -57,7 +66,7 @@ parseCommandLine(int argc, const char ** const argv,
 
     if (maxvalSpec) {
         if (cmdlineP->maxval > PGM_OVERALLMAXVAL)
-            pm_error("Maxval too large: %u.  Maximu is %u",
+            pm_error("Maxval too large: %u.  Maximum is %u",
                      cmdlineP->maxval, PGM_OVERALLMAXVAL);
         else if (cmdlineP->maxval == 0)
             pm_error("Maxval must not be zero");
@@ -95,38 +104,43 @@ parseCommandLine(int argc, const char ** const argv,
 
 
 static unsigned int
-randPool(unsigned int const digits) {
+randPool(unsigned int       const nDigits,
+         struct pm_randSt * const randStP) {
 /*----------------------------------------------------------------------------
-  Draw 'digits' bits from pool of random bits.  If the number of random bits
-  in pool is insufficient, call rand() and add 31 bits to it.
+  Draw 'nDigits' bits from pool of random bits.  If the number of random bits
+  in pool is insufficient, call pm_rand() and add N bits to it.
+
+  N is 31 or 32.  In raw mode we use N = 32 regardless of the actual number of
+  available bits.  If there are only 31 available, we use zero for the MSB.
 
-  'digits' must be at most 16.
+  'nDigits' must be at most 16.
 
-  We assume that each call to rand() generates 31 bits, or RAND_MAX ==
-  2147483647.
+  We assume that each call to pm_rand() generates 31 or 32 bits, or
+  randStP->max == 2147483647 or 4294967295.
 
-  The underlying logic is flexible and endian-free.  The above conditions
-  can be relaxed.
+  The underlying logic is flexible and endian-free.  The above conditions can
+  be relaxed.
 -----------------------------------------------------------------------------*/
     static unsigned long int hold=0;  /* entropy pool */
     static unsigned int len=0;        /* number of valid bits in pool */
 
-    unsigned int const mask = (1 << digits) - 1;
-
+    unsigned int const mask = (1 << nDigits) - 1;
+    unsigned int const randbits = (randStP->max == ceil31bits) ? 31 : 32;
     unsigned int retval;
 
-    assert(RAND_MAX == 2147483647 && digits <= 16);
+    assert(randStP->max == ceil31bits || randStP->max == ceil32bits);
+    assert(nDigits <= 16);
 
     retval = hold;  /* initial value */
 
-    if (len > digits) { /* Enough bits in hold to satisfy request */
-        hold >>= digits;
-        len   -= digits;
-    } else {              /* Load another 31 bits into hold */
-        hold    = rand();
+    if (len > nDigits) { /* Enough bits in hold to satisfy request */
+        hold >>= nDigits;
+        len   -= nDigits;
+    } else {            /* Load another 31 or 32 bits into hold */
+        hold    = pm_rand(randStP);
         retval |= (hold << len);
-        hold >>=  (digits - len);
-        len = 31 - digits + len;
+        hold >>=  (nDigits - len);
+        len = randbits - nDigits + len;
     }
     return (retval & mask);
 }
@@ -134,36 +148,61 @@ randPool(unsigned int const digits) {
 
 
 static void
-pgmnoise(FILE *       const ofP,
-         unsigned int const cols,
-         unsigned int const rows,
-         gray         const maxval) {
+reportVerbose(struct pm_randSt * const randStP,
+              gray               const maxval,
+              bool               const usingPool)  {
+
+    pm_message("random seed: %u", randStP->seed);
+    pm_message("random max: %u maxval: %u", randStP->max, maxval);
+    pm_message("method: %s", usingPool ? "pool" : "modulo");
+}
+
+
 
-    bool const usingPool = !(RAND_MAX==2147483647 && (maxval & (maxval+1)));
+static void
+pgmnoise(FILE *             const ofP,
+         unsigned int       const cols,
+         unsigned int       const rows,
+         gray               const maxval,
+         bool               const verbose,
+         struct pm_randSt * const randStP) {
+
+    bool const usingPool =
+        (randStP->max==ceil31bits || randStP->max==ceil32bits) &&
+        !(maxval & (maxval+1));
     unsigned int const bitLen = pm_maxvaltobits(maxval);
 
     unsigned int row;
     gray * destrow;
 
     /* If maxval is 2^n-1, we draw exactly n bits from the pool.
-       Otherwise call rand() and determine gray value by modulo.
+       Otherwise call pm_rand() and determine gray value by modulo.
 
        In the latter case, there is a minuscule skew toward 0 (=black)
        because smaller numbers are produced more frequently by modulo.
        Thus we employ the pool method only when it is certain that no
-       skew will ensue.
+       skew will result.
 
        To illustrate the point, consider converting the outcome of one
        roll of a fair, six-sided die to 5 values (0 to 4) by N % 5.  The
-       probability for values 1, 2, 3, 4 are 1/6, but 0 alone is 2/6.
+       probability for values 1, 2, 3, 4 is 1/6, but 0 alone is 2/6.
        Average is 10/6 or 1.6667, compared to 2.0 from an ideal
        generator which produces exactly 5 values.  With two dice
        average improves to 70/36 or 1.9444.
 
        The more (distinct) dice we roll, or the more binary digits we
        draw, the smaller the skew.
+
+       The pool method is economical.  But there is an additional merit:
+       No bits are lost this way.  This gives us a means to check the
+       integrity of the random number generator.
+
+       - Akira Urushibata, March 2021
     */
 
+    if (verbose)
+        reportVerbose(randStP, maxval, usingPool);
+
     destrow = pgm_allocrow(cols);
 
     pgm_writepgminit(ofP, cols, rows, maxval, 0);
@@ -172,12 +211,11 @@ pgmnoise(FILE *       const ofP,
         if (usingPool) {
             unsigned int col;
             for (col = 0; col < cols; ++col)
-                destrow[col] = randPool(bitLen);
-        }
-        else {
+                destrow[col] = randPool(bitLen, randStP);
+        } else {
             unsigned int col;
             for (col = 0; col < cols; ++col)
-                destrow[col] = rand() % (maxval + 1);
+                destrow[col] = pm_rand(randStP) % (maxval + 1);
         }
         pgm_writepgmrow(ofP, destrow, cols, maxval, 0);
     }
@@ -191,18 +229,22 @@ int
 main(int          argc,
      const char * argv[]) {
 
-    struct cmdlineInfo cmdline;
+    struct CmdlineInfo cmdline;
+    struct pm_randSt randSt;
 
     pm_proginit(&argc, argv);
 
     parseCommandLine(argc, argv, &cmdline);
 
-    srand(cmdline.randomseedSpec ? cmdline.randomseed : pm_randseed());
+    pm_randinit(&randSt);
+    pm_srand2(&randSt, cmdline.randomseedSpec, cmdline.randomseed);
+
+    pgmnoise(stdout, cmdline.width, cmdline.height, cmdline.maxval,
+             cmdline.verbose, &randSt);
 
-    pgmnoise(stdout, cmdline.width, cmdline.height, cmdline.maxval);
+    pm_randterm(&randSt);
 
     return 0;
 }
 
 
-
diff --git a/generator/ppmcie.c b/generator/ppmcie.c
index 86325ba6..26289c3c 100644
--- a/generator/ppmcie.c
+++ b/generator/ppmcie.c
@@ -26,8 +26,11 @@
   Introduced option to plot 1976 u' v' chromaticities.
 */
 
+#define _C99_SOURCE  /* Make sure snprintf() is in stdio.h */
+
 #include <assert.h>
 #include <math.h>
+#include <stdio.h>
 
 #include "pm_c_util.h"
 #include "ppm.h"
@@ -1114,15 +1117,15 @@ writeLabel(pixel **                   const pixels,
 
     PPM_ASSIGN(rgbcolor, maxval, maxval, maxval);
 
-    pm_snprintf(sysdesc, sizeof(sysdesc),
-                "System: %s\n"
-                "Primary illuminants (X, Y)\n"
-                "     Red:  %0.4f, %0.4f\n"
-                "     Green: %0.4f, %0.4f\n"
-                "     Blue:  %0.4f, %0.4f\n"
-                "White point (X, Y): %0.4f, %0.4f",
-                cs->name, cs->xRed, cs->yRed, cs->xGreen, cs->yGreen,
-                cs->xBlue, cs->yBlue, cs->xWhite, cs->yWhite);
+    snprintf(sysdesc, sizeof(sysdesc),
+             "System: %s\n"
+             "Primary illuminants (X, Y)\n"
+             "     Red:  %0.4f, %0.4f\n"
+             "     Green: %0.4f, %0.4f\n"
+             "     Blue:  %0.4f, %0.4f\n"
+             "White point (X, Y): %0.4f, %0.4f",
+             cs->name, cs->xRed, cs->yRed, cs->xGreen, cs->yGreen,
+             cs->xBlue, cs->yBlue, cs->xWhite, cs->yWhite);
     sysdesc[sizeof(sysdesc)-1] = '\0';  /* for robustness */
 
     ppmd_text(pixels, pixcols, pixrows, maxval,
diff --git a/generator/ppmcolors.c b/generator/ppmcolors.c
index 701812d1..ecbadf72 100644
--- a/generator/ppmcolors.c
+++ b/generator/ppmcolors.c
@@ -3,7 +3,7 @@
 *******************************************************************************
   Generate a color map containing all the colors representable with a certain
   maxval.
-  
+
   THIS PROGRAM IS OBSOLETE.  USE PAMSEQ INSTEAD.  WE KEEP THIS AROUND
   FOR BACKWARD COMPATIBILITY.
 
@@ -18,7 +18,7 @@ struct cmdlineInfo {
     /* All the information the user supplied in the command line,
        in a form easy for the program to use.
     */
-    unsigned int maxval;   
+    unsigned int maxval;
 };
 
 
@@ -34,9 +34,9 @@ parseCommandLine(int argc, char ** argv, struct cmdlineInfo *cmdlineP) {
     unsigned int option_def_index;
 
     option_def_index = 0;   /* incremented by OPTENT3 */
-    OPTENT3(0,   "maxval",       OPT_UINT,   
+    OPTENT3(0,   "maxval",       OPT_UINT,
             &cmdlineP->maxval,  NULL, 0);
-    
+
 
     opt.opt_table = option_def;
     opt.short_allowed = FALSE;  /* We have no short (old-fashioned) options */
@@ -76,7 +76,7 @@ main(int argc, char *argv[]) {
 
     rc = system(cmd);
 
-    if (rc != 0) 
+    if (rc != 0)
         pm_error("pamseq|pamtopnm pipeline failed.  system() rc = %d", rc);
 
     pm_strfree(cmd);
diff --git a/generator/ppmforge.c b/generator/ppmforge.c
index 114f7f18..e32d2b59 100644
--- a/generator/ppmforge.c
+++ b/generator/ppmforge.c
@@ -34,12 +34,14 @@
 #define _XOPEN_SOURCE 500  /* get M_PI in math.h */
 
 #include <math.h>
+#include <float.h>
 #include <assert.h>
 
 #include "pm_c_util.h"
-#include "ppm.h"
 #include "mallocvar.h"
+#include "rand.h"
 #include "shhopt.h"
+#include "ppm.h"
 
 static double const hugeVal = 1e50;
 
@@ -59,12 +61,8 @@ typedef struct {
 
 /* Definition for obtaining random numbers. */
 
-#define nrand 4               /* Gauss() sample count */
-#define Cast(low, high) ((low)+(((high)-(low)) * ((rand() & 0x7FFF) / arand)))
-
 /*  Local variables  */
 
-static double arand, gaussadd, gaussfac; /* Gaussian random parameters */
 static double fracdim;            /* Fractal dimension */
 static double powscale;           /* Power law scaling exponent */
 static int meshsize = 256;        /* FFT mesh size */
@@ -92,6 +90,7 @@ struct CmdlineInfo {
     float        ice;
     int          saturation;
     unsigned int seed;
+    unsigned int seedSpec;
     int          stars;
     unsigned int starsSpec;
     unsigned int width;
@@ -172,6 +171,11 @@ parseCommandLine(int argc, const char **argv,
         if (cmdlineP->dimension <= 0.0)
             pm_error("-dimension must be greater than zero.  "
                      "You specified %f", cmdlineP->dimension);
+        else if (cmdlineP->dimension > 5.0 + FLT_EPSILON)
+            pm_error("-dimension must not be greater than 5.  "
+                     "Results are not interesting with higher numbers, so "
+                     "we assume it is a mistake.  "
+                     "You specified %f", cmdlineP->dimension);
     } else
         cmdlineP->dimension = cmdlineP->clouds ? 2.15 : 2.4;
 
@@ -238,7 +242,7 @@ parseCommandLine(int argc, const char **argv,
 
 
 static void
-fourn(float *     const data,
+fourn(double *    const data,
       const int * const nn,
       int         const ndim,
       int         const isign) {
@@ -269,7 +273,7 @@ fourn(float *     const data,
     int i1, i2, i3;
     int i2rev, i3rev, ip1, ip2, ip3, ifp1, ifp2;
     int ibit, idim, k1, k2, n, nprev, nrem, ntot;
-    float tempi, tempr;
+    double tempi, tempr;
     double theta, wi, wpi, wpr, wr, wtemp;
 
 #define SWAP(a,b) tempr=(a); (a) = (b); (b) = tempr
@@ -332,72 +336,64 @@ fourn(float *     const data,
         nprev *= n;
     }
 }
-#undef SWAP
 
 
 
-static void
-initgauss(unsigned int const seed) {
+static double
+cast(double             const low,
+     double             const high,
+     struct pm_randSt * const randStP) {
 /*----------------------------------------------------------------------------
-  Initialize random number generators.
-
-  As given in Peitgen & Saupe, page 77.
+   A random number in the range ['low', 'high'].
 -----------------------------------------------------------------------------*/
-    /* Range of random generator */
-    arand = pow(2.0, 15.0) - 1.0;
-    gaussadd = sqrt(3.0 * nrand);
-    gaussfac = 2 * gaussadd / (nrand * arand);
-    srand(seed);
+    return (low + (high - low) * pm_drand(randStP));
+
 }
 
 
 
 static double
-gauss() {
+randPhase(struct pm_randSt * const randStP) {
 /*----------------------------------------------------------------------------
-  A Gaussian random number.
-
-  As given in Peitgen & Saupe, page 77.
+   A random number in the range [0, 2 * PI).
 -----------------------------------------------------------------------------*/
-    int i;
-    double sum;
+     return (2.0 * M_PI * pm_drand1(randStP));
 
-    for (i = 1, sum = 0.0; i <= nrand; ++i) {
-        sum += (rand() & 0x7FFF);
-    }
-    return gaussfac * sum - gaussadd;
 }
 
 
 
 static void
-spectralsynth(float **     const x,
-              unsigned int const n,
-              double        const h) {
+spectralsynth(double **          const aP,
+              unsigned int       const n,
+              double             const h,
+              struct pm_randSt * const randStP) {
 /*----------------------------------------------------------------------------
   Spectrally synthesized fractal motion in two dimensions.
 
   This algorithm is given under the name SpectralSynthesisFM2D on page 108 of
   Peitgen & Saupe.
 -----------------------------------------------------------------------------*/
-    unsigned bl;
+    unsigned int const bl = ((((unsigned long) n) * n) + 1) * 2;
+
     int i, j, i0, j0, nsize[3];
     double rad, phase, rcos, rsin;
-    float *a;
+    double * a;
+
+    MALLOCARRAY(a, bl);
 
-    bl = ((((unsigned long) n) * n) + 1) * 2 * sizeof(float);
-    a = (float *) calloc(bl, 1);
-    if (a == (float *) 0) {
-        pm_error("Cannot allocate %d x %d result array (% d bytes).",
+    if (!a) {
+        pm_error("Cannot allocate %u x %u result array (%u doubles).",
                  n, n, bl);
     }
-    *x = a;
+    for (i = 0; i < bl; ++i)
+        a[i] = 0.0;  /* initial value */
 
     for (i = 0; i <= n / 2; i++) {
         for (j = 0; j <= n / 2; j++) {
-            phase = 2 * M_PI * ((rand() & 0x7FFF) / arand);
+            phase = randPhase(randStP);
             if (i != 0 || j != 0) {
-                rad = pow((double) (i * i + j * j), -(h + 1) / 2) * gauss();
+                rad = pow( (double) (i*i + j*j), - (h + 1) / 2) * pm_gaussrand(randStP);
             } else {
                 rad = 0;
             }
@@ -416,8 +412,8 @@ spectralsynth(float **     const x,
     Imag(a, n / 2, n / 2) = 0;
     for (i = 1; i <= n / 2 - 1; i++) {
         for (j = 1; j <= n / 2 - 1; j++) {
-            phase = 2 * M_PI * ((rand() & 0x7FFF) / arand);
-            rad = pow((double) (i * i + j * j), -(h + 1) / 2) * gauss();
+            phase = randPhase(randStP);
+            rad = pow((double) (i * i + j * j), -(h + 1) / 2) * pm_gaussrand(randStP);
             rcos = rad * cos(phase);
             rsin = rad * sin(phase);
             Real(a, i, n - j) = rcos;
@@ -430,6 +426,8 @@ spectralsynth(float **     const x,
     nsize[0] = 0;
     nsize[1] = nsize[2] = n;          /* Dimension of frequency domain array */
     fourn(a, nsize, 2, -1);       /* Take inverse 2D Fourier transform */
+
+    *aP = a;
 }
 
 
@@ -469,22 +467,22 @@ temprgb(double   const temp,
 
 
 static void
-etoile(pixel * const pix) {
+etoile(pixel *        const pix,
+       struct pm_randSt * const randStP) {
 /*----------------------------------------------------------------------------
     Set a pixel in the starry sky.
 -----------------------------------------------------------------------------*/
-    if ((rand() % 1000) < starfraction) {
-#define StarQuality 0.5       /* Brightness distribution exponent */
-#define StarIntensity   8         /* Brightness scale factor */
-#define StarTintExp 0.5       /* Tint distribution exponent */
-        double v = StarIntensity * pow(1 / (1 - Cast(0, 0.9999)),
-                                       (double) StarQuality),
-            temp,
-            r, g, b;
-
-        if (v > 255) {
-            v = 255;
-        }
+    if ((pm_rand(randStP) % 1000) < starfraction) {
+        double const starQuality   = 0.5;
+            /* Brightness distribution exponent */
+        double const starIntensity = 8;
+            /* Brightness scale factor */
+        double const starTintExp = 0.5;
+            /* Tint distribution exponent */
+        double const v =
+            MIN(255.0,
+                starIntensity * pow(1 / (1 - cast(0, 0.9999, randStP)),
+                                    (double) starQuality));
 
         /* We make a special case for star color  of zero in order to
            prevent  floating  point  roundoff  which  would  otherwise
@@ -493,13 +491,17 @@ etoile(pixel * const pix) {
            256 shades in the image. */
 
         if (starcolor == 0) {
-            int vi = v;
+            pixval const vi = v;
 
             PPM_ASSIGN(*pix, vi, vi, vi);
         } else {
+            double temp;
+            double r, g, b;
+
             temp = 5500 + starcolor *
-                pow(1 / (1 - Cast(0, 0.9999)), StarTintExp) *
-                ((rand() & 7) ? -1 : 1);
+                pow(1 / (1 - cast(0, 0.9999, randStP)), starTintExp) *
+                ((pm_rand(randStP) & 7) ? -1 : 1);
+
             /* Constrain temperature to a reasonable value: >= 2600K
                (S Cephei/R Andromedae), <= 28,000 (Spica). */
             temp = MAX(2600, MIN(28000, temp));
@@ -535,9 +537,9 @@ atSat(double const x,
 
 
 static unsigned char *
-makeCp(float *      const a,
-       unsigned int const n,
-       pixval       const maxval) {
+makeCp(const double * const a,
+       unsigned int   const n,
+       pixval         const maxval) {
 
     /* Prescale the grid points into intensities. */
 
@@ -550,7 +552,7 @@ makeCp(float *      const a,
     if (cp == NULL)
         pm_error("Unable to allocate %u bytes for cp array", n);
 
-    ap = cp;
+    ap = cp;  /* initial value */
     {
         unsigned int i;
         for (i = 0; i < n; i++) {
@@ -565,17 +567,18 @@ makeCp(float *      const a,
 
 
 static void
-createPlanetStuff(bool             const clouds,
-                  float *          const a,
-                  unsigned int     const n,
-                  double **        const uP,
-                  double **        const u1P,
-                  unsigned int **  const bxfP,
-                  unsigned int **  const bxcP,
-                  unsigned char ** const cpP,
-                  Vector *         const sunvecP,
-                  unsigned int     const cols,
-                  pixval           const maxval) {
+createPlanetStuff(bool               const clouds,
+                  const double *     const a,
+                  unsigned int       const n,
+                  double **          const uP,
+                  double **          const u1P,
+                  unsigned int **    const bxfP,
+                  unsigned int **    const bxcP,
+                  unsigned char **   const cpP,
+                  Vector *           const sunvecP,
+                  unsigned int       const cols,
+                  pixval             const maxval,
+                  struct pm_randSt * const randStP) {
 
     double *u, *u1;
     unsigned int *bxf, *bxc;
@@ -585,8 +588,8 @@ createPlanetStuff(bool             const clouds,
 
     /* Compute incident light direction vector. */
 
-    shang = hourspec ? hourangle : Cast(0, 2 * M_PI);
-    siang = inclspec ? inclangle : Cast(-M_PI * 0.12, M_PI * 0.12);
+    shang = hourspec ? hourangle : randPhase(randStP);
+    siang = inclspec ? inclangle : cast(-M_PI * 0.12, M_PI * 0.12, randStP);
 
     sunvecP->x = sin(shang) * cos(siang);
     sunvecP->y = sin(siang);
@@ -594,7 +597,7 @@ createPlanetStuff(bool             const clouds,
 
     /* Allow only 25% of random pictures to be crescents */
 
-    if (!hourspec && ((rand() % 100) < 75)) {
+    if (!hourspec && ((pm_rand(randStP) % 100) < 75)) {
         flipped = (sunvecP->z < 0);
         sunvecP->z = fabs(sunvecP->z);
     } else
@@ -634,7 +637,7 @@ createPlanetStuff(bool             const clouds,
         pm_error("Cannot allocate %u element interpolation tables.", cols);
     {
         unsigned int j;
-        for (j = 0; j < cols; j++) {
+        for (j = 0; j < cols; ++j) {
             double const bx = (n - 1) * uprj(j, cols);
 
             bxf[j] = floor(bx);
@@ -651,8 +654,9 @@ createPlanetStuff(bool             const clouds,
 
 
 static void
-generateStarrySkyRow(pixel *      const pixels,
-                     unsigned int const cols) {
+generateStarrySkyRow(pixel *            const pixels,
+                     unsigned int       const cols,
+                     struct pm_randSt * const randStP) {
 /*----------------------------------------------------------------------------
   Generate a starry sky.  Note  that no FFT is performed;
   the output is  generated  directly  from  a  power  law
@@ -660,8 +664,8 @@ generateStarrySkyRow(pixel *      const pixels,
 -----------------------------------------------------------------------------*/
     unsigned int j;
 
-    for (j = 0; j < cols; j++)
-        etoile(pixels + j);
+    for (j = 0; j < cols; ++j)
+        etoile(pixels + j, randStP);
 }
 
 
@@ -867,21 +871,22 @@ limbDarken(int *          const irP,
 
 
 static void
-generatePlanetRow(pixel *         const pixelrow,
-                  unsigned int    const row,
-                  unsigned int    const rows,
-                  unsigned int    const cols,
-                  double          const t,
-                  double          const t1,
-                  double *        const u,
-                  double *        const u1,
-                  unsigned char * const cp,
-                  unsigned int *  const bxc,
-                  unsigned int *  const bxf,
-                  int             const byc,
-                  int             const byf,
-                  Vector          const sunvec,
-                  pixval          const maxval) {
+generatePlanetRow(pixel *            const pixelrow,
+                  unsigned int       const row,
+                  unsigned int       const rows,
+                  unsigned int       const cols,
+                  double             const t,
+                  double             const t1,
+                  double *           const u,
+                  double *           const u1,
+                  unsigned char *    const cp,
+                  unsigned int *     const bxc,
+                  unsigned int *     const bxf,
+                  int                const byc,
+                  int                const byf,
+                  Vector             const sunvec,
+                  pixval             const maxval,
+                  struct pm_randSt * const randStP) {
 
     unsigned int const StarClose = 2;
 
@@ -924,24 +929,25 @@ generatePlanetRow(pixel *         const pixelrow,
     /* Left stars */
 
     for (col = 0; (int)col < (int)(cols/2 - (lcos + StarClose)); ++col)
-        etoile(&pixelrow[col]);
+        etoile(&pixelrow[col], randStP);
 
     /* Right stars */
 
     for (col = cols/2 + (lcos + StarClose); col < cols; ++col)
-        etoile(&pixelrow[col]);
+        etoile(&pixelrow[col], randStP);
 }
 
 
 
 static void
-genplanet(bool         const stars,
-          bool         const clouds,
-          float *      const a,
-          unsigned int const cols,
-          unsigned int const rows,
-          unsigned int const n,
-          unsigned int const rseed) {
+genplanet(bool                const stars,
+          bool                const clouds,
+          const double *      const a,
+          unsigned int        const cols,
+          unsigned int        const rows,
+          unsigned int        const n,
+          unsigned int        const rseed,
+          struct pm_randSt *  const randStP) {
 /*----------------------------------------------------------------------------
   Generate planet from elevation array.
 
@@ -971,13 +977,13 @@ genplanet(bool         const stars,
                    clouds ? "clouds" : "planet",
                    rseed, fracdim, powscale, meshsize);
         createPlanetStuff(clouds, a, n, &u, &u1, &bxf, &bxc, &cp, &sunvec,
-                          cols, maxval);
+                          cols, maxval, randStP);
     }
 
     pixelrow = ppm_allocrow(cols);
     for (row = 0; row < rows; ++row) {
         if (stars)
-            generateStarrySkyRow(pixelrow, cols);
+            generateStarrySkyRow(pixelrow, cols, randStP);
         else {
             double const by = (n - 1) * uprj(row, rows);
             int    const byf = floor(by) * n;
@@ -993,7 +999,8 @@ genplanet(bool         const stars,
                 generatePlanetRow(pixelrow, row, rows, cols,
                                   t, t1, u, u1, cp, bxc, bxf, byc, byf,
                                   sunvec,
-                                  maxval);
+                                  maxval,
+                                  randStP);
         }
         ppm_writeppmrow(stdout, pixelrow, cols, maxval, FALSE);
     }
@@ -1010,9 +1017,9 @@ genplanet(bool         const stars,
 
 
 static void
-applyPowerLawScaling(float * const a,
-                     int     const meshsize,
-                     double  const powscale) {
+applyPowerLawScaling(double * const a,
+                     int      const meshsize,
+                     double   const powscale) {
 
     /* Apply power law scaling if non-unity scale is requested. */
 
@@ -1023,7 +1030,7 @@ applyPowerLawScaling(float * const a,
             for (j = 0; j < meshsize; j++) {
                 double const r = Real(a, i, j);
                 if (r > 0)
-                    Real(a, i, j) = pow(r, powscale);
+                    Real(a, i, j) = MIN(hugeVal, pow(r, powscale));
             }
         }
     }
@@ -1032,10 +1039,10 @@ applyPowerLawScaling(float * const a,
 
 
 static void
-computeExtremeReal(const float * const a,
-                   int           const meshsize,
-                   double *      const rminP,
-                   double *      const rmaxP) {
+computeExtremeReal(const double * const a,
+                   int            const meshsize,
+                   double *       const rminP,
+                   double *       const rmaxP) {
 
     /* Compute extrema for autoscaling. */
 
@@ -1061,8 +1068,8 @@ computeExtremeReal(const float * const a,
 
 
 static void
-replaceWithSpread(float * const a,
-                  int     const meshsize) {
+replaceWithSpread(double * const a,
+                  int      const meshsize) {
 /*----------------------------------------------------------------------------
   Replace the real part of each element of the 'a' array with a
   measure of how far the real is from the middle; sort of a standard
@@ -1096,17 +1103,19 @@ planet(unsigned int const cols,
 /*----------------------------------------------------------------------------
    Make a planet.
 -----------------------------------------------------------------------------*/
-    float * a;
+    double * a;
     bool error;
+    struct pm_randSt randSt;
 
-    initgauss(rseed);
+    pm_randinit(&randSt);
+    pm_srand(&randSt, rseed);
 
     if (stars) {
         a = NULL;
         error = FALSE;
     } else {
-        spectralsynth(&a, meshsize, 3.0 - fracdim);
-        if (a == NULL) {
+        spectralsynth(&a, meshsize, 3.0 - fracdim, &randSt);
+        if (!a) {
             error = TRUE;
         } else {
             applyPowerLawScaling(a, meshsize, powscale);
@@ -1117,7 +1126,7 @@ planet(unsigned int const cols,
         }
     }
     if (!error)
-        genplanet(stars, clouds, a, cols, rows, meshsize, rseed);
+        genplanet(stars, clouds, a, cols, rows, meshsize, rseed, &randSt);
 
     if (a != NULL)
         free(a);
@@ -1159,7 +1168,8 @@ main(int argc, const char ** argv) {
     cols = (MAX(cmdline.height, cmdline.width) + 1) & (~1);
     rows = cmdline.height;
 
-    success = planet(cols, rows, cmdline.night, cmdline.clouds, cmdline.seed);
+    success = planet(cols, rows, cmdline.night,
+                     cmdline.clouds, cmdline.seed);
 
     exit(success ? 0 : 1);
 }
diff --git a/generator/ppmmake.c b/generator/ppmmake.c
index 7bac9601..b1de7b52 100644
--- a/generator/ppmmake.c
+++ b/generator/ppmmake.c
@@ -31,11 +31,11 @@ static void
 parseCommandLine(int argc, char ** argv,
                  struct cmdlineInfo * const cmdlineP) {
 /*----------------------------------------------------------------------------
-  Convert program invocation arguments (argc,argv) into a format the 
+  Convert program invocation arguments (argc,argv) into a format the
   program can use easily, struct cmdlineInfo.  Validate arguments along
   the way and exit program with message if invalid.
 
-  Note that some string information we return as *cmdlineP is in the storage 
+  Note that some string information we return as *cmdlineP is in the storage
   argv[] points to.
 -----------------------------------------------------------------------------*/
     optEntry * option_def;
diff --git a/generator/ppmpat.c b/generator/ppmpat.c
index 908c200f..2f71914e 100644
--- a/generator/ppmpat.c
+++ b/generator/ppmpat.c
@@ -14,7 +14,6 @@
 #define _XOPEN_SOURCE 500  /* Make sure strdup() is in string.h */
                            /* get M_PI in math.h */
 #define _BSD_SOURCE  /* Make sure strdup() is in <string.h> */
-#define SPIROGRAPHS 0   /* Spirograph to be added soon */
 
 #include <assert.h>
 #include <math.h>
@@ -23,8 +22,9 @@
 
 #include "pm_c_util.h"
 #include "mallocvar.h"
-#include "shhopt.h"
 #include "nstring.h"
+#include "rand.h"
+#include "shhopt.h"
 #include "ppm.h"
 #include "ppmdraw.h"
 
@@ -39,12 +39,10 @@ typedef enum {
     PAT_POLES,
     PAT_SQUIG,
     PAT_CAMO,
-    PAT_ANTICAMO,
-    PAT_SPIRO1,
-    PAT_SPIRO2,
-    PAT_SPIRO3
+    PAT_ANTICAMO
 } Pattern;
 
+
 typedef struct {
 /*----------------------------------------------------------------------------
    An ordered list of colors with a cursor.
@@ -60,16 +58,79 @@ struct CmdlineInfo {
     /* All the information the user supplied in the command line,
        in a form easy for the program to use.
     */
-    Pattern      basePattern;
-    unsigned int width;
-    unsigned int height;
-    unsigned int colorSpec;
-    ColorTable   colorTable;
-    unsigned int randomseed;
-    unsigned int randomseedSpec;
+    Pattern         basePattern;
+    unsigned int    width;
+    unsigned int    height;
+    unsigned int    colorSpec;
+    ColorTable      colorTable;
+    unsigned int    randomseed;
+    unsigned int    randomseedSpec;
+    ppmd_drawproc * drawProc;
 };
 
 
+
+static pixel
+averageTwoColors(pixel const p1,
+                 pixel const p2) {
+
+    pixel p;
+
+    PPM_ASSIGN(p,
+               (PPM_GETR(p1) + PPM_GETR(p2)) / 2,
+               (PPM_GETG(p1) + PPM_GETG(p2)) / 2,
+               (PPM_GETB(p1) + PPM_GETB(p2)) / 2);
+
+    return p;
+}
+
+
+
+static ppmd_drawproc average_drawproc;
+
+static void
+average_drawproc(pixel **     const pixels,
+                 int          const cols,
+                 int          const rows,
+                 pixval       const maxval,
+                 int          const col,
+                 int          const row,
+                 const void * const clientdata) {
+/*----------------------------------------------------------------------------
+  Reset the pixel's color to the average of the original color and the color
+  indicated by * clientdata.
+-----------------------------------------------------------------------------*/
+
+    if (col >= 0 && col < cols && row >= 0 && row < rows)
+        pixels[row][col] =
+            averageTwoColors(pixels[row][col], *((const pixel*) clientdata));
+}
+
+
+
+static ppmd_drawproc checkerboard_drawproc;
+
+static void
+checkerboard_drawproc(pixel **     const pixels,
+                      int          const cols,
+                      int          const rows,
+                      pixval       const maxval,
+                      int          const col,
+                      int          const row,
+                      const void * const clientdata) {
+/*----------------------------------------------------------------------------
+  If col and row are both even or both odd, do nothing.
+  If one is even and the other is odd, set the pixel's color to that indicated
+  by * clientdata.
+-----------------------------------------------------------------------------*/
+    if (col >= 0 && col < cols && row >= 0 &&
+        row < rows && row % 2 != col % 2)
+
+        pixels[row][col] = *((const pixel*) clientdata);
+}
+
+
+
 static void
 validateColorCount(Pattern      const basePattern,
                    unsigned int const colorCount) {
@@ -80,7 +141,6 @@ validateColorCount(Pattern      const basePattern,
     switch (basePattern) {
     case PAT_GINGHAM2:
     case PAT_ARGYLE1:
-    case PAT_SPIRO1:
         if (colorCount != 2)
             pm_error("Wrong number of colors: %u. "
                      "2 colors are required for the specified pattern.",
@@ -112,8 +172,6 @@ validateColorCount(Pattern      const basePattern,
                      colorCount);
         break;
 
-    case PAT_SPIRO2:
-    case PAT_SPIRO3:
     default:
         pm_error("INTERNAL ERROR.");
     }
@@ -188,9 +246,7 @@ parseCommandLine(int argc, const char ** argv,
     unsigned int squig;
     unsigned int camo;
     unsigned int anticamo;
-    unsigned int spiro1;
-    unsigned int spiro2;
-    unsigned int spiro3;
+    unsigned int meshSpec;
 
     MALLOCARRAY_NOFAIL(option_def, 100);
 
@@ -208,9 +264,9 @@ parseCommandLine(int argc, const char ** argv,
     OPTENT3(0, "tartan",        OPT_FLAG,   NULL,
             &tartan,     0);
     OPTENT3(0, "argyle1",       OPT_FLAG,   NULL,
-            &argyle1,     0);
+            &argyle1,    0);
     OPTENT3(0, "argyle2",       OPT_FLAG,   NULL,
-            &argyle2,     0);
+            &argyle2,    0);
     OPTENT3(0, "poles",         OPT_FLAG,   NULL,
             &poles,      0);
     OPTENT3(0, "squig",         OPT_FLAG,   NULL,
@@ -219,20 +275,12 @@ parseCommandLine(int argc, const char ** argv,
             &camo,       0);
     OPTENT3(0, "anticamo",      OPT_FLAG,   NULL,
             &anticamo,   0);
-#if SPIROGRAPHS != 0
-    OPTENT3(0, "spiro1",        OPT_FLAG,   NULL,
-            &spiro1,     0);
-    OPTENT3(0, "spiro2",        OPT_FLAG,   NULL,
-            &spiro1,     0);
-    OPTENT3(0, "spiro3",        OPT_FLAG,   NULL,
-            &spiro1,     0);
-#else
-    spiro1 = spiro2 = spiro3 = 0;
-#endif
     OPTENT3(0, "color",         OPT_STRINGLIST, &colorText,
             &cmdlineP->colorSpec,           0);
     OPTENT3(0, "randomseed",    OPT_UINT,       &cmdlineP->randomseed,
             &cmdlineP->randomseedSpec,      0);
+    OPTENT3(0, "mesh",          OPT_FLAG,   NULL,
+            &meshSpec,   0);
 
     opt.opt_table = option_def;
     opt.short_allowed = FALSE;  /* We have no short (old-fashioned) options */
@@ -246,8 +294,7 @@ parseCommandLine(int argc, const char ** argv,
         gingham2 + gingham3 + madras + tartan + argyle1 + argyle2 +
         poles +
         squig +
-        camo + anticamo +
-        spiro1 + spiro2 + spiro3;
+        camo + anticamo;
 
     if (basePatternCount < 1)
         pm_error("You must specify a base pattern option such as -gingham2");
@@ -275,12 +322,6 @@ parseCommandLine(int argc, const char ** argv,
             cmdlineP->basePattern = PAT_CAMO;
         else if (anticamo)
             cmdlineP->basePattern = PAT_ANTICAMO;
-        else if (spiro1)
-            cmdlineP->basePattern = PAT_SPIRO1;
-        else if (spiro2)
-            cmdlineP->basePattern = PAT_SPIRO2;
-        else if (spiro3)
-            cmdlineP->basePattern = PAT_SPIRO3;
         else
             assert(false);  /* Every possibility is accounted for */
     }
@@ -291,12 +332,20 @@ parseCommandLine(int argc, const char ** argv,
     } else
         cmdlineP->colorTable.count = 0;
 
+    if (meshSpec) {
+        if (gingham2 + gingham3 + madras + tartan > 0)
+            cmdlineP->drawProc = &checkerboard_drawproc;
+        else
+            pm_message("-mesh ignored (no effect with specified pattern)");
+    } else
+        cmdlineP->drawProc = &average_drawproc;
+
     if (argc-1 != 2)
         pm_error("You must specify 2 non-option arguments: width and height "
                  "in pixels.  You specified %u", argc-1);
     else {
-        cmdlineP->width  = atoi(argv[1]);
-        cmdlineP->height = atoi(argv[2]);
+        cmdlineP->width  = pm_parse_width(argv[1]);
+        cmdlineP->height = pm_parse_height(argv[2]);
 
         if (cmdlineP->width < 1)
             pm_error("Width must be at least 1 pixel");
@@ -339,14 +388,15 @@ validateComputableDimensions(unsigned int const cols,
 
 
 static pixel
-randomColor(pixval const maxval) {
+randomColor(struct pm_randSt * const randStP,
+            pixval             const maxval) {
 
     pixel p;
 
     PPM_ASSIGN(p,
-               rand() % (maxval + 1),
-               rand() % (maxval + 1),
-               rand() % (maxval + 1)
+               pm_rand(randStP) % (maxval + 1),
+               pm_rand(randStP) % (maxval + 1),
+               pm_rand(randStP) % (maxval + 1)
         );
 
     return p;
@@ -354,15 +404,18 @@ randomColor(pixval const maxval) {
 
 
 
-#define DARK_THRESH 0.25
+static double const DARK_THRESH = 0.25;
+
+
 
 static pixel
-randomBrightColor(pixval const maxval) {
+randomBrightColor(struct pm_randSt * const randStP,
+                  pixval             const maxval) {
 
     pixel p;
 
     do {
-        p = randomColor(maxval);
+        p = randomColor(randStP, maxval);
     } while (PPM_LUMIN(p) <= maxval * DARK_THRESH);
 
     return p;
@@ -371,12 +424,13 @@ randomBrightColor(pixval const maxval) {
 
 
 static pixel
-randomDarkColor(pixval const maxval) {
+randomDarkColor(struct pm_randSt * const randStP,
+                pixval             const maxval) {
 
     pixel p;
 
     do {
-        p = randomColor(maxval);
+        p = randomColor(randStP, maxval);
     } while (PPM_LUMIN(p) > maxval * DARK_THRESH);
 
     return p;
@@ -384,40 +438,6 @@ randomDarkColor(pixval const maxval) {
 
 
 
-static pixel
-averageTwoColors(pixel const p1,
-                 pixel const p2) {
-
-    pixel p;
-
-    PPM_ASSIGN(p,
-               (PPM_GETR(p1) + PPM_GETR(p2)) / 2,
-               (PPM_GETG(p1) + PPM_GETG(p2)) / 2,
-               (PPM_GETB(p1) + PPM_GETB(p2)) / 2);
-
-    return p;
-}
-
-
-
-static ppmd_drawproc average_drawproc;
-
-static void
-average_drawproc(pixel **     const pixels,
-                 int          const cols,
-                 int          const rows,
-                 pixval       const maxval,
-                 int          const col,
-                 int          const row,
-                 const void * const clientdata) {
-
-    if (col >= 0 && col < cols && row >= 0 && row < rows)
-        pixels[row][col] =
-            averageTwoColors(pixels[row][col], *((const pixel*) clientdata));
-}
-
-
-
 static void
 nextColor(ColorTable * const colorTableP) {
 /*----------------------------------------------------------------------------
@@ -448,7 +468,8 @@ nextColorBg(ColorTable * const colorTableP) {
 
 
 static pixel
-randomAnticamoColor(pixval const maxval) {
+randomAnticamoColor(struct pm_randSt * const randStP,
+                    pixval             const maxval) {
 
     int v1, v2, v3;
     pixel p;
@@ -457,37 +478,49 @@ randomAnticamoColor(pixval const maxval) {
     v2 = (maxval + 1) / 2;
     v3 = 3 * v1;
 
-    switch (rand() % 15) {
+    switch (pm_rand(randStP) % 15) {
     case 0: case 1:
-        PPM_ASSIGN(p, rand() % v1 + v3, rand() % v2, rand() % v2);
+        PPM_ASSIGN(p, pm_rand(randStP) % v1 + v3,
+                      pm_rand(randStP) % v2,
+                      pm_rand(randStP) % v2);
         break;
 
     case 2:
     case 3:
-        PPM_ASSIGN(p, rand() % v2, rand() % v1 + v3, rand() % v2);
+        PPM_ASSIGN(p, pm_rand(randStP) % v2,
+                      pm_rand(randStP) % v1 + v3,
+                      pm_rand(randStP) % v2);
         break;
 
     case 4:
     case 5:
-        PPM_ASSIGN(p, rand() % v2, rand() % v2, rand() % v1 + v3);
+        PPM_ASSIGN(p, pm_rand(randStP) % v2,
+                      pm_rand(randStP) % v2,
+                      pm_rand(randStP) % v1 + v3);
         break;
 
     case 6:
     case 7:
     case 8:
-        PPM_ASSIGN(p, rand() % v2, rand() % v1 + v3, rand() % v1 + v3);
+        PPM_ASSIGN(p, pm_rand(randStP) % v2,
+                      pm_rand(randStP) % v1 + v3,
+                      pm_rand(randStP) % v1 + v3);
         break;
 
     case 9:
     case 10:
     case 11:
-        PPM_ASSIGN(p, rand() % v1 + v3, rand() % v2, rand() % v1 + v3);
+        PPM_ASSIGN(p, pm_rand(randStP) % v1 + v3,
+                      pm_rand(randStP) % v2,
+                      pm_rand(randStP) % v1 + v3);
         break;
 
     case 12:
     case 13:
     case 14:
-        PPM_ASSIGN(p, rand() % v1 + v3, rand() % v1 + v3, rand() % v2);
+        PPM_ASSIGN(p, pm_rand(randStP) % v1 + v3,
+                      pm_rand(randStP) % v1 + v3,
+                      pm_rand(randStP) % v2);
         break;
     }
 
@@ -497,7 +530,8 @@ randomAnticamoColor(pixval const maxval) {
 
 
 static pixel
-randomCamoColor(pixval const maxval) {
+randomCamoColor(struct pm_randSt * const randStP,
+                pixval             const maxval) {
 
     int const v1 = (maxval + 1 ) / 8;
     int const v2 = (maxval + 1 ) / 4;
@@ -505,31 +539,39 @@ randomCamoColor(pixval const maxval) {
 
     pixel p;
 
-    switch (rand() % 10) {
+    switch (pm_rand(randStP) % 10) {
     case 0:
     case 1:
     case 2:
         /* light brown */
-        PPM_ASSIGN(p, rand() % v3 + v3, rand() % v3 + v2, rand() % v3 + v2);
+        PPM_ASSIGN(p, pm_rand(randStP) % v3 + v3,
+                      pm_rand(randStP) % v3 + v2,
+                      pm_rand(randStP) % v3 + v2);
         break;
 
     case 3:
     case 4:
     case 5:
         /* dark green */
-        PPM_ASSIGN(p, rand() % v2, rand() % v2 + 3 * v1, rand() % v2);
+        PPM_ASSIGN(p, pm_rand(randStP) % v2,
+                      pm_rand(randStP) % v2 + 3 * v1,
+                      pm_rand(randStP) % v2);
         break;
 
     case 6:
     case 7:
         /* brown */
-        PPM_ASSIGN(p, rand() % v2 + v2, rand() % v2, rand() % v2);
+        PPM_ASSIGN(p, pm_rand(randStP) % v2 + v2,
+                      pm_rand(randStP) % v2,
+                      pm_rand(randStP) % v2);
         break;
 
     case 8:
     case 9:
         /* dark brown */
-        PPM_ASSIGN(p, rand() % v1 + v1, rand() % v1, rand() % v1);
+        PPM_ASSIGN(p, pm_rand(randStP) % v1 + v1,
+                      pm_rand(randStP) % v1,
+                      pm_rand(randStP) % v1);
         break;
     }
 
@@ -539,28 +581,30 @@ randomCamoColor(pixval const maxval) {
 
 
 static float
-rnduni(void) {
-    return rand() % 32767 / 32767.0;
+rnduni(struct pm_randSt * const randStP) {
+
+    return pm_rand(randStP) % 32767 / 32767.0;
 }
 
 
 
 static void
-clearBackgroundCamo(pixel **     const pixels,
-                    unsigned int const cols,
-                    unsigned int const rows,
-                    pixval       const maxval,
-                    ColorTable * const colorTableP,
-                    bool         const antiflag) {
+clearBackgroundCamo(pixel **           const pixels,
+                    unsigned int       const cols,
+                    unsigned int       const rows,
+                    pixval             const maxval,
+                    ColorTable *       const colorTableP,
+                    struct pm_randSt * const randStP,
+                    bool               const antiflag) {
 
     pixel color;
 
     if (colorTableP->count > 0) {
         color = colorTableP->color[0];
     } else if (antiflag)
-        color = randomAnticamoColor(maxval);
+        color = randomAnticamoColor(randStP, maxval);
     else
-        color = randomCamoColor(maxval);
+        color = randomCamoColor(randStP, maxval);
 
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, 0, cols, rows, PPMD_NULLDRAWPROC,
@@ -570,13 +614,14 @@ clearBackgroundCamo(pixel **     const pixels,
 
 
 static void
-camoFill(pixel **         const pixels,
-         unsigned int     const cols,
-         unsigned int     const rows,
-         pixval           const maxval,
-         struct fillobj * const fh,
-         ColorTable     * const colorTableP,
-         bool             const antiflag) {
+camoFill(pixel **           const pixels,
+         unsigned int       const cols,
+         unsigned int       const rows,
+         pixval             const maxval,
+         struct fillobj *   const fh,
+         ColorTable     *   const colorTableP,
+         struct pm_randSt * const randStP,
+         bool               const antiflag) {
 
     pixel color;
 
@@ -585,9 +630,9 @@ camoFill(pixel **         const pixels,
         color = colorTableP->color[colorTableP->index];
         nextColorBg(colorTableP);
     } else if (antiflag)
-        color = randomAnticamoColor(maxval);
+        color = randomAnticamoColor(randStP, maxval);
     else
-        color = randomCamoColor(maxval);
+        color = randomCamoColor(randStP, maxval);
 
     ppmd_fill(pixels, cols, rows, maxval, fh, PPMD_NULLDRAWPROC, &color);
 }
@@ -608,25 +653,29 @@ camoFill(pixel **         const pixels,
 
 
 static void
-computeXsYs(int *        const xs,
-            int *        const ys,
-            unsigned int const cols,
-            unsigned int const rows,
-            unsigned int const pointCt) {
-
-    unsigned int const cx = rand() % cols;
-    unsigned int const cy = rand() % rows;
-    double const a = rnduni() * (MAX_ELLIPSE_FACTOR - MIN_ELLIPSE_FACTOR) +
-        MIN_ELLIPSE_FACTOR;
-    double const b = rnduni() * (MAX_ELLIPSE_FACTOR - MIN_ELLIPSE_FACTOR) +
-        MIN_ELLIPSE_FACTOR;
-    double const theta = rnduni() * 2.0 * M_PI;
+computeXsYs(int *              const xs,
+            int *              const ys,
+            unsigned int       const cols,
+            unsigned int       const rows,
+            unsigned int       const pointCt,
+            struct pm_randSt * const randStP) {
+
+    unsigned int const cx = pm_rand(randStP) % cols;
+    unsigned int const cy = pm_rand(randStP) % rows;
+    double const a = rnduni(randStP) *
+                       (MAX_ELLIPSE_FACTOR - MIN_ELLIPSE_FACTOR) +
+                        MIN_ELLIPSE_FACTOR;
+    double const b = rnduni(randStP) *
+                       (MAX_ELLIPSE_FACTOR - MIN_ELLIPSE_FACTOR) +
+                        MIN_ELLIPSE_FACTOR;
+    double const theta = rnduni(randStP) * 2.0 * M_PI;
 
     unsigned int p;
 
     for (p = 0; p < pointCt; ++p) {
-        double const c = rnduni() * (MAX_POINT_FACTOR - MIN_POINT_FACTOR) +
-            MIN_POINT_FACTOR;
+        double const c = rnduni(randStP) *
+                           (MAX_POINT_FACTOR - MIN_POINT_FACTOR) +
+                            MIN_POINT_FACTOR;
         double const tx   = a * sin(p * 2.0 * M_PI / pointCt);
         double const ty   = b * cos(p * 2.0 * M_PI / pointCt);
         double const tang = atan2(ty, tx) + theta;
@@ -638,18 +687,20 @@ computeXsYs(int *        const xs,
 
 
 static void
-camo(pixel **     const pixels,
-     unsigned int const cols,
-     unsigned int const rows,
-     ColorTable * const colorTableP,
-     pixval       const maxval,
-     bool         const antiflag) {
+camo(pixel **           const pixels,
+     unsigned int       const cols,
+     unsigned int       const rows,
+     ColorTable *       const colorTableP,
+     struct pm_randSt * const randStP,
+     pixval             const maxval,
+     bool               const antiflag) {
 
     unsigned int const n = (rows * cols) / SQR(BLOBRAD) * 5;
 
     unsigned int i;
 
-    clearBackgroundCamo(pixels, cols, rows, maxval, colorTableP, antiflag);
+    clearBackgroundCamo(pixels, cols, rows, maxval,
+                        colorTableP, randStP, antiflag);
 
     if (colorTableP->count > 0) {
         assert(colorTableP->count > 1);
@@ -658,13 +709,13 @@ camo(pixel **     const pixels,
 
     for (i = 0; i < n; ++i) {
         unsigned int const pointCt =
-            rand() % (MAX_POINTS - MIN_POINTS + 1) + MIN_POINTS;
+            pm_rand(randStP) % (MAX_POINTS - MIN_POINTS + 1) + MIN_POINTS;
 
         int xs[MAX_POINTS], ys[MAX_POINTS];
         int x0, y0;
         struct fillobj * fh;
 
-        computeXsYs(xs, ys, cols, rows, pointCt);
+        computeXsYs(xs, ys, cols, rows, pointCt, randStP);
 
         x0 = (xs[0] + xs[pointCt - 1]) / 2;
         y0 = (ys[0] + ys[pointCt - 1]) / 2;
@@ -675,7 +726,8 @@ camo(pixel **     const pixels,
             pixels, cols, rows, maxval, x0, y0, pointCt, xs, ys, x0, y0,
             ppmd_fill_drawproc, fh);
 
-        camoFill(pixels, cols, rows, maxval, fh, colorTableP, antiflag);
+        camoFill(pixels, cols, rows, maxval, fh,
+                 colorTableP, randStP, antiflag);
 
         ppmd_fill_destroy(fh);
     }
@@ -688,17 +740,21 @@ camo(pixel **     const pixels,
 -----------------------------------------------------------------------------*/
 
 static void
-gingham2(pixel **     const pixels,
-         unsigned int const cols,
-         unsigned int const rows,
-         ColorTable   const colorTable,
-         pixval       const maxval) {
+gingham2(pixel **           const pixels,
+         unsigned int       const cols,
+         unsigned int       const rows,
+         ColorTable         const colorTable,
+         struct pm_randSt * const randStP,
+         ppmd_drawproc    * const drawproc,
+         pixval             const maxval) {
 
     bool  const colorSpec = (colorTable.count > 0);
     pixel const backcolor = colorSpec ?
-                            colorTable.color[0] : randomDarkColor(maxval);
+                              colorTable.color[0] :
+                              randomDarkColor(randStP, maxval);
     pixel const forecolor = colorSpec ?
-                            colorTable.color[1] : randomBrightColor(maxval);
+                              colorTable.color[1] :
+                              randomBrightColor(randStP, maxval);
     unsigned int const colso2 = cols / 2;
     unsigned int const rowso2 = rows / 2;
 
@@ -712,29 +768,34 @@ gingham2(pixel **     const pixels,
 
     /* Woof. */
     ppmd_filledrectangle(
-        pixels, cols, rows, maxval, 0, 0, cols, rowso2, average_drawproc,
+        pixels, cols, rows, maxval, 0, 0, cols, rowso2, drawproc,
         &backcolor);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, rowso2, cols, rows - rowso2,
-        average_drawproc, &forecolor);
+        drawproc, &forecolor);
 }
 
 
 
 static void
-gingham3(pixel **     const pixels,
-         unsigned int const cols,
-         unsigned int const rows,
-         ColorTable   const colorTable,
-         pixval       const maxval) {
+gingham3(pixel **           const pixels,
+         unsigned int       const cols,
+         unsigned int       const rows,
+         ColorTable         const colorTable,
+         struct pm_randSt * const randStP,
+         ppmd_drawproc    * const drawproc,
+         pixval             const maxval) {
 
     bool  const colorSpec = (colorTable.count > 0);
-    pixel const backcolor = colorSpec ?
-                            colorTable.color[0] : randomDarkColor(maxval);
+    pixel const backcolor  = colorSpec ?
+                               colorTable.color[0] :
+                               randomDarkColor(randStP, maxval);
     pixel const fore1color = colorSpec ?
-                            colorTable.color[1] : randomBrightColor(maxval);
+                               colorTable.color[1] :
+                               randomBrightColor(randStP, maxval);
     pixel const fore2color = colorSpec ?
-                            colorTable.color[2] : randomBrightColor(maxval);
+                               colorTable.color[2] :
+                               randomBrightColor(randStP, maxval);
     unsigned int const colso4 = cols / 4;
     unsigned int const rowso4 = rows / 4;
 
@@ -754,35 +815,40 @@ gingham3(pixel **     const pixels,
 
     /* Woof. */
     ppmd_filledrectangle(
-        pixels, cols, rows, maxval, 0, 0, cols, rowso4, average_drawproc,
+        pixels, cols, rows, maxval, 0, 0, cols, rowso4, drawproc,
         &backcolor);
     ppmd_filledrectangle(
-        pixels, cols, rows, maxval, 0, rowso4, cols, rowso4, average_drawproc,
+        pixels, cols, rows, maxval, 0, rowso4, cols, rowso4, drawproc,
         &fore1color);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, 2 * rowso4, cols, rowso4,
-        average_drawproc, &fore2color);
+        drawproc, &fore2color);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, 3 * rowso4, cols, rows - rowso4,
-        average_drawproc, &fore1color);
+        drawproc, &fore1color);
 }
 
 
 
 static void
-madras(pixel **     const pixels,
-       unsigned int const cols,
-       unsigned int const rows,
-       ColorTable   const colorTable,
-       pixval       const maxval) {
+madras(pixel **           const pixels,
+       unsigned int       const cols,
+       unsigned int       const rows,
+       ColorTable         const colorTable,
+       struct pm_randSt * const randStP,
+       ppmd_drawproc    * const drawproc,
+       pixval             const maxval) {
 
     bool  const colorSpec = (colorTable.count > 0);
-    pixel const backcolor = colorSpec ?
-                            colorTable.color[0] : randomDarkColor(maxval);
+    pixel const backcolor  = colorSpec ?
+                               colorTable.color[0] :
+                               randomDarkColor(randStP, maxval);
     pixel const fore1color = colorSpec ?
-                            colorTable.color[1] : randomBrightColor(maxval);
+                               colorTable.color[1] :
+                               randomBrightColor(randStP, maxval);
     pixel const fore2color = colorSpec ?
-                            colorTable.color[2] : randomBrightColor(maxval);
+                               colorTable.color[2] :
+                               randomBrightColor(randStP, maxval);
     unsigned int const cols2  = cols * 2 / 44;
     unsigned int const rows2  = rows * 2 / 44;
     unsigned int const cols3  = cols * 3 / 44;
@@ -846,72 +912,77 @@ madras(pixel **     const pixels,
 
     /* Woof. */
     ppmd_filledrectangle(
-        pixels, cols, rows, maxval, 0, 0, cols, rows2, average_drawproc,
+        pixels, cols, rows, maxval, 0, 0, cols, rows2, drawproc,
         &backcolor);
     ppmd_filledrectangle(
-        pixels, cols, rows, maxval, 0, rows2, cols, rows3, average_drawproc,
+        pixels, cols, rows, maxval, 0, rows2, cols, rows3, drawproc,
         &fore2color);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, rows2 + rows3, cols, rows2,
-        average_drawproc, &backcolor);
+        drawproc, &backcolor);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, 2 * rows2 + rows3, cols, rows2,
-        average_drawproc, &fore1color);
+        drawproc, &fore1color);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, 3 * rows2 + rows3, cols, rows2,
-        average_drawproc, &backcolor);
+        drawproc, &backcolor);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, 4 * rows2 + rows3, cols, rows6a,
-        average_drawproc, &fore2color);
+        drawproc, &fore2color);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, 4 * rows2 + rows3 + rows6a, cols, rows2,
-        average_drawproc, &backcolor);
+        drawproc, &backcolor);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, 5 * rows2 + rows3 + rows6a, cols, rows3,
-        average_drawproc, &fore1color);
+        drawproc, &fore1color);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, 5 * rows2 + 2 * rows3 + rows6a, cols,
-        rows2, average_drawproc, &backcolor);
+        rows2, drawproc, &backcolor);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, 6 * rows2 + 2 * rows3 + rows6a, cols,
-        rows3, average_drawproc, &fore1color);
+        rows3, drawproc, &fore1color);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, 6 * rows2 + 3 * rows3 + rows6a, cols,
-        rows2, average_drawproc, &backcolor);
+        rows2, drawproc, &backcolor);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, 7 * rows2 + 3 * rows3 + rows6a, cols,
-        rows6b, average_drawproc, &fore2color);
+        rows6b, drawproc, &fore2color);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, 7 * rows2 + 3 * rows3 + rows6a + rows6b,
-        cols, rows2, average_drawproc, &backcolor);
+        cols, rows2, drawproc, &backcolor);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, 8 * rows2 + 3 * rows3 + rows6a + rows6b,
-        cols, rows2, average_drawproc, &fore1color);
+        cols, rows2, drawproc, &fore1color);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, 9 * rows2 + 3 * rows3 + rows6a + rows6b,
-        cols, rows2, average_drawproc, &backcolor);
+        cols, rows2, drawproc, &backcolor);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0,
         10 * rows2 + 3 * rows3 + rows6a + rows6b,
-        cols, rows3, average_drawproc, &fore2color);
+        cols, rows3, drawproc, &fore2color);
 }
 
 
 
 static void
-tartan(pixel **     const pixels,
-       unsigned int const cols,
-       unsigned int const rows,
-       ColorTable   const colorTable,
-       pixval       const maxval) {
+tartan(pixel **           const pixels,
+       unsigned int       const cols,
+       unsigned int       const rows,
+       ColorTable         const colorTable,
+       struct pm_randSt * const randStP,
+       ppmd_drawproc    * const drawproc,
+       pixval             const maxval) {
 
     bool  const colorSpec = (colorTable.count > 0);
-    pixel const backcolor = colorSpec ?
-                            colorTable.color[0] : randomDarkColor(maxval);
+    pixel const backcolor  = colorSpec ?
+                               colorTable.color[0] :
+                               randomDarkColor(randStP, maxval);
     pixel const fore1color = colorSpec ?
-                            colorTable.color[1] : randomBrightColor(maxval);
+                               colorTable.color[1] :
+                               randomBrightColor(randStP, maxval);
     pixel const fore2color = colorSpec ?
-                            colorTable.color[2] : randomBrightColor(maxval);
+                               colorTable.color[2] :
+                               randomBrightColor(randStP, maxval);
     unsigned int const cols1  = cols / 22;
     unsigned int const rows1  = rows / 22;
     unsigned int const cols3  = cols * 3 / 22;
@@ -951,29 +1022,29 @@ tartan(pixel **     const pixels,
 
     /* Woof. */
     ppmd_filledrectangle(
-        pixels, cols, rows, maxval, 0, 0, cols, rows5a, average_drawproc,
+        pixels, cols, rows, maxval, 0, 0, cols, rows5a, drawproc,
         &backcolor);
     ppmd_filledrectangle(
-        pixels, cols, rows, maxval, 0, rows5a, cols, rows1, average_drawproc,
+        pixels, cols, rows, maxval, 0, rows5a, cols, rows1, drawproc,
         &fore1color);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, rows5a + rows1, cols, rows5b,
-        average_drawproc, &backcolor);
+        drawproc, &backcolor);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, rows10 + rows1, cols, rows3,
-        average_drawproc, &fore2color);
+        drawproc, &fore2color);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, rows10 + rows1 + rows3, cols, rows1,
-        average_drawproc, &backcolor);
+        drawproc, &backcolor);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, rows10 + 2 * rows1 + rows3, cols, rows3,
-        average_drawproc, &fore2color);
+        drawproc, &fore2color);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, rows10 + 2 * rows1 + 2 * rows3, cols,
-        rows1, average_drawproc, &backcolor);
+        rows1, drawproc, &backcolor);
     ppmd_filledrectangle(
         pixels, cols, rows, maxval, 0, rows10 + 3 * rows1 + 2 * rows3, cols,
-        rows3, average_drawproc, &fore2color);
+        rows3, drawproc, &fore2color);
 }
 
 
@@ -1013,14 +1084,15 @@ argyle(pixel **     const pixels,
        unsigned int const cols,
        unsigned int const rows,
        ColorTable   const colorTable,
+       struct pm_randSt * const randStP,
        pixval       const maxval,
        bool         const stripes) {
 
     bool  const colorSpec = (colorTable.count > 0);
     pixel const backcolor = colorSpec ?
-        colorTable.color[0] : randomDarkColor(maxval);
+        colorTable.color[0] : randomDarkColor(randStP, maxval);
     pixel const forecolor = colorSpec ?
-        colorTable.color[1] : randomBrightColor(maxval);
+        colorTable.color[1] : randomBrightColor(randStP, maxval);
 
     /* Fill canvas with background to start */
     ppmd_filledrectangle(
@@ -1032,7 +1104,8 @@ argyle(pixel **     const pixels,
     if (stripes) {
          /* Connect corners with thin stripes */
          pixel const stripecolor =
-             colorSpec ? colorTable.color[2] : randomBrightColor(maxval);
+             colorSpec ? colorTable.color[2] :
+                         randomBrightColor(randStP, maxval);
 
          ppmd_line(pixels, cols, rows, maxval, 0, 0, cols-1, rows-1,
               PPMD_NULLDRAWPROC, (char *) &stripecolor);
@@ -1054,27 +1127,28 @@ argyle(pixel **     const pixels,
 
 
 static void
-placeAndColorPolesRandomly(int *        const xs,
-                           int *        const ys,
-                           pixel *      const colors,
-                           unsigned int const cols,
-                           unsigned int const rows,
-                           pixval       const maxval,
-                           ColorTable * const colorTableP,
-                           unsigned int const poleCt) {
+placeAndColorPolesRandomly(int *              const xs,
+                           int *              const ys,
+                           pixel *            const colors,
+                           unsigned int       const cols,
+                           unsigned int       const rows,
+                           pixval             const maxval,
+                           ColorTable *       const colorTableP,
+                           struct pm_randSt * const randStP,
+                           unsigned int       const poleCt) {
 
     unsigned int i;
 
     for (i = 0; i < poleCt; ++i) {
 
-        xs[i] = rand() % cols;
-        ys[i] = rand() % rows;
+        xs[i] = pm_rand(randStP) % cols;
+        ys[i] = pm_rand(randStP) % rows;
 
         if (colorTableP->count > 0) {
             colors[i] = colorTableP->color[colorTableP->index];
             nextColor(colorTableP);
         } else
-            colors[i] = randomBrightColor(maxval);
+            colors[i] = randomBrightColor(randStP, maxval);
     }
 }
 
@@ -1104,11 +1178,12 @@ assignInterpolatedColor(pixel * const resultP,
 
 
 static void
-poles(pixel **     const pixels,
-      unsigned int const cols,
-      unsigned int const rows,
-      ColorTable * const colorTableP,
-      pixval       const maxval) {
+poles(pixel **           const pixels,
+      unsigned int       const cols,
+      unsigned int       const rows,
+      ColorTable *       const colorTableP,
+      struct pm_randSt * const randStP,
+      pixval             const maxval) {
 
     unsigned int const poleCt = MAX(2, MIN(MAXPOLES, cols * rows / 30000));
 
@@ -1117,7 +1192,7 @@ poles(pixel **     const pixels,
     unsigned int row;
 
     placeAndColorPolesRandomly(xs, ys, colors, cols, rows, maxval,
-                               colorTableP, poleCt);
+                               colorTableP, randStP, poleCt);
 
     /* Interpolate points */
 
@@ -1236,11 +1311,12 @@ sqRainbowCircleDrawproc(pixel **     const pixels,
 
 
 static void
-chooseSqPoleColors(ColorTable * const colorTableP,
-                   pixval       const maxval,
-                   pixel *      const color1P,
-                   pixel *      const color2P,
-                   pixel *      const color3P) {
+chooseSqPoleColors(ColorTable *       const colorTableP,
+                   pixval             const maxval,
+                   pixel *            const color1P,
+                   pixel *            const color2P,
+                   pixel *            const color3P,
+                   struct pm_randSt * const randStP) {
 
     if (colorTableP->count > 0) {
         *color1P = colorTableP->color[colorTableP->index];
@@ -1250,19 +1326,20 @@ chooseSqPoleColors(ColorTable * const colorTableP,
         *color3P = colorTableP->color[colorTableP->index];
         nextColor(colorTableP);
     } else {
-        *color1P = randomBrightColor(maxval);
-        *color2P = randomBrightColor(maxval);
-        *color3P = randomBrightColor(maxval);
+        *color1P = randomBrightColor(randStP, maxval);
+        *color2P = randomBrightColor(randStP, maxval);
+        *color3P = randomBrightColor(randStP, maxval);
     }
 }
 
 
 
 static void
-sqAssignColors(unsigned int const circlecount,
-               pixval       const maxval,
-               ColorTable * const colorTableP,
-               pixel *      const colors) {
+sqAssignColors(unsigned int       const circlecount,
+               pixval             const maxval,
+               ColorTable *       const colorTableP,
+               pixel *            const colors,
+               struct pm_randSt * const randStP) {
 
     float const cco3 = (circlecount - 1) / 3.0;
 
@@ -1271,7 +1348,7 @@ sqAssignColors(unsigned int const circlecount,
     pixel rc3;
     unsigned int i;
 
-    chooseSqPoleColors(colorTableP, maxval, &rc1, &rc2, &rc3);
+    chooseSqPoleColors(colorTableP, maxval, &rc1, &rc2, &rc3, randStP);
 
     for (i = 0; i < circlecount; ++i) {
         if (i < cco3) {
@@ -1333,24 +1410,25 @@ clearBackgroundSquig(pixel **     const pixels,
 
 
 static void
-chooseWrapAroundPoint(unsigned int const cols,
-                      unsigned int const rows,
-                      ppmd_point * const pFirstP,
-                      ppmd_point * const pLastP,
-                      ppmd_point * const p0P,
-                      ppmd_point * const p1P,
-                      ppmd_point * const p2P,
-                      ppmd_point * const p3P) {
-
-    switch (rand() % 4) {
+chooseWrapAroundPoint(unsigned int       const cols,
+                      unsigned int       const rows,
+                      ppmd_point *       const pFirstP,
+                      ppmd_point *       const pLastP,
+                      ppmd_point *       const p0P,
+                      ppmd_point *       const p1P,
+                      ppmd_point *       const p2P,
+                      ppmd_point *       const p3P,
+                      struct pm_randSt * const randStP) {
+
+    switch (pm_rand(randStP) % 4) {
     case 0:
-        p1P->x = rand() % cols;
+        p1P->x = pm_rand(randStP) % cols;
         p1P->y = 0;
         if (p1P->x < cols / 2)
-            pFirstP->x = rand() % (p1P->x * 2 + 1);
+            pFirstP->x = pm_rand(randStP) % (p1P->x * 2 + 1);
         else
-            pFirstP->x = cols - 1 - rand() % ((cols - p1P->x) * 2);
-        pFirstP->y = rand() % rows;
+            pFirstP->x = cols - 1 - pm_rand(randStP) % ((cols - p1P->x) * 2);
+        pFirstP->y = pm_rand(randStP) % rows;
         p2P->x = p1P->x;
         p2P->y = rows - 1;
         pLastP->x = 2 * p2P->x - pFirstP->x;
@@ -1362,13 +1440,13 @@ chooseWrapAroundPoint(unsigned int const cols,
         break;
 
     case 1:
-        p2P->x = rand() % cols;
+        p2P->x = pm_rand(randStP) % cols;
         p2P->y = 0;
         if (p2P->x < cols / 2)
-            pLastP->x = rand() % (p2P->x * 2 + 1);
+            pLastP->x = pm_rand(randStP) % (p2P->x * 2 + 1);
         else
-            pLastP->x = cols - 1 - rand() % ((cols - p2P->x) * 2);
-        pLastP->y = rand() % rows;
+            pLastP->x = cols - 1 - pm_rand(randStP) % ((cols - p2P->x) * 2);
+        pLastP->y = pm_rand(randStP) % rows;
         p1P->x = p2P->x;
         p1P->y = rows - 1;
         pFirstP->x = 2 * p1P->x - pLastP->x;
@@ -1381,12 +1459,12 @@ chooseWrapAroundPoint(unsigned int const cols,
 
     case 2:
         p1P->x = 0;
-        p1P->y = rand() % rows;
-        pFirstP->x = rand() % cols;
+        p1P->y = pm_rand(randStP) % rows;
+        pFirstP->x = pm_rand(randStP) % cols;
         if (p1P->y < rows / 2)
-            pFirstP->y = rand() % (p1P->y * 2 + 1);
+            pFirstP->y = pm_rand(randStP) % (p1P->y * 2 + 1);
         else
-            pFirstP->y = rows - 1 - rand() % ((rows - p1P->y) * 2);
+            pFirstP->y = rows - 1 - pm_rand(randStP) % ((rows - p1P->y) * 2);
         p2P->x = cols - 1;
         p2P->y = p1P->y;
         pLastP->x = p2P->x - pFirstP->x;
@@ -1399,12 +1477,12 @@ chooseWrapAroundPoint(unsigned int const cols,
 
     case 3:
         p2P->x = 0;
-        p2P->y = rand() % rows;
-        pLastP->x = rand() % cols;
+        p2P->y = pm_rand(randStP) % rows;
+        pLastP->x = pm_rand(randStP) % cols;
         if (p2P->y < rows / 2)
-            pLastP->y = rand() % (p2P->y * 2 + 1);
+            pLastP->y = pm_rand(randStP) % (p2P->y * 2 + 1);
         else
-            pLastP->y = rows - 1 - rand() % ((rows - p2P->y) * 2);
+            pLastP->y = rows - 1 - pm_rand(randStP) % ((rows - p2P->y) * 2);
         p1P->x = cols - 1;
         p1P->y = p2P->y;
         pFirstP->x = p1P->x - pLastP->x;
@@ -1420,11 +1498,12 @@ chooseWrapAroundPoint(unsigned int const cols,
 
 
 static void
-squig(pixel **     const pixels,
-      unsigned int const cols,
-      unsigned int const rows,
-      ColorTable * const colorTableP,
-      pixval       const maxval) {
+squig(pixel **           const pixels,
+      unsigned int       const cols,
+      unsigned int       const rows,
+      ColorTable *       const colorTableP,
+      struct pm_randSt * const randStP,
+      pixval             const maxval) {
 
     int i;
 
@@ -1453,10 +1532,11 @@ squig(pixel **     const pixels,
         ppmd_circlep(pixels, cols, rows, maxval,
                      ppmd_makePoint(0, 0), radius,
                      sqMeasureCircleDrawproc, &sqClientData);
-        sqAssignColors(squig.circleCt, maxval, colorTableP, squig.color);
+        sqAssignColors(squig.circleCt, maxval, colorTableP, squig.color,
+                       randStP);
 
         chooseWrapAroundPoint(cols, rows, &c[0], &c[SQ_POINTS-1],
-                              &p0, &p1, &p2, &p3);
+                              &p0, &p1, &p2, &p3, randStP);
 
         {
             /* Do the middle points */
@@ -1466,8 +1546,8 @@ squig(pixel **     const pixels,
               /* validateSquigAspect() assures that
                  cols - 2 * radius, rows -2 * radius are positive
               */
-                c[j].x = (rand() % (cols - 2 * radius)) + radius;
-                c[j].y = (rand() % (rows - 2 * radius)) + radius;
+                c[j].x = (pm_rand(randStP) % (cols - 2 * radius)) + radius;
+                c[j].y = (pm_rand(randStP) % (rows - 2 * radius)) + radius;
             }
         }
 
@@ -1490,6 +1570,7 @@ main(int argc, const char ** argv) {
 
     struct CmdlineInfo cmdline;
     pixel ** pixels;
+    struct pm_randSt randSt;
 
     pm_proginit(&argc, argv);
 
@@ -1497,65 +1578,68 @@ main(int argc, const char ** argv) {
 
     validateComputableDimensions(cmdline.width, cmdline.height);
 
-    srand(cmdline.randomseedSpec ? cmdline.randomseed : pm_randseed());
+    pm_randinit(&randSt);
+    pm_srand2(&randSt, cmdline.randomseedSpec, cmdline.randomseed);
 
     pixels = ppm_allocarray(cmdline.width, cmdline.height);
 
     switch (cmdline.basePattern) {
     case PAT_GINGHAM2:
         gingham2(pixels, cmdline.width, cmdline.height,
-                 cmdline.colorTable, PPM_MAXMAXVAL);
+                 cmdline.colorTable, &randSt, cmdline.drawProc, PPM_MAXMAXVAL);
         break;
 
     case PAT_GINGHAM3:
         gingham3(pixels, cmdline.width, cmdline.height,
-                 cmdline.colorTable, PPM_MAXMAXVAL);
+                 cmdline.colorTable, &randSt, cmdline.drawProc, PPM_MAXMAXVAL);
         break;
 
     case PAT_MADRAS:
         madras(pixels, cmdline.width, cmdline.height,
-               cmdline.colorTable, PPM_MAXMAXVAL);
+               cmdline.colorTable, &randSt, cmdline.drawProc, PPM_MAXMAXVAL);
         break;
 
     case PAT_TARTAN:
         tartan(pixels, cmdline.width, cmdline.height,
-               cmdline.colorTable, PPM_MAXMAXVAL);
+               cmdline.colorTable, &randSt, cmdline.drawProc, PPM_MAXMAXVAL);
         break;
 
     case PAT_ARGYLE1:
         argyle(pixels, cmdline.width, cmdline.height,
-               cmdline.colorTable, PPM_MAXMAXVAL, FALSE);
+               cmdline.colorTable, &randSt, PPM_MAXMAXVAL, FALSE);
         break;
 
     case PAT_ARGYLE2:
         argyle(pixels, cmdline.width, cmdline.height,
-               cmdline.colorTable, PPM_MAXMAXVAL, TRUE);
+               cmdline.colorTable, &randSt, PPM_MAXMAXVAL, TRUE);
         break;
 
     case PAT_POLES:
         poles(pixels, cmdline.width, cmdline.height,
-              &cmdline.colorTable, PPM_MAXMAXVAL);
+              &cmdline.colorTable, &randSt, PPM_MAXMAXVAL);
         break;
 
     case PAT_SQUIG:
         squig(pixels, cmdline.width, cmdline.height,
-              &cmdline.colorTable, PPM_MAXMAXVAL);
+              &cmdline.colorTable, &randSt, PPM_MAXMAXVAL);
         break;
 
     case PAT_CAMO:
         camo(pixels, cmdline.width, cmdline.height,
-             &cmdline.colorTable, PPM_MAXMAXVAL, 0);
+             &cmdline.colorTable, &randSt, PPM_MAXMAXVAL, 0);
         break;
 
     case PAT_ANTICAMO:
         camo(pixels, cmdline.width, cmdline.height,
-             &cmdline.colorTable, PPM_MAXMAXVAL, 1);
+             &cmdline.colorTable, &randSt, PPM_MAXMAXVAL, 1);
         break;
 
     default:
         pm_error("can't happen!");
     }
 
+    pm_randterm(&randSt);
+
     ppm_writeppm(stdout, pixels, cmdline.width, cmdline.height,
                  PPM_MAXMAXVAL, 0);
 
diff --git a/generator/ppmrainbow b/generator/ppmrainbow
index e8a329ff..68d519a0 100755
--- a/generator/ppmrainbow
+++ b/generator/ppmrainbow
@@ -26,6 +26,19 @@ exec perl -w -x -S -- "$0" "$@"
 use strict;
 use Getopt::Long;
 use File::Temp;
+use IO::Handle;
+
+
+
+sub pm_message($) {
+    STDERR->print("ppmrainbow: $_[0]\n");
+}
+
+sub pm_error($) {
+    pm_message($_[0]);
+    exit(1);
+}
+
 
 my ($FALSE, $TRUE) = (0,1);
 
@@ -46,15 +59,6 @@ sub doVersionHack($) {
 
 
 
-sub fatal($) {
-    my ($msg) = @_;
-
-    print(STDERR "ppmrainbow: $msg\n");
-    exit(1);
-}
-
-
-
 ##############################################################################
 #
 #                                 MAINLINE
@@ -79,14 +83,14 @@ GetOptions("width=i"   => \$Twid,
            "verbose!"  => \$verbose);
 
 if ($Twid < 1 || $Thgt < 1) {
-    fatal("invalid width and/or height");
+    pm_error("invalid width and/or height");
 }
 my $verboseCommand = $verbose ? "set -x;" : "";
 
 if (@ARGV < 1) {
-    fatal("You must specify at least one color as an argument");
+    pm_error("You must specify at least one color as an argument");
 } elsif (@ARGV < 2 && ! $repeat) {
-    fatal("With the -norepeat option, you must specify at least two colors " .
+    pm_error("With the -norepeat option, you must specify at least two colors " .
           "as arguments.");
 }
 
@@ -115,15 +119,18 @@ while (@colorlist >= 2) {
     my $rc = system("$verboseCommand pgmramp -lr $w $Thgt | " .
                     "pgmtoppm \"$colorlist[0]-$colorlist[1]\" >$outfile");
     if ($rc != 0) {
-        fatal("pgmramp|pgmtoppm pipe failed.");
+        pm_error("pgmramp|pgmtoppm pipe failed.");
     }
     $widthRemaining -= $w;
     $n++;
     shift @colorlist;
 }
 
-0 == system qq{$verboseCommand pnmcat -lr @outlist}
-    or exit 1;
+my $termStat =
+    system("$verboseCommand pamcat -leftright @outlist");
+if ($termStat != 0) {
+    exit 1;
+}
 
 exit 0;
 
diff --git a/generator/ppmrough.c b/generator/ppmrough.c
index c87a0364..d9f1bb65 100644
--- a/generator/ppmrough.c
+++ b/generator/ppmrough.c
@@ -16,11 +16,10 @@
 
 #include "pm_c_util.h"
 #include "mallocvar.h"
+#include "rand.h"
 #include "shhopt.h"
 #include "ppm.h"
 
-static pixel** PIX;
-static pixval BG_RED, BG_GREEN, BG_BLUE;
 
 
 struct CmdlineInfo {
@@ -28,9 +27,12 @@ struct CmdlineInfo {
      in a form easy for the program to use.
   */
     unsigned int left, right, top, bottom;
-    unsigned int width, height, var;
-    const char * bg_rgb;
-    const char * fg_rgb;
+    unsigned int leftSpec, rightSpec, topSpec, bottomSpec;
+    unsigned int width;
+    unsigned int height;
+    unsigned int var;
+    const char * bg;  /* Null if not specified */
+    const char * fg;  /* Null if not specified */
     unsigned int randomseed;
     unsigned int randomseedSpec;
     unsigned int verbose;
@@ -42,6 +44,8 @@ static void
 parseCommandLine(int argc, const char ** argv,
                  struct CmdlineInfo * const cmdlineP) {
 
+    unsigned int widthSpec, heightSpec, bgSpec, fgSpec, varSpec;
+
     optEntry * option_def;
         /* Instructions to OptParseOptions2 on how to parse our options.    */
     optStruct3 opt;
@@ -51,28 +55,30 @@ parseCommandLine(int argc, const char ** argv,
     MALLOCARRAY_NOFAIL(option_def, 100);
 
     option_def_index = 0;   /* incremented by OPTENTRY */
-    OPTENT3(0, "width",       OPT_UINT,   &cmdlineP->width,   NULL, 0);
-    OPTENT3(0, "height",      OPT_UINT,   &cmdlineP->height,  NULL, 0);
-    OPTENT3(0, "left",        OPT_UINT,   &cmdlineP->left,    NULL, 0);
-    OPTENT3(0, "right",       OPT_UINT,   &cmdlineP->right,   NULL, 0);
-    OPTENT3(0, "top",         OPT_UINT,   &cmdlineP->top,     NULL, 0);
-    OPTENT3(0, "bottom",      OPT_UINT,   &cmdlineP->bottom,  NULL, 0);
-    OPTENT3(0, "bg",          OPT_STRING, &cmdlineP->bg_rgb,  NULL, 0);
-    OPTENT3(0, "fg",          OPT_STRING, &cmdlineP->fg_rgb,  NULL, 0);
-    OPTENT3(0, "var",         OPT_UINT,   &cmdlineP->var,     NULL, 0);
+    OPTENT3(0, "width",       OPT_UINT,   &cmdlineP->width,
+            &widthSpec, 0);
+    OPTENT3(0, "height",      OPT_UINT,   &cmdlineP->height,
+            &heightSpec, 0);
+    OPTENT3(0, "left",        OPT_UINT,   &cmdlineP->left,
+            &cmdlineP->leftSpec, 0);
+    OPTENT3(0, "right",       OPT_UINT,   &cmdlineP->right,
+            &cmdlineP->rightSpec, 0);
+    OPTENT3(0, "top",         OPT_UINT,   &cmdlineP->top,
+            &cmdlineP->topSpec, 0);
+    OPTENT3(0, "bottom",      OPT_UINT,   &cmdlineP->bottom,
+            &cmdlineP->bottomSpec, 0);
+    OPTENT3(0, "bg",          OPT_STRING, &cmdlineP->bg,
+            &bgSpec, 0);
+    OPTENT3(0, "fg",          OPT_STRING, &cmdlineP->fg,
+            &fgSpec, 0);
+    OPTENT3(0, "var",         OPT_UINT,   &cmdlineP->var,
+            &varSpec, 0);
     OPTENT3(0, "randomseed",  OPT_UINT,   &cmdlineP->randomseed,
             &cmdlineP->randomseedSpec, 0);
     OPTENT3(0, "init",        OPT_UINT,   &cmdlineP->randomseed,
             &cmdlineP->randomseedSpec, 0);
-    OPTENT3(0, "verbose",     OPT_FLAG,   NULL, &cmdlineP->verbose, 0);
-
-    /* Set the defaults */
-    cmdlineP->width = 100;
-    cmdlineP->height = 100;
-    cmdlineP->left = cmdlineP->right = cmdlineP->top = cmdlineP->bottom = -1;
-    cmdlineP->bg_rgb = NULL;
-    cmdlineP->fg_rgb = NULL;
-    cmdlineP->var = 10;
+    OPTENT3(0, "verbose",     OPT_FLAG,   NULL,
+            &cmdlineP->verbose, 0);
 
     opt.opt_table = option_def;
     opt.short_allowed = FALSE;  /* We have no short (old-fashioned) options */
@@ -80,6 +86,26 @@ parseCommandLine(int argc, const char ** argv,
 
     pm_optParseOptions3(&argc, (char **)argv, opt, sizeof(opt), 0);
 
+    if (!widthSpec)
+        cmdlineP->width = 100;
+    if (!heightSpec)
+        cmdlineP->height = 100;
+    if (!bgSpec)
+        cmdlineP->bg = NULL;
+    if (!fgSpec)
+        cmdlineP->fg = NULL;
+    if (!varSpec)
+        cmdlineP->var = 10;
+
+    if (cmdlineP->topSpec && cmdlineP->top > cmdlineP->height)
+        pm_error("-top value too large.  Max is %u", cmdlineP->height);
+    if (cmdlineP->bottomSpec && cmdlineP->bottom > cmdlineP->height)
+        pm_error("-bottom value too large.  Max is %u", cmdlineP->height);
+    if (cmdlineP->leftSpec && cmdlineP->left > cmdlineP->width)
+        pm_error("-left value too large.  Max is %u", cmdlineP->width);
+    if (cmdlineP->rightSpec && cmdlineP->right > cmdlineP->width)
+        pm_error("-right value too large.  Max is %u", cmdlineP->width);
+
     if (argc-1 != 0)
         pm_error("There are no arguments.  You specified %d.", argc-1);
 
@@ -88,240 +114,353 @@ parseCommandLine(int argc, const char ** argv,
 
 
 
-static void
-procLeft(int          const r1,
-         int          const r2,
-         int          const c1,
-         int          const c2,
-         unsigned int const var) {
+static int
+mean(int const a,
+     int const b) {
 
-    int cm, rm, c;
+    return (a + b) / 2;
+}
 
-    if (r1 + 1 == r2)  return;
-    rm = (r1 + r2) >> 1;
-    cm = (c1 + c2) >> 1;
-    cm += (int)floor(((float)rand() / RAND_MAX - 0.5) * var + 0.5);
 
-    for (c = 0; c < cm; c++)
-        PPM_ASSIGN(PIX[rm][c], BG_RED, BG_GREEN, BG_BLUE);
 
-    procLeft(r1, rm, c1, cm, var);
-    procLeft(rm, r2, cm, c2, var);
+static void
+reportParameters(struct CmdlineInfo const cmdline,
+                 pixel              const bgcolor,
+                 pixel              const fgcolor) {
+
+    pm_message("width is %d, height is %d, variance is %d.",
+               cmdline.width, cmdline.height, cmdline.var);
+    if (cmdline.leftSpec)
+        pm_message("ragged left border is required");
+    if (cmdline.rightSpec)
+        pm_message("ragged right border is required");
+    if (cmdline.topSpec)
+        pm_message("ragged top border is required");
+    if (cmdline.bottomSpec)
+        pm_message("ragged bottom border is required");
+    pm_message("background is %s",
+               ppm_colorname(&bgcolor, PPM_MAXMAXVAL, 1));
+    pm_message("foreground is %s",
+               ppm_colorname(&fgcolor, PPM_MAXMAXVAL, 1));
+    if (cmdline.randomseedSpec)
+        pm_message("pm_rand() initialized with seed %u",
+                   cmdline.randomseed);
 }
 
 
 
 static void
-procRight(int          const r1,
-          int          const r2,
-          int          const c1,
-          int          const c2,
-          unsigned int const width,
-          unsigned int const var) {
+makeAllForegroundColor(pixel **     const pixels,
+                       unsigned int const rows,
+                       unsigned int const cols,
+                       pixel        const fgcolor) {
 
-    int cm, rm, c;
+    pixval const r = PPM_GETR(fgcolor);
+    pixval const g = PPM_GETG(fgcolor);
+    pixval const b = PPM_GETB(fgcolor);
 
-    if (r1 + 1 == r2)  return;
-    rm = (r1 + r2) >> 1;
-    cm = (c1 + c2) >> 1;
-    cm += (int)floor(((float)rand() / RAND_MAX - 0.5) * var + 0.5);
+    unsigned int row;
 
-    for (c = cm; c < width; c++)
-        PPM_ASSIGN(PIX[rm][c], BG_RED, BG_GREEN, BG_BLUE);
+    for (row = 0; row < rows; ++row) {
+        unsigned int col;
 
-    procRight(r1, rm, c1, cm, width, var);
-    procRight(rm, r2, cm, c2, width, var);
+        for (col = 0; col < cols; ++col)
+            PPM_ASSIGN(pixels[row][col], r, g, b);
+    }
 }
 
 
 
 static void
-procTop(int          const c1,
-        int          const c2,
-        int          const r1,
-        int          const r2,
-        unsigned int const var) {
+procLeft(pixel **           const pixels,
+         int                const r1,
+         int                const r2,
+         int                const c1,
+         int                const c2,
+         unsigned int       const width,
+         unsigned int       const var,
+         pixel              const bgcolor,
+         struct pm_randSt * const randStP) {
+
+    if (r1 + 1 != r2) {
+        int const rm = mean(r1, r2);
+        int const cm = mean(c1, c2) +
+            (int)floor(((float)pm_drand(randStP) - 0.5) * var + 0.5);
+
+        unsigned int c;
+
+        for (c = 0; c < MIN(width, MAX(0, cm)); ++c)
+            pixels[rm][c] = bgcolor;
+
+        procLeft(pixels, r1, rm, c1, cm, width, var, bgcolor, randStP);
+        procLeft(pixels, rm, r2, cm, c2, width, var, bgcolor, randStP);
+    }
+}
 
-    int rm, cm, r;
 
-    if (c1 + 1 == c2)  return;
-    cm = (c1 + c2) >> 1;
-    rm = (r1 + r2) >> 1;
-    rm += (int)floor(((float)rand() / RAND_MAX - 0.5) * var + 0.5);
 
-    for (r = 0; r < rm; r++)
-        PPM_ASSIGN(PIX[r][cm], BG_RED, BG_GREEN, BG_BLUE);
+static void
+procRight(pixel **           const pixels,
+          int                const r1,
+          int                const r2,
+          int                const c1,
+          int                const c2,
+          unsigned int       const width,
+          unsigned int       const var,
+          pixel              const bgcolor,
+          struct pm_randSt * const randStP) {
+
+    if (r1 + 1 != r2) {
+        int const rm = mean(r1, r2);
+        int const cm = mean(c1, c2) +
+            (int)floor(((float)pm_drand(randStP) - 0.5) * var + 0.5);
+
+        unsigned int c;
+
+        for (c = MAX(0, cm); c < width; ++c)
+            pixels[rm][c] = bgcolor;
+
+        procRight(pixels, r1, rm, c1, cm, width, var, bgcolor, randStP);
+        procRight(pixels, rm, r2, cm, c2, width, var, bgcolor, randStP);
+    }
+}
 
-    procTop(c1, cm, r1, rm, var);
-    procTop(cm, c2, rm, r2, var);
+
+
+static void
+procTop(pixel **           const pixels,
+        int                const c1,
+        int                const c2,
+        int                const r1,
+        int                const r2,
+        unsigned int       const height,
+        unsigned int       const var,
+        pixel              const bgcolor,
+        struct pm_randSt * const randStP) {
+
+    if (c1 + 1 != c2) {
+        int const cm = mean(c1, c2);
+        int const rm = mean(r1, r2) +
+            (int)floor(((float)pm_drand(randStP) - 0.5) * var + 0.5);
+
+        unsigned int r;
+
+        for (r = 0; r < MIN(height, MAX(0, rm)); ++r)
+            pixels[r][cm] = bgcolor;
+
+        procTop(pixels, c1, cm, r1, rm, height, var, bgcolor, randStP);
+        procTop(pixels, cm, c2, rm, r2, height, var, bgcolor, randStP);
+    }
 }
 
 
 
 static void
-procBottom(int          const c1,
-           int          const c2,
-           int          const r1,
-           int          const r2,
-           unsigned int const height,
-           unsigned int const var) {
+procBottom(pixel **           const pixels,
+           int                const c1,
+           int                const c2,
+           int                const r1,
+           int                const r2,
+           unsigned int       const height,
+           unsigned int       const var,
+           pixel              const bgcolor,
+           struct pm_randSt * const randStP) {
+
+    if (c1 + 1 != c2) {
+        int const cm = mean(c1, c2);
+        int const rm = mean(r1, r2) +
+            (int)floor(((float)pm_drand(randStP) - 0.5) * var + 0.5);
+
+        unsigned int r;
+
+        for (r = MAX(0, rm); r < height; ++r)
+            pixels[r][cm] = bgcolor;
+
+        procBottom(pixels, c1, cm, r1, rm, height, var, bgcolor, randStP);
+        procBottom(pixels, cm, c2, rm, r2, height, var, bgcolor, randStP);
+    }
+}
 
-    int rm, cm, r;
 
-    if (c1 + 1 == c2)  return;
-    cm = (c1 + c2) >> 1;
-    rm = (r1 + r2) >> 1;
-    rm += (int)floor(((float)rand() / RAND_MAX - 0.5) * var + 0.5);
 
-    for (r = rm; r < height; r++)
-        PPM_ASSIGN(PIX[r][cm], BG_RED, BG_GREEN, BG_BLUE);
+static void
+makeRaggedLeftBorder(pixel **           const pixels,
+                     unsigned int       const rows,
+                     unsigned int       const cols,
+                     bool               const leftSpec,
+                     unsigned int       const left,
+                     unsigned int       const var,
+                     pixel              const bgcolor,
+                     struct pm_randSt * const randStP) {
+
+    if (leftSpec) {
+        int const leftC1 = left;
+        int const leftC2 = left;
+        int const leftR1 = 0;
+        int const leftR2 = rows - 1;
+
+        unsigned int col;
 
-    procBottom(c1, cm, r1, rm, height, var);
-    procBottom(cm, c2, rm, r2, height, var);
+        for (col = 0; col < leftC1; ++col)
+            pixels[leftR1][col] = bgcolor;
+        for (col = 0; col < leftC2; ++col)
+            pixels[leftR2][col] = bgcolor;
+
+        procLeft(pixels, leftR1, leftR2, leftC1, leftC2, cols, var,
+                 bgcolor, randStP);
+    }
 }
 
 
 
-int
-main(int argc, const char * argv[]) {
+static void
+makeRaggedRightBorder(pixel **           const pixels,
+                      unsigned int       const rows,
+                      unsigned int       const cols,
+                      bool               const rightSpec,
+                      unsigned int       const right,
+                      unsigned int       const width,
+                      unsigned int       const var,
+                      pixel              const bgcolor,
+                      struct pm_randSt * const randStP) {
+
+    if (rightSpec) {
+        int const rightC1 = cols - right - 1;
+        int const rightC2 = cols - right - 1;
+        int const rightR1 = 0;
+        int const rightR2 = rows - 1;
 
-    struct CmdlineInfo cmdline;
-    pixel bgcolor, fgcolor;
-    pixval fg_red, fg_green, fg_blue;
-    int rows, cols, row;
-    int left, right, top, bottom;
+        unsigned int col;
 
-    pm_proginit(&argc, argv);
+        for (col = rightC1; col < cols; ++col)
+            pixels[rightR1][col] = bgcolor;
+        for (col = rightC2; col < cols; ++col)
+            pixels[rightR2][col] = bgcolor;
 
-    parseCommandLine(argc, argv, &cmdline);
+        procRight(pixels, rightR1, rightR2, rightC1, rightC2, width, var,
+                  bgcolor, randStP);
+    }
+}
 
-    srand(cmdline.randomseedSpec ? cmdline.randomseed : pm_randseed());
 
-    cols = cmdline.width;
-    rows = cmdline.height;
-    left = cmdline.left;
-    right = cmdline.right;
-    top = cmdline.top;
-    bottom = cmdline.bottom;
 
-    if (cmdline.bg_rgb)
-        bgcolor = ppm_parsecolor(cmdline.bg_rgb, PPM_MAXMAXVAL);
-    else
-        PPM_ASSIGN(bgcolor, 0, 0, 0);
-    BG_RED = PPM_GETR(bgcolor);
-    BG_GREEN = PPM_GETG(bgcolor);
-    BG_BLUE = PPM_GETB(bgcolor);
+static void
+makeRaggedTopBorder(pixel **           const pixels,
+                    unsigned int       const rows,
+                    unsigned int       const cols,
+                    bool               const topSpec,
+                    unsigned int       const top,
+                    unsigned int       const var,
+                    pixel              const bgcolor,
+                    struct pm_randSt * const randStP) {
+
+    if (topSpec) {
+        unsigned int const topR1 = top;
+        unsigned int const topR2 = top;
+        unsigned int const topC1 = 0;
+        unsigned int const topC2 = cols - 1;
 
-    if (cmdline.fg_rgb)
-        fgcolor = ppm_parsecolor(cmdline.fg_rgb, PPM_MAXMAXVAL);
-    else
-        PPM_ASSIGN(fgcolor, PPM_MAXMAXVAL, PPM_MAXMAXVAL, PPM_MAXMAXVAL);
-    fg_red = PPM_GETR(fgcolor);
-    fg_green = PPM_GETG(fgcolor);
-    fg_blue = PPM_GETB(fgcolor);
-
-    if (cmdline.verbose) {
-        pm_message("width is %d, height is %d, variance is %d.",
-                   cols, rows, cmdline.var);
-        if (left >= 0)
-            pm_message("ragged left border is required");
-        if (right >= 0)
-            pm_message("ragged right border is required");
-        if (top >= 0)
-            pm_message("ragged top border is required");
-        if (bottom >= 0)
-            pm_message("ragged bottom border is required");
-        pm_message("background is %s",
-                   ppm_colorname(&bgcolor, PPM_MAXMAXVAL, 1));
-        pm_message("foreground is %s",
-                   ppm_colorname(&fgcolor, PPM_MAXMAXVAL, 1));
-        if (cmdline.randomseedSpec)
-            pm_message("srand() initialized with seed %u", cmdline.randomseed);
-    }
+        unsigned int row;
 
-    /* Allocate memory for the whole pixmap */
-    PIX = ppm_allocarray(cols, rows);
+        for (row = 0; row < topR1; ++row)
+            pixels[row][topC1] = bgcolor;
+        for (row = 0; row < topR2; ++row)
+            pixels[row][topC2] = bgcolor;
 
-    /* First, set all pixel to foreground color */
-    for (row = 0; row < rows; row++) {
-        unsigned int col;
-        for (col = 0; col < cols; ++col)
-            PPM_ASSIGN(PIX[row][col], fg_red, fg_green, fg_blue);
+        procTop(pixels, topC1, topC2, topR1, topR2, rows,
+                var, bgcolor, randStP);
     }
-    /* Make a ragged left border */
-    if (left >= 0) {
-        int const left_c1 = left;
-        int const left_c2 = left;
-        int const left_r1 = 0;
-        int const left_r2 = rows - 1;
+}
 
-        unsigned int col;
 
-        for (col = 0; col < left_c1; ++col)
-            PPM_ASSIGN(PIX[left_r1][col], BG_RED, BG_GREEN, BG_BLUE);
-        for (col = 0; col < left_c2; ++col)
-            PPM_ASSIGN(PIX[left_r2][col], BG_RED, BG_GREEN, BG_BLUE);
 
-        procLeft(left_r1, left_r2, left_c1, left_c2, cmdline.var);
+static void
+makeRaggedBottomBorder(pixel **            const pixels,
+                       unsigned int        const rows,
+                       unsigned int        const cols,
+                       bool                const bottomSpec,
+                       unsigned int        const bottom,
+                       unsigned int        const height,
+                       unsigned int        const var,
+                       pixel               const bgcolor,
+                       struct pm_randSt *  const randStP) {
+
+    if (bottomSpec) {
+        unsigned int const bottomR1 = rows - bottom - 1;
+        unsigned int const bottomR2 = rows - bottom - 1;
+        unsigned int const bottomC1 = 0;
+        unsigned int const bottomC2 = cols - 1;
+
+        unsigned int row;
+
+        for (row = bottomR1; row < rows; ++row)
+            pixels[row][bottomC1] = bgcolor;
+        for (row = bottomR2; row < rows; ++row)
+            pixels[row][bottomC2] = bgcolor;
+
+        procBottom(pixels, bottomC1, bottomC2, bottomR1, bottomR2,
+                   height, var, bgcolor, randStP);
     }
+}
 
-    /* Make a ragged right border */
-    if (right >= 0) {
-        int const right_c1 = cols - right - 1;
-        int const right_c2 = cols - right - 1;
-        int const right_r1 = 0;
-        int const right_r2 = rows - 1;
 
-        unsigned int col;
 
-        for (col = right_c1; col < cols; col++)
-            PPM_ASSIGN(PIX[right_r1][col], BG_RED, BG_GREEN, BG_BLUE);
-        for (col = right_c2; col < cols; col++)
-            PPM_ASSIGN(PIX[right_r2][col], BG_RED, BG_GREEN, BG_BLUE);
+int
+main(int argc, const char ** const argv) {
 
-        procRight(right_r1, right_r2, right_c1, right_c2,
-                   cmdline.width, cmdline.var);
-    }
+    struct CmdlineInfo cmdline;
+    pixel bgcolor, fgcolor;
+    struct pm_randSt randSt;
+    static pixel** pixels;
 
-    /* Make a ragged top border */
-    if (top >= 0) {
-        int const top_r1 = top;
-        int const top_r2 = top;
-        int const top_c1 = 0;
-        int const top_c2 = cols - 1;
+    pm_proginit(&argc, argv);
 
-        unsigned int row;
+    parseCommandLine(argc, argv, &cmdline);
 
-        for (row = 0; row < top_r1; ++row)
-            PPM_ASSIGN(PIX[row][top_c1], BG_RED, BG_GREEN, BG_BLUE);
-        for (row = 0; row < top_r2; ++row)
-            PPM_ASSIGN(PIX[row][top_c2], BG_RED, BG_GREEN, BG_BLUE);
+    pm_randinit(&randSt);
+    pm_srand(&randSt,
+             cmdline.randomseedSpec ? cmdline.randomseed : pm_randseed());
 
-        procTop(top_c1, top_c2, top_r1, top_r2, cmdline.var);
-    }
+    if (cmdline.bg)
+        bgcolor = ppm_parsecolor(cmdline.bg, PPM_MAXMAXVAL);
+    else
+        PPM_ASSIGN(bgcolor, 0, 0, 0);
+
+    if (cmdline.fg)
+        fgcolor = ppm_parsecolor(cmdline.fg, PPM_MAXMAXVAL);
+    else
+        PPM_ASSIGN(fgcolor, PPM_MAXMAXVAL, PPM_MAXMAXVAL, PPM_MAXMAXVAL);
 
-    /* Make a ragged bottom border */
-    if (bottom >= 0) {
-        int const bottom_r1 = rows - bottom - 1;
-        int const bottom_r2 = rows - bottom - 1;
-        int const bottom_c1 = 0;
-        int const bottom_c2 = cols - 1;
+    if (cmdline.verbose)
+        reportParameters(cmdline, bgcolor, fgcolor);
 
-        unsigned int row;
+    pixels = ppm_allocarray(cmdline.width, cmdline.height);
 
-        for (row = bottom_r1; row < rows; ++row)
-            PPM_ASSIGN(PIX[row][bottom_c1], BG_RED, BG_GREEN, BG_BLUE);
-        for (row = bottom_r2; row < rows; ++row)
-            PPM_ASSIGN(PIX[row][bottom_c2], BG_RED, BG_GREEN, BG_BLUE);
+    makeAllForegroundColor(pixels, cmdline.height, cmdline.width, fgcolor);
 
-        procBottom(bottom_c1, bottom_c2, bottom_r1, bottom_r2,
-                   cmdline.height, cmdline.var);
-    }
+    makeRaggedLeftBorder(pixels, cmdline.height, cmdline.width,
+                         cmdline.leftSpec, cmdline.left,
+                         cmdline.var, bgcolor, &randSt);
+
+    makeRaggedRightBorder(pixels, cmdline.height, cmdline.width,
+                          cmdline.rightSpec, cmdline.right,
+                          cmdline.width, cmdline.var, bgcolor, &randSt);
+
+    makeRaggedTopBorder(pixels, cmdline.height, cmdline.width,
+                        cmdline.topSpec, cmdline.top,
+                        cmdline.var, bgcolor, &randSt);
+
+    makeRaggedBottomBorder(pixels, cmdline.height, cmdline.width,
+                           cmdline.bottomSpec, cmdline.bottom,
+                           cmdline.height, cmdline.var, bgcolor, &randSt);
+
+    pm_randterm(&randSt);
 
     /* Write pixmap */
-    ppm_writeppm(stdout, PIX, cols, rows, PPM_MAXMAXVAL, 0);
+    ppm_writeppm(stdout, pixels, cmdline.width, cmdline.height,
+                 PPM_MAXMAXVAL, 0);
 
-    ppm_freearray(PIX, rows);
+    ppm_freearray(pixels, cmdline.height);
 
     pm_close(stdout);
 
@@ -329,4 +468,3 @@ main(int argc, const char * argv[]) {
 }
 
 
-