]> git.imager.perl.org - imager.git/blobdiff - filters.im
he unpack code for ICO/CUR file handling could extend 32-bit unsigned values to 64...
[imager.git] / filters.im
index 9ebb0af74795acd51b7b1a2b95991a08f6ffcc95..e247c9ccc4c64df20ef37251655eac3f59d00e13 100644 (file)
@@ -1,3 +1,4 @@
+#define IMAGER_NO_CONTEXT
 #include "imager.h"
 #include "imageri.h"
 #include <stdlib.h>
@@ -69,12 +70,13 @@ Scales the pixel values by the amount specified.
 
 void
 i_contrast(i_img *im, float intensity) {
-  int x, y;
+  i_img_dim x, y;
   unsigned char ch;
   unsigned int new_color;
   i_color rcolor;
+  dIMCTXim(im);
   
-  mm_log((1,"i_contrast(im %p, intensity %f)\n", im, intensity));
+  im_log((aIMCTX, 1,"i_contrast(im %p, intensity %f)\n", im, intensity));
   
   if(intensity < 0) return;
   
@@ -97,11 +99,12 @@ i_contrast(i_img *im, float intensity) {
 
 static int
 s_hardinvert_low(i_img *im, int all) {
-  int x, y;
+  i_img_dim x, y;
   int ch;
   int invert_channels = all ? im->channels : i_img_color_channels(im);
+  dIMCTXim(im);
 
-  mm_log((1,"i_hardinvert(im %p)\n", im));
+  im_log((aIMCTX,1,"i_hardinvert)low(im %p, all %d)\n", im, all));
   
 #code im->bits <= 8  
   IM_COLOR *row, *entry;
@@ -159,7 +162,10 @@ i_hardinvertall(i_img *im) {
 /*
 =item i_noise(im, amount, type)
 
-Inverts the pixel values by the amount specified.
+Adjusts the sample values randomly by the amount specified.
+
+If type is 0, adjust all channels in a pixel by the same (random)
+amount amount, if non-zero adjust each sample independently.
 
   im     - image object
   amount - deviation in pixel values
@@ -175,14 +181,15 @@ Inverts the pixel values by the amount specified.
 
 void
 i_noise(i_img *im, float amount, unsigned char type) {
-  int x, y;
+  i_img_dim x, y;
   unsigned char ch;
   int new_color;
   float damount = amount * 2;
   i_color rcolor;
   int color_inc = 0;
+  dIMCTXim(im);
   
-  mm_log((1,"i_noise(im %p, intensity %.2f\n", im, amount));
+  im_log((aIMCTX, 1,"i_noise(im %p, intensity %.2f\n", im, amount));
   
   if(amount < 0) return;
   
@@ -216,63 +223,6 @@ i_noise(i_img *im, float amount, unsigned char type) {
   }
 }
 
-
-/* 
-=item i_noise(im, amount, type)
-
-Inverts the pixel values by the amount specified.
-
-  im     - image object
-  amount - deviation in pixel values
-  type   - noise individual for each channel if true
-
-=cut
-*/
-
-
-/*
-=item i_applyimage(im, add_im, mode)
-
-Apply's an image to another image
-
-  im     - target image
-  add_im - image that is applied to target
-  mode   - what method is used in applying:
-
-  0   Normal   
-  1   Multiply
-  2   Screen
-  3   Overlay
-  4   Soft Light
-  5   Hard Light
-  6   Color dodge
-  7   Color Burn
-  8   Darker
-  9   Lighter
-  10  Add
-  11  Subtract
-  12  Difference
-  13  Exclusion
-       
-=cut
-*/
-
-void i_applyimage(i_img *im, i_img *add_im, unsigned char mode) {
-  int x, y;
-  int mx, my;
-
-  mm_log((1, "i_applyimage(im %p, add_im %p, mode %d", im, add_im, mode));
-  
-  mx = (add_im->xsize <= im->xsize) ? add_im->xsize : add_im->xsize;
-  my = (add_im->ysize <= im->ysize) ? add_im->ysize : add_im->ysize;
-  
-  for(x = 0; x < mx; x++) {            
-    for(y = 0; y < my; y++) {
-    }
-  }
-}
-
-
 /* 
 =item i_bumpmap(im, bump, channel, light_x, light_y, st)
 
@@ -289,24 +239,25 @@ Makes a bumpmap on image im using the bump image as the elevation map.
 */
 
 void
-i_bumpmap(i_img *im, i_img *bump, int channel, int light_x, int light_y, int st) {
-  int x, y, ch;
-  int mx, my;
+i_bumpmap(i_img *im, i_img *bump, int channel, i_img_dim light_x, i_img_dim light_y, i_img_dim st) {
+  i_img_dim x, y;
+  int ch;
+  i_img_dim mx, my;
   i_color x1_color, y1_color, x2_color, y2_color, dst_color;
   double nX, nY;
   double tX, tY, tZ;
   double aX, aY, aL;
   double fZ;
   unsigned char px1, px2, py1, py2;
-  
+  dIMCTXim(im);
   i_img new_im;
 
-  mm_log((1, "i_bumpmap(im %p, add_im %p, channel %d, light_x %d, light_y %d, st %d)\n",
-         im, bump, channel, light_x, light_y, st));
+  im_log((aIMCTX, 1, "i_bumpmap(im %p, add_im %p, channel %d, light(" i_DFp "), st %" i_DF ")\n",
+         im, bump, channel, i_DFcp(light_x, light_y), i_DFc(st)));
 
 
   if(channel >= bump->channels) {
-    mm_log((1, "i_bumpmap: channel = %d while bump image only has %d channels\n", channel, bump->channels));
+    im_log((aIMCTX, 1, "i_bumpmap: channel = %d while bump image only has %d channels\n", channel, bump->channels));
     return;
   }
 
@@ -342,8 +293,8 @@ i_bumpmap(i_img *im, i_img *bump, int channel, int light_x, int light_y, int st)
 
       fZ = (sqrt((nX * nX) + (nY * nY)) / aL);
  
-      tX = abs(x - light_x) / aL;
-      tY = abs(y - light_y) / aL;
+      tX = i_abs(x - light_x) / aL;
+      tY = i_abs(y - light_y) / aL;
 
       tZ = 1 - (sqrt((tX * tX) + (tY * tY)) * fZ);
       
@@ -357,7 +308,7 @@ i_bumpmap(i_img *im, i_img *bump, int channel, int light_x, int light_y, int st)
     }
   }
 
-  i_copyto(im, &new_im, 0, 0, (int)im->xsize, (int)im->ysize, 0, 0);
+  i_copyto(im, &new_im, 0, 0, im->xsize, im->ysize, 0, 0);
   
   i_img_exorcise(&new_im);
 }
@@ -366,7 +317,7 @@ i_bumpmap(i_img *im, i_img *bump, int channel, int light_x, int light_y, int st)
 
 
 typedef struct {
-  float x,y,z;
+  double x,y,z;
 } fvec;
 
 
@@ -450,11 +401,11 @@ void
 i_bumpmap_complex(i_img *im,
                  i_img *bump,
                  int channel,
-                 int tx,
-                 int ty,
-                 float Lx,
-                 float Ly,
-                 float Lz,
+                 i_img_dim tx,
+                 i_img_dim ty,
+                 double Lx,
+                 double Ly,
+                 double Lz,
                  float cd,
                  float cs,
                  float n,
@@ -463,9 +414,9 @@ i_bumpmap_complex(i_img *im,
                  i_color *Is) {
   i_img new_im;
   
-  int inflight;
-  int x, y, ch;
-  int mx, Mx, my, My;
+  i_img_dim x, y;
+  int ch;
+  i_img_dim mx, Mx, my, My;
   
   float cdc[MAXCHANNELS];
   float csc[MAXCHANNELS];
@@ -479,11 +430,14 @@ i_bumpmap_complex(i_img *im,
   fvec R;         /* Reflection vector    */
   fvec V;         /* Vision vector        */
 
-  mm_log((1, "i_bumpmap_complex(im %p, bump %p, channel %d, tx %d, ty %d, Lx %.2f, Ly %.2f, Lz %.2f, cd %.2f, cs %.2f, n %.2f, Ia %p, Il %p, Is %p)\n",
-         im, bump, channel, tx, ty, Lx, Ly, Lz, cd, cs, n, Ia, Il, Is));
+  dIMCTXim(im);
+
+  im_log((aIMCTX, 1, "i_bumpmap_complex(im %p, bump %p, channel %d, t(" i_DFp 
+         "), Lx %.2f, Ly %.2f, Lz %.2f, cd %.2f, cs %.2f, n %.2f, Ia %p, Il %p, Is %p)\n",
+         im, bump, channel, i_DFcp(tx, ty), Lx, Ly, Lz, cd, cs, n, Ia, Il, Is));
   
   if (channel >= bump->channels) {
-    mm_log((1, "i_bumpmap_complex: channel = %d while bump image only has %d channels\n", channel, bump->channels));
+    im_log((aIMCTX, 1, "i_bumpmap_complex: channel = %d while bump image only has %d channels\n", channel, bump->channels));
     return;
   }
 
@@ -507,7 +461,6 @@ i_bumpmap_complex(i_img *im,
     L.z = -Lz;
     normalize(&L);
   } else {      /* Light is the position of the light source */
-    inflight = 0;
     L.x = -0.2;
     L.y = -0.4;
     L.z =  1;
@@ -566,7 +519,7 @@ i_bumpmap_complex(i_img *im,
     }
   }
   
-  i_copyto(im, &new_im, 0, 0, (int)im->xsize, (int)im->ysize, 0, 0);
+  i_copyto(im, &new_im, 0, 0, im->xsize, im->ysize, 0, 0);
   i_img_exorcise(&new_im);
 }
 
@@ -584,7 +537,8 @@ Quantizes Images to fewer levels.
 
 void
 i_postlevels(i_img *im, int levels) {
-  int x, y, ch;
+  i_img_dim x, y;
+  int ch;
   float pv;
   int rv;
   float av;
@@ -623,9 +577,10 @@ Makes an image looks like a mosaic with tilesize of size
 */
 
 void
-i_mosaic(i_img *im, int size) {
-  int x, y, ch;
-  int lx, ly, z;
+i_mosaic(i_img *im, i_img_dim size) {
+  i_img_dim x, y;
+  int ch, z;
+  i_img_dim lx, ly;
   long sqrsize;
 
   i_color rcolor;
@@ -673,12 +628,13 @@ Applies a watermark to the target image
 */
 
 void
-i_watermark(i_img *im, i_img *wmark, int tx, int ty, int pixdiff) {
-  int vx, vy, ch;
+i_watermark(i_img *im, i_img *wmark, i_img_dim tx, i_img_dim ty, int pixdiff) {
+  i_img_dim vx, vy;
+  int ch;
   i_color val, wval;
 
-       int mx = wmark->xsize;
-       int my = wmark->ysize;
+       i_img_dim mx = wmark->xsize;
+       i_img_dim my = wmark->ysize;
 
   for(vx=0;vx<mx;vx++) for(vy=0;vy<my;vy++) {
     
@@ -692,6 +648,102 @@ i_watermark(i_img *im, i_img *wmark, int tx, int ty, int pixdiff) {
   }
 }
 
+/*
+=item i_autolevels_mono(im, lsat, usat)
+
+Do autolevels, but monochromatically.
+
+=cut
+*/
+
+void
+i_autolevels_mono(i_img *im, float lsat, float usat) {
+  i_color val;
+  i_img_dim i, x, y, hist[256];
+  i_img_dim sum_lum, min_lum, max_lum;
+  i_img_dim upper_accum, lower_accum;
+  i_color *row;
+  dIMCTXim(im);
+  int adapt_channels = im->channels == 4 ? 2 : 1;
+  int color_channels = i_img_color_channels(im);
+  i_img_dim color_samples = im->xsize * color_channels;
+  
+
+  im_log((aIMCTX, 1,"i_autolevels_mono(im %p, lsat %f,usat %f)\n", im, lsat,usat));
+
+  /* build the histogram in 8-bits, unless the image has a very small
+     range it should make little difference to the result */
+  sum_lum = 0;
+  for (i = 0; i < 256; i++)
+    hist[i] = 0;
+
+  row = mymalloc(im->xsize * sizeof(i_color));
+  /* create histogram for each channel */
+  for (y = 0; y < im->ysize; y++) {
+    i_color *p = row;
+    i_glin(im, 0, im->xsize, y, row);
+    if (im->channels > 2)
+      i_adapt_colors(adapt_channels, im->channels, row, im->xsize);
+    for (x = 0; x < im->xsize; x++) {
+      hist[p->channel[0]]++;
+      ++p;
+    }
+  }
+  myfree(row);
+
+  for(i = 0; i < 256; i++) {
+    sum_lum += hist[i];
+  }
+  
+  min_lum = 0;
+  lower_accum = 0;
+  for (i = 0; i < 256; ++i) {
+    if (lower_accum < sum_lum * lsat)
+      min_lum = i;
+    lower_accum += hist[i];
+  }
+
+  max_lum = 255;
+  upper_accum = 0;
+  for(i = 255; i >= 0; i--) {
+    if (upper_accum < sum_lum * usat)
+      max_lum = i;
+    upper_accum  += hist[i];
+  }
+
+#code im->bits <= 8
+  IM_SAMPLE_T *srow = mymalloc(color_samples * sizeof(IM_SAMPLE_T));
+#ifdef IM_EIGHT_BIT
+  IM_WORK_T low = min_lum;
+  i_sample_t lookup[256];
+#else
+  IM_WORK_T low = min_lum / 255.0 * IM_SAMPLE_MAX;
+#endif
+  double scale = 255.0 / (max_lum - min_lum);
+
+#ifdef IM_EIGHT_BIT
+  for (i = 0; i < 256; ++i) {
+    IM_WORK_T tmp = (i - low) * scale;
+    lookup[i] = IM_LIMIT(tmp);
+  }
+#endif
+
+  for(y = 0; y < im->ysize; y++) {
+    IM_GSAMP(im, 0, im->xsize, y, srow, NULL, color_channels);
+    for(i = 0; i < color_samples; ++i) {
+#ifdef IM_EIGHT_BIT
+      srow[i] = lookup[srow[i]];
+#else
+      IM_WORK_T tmp = (srow[i] - low) * scale;
+      srow[i] = IM_LIMIT(tmp);
+#endif
+    }
+    IM_PSAMP(im, 0, im->xsize, y, srow, NULL, color_channels);
+  }
+  myfree(srow);
+#/code
+}
+
 
 /*
 =item i_autolevels(im, lsat, usat, skew)
@@ -705,19 +757,23 @@ occur when changing the contrast.
   usat - fraction of pixels that will be truncated at the higher end of the spectrum
   skew - not used yet
 
+Note: this code calculates levels and adjusts each channel separately,
+which will typically cause a color shift.
+
 =cut
 */
 
 void
 i_autolevels(i_img *im, float lsat, float usat, float skew) {
   i_color val;
-  int i, x, y, rhist[256], ghist[256], bhist[256];
-  int rsum, rmin, rmax;
-  int gsum, gmin, gmax;
-  int bsum, bmin, bmax;
-  int rcl, rcu, gcl, gcu, bcl, bcu;
+  i_img_dim i, x, y, rhist[256], ghist[256], bhist[256];
+  i_img_dim rsum, rmin, rmax;
+  i_img_dim gsum, gmin, gmax;
+  i_img_dim bsum, bmin, bmax;
+  i_img_dim rcl, rcu, gcl, gcu, bcl, bcu;
+  dIMCTXim(im);
 
-  mm_log((1,"i_autolevels(im %p, lsat %f,usat %f,skew %f)\n", im, lsat,usat,skew));
+  im_log((aIMCTX, 1,"i_autolevels(im %p, lsat %f,usat %f,skew %f)\n", im, lsat,usat,skew));
 
   rsum=gsum=bsum=0;
   for(i=0;i<256;i++) rhist[i]=ghist[i]=bhist[i] = 0;
@@ -772,9 +828,9 @@ Pseudo noise utility function used to generate perlin noise. (internal)
 */
 
 static
-float
-Noise(int x, int y) {
-  int n = x + y * 57; 
+double
+Noise(i_img_dim x, i_img_dim y) {
+  i_img_dim n = x + y * 57; 
   n = (n<<13) ^ n;
   return ( 1.0 - ( (n * (n * n * 15731 + 789221) + 1376312589) & 0x7fffffff) / 1073741824.0);
 }
@@ -791,11 +847,11 @@ Pseudo noise utility function used to generate perlin noise. (internal)
 */
 
 static
-float
-SmoothedNoise1(float x, float y) {
-  float corners = ( Noise(x-1, y-1)+Noise(x+1, y-1)+Noise(x-1, y+1)+Noise(x+1, y+1) ) / 16;
-  float sides   = ( Noise(x-1, y)  +Noise(x+1, y)  +Noise(x, y-1)  +Noise(x, y+1) ) /  8;
-  float center  =  Noise(x, y) / 4;
+double
+SmoothedNoise1(double x, double y) {
+  double corners = ( Noise(x-1, y-1)+Noise(x+1, y-1)+Noise(x-1, y+1)+Noise(x+1, y+1) ) / 16;
+  double sides   = ( Noise(x-1, y)  +Noise(x+1, y)  +Noise(x, y-1)  +Noise(x, y+1) ) /  8;
+  double center  =  Noise(x, y) / 4;
   return corners + sides + center;
 }
 
@@ -809,10 +865,11 @@ Utility function used to generate perlin noise. (internal)
 */
 
 static
-float C_Interpolate(float a, float b, float x) {
+double
+C_Interpolate(double a, double b, double x) {
   /*  float ft = x * 3.1415927; */
-  float ft = x * PI;
-  float f = (1 - cos(ft)) * .5;
+  double ft = x * PI;
+  double f = (1 - cos(ft)) * .5;
   return  a*(1-f) + b*f;
 }
 
@@ -826,21 +883,21 @@ Utility function used to generate perlin noise. (internal)
 */
 
 static
-float
-InterpolatedNoise(float x, float y) {
+double
+InterpolatedNoise(double x, double y) {
 
-  int integer_X = x;
-  float fractional_X = x - integer_X;
-  int integer_Y = y;
-  float fractional_Y = y - integer_Y;
+  i_img_dim integer_X = x;
+  double fractional_X = x - integer_X;
+  i_img_dim integer_Y = y;
+  double fractional_Y = y - integer_Y;
 
-  float v1 = SmoothedNoise1(integer_X,     integer_Y);
-  float v2 = SmoothedNoise1(integer_X + 1, integer_Y);
-  float v3 = SmoothedNoise1(integer_X,     integer_Y + 1);
-  float v4 = SmoothedNoise1(integer_X + 1, integer_Y + 1);
+  double v1 = SmoothedNoise1(integer_X,     integer_Y);
+  double v2 = SmoothedNoise1(integer_X + 1, integer_Y);
+  double v3 = SmoothedNoise1(integer_X,     integer_Y + 1);
+  double v4 = SmoothedNoise1(integer_X + 1, integer_Y + 1);
 
-  float i1 = C_Interpolate(v1 , v2 , fractional_X);
-  float i2 = C_Interpolate(v3 , v4 , fractional_X);
+  double i1 = C_Interpolate(v1 , v2 , fractional_X);
+  double i2 = C_Interpolate(v3 , v4 , fractional_X);
 
   return C_Interpolate(i1 , i2 , fractional_Y);
 }
@@ -859,8 +916,8 @@ static
 float
 PerlinNoise_2D(float x, float y) {
   int i,frequency;
-  float amplitude;
-  float total = 0;
+  double amplitude;
+  double total = 0;
   int Number_Of_Octaves=6;
   int n = Number_Of_Octaves - 1;
 
@@ -889,16 +946,17 @@ Perlin-like radial noise.
 */
 
 void
-i_radnoise(i_img *im, int xo, int yo, float rscale, float ascale) {
-  int x, y, ch;
+i_radnoise(i_img *im, i_img_dim xo, i_img_dim yo, double rscale, double ascale) {
+  i_img_dim x, y;
+  int ch;
   i_color val;
   unsigned char v;
-  float xc, yc, r;
+  double xc, yc, r;
   double a;
   
   for(y = 0; y < im->ysize; y++) for(x = 0; x < im->xsize; x++) {
-    xc = (float)x-xo+0.5;
-    yc = (float)y-yo+0.5;
+    xc = (double)x-xo+0.5;
+    yc = (double)y-yo+0.5;
     r = rscale*sqrt(xc*xc+yc*yc)+1.2;
     a = (PI+atan2(yc,xc))*ascale;
     v = saturate(128+100*(PerlinNoise_2D(a,r)));
@@ -923,14 +981,15 @@ Perlin-like 2d noise noise.
 */
 
 void
-i_turbnoise(i_img *im, float xo, float yo, float scale) {
-  int x,y,ch;
+i_turbnoise(i_img *im, double xo, double yo, double scale) {
+  i_img_dim x,y;
+  int ch;
   unsigned char v;
   i_color val;
 
   for(y = 0; y < im->ysize; y++) for(x = 0; x < im->xsize; x++) {
     /*    v=saturate(125*(1.0+PerlinNoise_2D(xo+(float)x/scale,yo+(float)y/scale))); */
-    v = saturate(120*(1.0+sin(xo+(float)x/scale+PerlinNoise_2D(xo+(float)x/scale,yo+(float)y/scale))));
+    v = saturate(120*(1.0+sin(xo+(double)x/scale+PerlinNoise_2D(xo+(double)x/scale,yo+(float)y/scale))));
     for(ch=0; ch<im->channels; ch++) val.channel[ch] = v;
     i_ppix(im, x, y, &val);
   }
@@ -958,21 +1017,23 @@ Gradient generating function.
 
 
 void
-i_gradgen(i_img *im, int num, int *xo, int *yo, i_color *ival, int dmeasure) {
+i_gradgen(i_img *im, int num, i_img_dim *xo, i_img_dim *yo, i_color *ival, int dmeasure) {
   
   i_color val;
-  int p, x, y, ch;
+  int p, ch;
+  i_img_dim x, y;
   int channels = im->channels;
-  int xsize    = im->xsize;
-  int ysize    = im->ysize;
-  int bytes;
+  i_img_dim xsize    = im->xsize;
+  i_img_dim ysize    = im->ysize;
+  size_t bytes;
 
-  float *fdist;
+  double *fdist;
+  dIMCTXim(im);
 
-  mm_log((1,"i_gradgen(im %p, num %d, xo %p, yo %p, ival %p, dmeasure %d)\n", im, num, xo, yo, ival, dmeasure));
+  im_log((aIMCTX, 1,"i_gradgen(im %p, num %d, xo %p, yo %p, ival %p, dmeasure %d)\n", im, num, xo, yo, ival, dmeasure));
   
   for(p = 0; p<num; p++) {
-    mm_log((1,"i_gradgen: (%d, %d)\n", xo[p], yo[p]));
+    im_log((aIMCTX,1,"i_gradgen: p%d(" i_DFp ")\n", p, i_DFcp(xo[p], yo[p])));
     ICL_info(&ival[p]);
   }
 
@@ -984,19 +1045,19 @@ i_gradgen(i_img *im, int num, int *xo, int *yo, i_color *ival, int dmeasure) {
      overflow is a programming error rather than an end-user error, so
      calling exit() is justified.
   */
-  bytes = sizeof(float) * num;
-  if (bytes / num != sizeof(float)) {
+  bytes = sizeof(double) * num;
+  if (bytes / num != sizeof(double)) {
     fprintf(stderr, "integer overflow calculating memory allocation");
     exit(1);
   }
   fdist = mymalloc( bytes ); /* checked 14jul05 tonyc */
   
   for(y = 0; y<ysize; y++) for(x = 0; x<xsize; x++) {
-    float cs = 0;
-    float csd = 0;
+    double cs = 0;
+    double csd = 0;
     for(p = 0; p<num; p++) {
-      int xd    = x-xo[p];
-      int yd    = y-yo[p];
+      i_img_dim xd    = x-xo[p];
+      i_img_dim yd    = y-yo[p];
       switch (dmeasure) {
       case 0: /* euclidean */
        fdist[p]  = sqrt(xd*xd + yd*yd); /* euclidean distance */
@@ -1008,7 +1069,7 @@ i_gradgen(i_img *im, int num, int *xo, int *yo, i_color *ival, int dmeasure) {
        fdist[p]  = i_max(xd*xd, yd*yd); /* manhattan distance */
        break;
       default:
-       i_fatal(3,"i_gradgen: Unknown distance measure\n");
+       im_fatal(aIMCTX, 3,"i_gradgen: Unknown distance measure\n");
       }
       cs += fdist[p];
     }
@@ -1029,26 +1090,28 @@ i_gradgen(i_img *im, int num, int *xo, int *yo, i_color *ival, int dmeasure) {
 }
 
 void
-i_nearest_color_foo(i_img *im, int num, int *xo, int *yo, i_color *ival, int dmeasure) {
+i_nearest_color_foo(i_img *im, int num, i_img_dim *xo, i_img_dim *yo, i_color *ival, int dmeasure) {
 
-  int p, x, y;
-  int xsize    = im->xsize;
-  int ysize    = im->ysize;
+  int p;
+  i_img_dim x, y;
+  i_img_dim xsize    = im->xsize;
+  i_img_dim ysize    = im->ysize;
+  dIMCTXim(im);
 
-  mm_log((1,"i_gradgen(im %p, num %d, xo %p, yo %p, ival %p, dmeasure %d)\n", im, num, xo, yo, ival, dmeasure));
+  im_log((aIMCTX,1,"i_gradgen(im %p, num %d, xo %p, yo %p, ival %p, dmeasure %d)\n", im, num, xo, yo, ival, dmeasure));
   
   for(p = 0; p<num; p++) {
-    mm_log((1,"i_gradgen: (%d, %d)\n", xo[p], yo[p]));
+    im_log((aIMCTX, 1,"i_gradgen: p%d(" i_DFp ")\n", p, i_DFcp(xo[p], yo[p])));
     ICL_info(&ival[p]);
   }
 
   for(y = 0; y<ysize; y++) for(x = 0; x<xsize; x++) {
-    int   midx    = 0;
-    float mindist = 0;
-    float curdist = 0;
+    int    midx    = 0;
+    double mindist = 0;
+    double curdist = 0;
 
-    int xd        = x-xo[0];
-    int yd        = y-yo[0];
+    i_img_dim xd        = x-xo[0];
+    i_img_dim yd        = y-yo[0];
 
     switch (dmeasure) {
     case 0: /* euclidean */
@@ -1061,7 +1124,7 @@ i_nearest_color_foo(i_img *im, int num, int *xo, int *yo, i_color *ival, int dme
       mindist = i_max(xd*xd, yd*yd); /* manhattan distance */
       break;
     default:
-      i_fatal(3,"i_nearest_color: Unknown distance measure\n");
+      im_fatal(aIMCTX, 3,"i_nearest_color: Unknown distance measure\n");
     }
 
     for(p = 1; p<num; p++) {
@@ -1078,7 +1141,7 @@ i_nearest_color_foo(i_img *im, int num, int *xo, int *yo, i_color *ival, int dme
        curdist = i_max(xd*xd, yd*yd); /* manhattan distance */
        break;
       default:
-       i_fatal(3,"i_nearest_color: Unknown distance measure\n");
+       im_fatal(aIMCTX, 3,"i_nearest_color: Unknown distance measure\n");
       }
       if (curdist < mindist) {
        mindist = curdist;
@@ -1112,11 +1175,11 @@ int num - number of points/colors in xo, yo, oval
 
 =item *
 
-int *xo - array of I<num> x positions
+i_img_dim *xo - array of I<num> x positions
 
 =item *
 
-int *yo - array of I<num> y positions
+i_img_dim *yo - array of I<num> y positions
 
 =item *
 
@@ -1157,18 +1220,20 @@ An invalid value causes an error exit (the program is aborted).
  */
 
 int
-i_nearest_color(i_img *im, int num, int *xo, int *yo, i_color *oval, int dmeasure) {
+i_nearest_color(i_img *im, int num, i_img_dim *xo, i_img_dim *yo, i_color *oval, int dmeasure) {
   i_color *ival;
   float *tval;
-  float c1, c2;
+  double c1, c2;
   i_color val;
-  int p, x, y, ch;
-  int xsize    = im->xsize;
-  int ysize    = im->ysize;
+  int p, ch;
+  i_img_dim x, y;
+  i_img_dim xsize    = im->xsize;
+  i_img_dim ysize    = im->ysize;
   int *cmatch;
-  int ival_bytes, tval_bytes;
+  size_t ival_bytes, tval_bytes;
+  dIMCTXim(im);
 
-  mm_log((1,"i_nearest_color(im %p, num %d, xo %p, yo %p, oval %p, dmeasure %d)\n", im, num, xo, yo, oval, dmeasure));
+  im_log((aIMCTX, 1,"i_nearest_color(im %p, num %d, xo %p, yo %p, oval %p, dmeasure %d)\n", im, num, xo, yo, oval, dmeasure));
 
   i_clear_error();
 
@@ -1204,11 +1269,11 @@ i_nearest_color(i_img *im, int num, int *xo, int *yo, i_color *oval, int dmeasur
   
   for(y = 0; y<ysize; y++) for(x = 0; x<xsize; x++) {
     int   midx    = 0;
-    float mindist = 0;
-    float curdist = 0;
+    double mindist = 0;
+    double curdist = 0;
     
-    int xd        = x-xo[0];
-    int yd        = y-yo[0];
+    i_img_dim xd        = x-xo[0];
+    i_img_dim yd        = y-yo[0];
 
     switch (dmeasure) {
     case 0: /* euclidean */
@@ -1221,7 +1286,7 @@ i_nearest_color(i_img *im, int num, int *xo, int *yo, i_color *oval, int dmeasur
       mindist = i_max(xd*xd, yd*yd); /* manhattan distance */
       break;
     default:
-      i_fatal(3,"i_nearest_color: Unknown distance measure\n");
+      im_fatal(aIMCTX, 3,"i_nearest_color: Unknown distance measure\n");
     }
     
     for(p = 1; p<num; p++) {
@@ -1238,7 +1303,7 @@ i_nearest_color(i_img *im, int num, int *xo, int *yo, i_color *oval, int dmeasur
        curdist = i_max(xd*xd, yd*yd); /* manhattan distance */
        break;
       default:
-       i_fatal(3,"i_nearest_color: Unknown distance measure\n");
+       im_fatal(aIMCTX, 3,"i_nearest_color: Unknown distance measure\n");
       }
       if (curdist < mindist) {
        mindist = curdist;
@@ -1257,11 +1322,21 @@ i_nearest_color(i_img *im, int num, int *xo, int *yo, i_color *oval, int dmeasur
   
   }
 
-  for(p = 0; p<num; p++) for(ch = 0; ch<im->channels; ch++)
-    ival[p].channel[ch] = tval[p*im->channels + ch];
+  for(p = 0; p<num; p++) {
+    for(ch = 0; ch<im->channels; ch++)
+      ival[p].channel[ch] = tval[p*im->channels + ch];
+
+    /* avoid uninitialized value messages from valgrind */
+    while (ch < MAXCHANNELS)
+      ival[p].channel[ch++] = 0;
+  }
 
   i_nearest_color_foo(im, num, xo, yo, ival, dmeasure);
 
+  myfree(cmatch);
+  myfree(ival);
+  myfree(tval);
+
   return 1;
 }
 
@@ -1277,7 +1352,8 @@ image from double the original.
 void
 i_unsharp_mask(i_img *im, double stddev, double scale) {
   i_img *copy;
-  int x, y, ch;
+  i_img_dim x, y;
+  int ch;
 
   if (scale < 0)
     return;
@@ -1354,7 +1430,8 @@ i_img *
 i_diff_image(i_img *im1, i_img *im2, double mindist) {
   i_img *out;
   int outchans, diffchans;
-  int xsize, ysize;
+  i_img_dim xsize, ysize;
+  dIMCTXim(im1);
 
   i_clear_error();
   if (im1->channels != im2->channels) {
@@ -1375,7 +1452,8 @@ i_diff_image(i_img *im1, i_img *im2, double mindist) {
     i_color *line1 = mymalloc(xsize * sizeof(*line1)); /* checked 17feb2005 tonyc */
     i_color *line2 = mymalloc(xsize * sizeof(*line1)); /* checked 17feb2005 tonyc */
     i_color empty;
-    int x, y, ch;
+    i_img_dim x, y;
+    int ch;
     int imindist = (int)mindist;
 
     for (ch = 0; ch < MAXCHANNELS; ++ch)
@@ -1410,7 +1488,8 @@ i_diff_image(i_img *im1, i_img *im2, double mindist) {
     i_fcolor *line1 = mymalloc(xsize * sizeof(*line1)); /* checked 17feb2005 tonyc */
     i_fcolor *line2 = mymalloc(xsize * sizeof(*line2)); /* checked 17feb2005 tonyc */
     i_fcolor empty;
-    int x, y, ch;
+    i_img_dim x, y;
+    int ch;
     double dist = mindist / 255.0;
 
     for (ch = 0; ch < MAXCHANNELS; ++ch)
@@ -1536,7 +1615,7 @@ struct fount_state {
   double cos;
   double sin;
   double theta;
-  int xa, ya;
+  i_img_dim xa, ya;
   void *ssample_data;
   fount_func ffunc;
   fount_repeat rpfunc;
@@ -1668,13 +1747,13 @@ i_fountain(i_img *im, double xa, double ya, double xb, double yb,
            int combine, int super_sample, double ssample_param, 
            int count, i_fountain_seg *segs) {
   struct fount_state state;
-  int x, y;
+  i_img_dim x, y;
   i_fcolor *line = NULL;
   i_fcolor *work = NULL;
-  int line_bytes;
-  i_fountain_seg *my_segs;
+  size_t line_bytes;
   i_fill_combine_f combine_func = NULL;
   i_fill_combinef_f combinef_func = NULL;
+  dIMCTXim(im);
 
   i_clear_error();
 
@@ -1694,7 +1773,6 @@ i_fountain(i_img *im, double xa, double ya, double xb, double yb,
 
   fount_init_state(&state, xa, ya, xb, yb, type, repeat, combine, 
                    super_sample, ssample_param, count, segs);
-  my_segs = state.segs;
 
   for (y = 0; y < im->ysize; ++y) {
     i_glinf(im, 0, im->xsize, y, line);
@@ -1798,7 +1876,7 @@ fount_init_state(struct fount_state *state, double xa, double ya,
                  i_fountain_repeat repeat, int combine, int super_sample, 
                  double ssample_param, int count, i_fountain_seg *segs) {
   int i, j;
-  int bytes;
+  size_t bytes;
   i_fountain_seg *my_segs = mymalloc(sizeof(i_fountain_seg) * count); /* checked 2jul06 - duplicating original */
   /*int have_alpha = im->channels == 2 || im->channels == 4;*/
   
@@ -2005,7 +2083,7 @@ square.
 */
 static double
 square_fount_f(double x, double y, struct fount_state *state) {
-  int xc, yc; /* centred on A */
+  i_img_dim xc, yc; /* centred on A */
   double xt, yt; /* rotated by theta */
   xc = x - state->xa;
   yc = y - state->ya;
@@ -2187,7 +2265,7 @@ Simple grid-based super-sampling.
 static int
 simple_ssample(i_fcolor *out, double x, double y, struct fount_state *state) {
   i_fcolor *work = state->ssample_data;
-  int dx, dy;
+  i_img_dim dx, dy;
   int grid = state->parm;
   double base = -0.5 + 0.5 / grid;
   double step = 1.0 / grid;
@@ -2379,8 +2457,9 @@ fill_fountf(i_fill_t *fill, i_img_dim x, i_img_dim y, i_img_dim width,
       got_one = f->state.ssfunc(&c, x, y, &f->state);
     else
       got_one = fount_getat(&c, x, y, &f->state);
-    
-    *data++ = c;
+
+    if (got_one)
+      *data++ = c;
     
     ++x;
   }