/* ichika0.c - Don Yang (uguu.org)

   11/17/04
                                                                           */
/*@ -realcompare @*/

#include<assert.h>
#include<stdio.h>
#include<stdlib.h>


#define SHIFT_HUE       2.0f
#define SHIFT_SAT       0.70f
#define SHIFT_VAL       0.09f

#define EYE_HUE_MIN     3.60f
#define EYE_HUE_MAX     4.24f
#define EYE_SAT1        0.65f
#define EYE_VAL1        0.42f
#define EYE_SAT2        0.36f
#define EYE_VAL2        0.60f
#define EYE_THRESHOLD   0.6f
#define EYE_SCALE_VAL   1.29f

#define FACE_HUE1_MIN   0.30f
#define FACE_HUE1_MAX   0.45f
#define FACE_HUE2_MIN   4.28f
#define FACE_HUE2_MAX   6.00f
#define HAIR_HUE_MIN    0.0f
#define HAIR_HUE_MAX    0.1f

typedef struct
{
   /*@only@*/unsigned char *image;
   int width, height, size, align, headersize;
   char header[64];
} Image;

#if 0
   typedef double flt;
#else
   typedef float flt;
#endif

static void rgb2hsv(flt r, flt g, flt b,
                    /*@out@*/flt *h, /*@out@*/flt *s, /*@out@*/flt *v);
static void hsv2rgb(flt h, flt s, flt v,
                    /*@out@*/flt *r, /*@out@*/flt *g, /*@out@*/flt *b);

static int CheckNeighbors(Image *img, int offset);
static int CheckNeighbors0(Image *img, int offset, int delta);
static void Erase(Image *img, unsigned char *mask, int offset);
static /*@only@*//*@null@*/Image *LoadImage(char *name);
static void ProcessImage(Image *img);
static void ReorderBytes(Image *img);
static void SaveImage(Image *img, char *name);

int main(int argc, char **argv)
{
   Image *img;

   if( argc < 3 )
      return printf("%s <infile.bmp|infile.ppm> <outfile>\n", *argv);

   if( (img = LoadImage(argv[1])) == NULL )
      return -1;
   if( img->header[0] == 'B' )
      ReorderBytes(img);
   ProcessImage(img);
   if( img->header[0] == 'B' )
      ReorderBytes(img);
   SaveImage(img, argv[2]);
   free(img->image);
   free(img);
   return 0;
}


static void rgb2hsv(flt r, flt g, flt b,
                    /*@out@*/flt *h, /*@out@*/flt *s, /*@out@*/flt *v)
{
   flt p, q;

   assert(r >= 0.0f && r <= 1.0f);
   assert(g >= 0.0f && g <= 1.0f);
   assert(b >= 0.0f && b <= 1.0f);
   p = q = r;
   if( p > g ) p = g; if( q < g ) q = g;
   if( p > b ) p = b; if( q < b ) q = b;
   if( p == q )
   {
      *h = *s = *v = 0.0f;
   }
   else
   {
      assert(q > 0.0f);
      p = q - p;
      if( r == q )      { *h = (g - b) / p; }
      else if( g == q ) { *h = (b - r) / p + 2.0f; }
      else /* b == q */ { *h = (r - g) / p + 4.0f; }
      *s = p / (*v = q);
   }
}

static void hsv2rgb(flt h, flt s, flt v,
                    /*@out@*/flt *r, /*@out@*/flt *g, /*@out@*/flt *b)
{
   flt f, p, q, t;
   int i;

   assert(h >= 0.0f && h < 6.0f);
   assert(s >= 0.0f && s <= 1.0f);
   assert(v >= 0.0f && v <= 1.0f);
   i = (int)h;
   f = h - (flt)i;
   p = v * (1.0f - s);
   q = v * (1.0f - s * f);
   t = v * (1.0f - s * (1.0f - f));
   switch( i )
   {
      case 0:  *r = v; *g = t; *b = p; break;
      case 1:  *r = q; *g = v; *b = p; break;
      case 2:  *r = p; *g = v; *b = t; break;
      case 3:  *r = p; *g = q; *b = v; break;
      case 4:  *r = t; *g = p; *b = v; break;
      default: *r = v; *g = p; *b = q; break;
   }
}


static int CheckNeighbors(Image *img, int offset)
{
   return
      ((CheckNeighbors0(img, offset, 1) +
        CheckNeighbors0(img, offset, -1) +
        CheckNeighbors0(img, offset, img->width) +
        CheckNeighbors0(img, offset, -img->width)) == 4) ? 0 : 1;
}
static int CheckNeighbors0(Image *img, int offset, int delta)
{
   flt h, s, v;

   if( offset < img->width || offset > img->size - img->width ||
       (offset % img->width) == 0 || ((offset + 1) % img->width) == 0 )
      return 0;
   rgb2hsv((flt)(img->image[offset * 3]) / 255.0f,
           (flt)(img->image[offset * 3 + 1]) / 255.0f,
           (flt)(img->image[offset * 3 + 2]) / 255.0f,
           &h, &s, &v);
   if( (FACE_HUE1_MIN < h && h < FACE_HUE1_MAX) ||
       (FACE_HUE2_MIN < h && h < FACE_HUE2_MAX) ||
       (HAIR_HUE_MIN < h && h < HAIR_HUE_MAX) )
      return 1;
   return CheckNeighbors0(img, offset + delta, delta);
}

static void Erase(Image *img, unsigned char *mask, int offset)
{
   int x, y;

   if( mask[offset] == (unsigned char)0xff )
      return;
   mask[offset] = (unsigned char)0xff;
   x = offset % img->width;
   y = offset / img->width;
   if( x < img->width - 1 )   Erase(img, mask, offset + 1);
   if( y < img->height - 1 )  Erase(img, mask, offset + img->width);
   if( x > 0 )                Erase(img, mask, offset - 1);
   if( y > 0 )                Erase(img, mask, offset - img->width);
}

static /*@only@*//*@null@*/Image *LoadImage(char *name)
{
   Image *img;
   FILE *infile;
   unsigned char *u;
   int i, j;

   if( (infile = fopen(name, "rb")) == NULL )
   {
      printf("can not open %s\n", name);
      return NULL;
   }
   do
   {
      if( (img = (Image*)malloc(sizeof(Image))) == NULL )
         break;
      img->image = NULL;

      /* Parse header */
      if( fread(img->header, 64, 1, infile) != 1 )
      {
         (void)puts("read error");
         break;
      }
      if( img->header[0] == 'B' && img->header[1] == 'M' )
      {
         if( *((int*)&(img->header[28])) != 24 )
         {
            (void)puts("invalid BMP");
            break;
         }
         img->width = *((int*)&(img->header[18]));    /* endian-specific */
         img->height = *((int*)&(img->header[22]));
         img->align = (4 - ((img->width * 3) & 3)) & 3;
         img->headersize = 54;
      }
      else if( img->header[0] == 'P' && img->header[1] == '6' )
      {
         for(i = j = img->headersize = 0; i < 63; i++)
         {
            if( img->header[i] == '\n' )
            {
               if( ++j == 3 )
               {
                  img->header[img->headersize = ++i] = '\0';
                  break;
               }
            }
         }
         if( img->headersize == 0 )
         {
            (void)puts("invalid PPM");
            break;
         }
         if( sscanf(img->header, "P6\n%d %d\n",
               &(img->width), &(img->height)) != 2 )
         {
            (void)puts("invalid PPM");
            break;
         }
         img->align = 0;
      }
      else
      {
         (void)puts("unrecognized input format");
         break;
      }
      if( img->width < 16 || img->height < 16 )
      {
         (void)puts("invalid image size");
         break;
      }

      /* Load image
         Note that BMPs will be loaded upside down, but we don't care. */
      if( (img->image = (unsigned char*)malloc(
            (size_t)(img->size = (img->width * img->height)) * 3)) == NULL )
      {
         (void)puts("out of memory");
         break;
      }
      u = img->image;
      (void)fseek(infile, (long)img->headersize, SEEK_SET);
      for(i = 0; i < img->height; i++)
      {
         if( fread(u, (size_t)img->width * 3, 1, infile) < 1 )
         {
            (void)puts("read error");
            break;
         }
         for(j = 0; j < img->align; j++)
            (void)fgetc(infile);
         u += img->width * 3;
      }
      if( i < img->height )
         break;

      (void)fclose(infile);
      return /*@i1@*/img;

   } while(0);
   (void)fclose(infile);
   if( /*@i1@*/img != NULL )
   {
      if( img->image != NULL )
         free(img->image);
      free(img);
   }
   return NULL;
}


static void ProcessImage(Image *img)
{
   flt h, s, v, r, g, b, r0, g0, b0, d1, d2;
   unsigned char *u, *m, *mask;
   int i;

   /* Create mask */
   if( (mask = (unsigned char*)malloc((size_t)img->size)) == NULL )
      return;
   u = img->image;
   m = mask;
   for(i = 0; i < img->size; i++)
   {
      rgb2hsv((flt)(u[0]) / 255.0f,
              (flt)(u[1]) / 255.0f,
              (flt)(u[2]) / 255.0f,
              &h, &s, &v);
      *m = (unsigned char)0xff;
      if( EYE_HUE_MIN < h && h < EYE_HUE_MAX )
      {
         d1 = (s-EYE_SAT1)*(s-EYE_SAT1)
            + (v-EYE_VAL1)*(v-EYE_VAL1) * EYE_SCALE_VAL;
         d2 = (s-EYE_SAT2)*(s-EYE_SAT2)
            + (v-EYE_VAL2)*(v-EYE_VAL2) * EYE_SCALE_VAL;
         if( d1 < EYE_THRESHOLD || d2 < EYE_THRESHOLD )
         {
            d1 = (d1 < EYE_THRESHOLD) ? (d1 / EYE_THRESHOLD)
                                      : (d2 / EYE_THRESHOLD);
            *m = (unsigned char)(255.0f * d1);
         }
      }
      u += 3;
      m++;
   }

   /* Remove mask areas at edges */
   for(i = 0; i < img->width; i++)
   {
      Erase(img, mask, i);
      Erase(img, mask, i + img->width);
      Erase(img, mask, i + img->size - img->width * 2);
      Erase(img, mask, i + img->size - img->width);
   }
   for(i = 0; i < img->height; i++)
   {
      Erase(img, mask, i * img->width);
      Erase(img, mask, i * img->width + 1);
      Erase(img, mask, (i + 1) * img->width - 2);
      Erase(img, mask, (i + 1) * img->width - 1);
   }

   /* Remove areas not surrounded by skin/hair */
   m = mask;
   for(i = 0; i < img->size; i++)
   {
      if( *m != (unsigned char)0xff )
      {
         if( CheckNeighbors(img, i) != 0 )
            Erase(img, mask, i);
      }
      m++;
   }

   /* Shift hue */
   u = img->image;
   m = mask;
   for(i = 0; i < img->size; i++)
   {
      if( *m != (unsigned char)0xff )
      {
         assert(i > img->width * 2 && i < img->size - img->width * 2);
         d1 = (flt)*m / 255.0f;
         d2 = 1.0f - d1;

         rgb2hsv(r0 = (flt)(u[0]) / 255.0f,
                 g0 = (flt)(u[1]) / 255.0f,
                 b0 = (flt)(u[2]) / 255.0f,
                 &h, &s, &v);
         if( (h += SHIFT_HUE) >= 6.0f )   h -= 6.0f;
         if( (s += SHIFT_SAT) > 1.0f )    s = 1.0f;
         if( (v += SHIFT_VAL) > 1.0f )    v = 1.0f;
         hsv2rgb(h, s, v, &r, &g, &b);
         r = r0 * d1 + r * d2;
         g = r0 * d1 + g * d2;
         b = r0 * d1 + b * d2;
         u[0] = (unsigned char)(r * 255.0f);
         u[1] = (unsigned char)(g * 255.0f);
         u[2] = (unsigned char)(b * 255.0f);
      }
      u += 3;
      m++;
   }

   free(mask);
}

static void ReorderBytes(Image *img)
{
   unsigned char *u, c;
   int i;

   /* Convert BMP's BGR to RGB */
   u = img->image - 1;
   for(i = 0; i < img->size; i++)
   {
      c = *++u;
      *u = *(u + 2);
      *(u += 2) = c;
   }
}

static void SaveImage(Image *img, char *name)
{
   FILE *outfile;
   int i, j;
   unsigned char *u;

   if( (outfile = fopen(name, "wb+")) == NULL )
   {
      printf("can not create %s\n", name);
      return;
   }
   (void)fwrite(img->header, (size_t)img->headersize, 1, outfile);
   u = img->image;
   for(i = 0; i < img->height; i++)
   {
      (void)fwrite(u, (size_t)img->width * 3, 1, outfile);
      for(j = 0; j < img->align; j++)
         (void)fputc(0, outfile);
      u += img->width * 3;
   }
   (void)fclose(outfile);
}