home *** CD-ROM | disk | FTP | other *** search
/ PC Pro 2002 April / pcpro0402.iso / essentials / graphics / Gimp / gimp-src-20001226.exe / src / gimp / plug-ins / common / nlfilt.c < prev    next >
Encoding:
C/C++ Source or Header  |  2000-10-22  |  40.8 KB  |  1,286 lines

  1. /**************************************************
  2.  * file: nlfilt/nlfilt.c
  3.  *
  4.  * Copyright (c) 1997 Eric L. Hernes (erich@rrnet.com)
  5.  * All rights reserved.
  6.  *
  7.  * Redistribution and use in source and binary forms, with or without
  8.  * modification, are permitted provided that the following conditions
  9.  * are met:
  10.  * 1. Redistributions of source code must retain the above copyright
  11.  *    notice, this list of conditions and the following disclaimer.
  12.  * 2. The name of the author may not be used to endorse or promote products
  13.  *    derived from this software withough specific prior written permission
  14.  *
  15.  * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR
  16.  * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
  17.  * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
  18.  * IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT,
  19.  * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT
  20.  * NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
  21.  * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
  22.  * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
  23.  * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF
  24.  * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
  25.  *
  26.  * $Id: nlfilt.c,v 1.25 2000/10/21 13:15:55 austin Exp $
  27.  */
  28.  
  29. /*
  30.  * Algorithm fixes, V2.0 compatibility by David Hodson  hodsond@ozemail.com.au
  31.  */
  32.  
  33. /* add any necessary includes  */
  34.  
  35. #ifdef HAVE_CONFIG_H
  36. #include "config.h"
  37. #endif
  38.  
  39. #include <string.h>
  40. #include <stdio.h>
  41. #include <stdlib.h>
  42. #include <math.h>
  43.  
  44. #include <gtk/gtk.h>
  45.  
  46. #include <libgimp/gimp.h>
  47. #include <libgimp/gimpui.h>
  48.  
  49. #include "libgimp/stdplugins-intl.h"
  50.  
  51.  
  52. struct Grgb
  53. {
  54.   guint8 red;
  55.   guint8 green;
  56.   guint8 blue;
  57. };
  58.  
  59. struct GRegion
  60. {
  61.   gint32 x;
  62.   gint32 y;
  63.   gint32 width;
  64.   gint32 height;
  65. };
  66.  
  67. struct piArgs
  68. {
  69.   gint32  img;
  70.   gint32  drw;
  71.   gdouble alpha;
  72.   gdouble radius;
  73.   gint    filter;
  74. };
  75.  
  76. typedef enum
  77. {
  78.   filter_alpha_trim,
  79.   filter_opt_est,
  80.   filter_edge_enhance
  81. } FilterType;
  82.  
  83. /*  preview stuff -- to be removed as soon as we have a real libgimp preview  */
  84.  
  85. struct mwPreview
  86. {
  87.   gint     width;
  88.   gint     height;
  89.   gint     bpp;
  90.   gdouble  scale;
  91.   guchar  *bits;
  92. };
  93.  
  94. #define PREVIEW_SIZE 100
  95.  
  96. static gint   do_preview = TRUE;
  97. static struct mwPreview *thePreview;
  98.  
  99. static GtkWidget        * mw_preview_new   (GtkWidget        *parent,
  100.                                             struct mwPreview *mwp);
  101. static struct mwPreview * mw_preview_build (GimpDrawable     *drw);
  102.  
  103.  
  104. /* function protos */
  105.  
  106. static void query (void);
  107. static void run   (gchar      *name,
  108.            gint        nparam,
  109.            GimpParam  *param,
  110.            gint       *nretvals,
  111.            GimpParam **retvals);
  112.  
  113. static gint pluginCore        (struct piArgs *argp);
  114. static gint pluginCoreIA      (struct piArgs *argp);
  115.  
  116. static void nlfilt_do_preview (GtkWidget  *preview);
  117.  
  118. static inline gint nlfiltInit (gdouble     alpha,
  119.                    gdouble     radius,
  120.                    FilterType  filter);
  121. static inline void nlfiltRow  (guchar     *srclast,
  122.                                guchar     *srcthis,
  123.                                guchar     *srcnext,
  124.                    guchar     *dst,
  125.                    gint        width,
  126.                    gint        Bpp,
  127.                    gint        filtno);
  128.  
  129. GimpPlugInInfo PLUG_IN_INFO =
  130. {
  131.   NULL,  /* init_proc  */
  132.   NULL,  /* quit_proc  */
  133.   query, /* query_proc */
  134.   run,   /* run_proc   */
  135. };
  136.  
  137. MAIN ()
  138.  
  139. static void
  140. query (void)
  141. {
  142.   static GimpParamDef args[] =
  143.   {
  144.     { GIMP_PDB_INT32, "run_mode", "Interactive, non-interactive" },
  145.     { GIMP_PDB_IMAGE, "img", "The Image to Filter" },
  146.     { GIMP_PDB_DRAWABLE, "drw", "The Drawable" },
  147.     { GIMP_PDB_FLOAT, "alpha", "The amount of the filter to apply" },
  148.     { GIMP_PDB_FLOAT, "radius", "The filter radius" },
  149.     { GIMP_PDB_INT32, "filter", "The Filter to Run, 0 - alpha trimmed mean; 1 - optimal estimation (alpha controls noise variance); 2 - edge enhancement" }
  150.   };
  151.   static gint nargs = sizeof (args) / sizeof (args[0]);
  152.  
  153.   gimp_install_procedure ("plug_in_nlfilt",
  154.               "Nonlinear swiss army knife filter",
  155.               "This is the pnmnlfilt, in gimp's clothing.  See the pnmnlfilt manpage for details.",
  156.               "Graeme W. Gill, gimp 0.99 plugin by Eric L. Hernes",
  157.               "Graeme W. Gill, Eric L. Hernes",
  158.               "1997",
  159.               N_("<Image>/Filters/Enhance/NL Filter..."),
  160.               "RGB,GRAY",
  161.               GIMP_PLUGIN,
  162.               nargs, 0,
  163.               args, NULL);
  164. }
  165.  
  166. static void
  167. run (gchar      *name,
  168.      gint        nparam,
  169.      GimpParam  *param,
  170.      gint       *nretvals,
  171.      GimpParam **retvals)
  172. {
  173.   static GimpParam rvals[1];
  174.  
  175.   struct piArgs args;
  176.  
  177.   *nretvals = 1;
  178.   *retvals = rvals;
  179.  
  180.   memset (&args, (int) 0, sizeof (struct piArgs));
  181.  
  182.   args.radius = -1.0;
  183.   gimp_get_data ("plug_in_nlfilt", &args);
  184.   args.img = param[1].data.d_image;
  185.   args.drw = param[2].data.d_drawable;
  186.  
  187.   rvals[0].type = GIMP_PDB_STATUS;
  188.   rvals[0].data.d_status = GIMP_PDB_SUCCESS;
  189.  
  190.   switch (param[0].data.d_int32)
  191.     {
  192.       GimpDrawable *drw;
  193.     case GIMP_RUN_INTERACTIVE:
  194.       INIT_I18N_UI();
  195.       /* XXX: add code here for interactive running */
  196.       if (args.radius == -1)
  197.     {
  198.       args.alpha  = (gdouble) 0.3;
  199.       args.radius = (gdouble) 0.3;
  200.       args.filter = 0;
  201.     }
  202.       drw = gimp_drawable_get (args.drw);
  203.       thePreview = mw_preview_build (drw);
  204.  
  205.       if (pluginCoreIA (&args) == -1)
  206.     {
  207.       rvals[0].data.d_status = GIMP_PDB_EXECUTION_ERROR;
  208.     }
  209.       else
  210.     {
  211.       gimp_set_data ("plug_in_nlfilt", &args, sizeof (struct piArgs));
  212.     }
  213.       break;
  214.  
  215.     case GIMP_RUN_NONINTERACTIVE:
  216.       INIT_I18N();
  217.       /* XXX: add code here for non-interactive running */
  218.       if (nparam != 6)
  219.     {
  220.       rvals[0].data.d_status = GIMP_PDB_CALLING_ERROR;
  221.       break;
  222.     }
  223.       args.alpha  = param[3].data.d_float;
  224.       args.radius = param[4].data.d_float;
  225.       args.filter = param[5].data.d_int32;
  226.  
  227.       if (pluginCore (&args) == -1)
  228.     {
  229.       rvals[0].data.d_status = GIMP_PDB_EXECUTION_ERROR;
  230.       break;
  231.     }
  232.       break;
  233.  
  234.     case GIMP_RUN_WITH_LAST_VALS:
  235.       INIT_I18N();
  236.       /* XXX: add code here for last-values running */
  237.       if (pluginCore (&args) == -1)
  238.     {
  239.       rvals[0].data.d_status = GIMP_PDB_EXECUTION_ERROR;
  240.     }
  241.       break;
  242.   }
  243. }
  244.  
  245. static gint
  246. pluginCore (struct piArgs *argp)
  247. {
  248.   GimpDrawable *drw;
  249.   GimpPixelRgn srcPr, dstPr;
  250.   guchar *srcbuf, *dstbuf;
  251.   guchar *lastrow, *thisrow, *nextrow, *temprow;
  252.   guint width, height, Bpp;
  253.   gint filtno, y, rowsize, exrowsize, p_update;
  254.  
  255.   drw = gimp_drawable_get (argp->drw);
  256.  
  257.   width = drw->width;
  258.   height = drw->height;
  259.   Bpp = drw->bpp;
  260.   rowsize = width * Bpp;
  261.   exrowsize = (width + 2) * Bpp;
  262.   p_update = width / 20 + 1;
  263.  
  264.   gimp_pixel_rgn_init (&srcPr, drw, 0, 0, width, height, FALSE, FALSE);
  265.   gimp_pixel_rgn_init (&dstPr, drw, 0, 0, width, height, TRUE, TRUE);
  266.  
  267.   /* source buffer gives one pixel margin all around destination buffer */
  268.   srcbuf = g_new0 (guchar, exrowsize * 3);
  269.   dstbuf = g_new0 (guchar, rowsize);
  270.  
  271.   /* pointers to second pixel in each source row */
  272.   lastrow = srcbuf + Bpp;
  273.   thisrow = lastrow + exrowsize;
  274.   nextrow = thisrow + exrowsize;
  275.  
  276.   filtno = nlfiltInit (argp->alpha, argp->radius, argp->filter);
  277.   gimp_progress_init (_("NL Filter"));
  278.  
  279.   /* first row */
  280.   gimp_pixel_rgn_get_row (&srcPr, thisrow, 0, 0, width);
  281.   /* copy thisrow[0] to thisrow[-1], thisrow[width-1] to thisrow[width] */
  282.   memcpy (thisrow - Bpp, thisrow, Bpp);
  283.   memcpy (thisrow + rowsize, thisrow + rowsize - Bpp, Bpp);
  284.   /* copy whole thisrow to lastrow */
  285.   memcpy (lastrow - Bpp, thisrow - Bpp, exrowsize);
  286.  
  287.   for (y = 0; y < height - 1; y++)
  288.     {
  289.       if ((y % p_update) == 0)
  290.     gimp_progress_update ((gdouble) y / (gdouble) height);
  291.  
  292.       gimp_pixel_rgn_get_row (&srcPr, nextrow, 0, y + 1, width);
  293.       memcpy (nextrow - Bpp, nextrow, Bpp);
  294.       memcpy (nextrow + rowsize, nextrow + rowsize - Bpp, Bpp);
  295.       nlfiltRow (lastrow, thisrow, nextrow, dstbuf, width, Bpp, filtno);
  296.       gimp_pixel_rgn_set_row (&dstPr, dstbuf, 0, y, width);
  297.       /* rotate row buffers */
  298.       temprow = lastrow; lastrow = thisrow;
  299.       thisrow = nextrow; nextrow = temprow;
  300.     }
  301.  
  302.   /* last row */
  303.   memcpy (nextrow - Bpp, thisrow - Bpp, exrowsize);
  304.   nlfiltRow (lastrow, thisrow, nextrow, dstbuf, width, Bpp, filtno);
  305.   gimp_pixel_rgn_set_row (&dstPr, dstbuf, 0, height - 1, width);
  306.  
  307.   g_free (srcbuf);
  308.   g_free (dstbuf);
  309.  
  310.   gimp_drawable_flush (drw);
  311.   gimp_drawable_merge_shadow (drw->id, TRUE);
  312.   gimp_drawable_update (drw->id, 0, 0, width, height);
  313.   gimp_displays_flush ();
  314.  
  315.   return 0;
  316. }
  317.  
  318. gboolean run_flag = FALSE;
  319.  
  320. static void
  321. nlfilt_ok_callback (GtkWidget *widget,
  322.             gpointer   data)
  323. {
  324.   run_flag = TRUE;
  325.  
  326.   gtk_widget_destroy (GTK_WIDGET (data));
  327. }
  328.  
  329. static void
  330. nlfilt_radio_button_update (GtkWidget *widget,
  331.                 gpointer   data)
  332. {
  333.   gimp_radio_button_update (widget, data);
  334.  
  335.   if (do_preview && GTK_TOGGLE_BUTTON (widget)->active)
  336.     nlfilt_do_preview (NULL);
  337. }
  338.  
  339. static void
  340. nlfilt_double_adjustment_update (GtkAdjustment *adjustment,
  341.                  gpointer       data)
  342. {
  343.   gimp_double_adjustment_update (adjustment, data);
  344.  
  345.   if (do_preview)
  346.     nlfilt_do_preview (NULL);
  347. }
  348.  
  349. static gint
  350. pluginCoreIA (struct piArgs *argp)
  351. {
  352.   gint retval = -1; /* default to error return */
  353.   GtkWidget *dlg;
  354.   GtkWidget *main_vbox;
  355.   GtkWidget *frame;
  356.   GtkWidget *hbox;
  357.   GtkWidget *table;
  358.   GtkWidget *preview;
  359.   GtkObject *adj;
  360.  
  361.   gimp_ui_init ("nlfilt", TRUE);
  362.  
  363.   dlg = gimp_dialog_new (_("NL Filter"), "nlfilt",
  364.              gimp_standard_help_func, "filters/nlfilt.html",
  365.              GTK_WIN_POS_MOUSE,
  366.              FALSE, TRUE, FALSE,
  367.  
  368.              _("OK"), nlfilt_ok_callback,
  369.              NULL, NULL, NULL, TRUE, FALSE,
  370.              _("Cancel"), gtk_widget_destroy,
  371.              NULL, 1, NULL, FALSE, TRUE,
  372.  
  373.              NULL);
  374.  
  375.   gtk_signal_connect (GTK_OBJECT (dlg), "destroy",
  376.               GTK_SIGNAL_FUNC (gtk_main_quit),
  377.               NULL);
  378.  
  379.   main_vbox = gtk_vbox_new (FALSE, 4);
  380.   gtk_container_set_border_width (GTK_CONTAINER (main_vbox), 6);
  381.   gtk_box_pack_start (GTK_BOX (GTK_DIALOG (dlg)->vbox), main_vbox,
  382.               TRUE, TRUE, 0);
  383.   gtk_widget_show (main_vbox);
  384.  
  385.   hbox = gtk_hbox_new (FALSE, 4);
  386.   gtk_box_pack_start (GTK_BOX (main_vbox), hbox, FALSE, FALSE, 0);
  387.   gtk_widget_show (hbox);
  388.  
  389.   preview = mw_preview_new (hbox, thePreview);
  390.   gtk_object_set_data (GTK_OBJECT (preview), "piArgs", argp);
  391.   nlfilt_do_preview (preview);
  392.  
  393.   frame = gimp_radio_group_new2 (TRUE, _("Filter"),
  394.                  nlfilt_radio_button_update,
  395.                  &argp->filter, (gpointer) argp->filter,
  396.  
  397.                  _("Alpha Trimmed Mean"),
  398.                  (gpointer) filter_alpha_trim, NULL,
  399.                  _("Optimal Estimation"),
  400.                  (gpointer) filter_opt_est, NULL,
  401.                  _("Edge Enhancement"),
  402.                  (gpointer) filter_edge_enhance, NULL,
  403.  
  404.                  NULL);
  405.  
  406.   gtk_box_pack_start (GTK_BOX (hbox), frame, FALSE, FALSE, 0);
  407.   gtk_widget_show (frame);
  408.  
  409.   frame = gtk_frame_new (_("Parameter Settings"));
  410.   gtk_frame_set_shadow_type(GTK_FRAME(frame), GTK_SHADOW_ETCHED_IN);
  411.   gtk_box_pack_start (GTK_BOX (main_vbox), frame, FALSE, FALSE, 0);
  412.   gtk_widget_show (frame);
  413.  
  414.   table = gtk_table_new (2, 3, FALSE);
  415.   gtk_table_set_col_spacings (GTK_TABLE (table), 4);
  416.   gtk_table_set_row_spacings (GTK_TABLE (table), 2);
  417.   gtk_container_set_border_width (GTK_CONTAINER (table), 4);
  418.   gtk_container_add (GTK_CONTAINER (frame), table);
  419.  
  420.   adj = gimp_scale_entry_new (GTK_TABLE (table), 0, 0,
  421.                   _("Alpha:"), 0, 0,
  422.                   argp->alpha, 0.0, 1.0, 0.05, 0.1, 2,
  423.                   TRUE, 0, 0,
  424.                   NULL, NULL);
  425.   gtk_signal_connect (GTK_OBJECT (adj), "value_changed",
  426.               GTK_SIGNAL_FUNC (nlfilt_double_adjustment_update),
  427.               &argp->alpha);
  428.  
  429.   adj = gimp_scale_entry_new (GTK_TABLE (table), 0, 1,
  430.                   _("Radius:"), 0, 0,
  431.                   argp->radius, 1.0 / 3.0, 1.0, 0.05, 0.1, 2,
  432.                   TRUE, 0, 0,
  433.                   NULL, NULL);
  434.   gtk_signal_connect (GTK_OBJECT (adj), "value_changed",
  435.               GTK_SIGNAL_FUNC (nlfilt_double_adjustment_update),
  436.               &argp->radius);
  437.  
  438.   gtk_widget_show (table);
  439.  
  440.   gtk_widget_show (dlg);
  441.  
  442.   gtk_main ();
  443.   gdk_flush ();
  444.  
  445.   if (run_flag)
  446.     {
  447. #if 0
  448.       fprintf (stderr, "running:\n");
  449.       fprintf (stderr, "\t(image %d)\n", argp->img);
  450.       fprintf (stderr, "\t(drawable %d)\n", argp->drw);
  451.       fprintf (stderr, "\t(alpha %f)\n", argp->alpha);
  452.       fprintf (stderr, "\t(radius %f)\n", argp->radius);
  453. #endif
  454.       return pluginCore (argp);
  455.     }
  456.   else
  457.     {
  458.       return retval;
  459.     }
  460. }
  461.  
  462. static void
  463. nlfilt_do_preview (GtkWidget *w)
  464. {
  465.   static GtkWidget *theWidget = NULL;
  466.   struct piArgs *ap;
  467.   guchar *dst, *src0, *src1, *src2;
  468.   gint y, rowsize, filtno;
  469.   
  470.   if (theWidget == NULL)
  471.     {
  472.       theWidget = w;
  473.     }
  474.  
  475.   ap = gtk_object_get_data (GTK_OBJECT (theWidget), "piArgs");
  476.  
  477.   rowsize = thePreview->width * thePreview->bpp;
  478.   filtno =  nlfiltInit (ap->alpha, ap->radius, ap->filter);
  479.  
  480.   src0 = thePreview->bits + thePreview->bpp;
  481.   src1 = src0 + rowsize;
  482.   src2 = src1 + rowsize;
  483.   dst = g_malloc (rowsize);
  484.  
  485.   /* for preview, don't worry about edge effects */
  486.   for (y = 1; y < thePreview->height - 1; y++)
  487.     {
  488.       nlfiltRow (src0, src1, src2, dst + thePreview->bpp,
  489.                  thePreview->width - 2, thePreview->bpp, filtno);
  490.       gtk_preview_draw_row (GTK_PREVIEW (theWidget),
  491.                 dst + thePreview->bpp, 1, y, thePreview->width - 2);
  492.       src0 = src1; src1 = src2; src2 += rowsize;
  493.     }
  494.  
  495.   gtk_widget_draw (theWidget, NULL);
  496.   gdk_flush ();
  497.   g_free (dst);
  498. }
  499.  
  500. static void
  501. mw_preview_toggle_callback (GtkWidget *widget,
  502.                             gpointer   data)
  503. {
  504.   gimp_toggle_button_update (widget, data);
  505.  
  506.   if (do_preview)
  507.     nlfilt_do_preview (NULL);
  508. }
  509.  
  510. static struct mwPreview *
  511. mw_preview_build_virgin (GimpDrawable *drw)
  512. {
  513.   struct mwPreview *mwp;
  514.  
  515.   mwp = g_new (struct mwPreview, 1);
  516.  
  517.   if (drw->width > drw->height)
  518.     {
  519.       mwp->scale  = (gdouble) drw->width / (gdouble) PREVIEW_SIZE;
  520.       mwp->width  = PREVIEW_SIZE;
  521.       mwp->height = drw->height / mwp->scale;
  522.     }
  523.   else
  524.     {
  525.       mwp->scale  = (gdouble) drw->height / (gdouble) PREVIEW_SIZE;
  526.       mwp->height = PREVIEW_SIZE;
  527.       mwp->width  = drw->width / mwp->scale;
  528.     }
  529.  
  530.   mwp->bpp  = 3;
  531.   mwp->bits = NULL;
  532.  
  533.   return mwp;
  534. }
  535.  
  536. static struct mwPreview *
  537. mw_preview_build (GimpDrawable *drw)
  538. {
  539.   struct mwPreview *mwp;
  540.   gint x, y, b;
  541.   guchar *bc;
  542.   guchar *drwBits;
  543.   GimpPixelRgn pr;
  544.  
  545.   mwp = mw_preview_build_virgin (drw);
  546.  
  547.   gimp_pixel_rgn_init (&pr, drw, 0, 0, drw->width, drw->height, FALSE, FALSE);
  548.   drwBits = g_new (guchar, drw->width * drw->bpp);
  549.  
  550.   bc = mwp->bits = g_new (guchar, mwp->width * mwp->height * mwp->bpp);
  551.   for (y = 0; y < mwp->height; y++)
  552.     {
  553.       gimp_pixel_rgn_get_row (&pr, drwBits, 0, (int)(y*mwp->scale), drw->width);
  554.  
  555.       for (x = 0; x < mwp->width; x++)
  556.         {
  557.           for (b = 0; b < mwp->bpp; b++)
  558.             *bc++ = *(drwBits +
  559.                       ((gint) (x * mwp->scale) * drw->bpp) + b % drw->bpp);
  560.         }
  561.     }
  562.   g_free (drwBits);
  563.  
  564.   return mwp;
  565. }
  566.  
  567. static GtkWidget *
  568. mw_preview_new (GtkWidget        *parent,
  569.                 struct mwPreview *mwp)
  570. {
  571.   GtkWidget *preview;
  572.   GtkWidget *frame;
  573.   GtkWidget *pframe;
  574.   GtkWidget *vbox;
  575.   GtkWidget *button;
  576.   guchar *color_cube;
  577.    
  578.   gtk_preview_set_gamma (gimp_gamma ());
  579.   gtk_preview_set_install_cmap (gimp_install_cmap ());
  580.   color_cube = gimp_color_cube ();
  581.   gtk_preview_set_color_cube (color_cube[0], color_cube[1],
  582.                               color_cube[2], color_cube[3]);
  583.  
  584.   gtk_widget_set_default_visual (gtk_preview_get_visual ());
  585.   gtk_widget_set_default_colormap (gtk_preview_get_cmap ());
  586.  
  587.   frame = gtk_frame_new (_("Preview"));
  588.   gtk_frame_set_shadow_type (GTK_FRAME (frame), GTK_SHADOW_ETCHED_IN);
  589.   gtk_box_pack_start (GTK_BOX (parent), frame, FALSE, FALSE, 0);
  590.   gtk_widget_show (frame);
  591.  
  592.   vbox = gtk_vbox_new (FALSE, 2);
  593.   gtk_container_set_border_width (GTK_CONTAINER (vbox), 4);
  594.   gtk_container_add (GTK_CONTAINER (frame), vbox);
  595.   gtk_widget_show (vbox);
  596.  
  597.   pframe = gtk_frame_new (NULL);
  598.   gtk_frame_set_shadow_type (GTK_FRAME(pframe), GTK_SHADOW_IN);
  599.   gtk_box_pack_start (GTK_BOX (vbox), pframe, FALSE, FALSE, 0);
  600.   gtk_widget_show (pframe);
  601.  
  602.   preview = gtk_preview_new (GTK_PREVIEW_COLOR);
  603.   gtk_preview_size (GTK_PREVIEW (preview), mwp->width, mwp->height);
  604.   gtk_container_add (GTK_CONTAINER (pframe), preview);
  605.   gtk_widget_show (preview);
  606.  
  607.   button = gtk_check_button_new_with_label (_("Do Preview"));
  608.   gtk_toggle_button_set_active (GTK_TOGGLE_BUTTON (button), do_preview);
  609.   gtk_signal_connect (GTK_OBJECT (button), "toggled",
  610.                       GTK_SIGNAL_FUNC (mw_preview_toggle_callback),
  611.                       &do_preview);
  612.   gtk_box_pack_start (GTK_BOX (vbox), button, FALSE, FALSE, 0);
  613.   gtk_widget_show (button);
  614.  
  615.   return preview;
  616. }
  617.  
  618. /* pnmnlfilt.c - 4 in 1 (2 non-linear) filter
  619. **             - smooth an anyimage
  620. **             - do alpha trimmed mean filtering on an anyimage
  621. **             - do optimal estimation smoothing on an anyimage
  622. **             - do edge enhancement on an anyimage
  623. **
  624. ** Version 1.0
  625. **
  626. ** The implementation of an alpha-trimmed mean filter
  627. ** is based on the description in IEEE CG&A May 1990
  628. ** Page 23 by Mark E. Lee and Richard A. Redner.
  629. **
  630. ** The paper recommends using a hexagon sampling region around each
  631. ** pixel being processed, allowing an effective sub pixel radius to be
  632. ** specified. The hexagon values are sythesised by area sampling the
  633. ** rectangular pixels with a hexagon grid. The seven hexagon values
  634. ** obtained from the 3x3 pixel grid are used to compute the alpha
  635. ** trimmed mean. Note that an alpha value of 0.0 gives a conventional
  636. ** mean filter (where the radius controls the contribution of
  637. ** surrounding pixels), while a value of 0.5 gives a median filter.
  638. ** Although there are only seven values to trim from before finding
  639. ** the mean, the algorithm has been extended from that described in
  640. ** CG&A by using interpolation, to allow a continuous selection of
  641. ** alpha value between and including 0.0 to 0.5  The useful values
  642. ** for radius are between 0.3333333 (where the filter will have no
  643. ** effect because only one pixel is sampled), to 1.0, where all
  644. ** pixels in the 3x3 grid are sampled.
  645. **
  646. ** The optimal estimation filter is taken from an article "Converting Dithered
  647. ** Images Back to Gray Scale" by Allen Stenger, Dr Dobb's Journal, November
  648. ** 1992, and this article references "Digital Image Enhancement andNoise Filtering by
  649. ** Use of Local Statistics", Jong-Sen Lee, IEEE Transactions on Pattern Analysis and
  650. ** Machine Intelligence, March 1980.
  651. **
  652. ** Also borrow the  technique used in pgmenhance(1) to allow edge
  653. ** enhancement if the alpha value is negative.
  654. **
  655. ** Author:
  656. **         Graeme W. Gill, 30th Jan 1993
  657. **         graeme@labtam.oz.au
  658. **
  659. ** Permission to use, copy, modify, and distribute this software and its
  660. ** documentation for any purpose and without fee is hereby granted, provided
  661. ** that the above copyright notice appear in all copies and that both that
  662. ** copyright notice and this permission notice appear in supporting
  663. ** documentation.  This software is provided "as is" without express or
  664. ** implied warranty.
  665. */
  666.  
  667. /* ************************************************** */
  668. /* Hexagon intersecting square area functions */
  669. /* Compute the area of the intersection of a triangle */
  670. /* and a rectangle */
  671.  
  672. gdouble triang_area(gdouble, gdouble, gdouble, gdouble, gdouble,
  673.                    gdouble, gdouble, gdouble, gint);
  674. gdouble rectang_area(gdouble, gdouble, gdouble, gdouble,
  675.                     gdouble, gdouble, gdouble, gdouble);
  676. gdouble hex_area(gdouble, gdouble, gdouble, gdouble, gdouble);
  677.  
  678. gint atfilt0(gint *p);
  679. gint atfilt1(gint *p);
  680. gint atfilt2(gint *p);
  681. gint atfilt3(gint *p);
  682. gint atfilt4(gint *p);
  683. gint atfilt5(gint *p);
  684. gint (*atfuncs[6])(gint *) = {
  685.   atfilt0,
  686.   atfilt1,
  687.   atfilt2,
  688.   atfilt3,
  689.   atfilt4,
  690.   atfilt5
  691. };
  692.  
  693. gint noisevariance;      /* global so that pixel processing code can get at it quickly */
  694.  
  695. #define MXIVAL 255    /* maximum input value */
  696. #define NOIVAL (MXIVAL + 1)             /* number of possible input values */
  697.  
  698. #define SCALEB 8                                /* scale bits */
  699. #define SCALE (1 << SCALEB)     /* scale factor */
  700. #define MXSVAL (MXIVAL * SCALE) /* maximum scaled values */
  701.  
  702. #define CSCALEB 2                               /* coarse scale bits */
  703. #define CSCALE (1 << CSCALEB)   /* coarse scale factor */
  704. #define MXCSVAL (MXIVAL * CSCALE)       /* maximum coarse scaled values */
  705. #define NOCSVAL (MXCSVAL + 1)   /* number of coarse scaled values */
  706. #define SCTOCSC(x) ((x) >> (SCALEB - CSCALEB))  /* convert from scaled to coarse scaled */
  707. #define CSCTOSC(x) ((x) << (SCALEB - CSCALEB))  /* convert from course scaled to scaled */
  708.  
  709. #ifndef MAXINT
  710. # define MAXINT 0x7fffffff      /* assume this is a 32 bit machine */
  711. #endif
  712.  
  713. /* round and scale floating point to scaled integer */
  714. #define SROUND(x) ((gint)(((x) * (gdouble)SCALE) + 0.5))
  715. /* round and un-scale scaled integer value */
  716. #define RUNSCALE(x) (((x) + (1 << (SCALEB-1))) >> SCALEB)       /* rounded un-scale */
  717. #define UNSCALE(x) ((x) >> SCALEB)
  718.  
  719. /* Note: modified by David Hodson, nlfiltRow now accesses
  720.  * srclast, srcthis, and srcnext from [-Bpp] to [width*Bpp-1].
  721.  * Beware if you use this code anywhere else!
  722.  */
  723. static void
  724. nlfiltRow(guchar *srclast, guchar *srcthis, guchar *srcnext, guchar *dst,
  725.           gint width, gint Bpp, gint filtno) {
  726.  
  727.    gint pf[9];
  728.    guchar *ip0, *ip1, *ip2, *or, *orend;
  729.  
  730.    or = dst;
  731.    orend = dst + width * Bpp;
  732.    ip0 = srclast;
  733.    ip1 = srcthis;
  734.    ip2 = srcnext;
  735.    for (or = dst; or < orend; ip0++, ip1++, ip2++, or++) {
  736.       pf[0] = *ip1;
  737.       pf[1] = *(ip1 - Bpp);
  738.       pf[2] = *(ip2 - Bpp);
  739.       pf[3] = *(ip2);
  740.       pf[4] = *(ip2 + Bpp);
  741.       pf[5] = *(ip1 + Bpp);
  742.       pf[6] = *(ip0 + Bpp);
  743.       pf[7] = *(ip0);
  744.       pf[8] = *(ip0 - Bpp);
  745.       *or=(atfuncs[filtno])(pf);
  746.    }
  747. }
  748.  
  749. /* We restrict radius to the values: 0.333333 <= radius <= 1.0 */
  750. /* so that no fewer and no more than a 3x3 grid of pixels around */
  751. /* the pixel in question needs to be read. Given this, we only */
  752. /* need 3 or 4 weightings per hexagon, as follows: */
  753. /*                  _ _                         */
  754. /* Virtical hex:   |_|_|  1 2                   */
  755. /*                 |X|_|  0 3                   */
  756. /*                                       _      */
  757. /*              _                      _|_|   1 */
  758. /* Middle hex: |_| 1  Horizontal hex: |X|_| 0 2 */
  759. /*             |X| 0                    |_|   3 */
  760. /*             |_| 2                            */
  761.  
  762. /* all filters */
  763. gint V0[NOIVAL],V1[NOIVAL],V2[NOIVAL],V3[NOIVAL];   /* vertical hex */
  764. gint M0[NOIVAL],M1[NOIVAL],M2[NOIVAL];              /* middle hex */
  765. gint H0[NOIVAL],H1[NOIVAL],H2[NOIVAL],H3[NOIVAL];   /* horizontal hex */
  766.  
  767. /* alpha trimmed and edge enhancement only */
  768. gint ALFRAC[NOIVAL * 8];               /* fractional alpha divider table */
  769.  
  770. /* optimal estimation only */
  771. gint AVEDIV[7 * NOCSVAL];              /* divide by 7 to give average value */
  772. gint SQUARE[2 * NOCSVAL];              /* scaled square lookup table */
  773.  
  774. /* Table initialisation function - return alpha range */
  775. static inline gint
  776. nlfiltInit(gdouble alpha, gdouble radius, FilterType filter) {
  777.    gint alpharange;                 /* alpha range value 0 - 3 */
  778.    gdouble meanscale;               /* scale for finding mean */
  779.    gdouble mmeanscale;              /* scale for finding mean - midle hex */
  780.    gdouble alphafraction;   /* fraction of next largest/smallest
  781.                              *  to subtract from sum
  782.                              */
  783.    switch (filter) {
  784.        case filter_alpha_trim: {
  785.           gdouble noinmean;
  786.           /* alpha only makes sense in range 0.0 - 0.5 */
  787.           alpha /= 2.0;
  788.               /* number of elements (out of a possible 7) used in the mean */
  789.           noinmean = ((0.5 - alpha) * 12.0) + 1.0;
  790.           mmeanscale = meanscale = 1.0/noinmean;
  791.           if (alpha == 0.0) {                    /* mean filter */
  792.              alpharange = 0;
  793.              alphafraction = 0.0;            /* not used */
  794.           } else if (alpha < (1.0/6.0)) {    /* mean of 5 to 7 middle values */
  795.              alpharange = 1;
  796.              alphafraction = (7.0 - noinmean)/2.0;
  797.           } else if (alpha < (1.0/3.0)) {    /* mean of 3 to 5 middle values */
  798.              alpharange = 2;
  799.              alphafraction = (5.0 - noinmean)/2.0;
  800.           } else {                           /* mean of 1 to 3 middle values */
  801.                                              /* alpha==0.5  => median filter */
  802.              alpharange = 3;
  803.              alphafraction = (3.0 - noinmean)/2.0;
  804.           }
  805.        }
  806.        break;
  807.        case filter_opt_est: {
  808.           gint i;
  809.           gdouble noinmean = 7.0;
  810.  
  811.               /* edge enhancement function */
  812.           alpharange = 5;
  813.  
  814.               /* compute scaled hex values */
  815.           mmeanscale=meanscale=1.0;
  816.  
  817.               /* Set up 1:1 division lookup - not used */
  818.           alphafraction=1.0/noinmean;
  819.  
  820.               /* estimate of noise variance */
  821.           noisevariance = alpha * (gdouble)255;
  822.           noisevariance = noisevariance * noisevariance / 8.0;
  823.  
  824.               /* set yp optimal estimation specific stuff */
  825.  
  826.           for (i=0;i<(7*NOCSVAL);i++) { /* divide scaled value by 7 lookup */
  827.              AVEDIV[i] = CSCTOSC(i)/7;       /* scaled divide by 7 */
  828.           }
  829.               /* compute square and rescale by
  830.                * (val >> (2 * SCALEB + 2)) table
  831.                */
  832.           for (i=0;i<(2*NOCSVAL);i++) {
  833.              gint val;
  834.                  /* NOCSVAL offset to cope with -ve input values */
  835.              val = CSCTOSC(i - NOCSVAL);
  836.              SQUARE[i] = (val * val) >> (2 * SCALEB + 2);
  837.           }
  838.        }
  839.        break;
  840.        case filter_edge_enhance: {
  841.           if (alpha == 1.0) alpha = 0.99;
  842.           alpharange = 4;
  843.               /* mean of 7 and scaled by -alpha/(1-alpha) */
  844.           meanscale = 1.0 * (-alpha/((1.0 - alpha) * 7.0));
  845.  
  846.               /* middle pixel has 1/(1-alpha) as well */
  847.           mmeanscale = 1.0 * (1.0/(1.0 - alpha) + meanscale);
  848.           alphafraction = 0.0;    /* not used */
  849.        }
  850.        break;
  851.        default:
  852.           fprintf(stderr, "unknown filter %d\n", filter);
  853.           return -1;
  854.    }
  855.        /*
  856.         * Setup pixel weighting tables -
  857.         * note we pre-compute mean division here too.
  858.         */
  859.    {
  860.       gint i;
  861.       gdouble hexhoff,hexvoff;
  862.       gdouble tabscale,mtabscale;
  863.       gdouble v0,v1,v2,v3,m0,m1,m2,h0,h1,h2,h3;
  864.  
  865.           /* horizontal offset of virtical hex centers */
  866.       hexhoff = radius/2;
  867.  
  868.           /* vertical offset of virtical hex centers */
  869.       hexvoff = 3.0 * radius/sqrt(12.0);
  870.  
  871.           /*
  872.            * scale tables to normalise by hexagon
  873.            * area, and number of hexes used in mean
  874.            */
  875.       tabscale = meanscale / (radius * hexvoff);
  876.       mtabscale = mmeanscale / (radius * hexvoff);
  877.       v0 = hex_area(0.0,0.0,hexhoff,hexvoff,radius) * tabscale;
  878.       v1 = hex_area(0.0,1.0,hexhoff,hexvoff,radius) * tabscale;
  879.       v2 = hex_area(1.0,1.0,hexhoff,hexvoff,radius) * tabscale;
  880.       v3 = hex_area(1.0,0.0,hexhoff,hexvoff,radius) * tabscale;
  881.       m0 = hex_area(0.0,0.0,0.0,0.0,radius) * mtabscale;
  882.       m1 = hex_area(0.0,1.0,0.0,0.0,radius) * mtabscale;
  883.       m2 = hex_area(0.0,-1.0,0.0,0.0,radius) * mtabscale;
  884.       h0 = hex_area(0.0,0.0,radius,0.0,radius) * tabscale;
  885.       h1 = hex_area(1.0,1.0,radius,0.0,radius) * tabscale;
  886.       h2 = hex_area(1.0,0.0,radius,0.0,radius) * tabscale;
  887.       h3 = hex_area(1.0,-1.0,radius,0.0,radius) * tabscale;
  888.  
  889.       for (i=0; i <= MXIVAL; i++) {
  890.          gdouble fi;
  891.          fi = (gdouble)i;
  892.          V0[i] = SROUND(fi * v0);
  893.          V1[i] = SROUND(fi * v1);
  894.          V2[i] = SROUND(fi * v2);
  895.          V3[i] = SROUND(fi * v3);
  896.          M0[i] = SROUND(fi * m0);
  897.          M1[i] = SROUND(fi * m1);
  898.          M2[i] = SROUND(fi * m2);
  899.          H0[i] = SROUND(fi * h0);
  900.          H1[i] = SROUND(fi * h1);
  901.          H2[i] = SROUND(fi * h2);
  902.          H3[i] = SROUND(fi * h3);
  903.       }
  904.           /* set up alpha fraction lookup table used on big/small */
  905.       for (i=0; i < (NOIVAL * 8); i++) {
  906.          ALFRAC[i] = SROUND((gdouble)i * alphafraction);
  907.       }
  908.    }
  909.    return alpharange;
  910. }
  911.  
  912. /* Core pixel processing function - hand it 3x3 pixels and return result. */
  913. /* Mean filter */
  914. gint
  915. atfilt0(gint32 *p) {
  916.    gint retv;
  917.        /* map to scaled hexagon values */
  918.    retv = M0[p[0]] + M1[p[3]] + M2[p[7]];
  919.    retv += H0[p[0]] + H1[p[2]] + H2[p[1]] + H3[p[8]];
  920.    retv += V0[p[0]] + V1[p[3]] + V2[p[2]] + V3[p[1]];
  921.    retv += V0[p[0]] + V1[p[3]] + V2[p[4]] + V3[p[5]];
  922.    retv += H0[p[0]] + H1[p[4]] + H2[p[5]] + H3[p[6]];
  923.    retv += V0[p[0]] + V1[p[7]] + V2[p[6]] + V3[p[5]];
  924.    retv += V0[p[0]] + V1[p[7]] + V2[p[8]] + V3[p[1]];
  925.    return UNSCALE(retv);
  926. }
  927.  
  928. /* Mean of 5 - 7 middle values */
  929. gint
  930. atfilt1(gint32 *p) {
  931.    gint h0,h1,h2,h3,h4,h5,h6;       /* hexagon values    2 3   */
  932.                                     /*                  1 0 4  */
  933.                                     /*                   6 5   */
  934.    gint big,small;
  935.        /* map to scaled hexagon values */
  936.    h0 = M0[p[0]] + M1[p[3]] + M2[p[7]];
  937.    h1 = H0[p[0]] + H1[p[2]] + H2[p[1]] + H3[p[8]];
  938.    h2 = V0[p[0]] + V1[p[3]] + V2[p[2]] + V3[p[1]];
  939.    h3 = V0[p[0]] + V1[p[3]] + V2[p[4]] + V3[p[5]];
  940.    h4 = H0[p[0]] + H1[p[4]] + H2[p[5]] + H3[p[6]];
  941.    h5 = V0[p[0]] + V1[p[7]] + V2[p[6]] + V3[p[5]];
  942.    h6 = V0[p[0]] + V1[p[7]] + V2[p[8]] + V3[p[1]];
  943.        /* sum values and also discover the largest and smallest */
  944.    big = small = h0;
  945. #define CHECK(xx) \
  946.         h0 += xx; \
  947.         if (xx > big) \
  948.                 big = xx; \
  949.         else if (xx < small) \
  950.                 small = xx;
  951.         CHECK(h1)
  952.         CHECK(h2)
  953.         CHECK(h3)
  954.         CHECK(h4)
  955.         CHECK(h5)
  956.         CHECK(h6)
  957. #undef CHECK
  958.        /* Compute mean of middle 5-7 values */
  959.    return UNSCALE(h0 -ALFRAC[(big + small)>>SCALEB]);
  960. }
  961.  
  962. /* Mean of 3 - 5 middle values */
  963. gint
  964. atfilt2(gint32 *p) {
  965.    gint h0,h1,h2,h3,h4,h5,h6;       /* hexagon values    2 3   */
  966.                                     /*                  1 0 4  */
  967.                                     /*                   6 5   */
  968.    gint big0,big1,small0,small1;
  969.        /* map to scaled hexagon values */
  970.    h0 = M0[p[0]] + M1[p[3]] + M2[p[7]];
  971.    h1 = H0[p[0]] + H1[p[2]] + H2[p[1]] + H3[p[8]];
  972.    h2 = V0[p[0]] + V1[p[3]] + V2[p[2]] + V3[p[1]];
  973.    h3 = V0[p[0]] + V1[p[3]] + V2[p[4]] + V3[p[5]];
  974.    h4 = H0[p[0]] + H1[p[4]] + H2[p[5]] + H3[p[6]];
  975.    h5 = V0[p[0]] + V1[p[7]] + V2[p[6]] + V3[p[5]];
  976.    h6 = V0[p[0]] + V1[p[7]] + V2[p[8]] + V3[p[1]];
  977.        /* sum values and also discover the 2 largest and 2 smallest */
  978.    big0 = small0 = h0;
  979.    small1 = MAXINT;
  980.    big1 = 0;
  981. #define CHECK(xx) \
  982.         h0 += xx; \
  983.         if (xx > big1) \
  984.                 { \
  985.                 if (xx > big0) \
  986.                         { \
  987.                         big1 = big0; \
  988.                         big0 = xx; \
  989.                         } \
  990.                 else \
  991.                         big1 = xx; \
  992.                 } \
  993.         if (xx < small1) \
  994.                 { \
  995.                 if (xx < small0) \
  996.                         { \
  997.                         small1 = small0; \
  998.                         small0 = xx; \
  999.                         } \
  1000.                 else \
  1001.                         small1 = xx; \
  1002.                 }
  1003.         CHECK(h1)
  1004.         CHECK(h2)
  1005.         CHECK(h3)
  1006.         CHECK(h4)
  1007.         CHECK(h5)
  1008.         CHECK(h6)
  1009. #undef CHECK
  1010.        /* Compute mean of middle 3-5 values */
  1011.   return UNSCALE(h0 -big0 -small0 -ALFRAC[(big1 + small1)>>SCALEB]);
  1012. }
  1013.  
  1014. /*
  1015.  * Mean of 1 - 3 middle values.
  1016.  * If only 1 value, then this is a median filter.
  1017.  */
  1018. gint32
  1019. atfilt3(gint32 *p) {
  1020.    gint h0,h1,h2,h3,h4,h5,h6;       /* hexagon values    2 3   */
  1021.                                    /*                  1 0 4  */
  1022.                                    /*                   6 5   */
  1023.    gint big0,big1,big2,small0,small1,small2;
  1024.        /* map to scaled hexagon values */
  1025.    h0 = M0[p[0]] + M1[p[3]] + M2[p[7]];
  1026.    h1 = H0[p[0]] + H1[p[2]] + H2[p[1]] + H3[p[8]];
  1027.    h2 = V0[p[0]] + V1[p[3]] + V2[p[2]] + V3[p[1]];
  1028.    h3 = V0[p[0]] + V1[p[3]] + V2[p[4]] + V3[p[5]];
  1029.    h4 = H0[p[0]] + H1[p[4]] + H2[p[5]] + H3[p[6]];
  1030.    h5 = V0[p[0]] + V1[p[7]] + V2[p[6]] + V3[p[5]];
  1031.    h6 = V0[p[0]] + V1[p[7]] + V2[p[8]] + V3[p[1]];
  1032.        /* sum values and also discover the 3 largest and 3 smallest */
  1033.    big0 = small0 = h0;
  1034.    small1 = small2 = MAXINT;
  1035.    big1 = big2 = 0;
  1036. #define CHECK(xx) \
  1037.         h0 += xx; \
  1038.         if (xx > big2) \
  1039.                 { \
  1040.                 if (xx > big1) \
  1041.                         { \
  1042.                         if (xx > big0) \
  1043.                                 { \
  1044.                                 big2 = big1; \
  1045.                                 big1 = big0; \
  1046.                                 big0 = xx; \
  1047.                                 } \
  1048.                         else \
  1049.                                 { \
  1050.                                 big2 = big1; \
  1051.                                 big1 = xx; \
  1052.                                 } \
  1053.                         } \
  1054.                 else \
  1055.                         big2 = xx; \
  1056.                 } \
  1057.         if (xx < small2) \
  1058.                 { \
  1059.                 if (xx < small1) \
  1060.                         { \
  1061.                         if (xx < small0) \
  1062.                                 { \
  1063.                                 small2 = small1; \
  1064.                                 small1 = small0; \
  1065.                                 small0 = xx; \
  1066.                                 } \
  1067.                         else \
  1068.                                 { \
  1069.                                 small2 = small1; \
  1070.                                 small1 = xx; \
  1071.                                 } \
  1072.                         } \
  1073.                 else \
  1074.                         small2 = xx; \
  1075.                 }
  1076.         CHECK(h1)
  1077.         CHECK(h2)
  1078.         CHECK(h3)
  1079.         CHECK(h4)
  1080.         CHECK(h5)
  1081.         CHECK(h6)
  1082. #undef CHECK
  1083.        /* Compute mean of middle 1-3 values */
  1084.    return  UNSCALE(h0-big0-big1-small0-small1-ALFRAC[(big2+small2)>>SCALEB]);
  1085. }
  1086.  
  1087. /* Edge enhancement */
  1088. gint
  1089. atfilt4(gint *p) {
  1090.    gint hav;
  1091.        /* map to scaled hexagon values and compute enhance value */
  1092.    hav = M0[p[0]] + M1[p[3]] + M2[p[7]];
  1093.    hav += H0[p[0]] + H1[p[2]] + H2[p[1]] + H3[p[8]];
  1094.    hav += V0[p[0]] + V1[p[3]] + V2[p[2]] + V3[p[1]];
  1095.    hav += V0[p[0]] + V1[p[3]] + V2[p[4]] + V3[p[5]];
  1096.    hav += H0[p[0]] + H1[p[4]] + H2[p[5]] + H3[p[6]];
  1097.    hav += V0[p[0]] + V1[p[7]] + V2[p[6]] + V3[p[5]];
  1098.    hav += V0[p[0]] + V1[p[7]] + V2[p[8]] + V3[p[1]];
  1099.    if (hav < 0)
  1100.       hav = 0;
  1101.    hav = UNSCALE(hav);
  1102.    if (hav > (gdouble)255)
  1103.       hav = (gdouble)255;
  1104.    return hav;
  1105. }
  1106.  
  1107. /* Optimal estimation - do smoothing in inverse proportion */
  1108. /* to the local variance. */
  1109. /* notice we use the globals noisevariance */
  1110. gint
  1111. atfilt5(gint *p) {
  1112.    gint mean,variance,temp;
  1113.    gint h0,h1,h2,h3,h4,h5,h6;       /* hexagon values    2 3   */
  1114.                                    /*                  1 0 4  */
  1115.                                    /*                   6 5   */
  1116.        /* map to scaled hexagon values */
  1117.    h0 = M0[p[0]] + M1[p[3]] + M2[p[7]];
  1118.    h1 = H0[p[0]] + H1[p[2]] + H2[p[1]] + H3[p[8]];
  1119.    h2 = V0[p[0]] + V1[p[3]] + V2[p[2]] + V3[p[1]];
  1120.    h3 = V0[p[0]] + V1[p[3]] + V2[p[4]] + V3[p[5]];
  1121.    h4 = H0[p[0]] + H1[p[4]] + H2[p[5]] + H3[p[6]];
  1122.    h5 = V0[p[0]] + V1[p[7]] + V2[p[6]] + V3[p[5]];
  1123.    h6 = V0[p[0]] + V1[p[7]] + V2[p[8]] + V3[p[1]];
  1124.    mean = h0 + h1 + h2 + h3 + h4 + h5 + h6;
  1125.        /* compute scaled mean by dividing by 7 */
  1126.    mean = AVEDIV[SCTOCSC(mean)];
  1127.  
  1128.        /* compute scaled variance */
  1129.    temp = (h1 - mean); variance = SQUARE[NOCSVAL + SCTOCSC(temp)];
  1130.  
  1131.        /* and rescale to keep */
  1132.    temp = (h2 - mean); variance += SQUARE[NOCSVAL + SCTOCSC(temp)];
  1133.  
  1134.  /* within 32 bit limits */
  1135.    temp = (h3 - mean); variance += SQUARE[NOCSVAL + SCTOCSC(temp)];
  1136.    temp = (h4 - mean); variance += SQUARE[NOCSVAL + SCTOCSC(temp)];
  1137.    temp = (h5 - mean); variance += SQUARE[NOCSVAL + SCTOCSC(temp)];
  1138.    temp = (h6 - mean); variance += SQUARE[NOCSVAL + SCTOCSC(temp)];
  1139.        /* (temp = h0 - mean) */
  1140.    temp = (h0 - mean); variance += SQUARE[NOCSVAL + SCTOCSC(temp)];
  1141.    if (variance != 0)      /* avoid possible divide by 0 */
  1142.           /* optimal estimate */
  1143.       temp = mean + (variance * temp) / (variance + noisevariance);
  1144.    else temp = h0;
  1145.    if (temp < 0)
  1146.       temp = 0;
  1147.    temp = RUNSCALE(temp);
  1148.    if (temp > (gdouble)255) temp = (gdouble)255;
  1149.    return temp;
  1150. }
  1151.  
  1152.  
  1153. /* Triangle orientation is per geometric axes (not graphical axies) */
  1154.  
  1155. #define NW 0    /* North west triangle /| */
  1156. #define NE 1    /* North east triangle |\ */
  1157. #define SW 2    /* South west triangle \| */
  1158. #define SE 3    /* South east triangle |/ */
  1159. #define STH 2
  1160. #define EST 1
  1161.  
  1162. #define SWAPI(a,b) (t = a, a = -b, b = -t)
  1163.  
  1164. /* compute the area of overlap of a hexagon diameter d, */
  1165. /* centered at hx,hy, with a unit square of center sx,sy. */
  1166. gdouble
  1167. hex_area(gdouble sx, gdouble sy, gdouble hx, gdouble hy, gdouble d) {
  1168.    gdouble hx0,hx1,hx2,hy0,hy1,hy2,hy3;
  1169.    gdouble sx0,sx1,sy0,sy1;
  1170.  
  1171.        /* compute square co-ordinates */
  1172.    sx0 = sx - 0.5;
  1173.    sy0 = sy - 0.5;
  1174.    sx1 = sx + 0.5;
  1175.    sy1 = sy + 0.5;
  1176.        /* compute hexagon co-ordinates */
  1177.    hx0 = hx - d/2.0;
  1178.    hx1 = hx;
  1179.    hx2 = hx + d/2.0;
  1180.    hy0 = hy - 0.5773502692 * d;    /* d / sqrt(3) */
  1181.    hy1 = hy - 0.2886751346 * d;    /* d / sqrt(12) */
  1182.    hy2 = hy + 0.2886751346 * d;    /* d / sqrt(12) */
  1183.    hy3 = hy + 0.5773502692 * d;    /* d / sqrt(3) */
  1184.  
  1185.    return
  1186.       triang_area(sx0,sy0,sx1,sy1,hx0,hy2,hx1,hy3,NW) +
  1187.       triang_area(sx0,sy0,sx1,sy1,hx1,hy2,hx2,hy3,NE) +
  1188.       rectang_area(sx0,sy0,sx1,sy1,hx0,hy1,hx2,hy2) +
  1189.       triang_area(sx0,sy0,sx1,sy1,hx0,hy0,hx1,hy1,SW) +
  1190.       triang_area(sx0,sy0,sx1,sy1,hx1,hy0,hx2,hy1,SE);
  1191. }
  1192.  
  1193. gdouble
  1194. triang_area(gdouble rx0, gdouble ry0, gdouble rx1, gdouble ry1, gdouble tx0,
  1195.             gdouble ty0, gdouble tx1, gdouble ty1, gint tt) {
  1196.    gdouble a,b,c,d;
  1197.    gdouble lx0,ly0,lx1,ly1;
  1198.        /* Convert everything to a NW triangle */
  1199.    if (tt & STH) {
  1200.       gdouble t;
  1201.       SWAPI(ry0,ry1);
  1202.       SWAPI(ty0,ty1);
  1203.    } if (tt & EST) {
  1204.       gdouble t;
  1205.       SWAPI(rx0,rx1);
  1206.       SWAPI(tx0,tx1);
  1207.    }
  1208.        /* Compute overlapping box */
  1209.    if (tx0 > rx0)
  1210.       rx0 = tx0;
  1211.    if (ty0 > ry0)
  1212.       ry0 = ty0;
  1213.    if (tx1 < rx1)
  1214.       rx1 = tx1;
  1215.    if (ty1 < ry1)
  1216.       ry1 = ty1;
  1217.    if (rx1 <= rx0 || ry1 <= ry0)
  1218.       return 0.0;
  1219.        /* Need to compute diagonal line intersection with the box */
  1220.        /* First compute co-efficients to formulas x = a + by and y = c + dx */
  1221.    b = (tx1 - tx0)/(ty1 - ty0);
  1222.    a = tx0 - b * ty0;
  1223.    d = (ty1 - ty0)/(tx1 - tx0);
  1224.    c = ty0 - d * tx0;
  1225.  
  1226.        /* compute top or right intersection */
  1227.    tt = 0;
  1228.    ly1 = ry1;
  1229.    lx1 = a + b * ly1;
  1230.    if (lx1 <= rx0)
  1231.       return (rx1 - rx0) * (ry1 - ry0);
  1232.    else if (lx1 > rx1) {     /* could be right hand side */
  1233.       lx1 = rx1;
  1234.       ly1 = c + d * lx1;
  1235.       if (ly1 <= ry0)
  1236.          return (rx1 - rx0) * (ry1 - ry0);
  1237.       tt = 1; /* right hand side intersection */
  1238.    }
  1239.        /* compute left or bottom intersection */
  1240.    lx0 = rx0;
  1241.    ly0 = c + d * lx0;
  1242.    if (ly0 >= ry1)
  1243.       return (rx1 - rx0) * (ry1 - ry0);
  1244.    else if (ly0 < ry0) {    /* could be right hand side */
  1245.       ly0 = ry0;
  1246.       lx0 = a + b * ly0;
  1247.       if (lx0 >= rx1)
  1248.          return (rx1 - rx0) * (ry1 - ry0);
  1249.       tt |= 2;        /* bottom intersection */
  1250.    }
  1251.  
  1252.    if (tt == 0) {    /* top and left intersection */
  1253.                        /* rectangle minus triangle */
  1254.       return ((rx1 - rx0) * (ry1 - ry0))
  1255.          - (0.5 * (lx1 - rx0) * (ry1 - ly0));
  1256.    }
  1257.    else if (tt == 1) {       /* right and left intersection */
  1258.       return ((rx1 - rx0) * (ly0 - ry0))
  1259.          + (0.5 * (rx1 - rx0) * (ly1 - ly0));
  1260.    } else if (tt == 2) {      /* top and bottom intersection */
  1261.       return ((rx1 - lx1) * (ry1 - ry0))
  1262.          + (0.5 * (lx1 - lx0) * (ry1 - ry0));
  1263.    } else { /* tt == 3 */      /* right and bottom intersection */
  1264.           /* triangle */
  1265.       return (0.5 * (rx1 - lx0) * (ly1 - ry0));
  1266.    }
  1267. }
  1268.  
  1269. /* Compute rectangle area */
  1270. gdouble
  1271. rectang_area(gdouble rx0, gdouble ry0, gdouble rx1, gdouble ry1, gdouble tx0,
  1272.              gdouble ty0, gdouble tx1, gdouble ty1) {
  1273.        /* Compute overlapping box */
  1274.    if (tx0 > rx0)
  1275.       rx0 = tx0;
  1276.    if (ty0 > ry0)
  1277.       ry0 = ty0;
  1278.    if (tx1 < rx1)
  1279.       rx1 = tx1;
  1280.    if (ty1 < ry1)
  1281.       ry1 = ty1;
  1282.    if (rx1 <= rx0 || ry1 <= ry0)
  1283.       return 0.0;
  1284.    return (rx1 - rx0) * (ry1 - ry0);
  1285. }
  1286.