#include<math.h>
#include<stdio.h>
#include<stdlib.h>
#include<string.h>
#include<time.h>
#include<png.h>

typedef int _;
typedef double f;

png_image N;
png_bytep Z[2];

_ C,H,I,S,E,
  e[50][2],
  s, i, j, x, y, p, q, *o,
  t, r,
  g[50][50],
  c[50];

f T[6], P, Q;

_ b(_ u, char *v)
{
   return printf(u ? u - 1 ? "Error writing %s\n"
                           : "Error reading %s\n"
                   : "%s <input.png> <output1.png> <output2.png>\n",
                 v);
}

f R(f u)
{
   return (f)rand() / (RAND_MAX + 1u) * u;
}

_ a()
{
   _ d[50],
     W = 5,
     z = -1,
     u, v, w;

   for(memset(d, u = 0, sizeof(d)); u < 50; u++)
      if( !c[u] )
      {
         w = 4;
         for(v = 0; v < 50; v++)
            if( u - v && g[u][v] && c[v] &&
                !(d[u] & (1 << c[v])) )
            {
               d[u] |= 1 << c[v];
               w--;
               if( w == 0 )
                  return 0;
            }
         if( W > w )
         {
            W = w;
            z = u;
         }
      }

   if( z < 0 )
      return 1;

   for(u = 1; u <= 4; u++)
      if( !(d[z] & (1 << u)) )
      {
         c[z] = u;
         if( a() )
            return 1;

         c[z] = 0;
      }

   return 0;
}

_ D(_ u)
{
   r = e[--u][0] - p;
   t = e[u][1] - q;
   return r * r + t * t;
}

_ M(_ u)
{
   return u > 0 ? u % 50
                : (50 - (-u % 50)) % 50;
}

_ main(_ u, char **v)
{
   if( u - 4 )
      return b(0, *v);

   N.version = PNG_IMAGE_VERSION;
   if( png_image_begin_read_from_file(&N, v[1]) )
   {
      N.format = PNG_FORMAT_RGBA;
      C = N.width;
      H = N.height;
      for(s = C * H; i < 2; i++)
         Z[i] = (png_bytep)calloc(s, 4);
      t = *Z && Z[1] &&
          png_image_finish_read(&N, NULL, *Z, 0, NULL);
   }
   if( !t )
      return b(1, v[1]);

   srand(time(NULL));

   o = (_*)Z[1];

   if( sizeof('c') > 1 )
   {
      P = R(40) + 10;
      P /= C > H ? C : H;
      T[0] = T[4] = cos(Q = R(22.0 / 7)) * P;
      T[3] = -(T[1] = sin(Q) * P);
      T[2] = R(4);
      T[5] = R(4);

      T[0] /= P = sqrt(3);
      T[1] /= P;

      for(y = 0; y < 50; y++)
         for(x = 0; x < 50;)
            g[y][x++] = R(4);

      for(y = 0; y < H; y++)
         for(x = 0; x < C;)
         {
            P = T[0] * x + T[1] * y + T[2];
            Q = T[3] * x + T[4] * y + T[5];

            t = floor(P + Q);
            r = floor((floor(Q - P) + t) / 3);

            o[y * C + x++] =
               g[M(floor((floor(2 * P + 1) + t) / 3) - r)][M(r)];
         }
   }
   else
   {
      for(i = 0; i < 50; o[y * C + x] = ++i)
      {
         e[i][0] = x = R(C);
         e[i][1] = y = R(H);
      }

      for(S = 1; S * 2 < C && S * 2 < H; S *= 2);
      for(; S > 0; S /= 2)
         for(y = 0; y < H; y++)
            for(x = 0; x < C; x++)
               if( (I = o[y * C + x]) )
                  for(i = -S; i <= S; i += S)
                  {
                     q = y + i;
                     if( q >= 0 && q < H )
                        for(j = -S; j <= S; j += S)
                        {
                           p = x + j;
                           if( p >= 0 && p < C )
                           {
                              E = o[q * C + p];

                              o[q * C + p] =
                                 !E || D(E) > D(I)
                                    ? I
                                    : E;
                           }
                        }
                  }

      for(i = 0; i < s; i++)
      {
         I = o[i] - 1;
         if( i % C )
         {
            E = o[i - 1] - 1;
            g[I][E] =
            g[E][I] = 1;
         }
         if( i >= C )
         {
            E = o[i - C] - 1;
            g[I][E] =
            g[E][I] = 1;
         }
      }

      a();

      for(i = 0; i < s; i++)
         o[i] = c[o[i] - 1] - 1;
   }

   r = R(8) + 2;
   t = R(8) + 2;
   p = 2 * t;
   q = 2 * r;
   for(y = i = 0; y < H; y++)
      for(x = 0; x < C; x++, i++)
         if( (s = o[i])
                ? s > 2
                   ? 1
                   : s < 2
                      ? R(2) < ((x % r * t > y % t * r) ^
                                ((x / r + y / t) % 2))
                      : R(2) > ((x % p * q > (q - 1 - y % q) * p) ^
                                ((x / p + y / q) % 2))
               : 0 )
         {
            memcpy(Z[1] + i * 4, *Z + i * 4, 4);
            memset(*Z + i * 4, 0, 4);
         }
         else
            o[i] = 0;

   for(i = 0; i < 2; i++)
      if( !png_image_write_to_file(&N, v[i + 2], 0, Z[i], 0, NULL) )
         return b(2, v[i + 2]);

   return 0;
}