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

   01/02/12
                                                                           */
#include<stdio.h>
#include<stdlib.h>

FILE *file[3];
int width[2], height[2], offset[2];
int i, c;
int x, y, x0, y0;
int step, ip0, ip1, iq0, iq1;
unsigned char *buffer[2];
unsigned char pixel[12];
char header[256];
char pixels[256][256];  /* [primary][secondary] */
int best[4], z[4], scale[4];
int best_area, z_area;
float p, q, r, s, fpixel[12];

/* Pixel schedules */
int index[2][3][12] =
{
   /* Rc Gp Bc Rc Gc Bp
      Rp Gc Bc Rq Gq Bq */
   {
      {6, 9, 0, 3},
      {1, 10, 4, 7},
      {5, 11, 2, 8}
   },
   /* Rq Gq Bq Rp Gc Bc
      Rc Gc Bp Rc Gp Bc */
   {
      {3, 0, 6, 9},
      {10, 1, 4, 7},
      {8, 2, 5, 11}
   }
};

int main(int argc, char **argv)
{
   /* Check command line arguments */
   if( argc < 4 )
   {
      printf("%s <primary.ppm> <secondary.ppm> <output.ppm> [odd]\n", *argv);
   }
   else
   {
      /* Open files */
      for(x = 0; x < 2; x++)
      {
         if( (file[x] = fopen(argv[x + 1], "rb")) == NULL )
            goto error;
         fseek(file[x], 0, SEEK_SET);
         if( (int)fread(header, 1, sizeof(header), file[x]) <= 0 ||
             sscanf(header, "P6\n%d %d\n%d", &width[x], &height[x], &i) != 3 ||
             i != 255 )
            goto error;
         for(i = c = offset[x] = 0; i < sizeof(header) && offset[x] == 0; i++)
            if( header[i] == '\n' )
               if( ++c == 3 )
                  offset[x] = i + 1;
         if( offset[x] == 0 ||
             (buffer[x] = (unsigned char*)malloc(width[x] * 3)) == NULL )
            goto error;
         fseek(file[x], offset[x], SEEK_SET);
      }

      /* Mark all pixel combinations as unused */
      for(x = 0; x < 256; x++)
         for(y = 0; y < 256; y++)
            pixels[x][y] = 0;

      /* Record all pixel combinations */
      for(y = 0; y < *height; y++)
      {
         for(i = 0; i < 2; i++)
            if( fread(buffer[i], *width * 3, 1, file[i]) != 1 )
               goto error;
         for(x = 0; x < *width * 3; x++)
            pixels[buffer[0][x]][buffer[1][x]] = 1;
      }

      /* Set initial range and try to expand outward */
      best[0] = best[2] = scale[0] = scale[2] = 64;
      best[1] = best[3] = scale[1] = scale[3] = 256 - 64;
      best_area = (best[1] - best[0]) * (best[3] - best[2]);
      for(step = 64; step > 0; step /= 2)
      {
         for(ip0 = -1; ip0 <= 0; ip0++)
         {
            z[0] = scale[0] + ip0 * step;
            if( z[0] >= 0 )
            {
               for(iq0 = -1; iq0 <= 0; iq0++)
               {
                  z[2] = scale[2] + iq0 * step;
                  if( z[2] >= 0 )
                  {
                     for(ip1 = 1; ip1 >= 0; ip1--)
                     {
                        z[1] = scale[1] + ip1 * step;
                        if( z[1] <= 256 )
                        {
                           for(iq1 = 1; iq1 >= 0; iq1--)
                           {
                              z[3] = scale[3] + iq1 * step;
                              if( z[3] <= 256 )
                              {
                                 z_area = (z[1] - z[0]) * (z[3] - z[2]);
                                 if( z_area > best_area )
                                 {
                                    for(x = c = 0; x < 256 && c >= 0 && c < 256; x++)
                                    {
                                       x0 = x * (z[1] - z[0]) / 256 + z[0];
                                       for(y = 0; y < 256 && c >= 0 && c < 256; y++)
                                       {
                                          y0 = y * (z[3] - z[2]) / 256 + z[2];
                                          if( pixels[x][y] != 0 )
                                             c = (3 * x0 - y0) / 2;
                                       }
                                    }
                                    if( c >= 0 && c < 256 )
                                    {
                                       for(i = 0; i < 4; i++)
                                          best[i] = z[i];
                                       best_area = z_area;
                                    }
                                 }
                              }
                           }
                        }
                     }
                  }
               }
            }
         }
         for(i = 0; i < 4; i++)
            scale[i] = best[i];
      }

      /* Open output */
      if( (file[2] = fopen(argv[3], "wb+")) == NULL )
         goto error;
      else
      {
         /* Write output header */
         fprintf(file[2], "P6\n%d %d\n255\n", *width * 2, *height * 2);

         /* Merge images */
         for(i = 0; i < 2; i++)
            fseek(file[i], offset[i], SEEK_SET);
         c = argc > 4 ? 1 : 0;
         for(y = 0; y < *height; y++)
         {
            for(i = 0; i < 2; i++)
               if( fread(buffer[i], *width * 3, 1, file[i]) != 1 )
                  goto error;

            for(step = 0; step <= 6; step += 6)
            {
               for(x = 0; x < *width; x++)
               {
                  for(x0 = 0; x0 < 3; x0++)
                  {
                     p = (float)buffer[0][x * 3 + x0] * (scale[1] - scale[0]) / 256.f + scale[0];
                     q = (float)buffer[1][x * 3 + x0] * (scale[3] - scale[2]) / 256.f + scale[2];
                     r = (3 * p - q) / 2.0f;
                     fpixel[index[c][x0][0]] = p;
                     fpixel[index[c][x0][1]] = q;
                     fpixel[index[c][x0][2]] = r;
                     fpixel[index[c][x0][3]] = r;
                  }
                  for(x0 = 0; x0 < 12; x0++)
                  {
                     ip0 = fpixel[x0];
                     q = fpixel[x0] - ip0;
                     if( (float)rand() / (float)RAND_MAX < q )
                        ip0++;
                     pixel[x0] = (unsigned char)(ip0 < 0 ? 0 : (ip0 > 255 ? 255 : ip0));
                  }
                  fwrite(pixel + step, 6, 1, file[2]);
               }
            }
         }
      }
   }

   for(i = 0; i < 2; i++)
   {
      free(buffer[i]);
      fclose(file[i]);
   }
   return fclose(file[2]);

error:
   return !puts("error");
}