Create Subversion repository

git-svn-id: http://svn.code.sf.net/p/netpbm/code/trunk@1 9d0c8265-081b-0410-96cb-a4ca84ce46f8

1 files changed, 437 insertions, 0 deletions

diff --git a/analyzer/pamtilt.c b/analyzer/pamtilt.c
new file mode 100644
index 00000000..3753955b
--- /dev/null
+++ b/analyzer/pamtilt.c
@@ -0,0 +1,437 @@
+/*=============================================================================
+                             pgmtilt
+===============================================================================
+  Print the tilt angle of a PGM file
+
+  Based on pgmskew by Gregg Townsend, August 2005.
+
+  Adapted to Netpbm by Bryan Henderson, August 2005.
+
+  All work has been contributed to the public domain by its authors.
+=============================================================================*/
+
+#include <math.h>
+
+#include "pam.h"
+#include "mallocvar.h"
+#include "shhopt.h"
+
+/* program constants */
+#define DIFFMAX 255                     /* maximum scaling for differences */
+#define BARLENGTH 60                    /* length of bars printed with -v */
+
+struct cmdlineInfo {
+    /* command parameters */
+    const char * inputFilename;
+    float maxangle;     /* maximum angle attempted */
+    float astep;        /* initial angle increment */
+    float qmin;         /* minimum quality of solution */
+    unsigned int hstep; /* horizontal step size */
+    unsigned int vstep; /* vertical step size */
+    unsigned int dstep; /* difference distance */
+    unsigned int fast;    /* skip third iteration */
+    unsigned int verbose; /* generate commentary */
+};
+
+static void
+abandon(void) {
+
+    printf("00.00\n");
+    exit(0);
+}
+
+
+
+static void
+parseCommandLine(int argc, char *argv[],
+                 struct cmdlineInfo * const cmdlineP) {
+
+    static optEntry option_def[50];
+    static optStruct3 opt;
+    unsigned int option_def_index;
+
+    /* set defaults */
+    cmdlineP->hstep = 11;        /* read only every 11th column */
+    cmdlineP->vstep = 5;         /* calc differences every 5th row */
+    cmdlineP->dstep = 2;         /* check for differences two rows down */
+    cmdlineP->maxangle = 10.0;   /* assume skew is less than +/- ten degrees */
+    cmdlineP->astep = 1.0;       /* initially check by one-degree increments */
+    cmdlineP->qmin = 1.0;        /* don't require S/N better than 1.0 */
+
+    /* initialize option table */
+    option_def_index = 0;       /* incremented by OPTENT3 */
+    OPTENT3(0, "fast",     OPT_FLAG,  NULL, &cmdlineP->fast,     0);
+    OPTENT3(0, "verbose",  OPT_FLAG,  NULL, &cmdlineP->verbose,  0);
+    OPTENT3(0, "angle",    OPT_FLOAT, &cmdlineP->maxangle, NULL, 0);
+    OPTENT3(0, "quality",  OPT_FLOAT, &cmdlineP->qmin,     NULL, 0);
+    OPTENT3(0, "astep",    OPT_FLOAT, &cmdlineP->astep,    NULL, 0);
+    OPTENT3(0, "hstep",    OPT_UINT,  &cmdlineP->hstep,    NULL, 0);
+    OPTENT3(0, "vstep",    OPT_UINT,  &cmdlineP->vstep,    NULL, 0);
+    OPTENT3(0, "dstep",    OPT_UINT,  &cmdlineP->dstep,    NULL, 0);
+
+    opt.opt_table = option_def;
+    opt.short_allowed = FALSE;          /* no short options used */
+    opt.allowNegNum = FALSE;            /* don't allow negative values */
+    optParseOptions3(&argc, argv, opt, sizeof(opt), 0);
+
+    if (cmdlineP->hstep < 1)
+        pm_error("-hstep must be at least 1 column.");
+    if (cmdlineP->vstep < 1)
+        pm_error("-vstep must be at least 1 row.");
+    if (cmdlineP->dstep < 1)
+        pm_error("-dstep must be at least 1 row.");
+    if (cmdlineP->maxangle < 1 || cmdlineP->maxangle > 45)
+        pm_error("-maxangle must be between 1 and 45 degrees.");
+
+    if (argc-1 < 1)                      /* if input file name given */
+        cmdlineP->inputFilename = "-";
+    else {
+        cmdlineP->inputFilename = argv[1];
+
+        if (argc-1 > 1)
+            pm_error("There is at most one argument.  You specified %d",
+                     argc-1);
+    }
+}
+
+
+
+static void
+computeSteps(const struct pam * const pamP,
+             unsigned int       const hstepReq,
+             unsigned int       const vstepReq,
+             unsigned int *     const hstepP,
+             unsigned int *     const vstepP) {
+/*----------------------------------------------------------------------------
+  Adjust parameters if necessary now that we have the image size
+-----------------------------------------------------------------------------*/
+    if (pamP->width < 10 || pamP->height < 10)
+        abandon();
+
+    if (pamP->width < 10 * hstepReq)
+        *hstepP = pamP->width / 10;
+    else
+        *hstepP = hstepReq;
+
+    if (pamP->height < 10 * vstepReq)
+        *vstepP = pamP->height / 10;
+    else
+        *vstepP = vstepReq;
+}
+
+
+
+static void
+load(const struct pam * const pamP,
+     unsigned int       const hstep,
+     sample ***         const pixelsP,
+     unsigned int *     const hsamplesP) {
+/*----------------------------------------------------------------------------
+  read file into memory, returning array of rows
+-----------------------------------------------------------------------------*/
+    unsigned int const hsamples = 1 + (pamP->width - 1) / hstep;
+        /* use only this many cols */
+
+    unsigned int row;
+    tuple *  tuplerow;
+    sample ** pixels;
+
+    tuplerow = pnm_allocpamrow(pamP);
+    
+    MALLOCARRAY(pixels, pamP->height);
+    if (pixels == NULL)
+        pm_error("Unable to allocate array of %u pixel rows",
+                 pamP->height);
+
+    if (pixels == NULL)
+        pm_error("Unable to allocate array of %u rows", pamP->height);
+
+    for (row = 0; row < pamP->height; ++row) {
+        unsigned int i;
+        unsigned int col;
+
+        pnm_readpamrow(pamP, tuplerow);
+
+        MALLOCARRAY(pixels[row], hsamples);
+        if (pixels[row] == NULL)
+            pm_error("Unable to allocate %u-sample array for Row %u",
+                     hsamples, row);
+
+        /* save every hstep'th column */
+        for (i = col = 0; i < hsamples; ++i, col += hstep)
+            pixels[row][i] = tuplerow[col][0];
+    }
+
+    pnm_freepamrow(tuplerow);
+
+    *hsamplesP = hsamples;
+    *pixelsP   = pixels;
+}
+
+
+
+static void
+freePixels(sample **    const samples,
+           unsigned int const rows) {
+
+    unsigned int row;
+
+    for (row = 0; row < rows; ++row)
+        free(samples[row]);
+
+    free(samples);
+}
+
+
+
+static void
+replacePixelValuesWithScaledDiffs(
+    const struct pam * const pamP,
+    sample **          const pixels,
+    unsigned int       const hsamples,
+    unsigned int       const dstep) {
+/*----------------------------------------------------------------------------
+  Replace pixel values with scaled differences used for all
+  calculations
+-----------------------------------------------------------------------------*/
+    float const m = DIFFMAX / (float) pamP->maxval;  /* scale multiplier */
+
+    unsigned int row;
+
+    for (row = dstep; row < pamP->height; ++row) {
+        unsigned int col;
+        for (col = 0; col < hsamples; ++col) {
+            int const d = pixels[row - dstep][col] - pixels[row][col];
+            unsigned int const absd = d < 0 ? -d : d;
+            pixels[row - dstep][col] = m * absd;  /* scale the difference */
+        }
+    }
+}
+
+
+
+static void
+scoreAngle(const struct pam * const pamP,
+           sample **          const pixels,
+           unsigned int       const hstep,
+           unsigned int       const vstep,
+           unsigned int       const hsamples,
+           float              const deg,
+           float *            const scoreP) {
+/*----------------------------------------------------------------------------
+  calculate score for a given angle
+-----------------------------------------------------------------------------*/
+    float  const radians = (float)deg/360 * 2 * M_PI;
+    float  const dy      = hstep * tan(radians);
+    int    const dtotal  = pamP->width * tan(radians);
+    double const tscale  = 1.0 / hsamples;
+
+    double total;
+    int first;
+    int last;
+    
+    unsigned int row;
+
+    total = 0.0; /* initial value */
+
+    if (dtotal > 0) {
+        first = 0;
+        last = pamP->height - 1 - (dtotal + 1);
+    } else {
+        first = -(dtotal - 1);
+        last = pamP->height - 1;
+    }
+    for (row = first; row < last; row += vstep) {
+        float o;
+        long t;
+        double dt;
+
+        unsigned int i;
+
+        for (i = 0, t = 0, o = 0.5;
+             i < hsamples;
+             ++i, t += pixels[(int)(row + o)][i], o += dy) {
+        }
+        dt = tscale * t;
+        total += dt * dt;
+    }
+    *scoreP = total / (last - first);
+}
+
+
+
+static void
+getBestAngleLocal(
+    const struct pam * const pamP,
+    sample **          const pixels,
+    unsigned int       const hstep,
+    unsigned int       const vstep,
+    unsigned int       const hsamples,
+    float              const minangle,
+    float              const maxangle,
+    float              const incr,
+    bool               const verbose,
+    float *            const bestAngleP,
+    float *            const qualityP) {
+/*----------------------------------------------------------------------------
+  find angle of highest score within a range
+-----------------------------------------------------------------------------*/
+    int const nsamples = ((maxangle - minangle) / incr + 1.5);
+
+    float score;
+    float quality;  /* signal/noise ratio of the best angle */
+    float angle;
+    float bestangle;
+    float bestscore;
+    float total;
+    float others;
+    float * results;
+    int i;
+
+    MALLOCARRAY_NOFAIL(results, nsamples);      /* allocate array of results */
+
+    /* try all angles within the given range, stepping by incr */
+    bestangle = minangle;
+    bestscore = 0;
+    total = 0;
+    for (i = 0; i < nsamples; i++) {
+        angle = minangle + i * incr;
+        scoreAngle(pamP, pixels, hstep, vstep, hsamples, angle, &score);
+        results[i] = score;
+        if (score > bestscore ||
+            (score == bestscore && fabs(angle) < fabs(bestangle))) {
+            bestscore = score;
+            bestangle = angle;
+        }
+        total += score;
+    }
+
+    others = (total-bestscore) / (nsamples-1);  /* get mean of other scores */
+    quality = bestscore / others;
+
+    if (verbose) {
+        fprintf(stderr,
+                "\n%2d angles from %6.2f to %5.2f by %4.2f:  "
+                "best = %6.2f, S:N = %4.2f\n",
+                nsamples, minangle, maxangle, incr, bestangle, quality);
+        for (i = 0; i < nsamples; ++i) {
+            float const angle = minangle + i * incr;
+            int const n = (int) (BARLENGTH * results[i] / bestscore + 0.5);
+            fprintf(stderr, "%6.2f: %8.2f %0*d\n", angle, results[i], n, 0);
+        }
+    }
+    free(results);
+
+    *qualityP   = quality;
+    *bestAngleP = bestangle;
+}
+
+
+
+static void
+readRelevantPixels(const char *   const inputFilename,
+                   unsigned int   const hstepReq,
+                   unsigned int   const vstepReq,
+                   unsigned int * const hstepP,
+                   unsigned int * const vstepP,
+                   sample ***     const pixelsP,
+                   struct pam *   const pamP,
+                   unsigned int * const hsamplesP) {
+/*----------------------------------------------------------------------------
+  load the image, saving only the pixels we might actually inspect
+-----------------------------------------------------------------------------*/
+    FILE * ifP;
+    unsigned int hstep;
+    unsigned int vstep;
+
+    ifP = pm_openr(inputFilename);
+    pnm_readpaminit(ifP, pamP, PAM_STRUCT_SIZE(tuple_type));
+    computeSteps(pamP, hstepReq, vstepReq, &hstep, &vstep);
+
+    load(pamP, hstep, pixelsP, hsamplesP);
+
+    *hstepP = hstep;
+    *vstepP = vstep;
+    
+    pm_close(ifP);
+}
+
+
+
+static void
+getAngle(const struct pam * const pamP,
+         sample **          const pixels,
+         unsigned int       const hstep,
+         unsigned int       const vstep,
+         unsigned int       const hsamples,
+         float              const maxangle,
+         float              const astep,
+         float              const qmin,
+         bool               const fast,
+         bool               const verbose,
+         float *            const angleP) {
+
+    float a;
+    float da;
+    float lastq;        /* quality (s/n ratio) of last measurement */
+    
+    getBestAngleLocal(pamP, pixels, hstep, vstep, hsamples,
+                      -maxangle, maxangle, astep, verbose,
+                      &a, &lastq);
+
+    if ((a < -maxangle + astep / 2) || (a > maxangle - astep / 2))
+        /* extreme val almost certainly wrong */
+        abandon();
+    if (lastq < qmin)
+        /* insufficient s/n ratio */
+        abandon();
+
+    /* make a finer search in the neighborhood */
+    da = astep / 10;
+    getBestAngleLocal(pamP, pixels, hstep, vstep, hsamples,
+                      a - 9 * da, a + 9 * da, da, verbose,
+                      &a, &lastq);
+
+    /* iterate once more unless we don't need that much accuracy */
+    if (!fast) {
+        da /= 10;
+        getBestAngleLocal(pamP, pixels, hstep, vstep, hsamples,
+                          a - 9 * da, a + 9 * da, da, verbose,
+                          &a, &lastq);
+    }
+    *angleP = a;
+}
+
+
+
+int
+main(int argc, char *argv[]) {
+
+    struct cmdlineInfo cmdline;
+    struct pam pam;
+    sample ** pixels;       /* pixel data */
+    unsigned int hsamples; /* horizontal samples used */
+    unsigned int hstep;    /* horizontal step size */
+    unsigned int vstep;    /* vertical step size */
+    float angle;
+
+    pgm_init(&argc, argv);              /* initialize netpbm system */
+
+    parseCommandLine(argc, argv, &cmdline);
+
+    readRelevantPixels(cmdline.inputFilename, cmdline.hstep, cmdline.vstep,
+                       &hstep, &vstep, &pixels, &pam, &hsamples);
+
+    replacePixelValuesWithScaledDiffs(&pam, pixels, hsamples, cmdline.dstep);
+
+    getAngle(&pam, pixels, hstep, vstep, hsamples,
+             cmdline.maxangle, cmdline.astep, cmdline.qmin,
+             cmdline.fast, cmdline.verbose, &angle);
+
+    /* report the result on stdout */
+    printf("%.2f\n", angle);
+
+    freePixels(pixels, pam.height);
+
+    return 0;
+}