cacalabs
/
libpipi
1
0
Fork 0
Code Issues 0 Pull Requests 0 Releases 0 Wiki Activity
Browse Source

* More img2rubik experimentation, with a lot of debugging messages. 

git-svn-id: file:///srv/caca.zoy.org/var/lib/svn/libpipi/trunk@2731 92316355-f0b4-4df1-b90c-862c8a59935f
remotes/tiles
sam 16 years ago
parent
e23268c286
commit
0151c722cd
1 changed files with 193 additions and 76 deletions
Split View
Diff Options
Show Stats Download Patch File Download Diff File
  1. +193
    -76
      examples/img2rubik.c

+ 193
- 76
examples/img2rubik.c View File

@@ -8,6 +8,9 @@

#include <pipi.h>

#define R 0
#define G 1
#define B 2
#define X 3
#define Y 4
#define A 5
@@ -15,7 +18,7 @@
#define BRIGHT(x) (0.299*(x)[0] + 0.587*(x)[1] + 0.114*(x)[2])

#define MAXCOLORS 16
#define STEPS 256
#define STEPS 32
#define EPSILON (0.000001)

typedef struct
@@ -37,19 +40,19 @@ static void init_uv(void)
{
double tmp;

ylen = sqrt(y[0] * y[0] + y[1] * y[1] + y[2] * y[2]);
ylen = sqrt(y[R] * y[R] + y[G] * y[G] + y[B] * y[B]);

u[0] = y[1];
u[1] = -y[0];
u[2] = 0;
tmp = sqrt(u[0] * u[0] + u[1] * u[1] + u[2] * u[2]);
u[0] /= tmp; u[1] /= tmp; u[2] /= tmp;
u[R] = y[1];
u[G] = -y[0];
u[B] = 0;
tmp = sqrt(u[R] * u[R] + u[G] * u[G] + u[B] * u[B]);
u[R] /= tmp; u[G] /= tmp; u[B] /= tmp;

v[0] = y[1] * u[2] - y[2] * u[1];
v[1] = y[2] * u[0] - y[0] * u[2];
v[2] = y[0] * u[1] - y[1] * u[0];
tmp = sqrt(v[0] * v[0] + v[1] * v[1] + v[2] * v[2]);
v[0] /= tmp; v[1] /= tmp; v[2] /= tmp;
v[R] = y[G] * u[B] - y[B] * u[G];
v[G] = y[B] * u[R] - y[R] * u[B];
v[B] = y[R] * u[G] - y[G] * u[R];
tmp = sqrt(v[R] * v[R] + v[G] * v[G] + v[B] * v[B]);
v[R] /= tmp; v[G] /= tmp; v[B] /= tmp;
}

/*
@@ -64,9 +67,9 @@ static hull_t *compute_hull(int ncolors, double const *palette)
double pal[ncolors][3];
for(i = 0; i < ncolors; i++)
{
pal[i][0] = palette[i * 3];
pal[i][1] = palette[i * 3 + 1];
pal[i][2] = palette[i * 3 + 2];
pal[i][R] = palette[i * 3];
pal[i][G] = palette[i * 3 + 1];
pal[i][B] = palette[i * 3 + 2];
}

/*
@@ -91,9 +94,9 @@ static hull_t *compute_hull(int ncolors, double const *palette)

double gray[3];

gray[0] = light[0] - dark[0];
gray[1] = light[1] - dark[1];
gray[2] = light[2] - dark[2];
gray[R] = light[R] - dark[R];
gray[G] = light[G] - dark[G];
gray[B] = light[B] - dark[B];

/*
* 3. Browse the grey axis and do stuff
@@ -101,8 +104,9 @@ static hull_t *compute_hull(int ncolors, double const *palette)
int n;
for(n = 0; n <= STEPS; n++)
{
double pts[ncolors * (ncolors - 1) / 2][6];
double ptmp[6];
printf("slice %i\n", n);
double pts[ncolors * (ncolors - 1) / 2][5];
double ptmp[5];
#define SWAP(p1,p2) do { memcpy(ptmp, p1, sizeof(ptmp)); \
memcpy(p1, p2, sizeof(ptmp)); \
memcpy(p2, ptmp, sizeof(ptmp)); } while(0)
@@ -110,9 +114,9 @@ static hull_t *compute_hull(int ncolors, double const *palette)
int npts = 0;

double p0[3];
p0[0] = dark[0] + t * gray[0];
p0[1] = dark[1] + t * gray[1];
p0[2] = dark[2] + t * gray[2];
p0[R] = dark[R] + t * gray[R];
p0[G] = dark[G] + t * gray[G];
p0[B] = dark[B] + t * gray[B];

/*
* 3.1. Find all edges that intersect the t.y + (u,v) plane
@@ -120,13 +124,13 @@ static hull_t *compute_hull(int ncolors, double const *palette)
for(i = 0; i < ncolors; i++)
{
double k1[3];
k1[0] = pal[i][0] - p0[0];
k1[1] = pal[i][1] - p0[1];
k1[2] = pal[i][2] - p0[2];
tmp = sqrt(k1[0] * k1[0] + k1[1] * k1[1] + k1[2] * k1[2]);
k1[R] = pal[i][R] - p0[R];
k1[G] = pal[i][G] - p0[G];
k1[B] = pal[i][B] - p0[B];
tmp = sqrt(k1[R] * k1[R] + k1[G] * k1[G] + k1[B] * k1[B]);

/* If k1.y > t.y.y, we don't want this point */
double yk1 = y[0] * k1[0] + y[1] * k1[1] + y[2] * k1[2];
double yk1 = y[R] * k1[R] + y[G] * k1[G] + y[B] * k1[B];
if(yk1 > t * ylen * ylen + EPSILON)
continue;

@@ -136,13 +140,13 @@ static hull_t *compute_hull(int ncolors, double const *palette)
continue;

double k2[3];
k2[0] = pal[j][0] - p0[0];
k2[1] = pal[j][1] - p0[1];
k2[2] = pal[j][2] - p0[2];
tmp = sqrt(k2[0] * k2[0] + k2[1] * k2[1] + k2[2] * k2[2]);
k2[R] = pal[j][R] - p0[R];
k2[G] = pal[j][G] - p0[G];
k2[B] = pal[j][B] - p0[B];
tmp = sqrt(k2[R] * k2[R] + k2[G] * k2[G] + k2[B] * k2[B]);

/* If k2.y < t.y.y, we don't want this point */
double yk2 = y[0] * k2[0] + y[1] * k2[1] + y[2] * k2[2];
double yk2 = y[R] * k2[R] + y[G] * k2[G] + y[B] * k2[B];
if(yk2 < t * ylen * ylen - EPSILON)
continue;

@@ -152,9 +156,9 @@ static hull_t *compute_hull(int ncolors, double const *palette)
double s = yk1 == yk2 ?
0.5 : (t * ylen * ylen - yk1) / (yk2 - yk1);

pts[npts][0] = p0[0] + k1[0] + s * (k2[0] - k1[0]);
pts[npts][1] = p0[1] + k1[1] + s * (k2[1] - k1[1]);
pts[npts][2] = p0[2] + k1[2] + s * (k2[2] - k1[2]);
pts[npts][R] = p0[R] + k1[R] + s * (k2[R] - k1[R]);
pts[npts][G] = p0[G] + k1[G] + s * (k2[G] - k1[G]);
pts[npts][B] = p0[B] + k1[B] + s * (k2[B] - k1[B]);
npts++;
}
}
@@ -167,12 +171,12 @@ static hull_t *compute_hull(int ncolors, double const *palette)
/* Make our problem a 2-D problem. */
for(i = 0; i < npts; i++)
{
pts[i][X] = (pts[i][0] - p0[0]) * u[0]
+ (pts[i][1] - p0[1]) * u[1]
+ (pts[i][2] - p0[2]) * u[2];
pts[i][Y] = (pts[i][0] - p0[0]) * v[0]
+ (pts[i][1] - p0[1]) * v[1]
+ (pts[i][2] - p0[2]) * v[2];
pts[i][X] = (pts[i][R] - p0[R]) * u[R]
+ (pts[i][G] - p0[G]) * u[G]
+ (pts[i][B] - p0[B]) * u[B];
pts[i][Y] = (pts[i][R] - p0[R]) * v[R]
+ (pts[i][G] - p0[G]) * v[G]
+ (pts[i][B] - p0[B]) * v[B];
}

/* Find the leftmost point */
@@ -186,32 +190,63 @@ static hull_t *compute_hull(int ncolors, double const *palette)
}
SWAP(pts[0], pts[left]);

/* Sort the remaining points radially around pts[0] */
/* Sort the remaining points radially around pts[0]. Bubble sort
* is okay for small sizes, I don't care. */
for(i = 1; i < npts; i++)
for(j = 1; j < npts - i; j++)
if((pts[j][X] - pts[0][X]) * (pts[j + 1][Y] - pts[0][Y])
< (pts[j + 1][X] - pts[0][X]) * (pts[j][Y] - pts[0][Y]))
{
double k1 = (pts[j][X] - pts[0][X])
* (pts[j + 1][Y] - pts[0][Y]);
double k2 = (pts[j + 1][X] - pts[0][X])
* (pts[j][Y] - pts[0][Y]);
if(k1 < k2)
SWAP(pts[j], pts[j + 1]);
else if(k1 == k2)
{
/* Aligned! keep the farthest point */
double ax = pts[j][X] - pts[0][X];
double ay = pts[j][Y] - pts[0][Y];
double bx = pts[j + 1][X] - pts[0][X];
double by = pts[j + 1][Y] - pts[0][Y];

if(ax * ax + ay * ay > bx * bx + by * by)
SWAP(pts[j], pts[j + 1]);
}
}

for(i = 0; i < npts; i++)
printf(" : %g %g\n", pts[i][X], pts[i][Y]);
/* Remove points not in the convex hull */
for(i = 2; i < npts; /* */)
{
printf(" # %g %g", pts[i][X], pts[i][Y]);
if(i < 2)
{
printf(" skipping\n");
i++;
continue;
}

if((pts[i - 1][X] - pts[i - 2][X]) * (pts[i][Y] - pts[i - 2][Y])
< (pts[i][X] - pts[i - 2][X]) * (pts[i - 1][Y] - pts[i - 2][Y]))
double k1 = (pts[i - 1][X] - pts[i - 2][X])
* (pts[i][Y] - pts[i - 2][Y]);
double k2 = (pts[i][X] - pts[i - 2][X])
* (pts[i - 1][Y] - pts[i - 2][Y]);
printf(" %g < %g ? ", k1, k2);
if(k1 <= k2 + EPSILON)
{
printf("yes, removing\n");
for(j = i - 1; j < npts - 1; j++)
SWAP(pts[j], pts[j + 1]);
npts--;
}
else
{
i++;
printf("no\n");
}
}
for(i = 0; i < npts; i++)
printf(" : %g %g\n", pts[i][X], pts[i][Y]);

/* Compute the barycentre coordinates */
double ctx = 0., cty = 0., weight = 0.;
@@ -255,17 +290,22 @@ static hull_t *compute_hull(int ncolors, double const *palette)
* 3.3. Store the barycentre and convex hull information.
*/

ret->bary[n][0] = p0[0] + ctx * u[0] + cty * v[0];
ret->bary[n][1] = p0[1] + ctx * u[1] + cty * v[1];
ret->bary[n][2] = p0[2] + ctx * u[2] + cty * v[2];
//if(ncolors == 6) printf("%i %g %g %g\n", n, ret->bary[n][0], ret->bary[n][1], ret->bary[n][2]);
ret->bary[n][R] = p0[R] + ctx * u[R] + cty * v[R];
ret->bary[n][G] = p0[G] + ctx * u[G] + cty * v[G];
ret->bary[n][B] = p0[B] + ctx * u[B] + cty * v[B];
//if(ncolors == 6) printf("%i %g %g %g\n", n, ret->bary[n][R], ret->bary[n][G], ret->bary[n][B]);

for(i = 0; i < npts; i++)
{
ret->pts[n][i][0] = pts[i][0];
ret->pts[n][i][1] = pts[i][1];
ret->pts[n][i][2] = pts[i][2];
ret->pts[n][i][R] = pts[i][R];
ret->pts[n][i][G] = pts[i][G];
ret->pts[n][i][B] = pts[i][B];
ret->pts[n][i][X] = pts[i][X] - ctx;
ret->pts[n][i][Y] = pts[i][Y] - cty;
ret->pts[n][i][A] = atan2(pts[i][Y] - cty, pts[i][X] - ctx);
printf("%g %g -> %g\n", ret->pts[n][i][X], ret->pts[n][i][Y], ret->pts[n][i][A]);
}
printf("\n");
ret->hullsize[n] = npts;
}

@@ -283,12 +323,12 @@ int main(int argc, char *argv[])

static double const mypal[] =
{
0.8, 0.0, 0.0, /* red */
0.0, 0.8, 0.0, /* green */
0.0, 0.0, 0.6, /* blue */
1.0, 1.0, 1.0, /* white */
0.9, 0.9, 0.0, /* yellow */
0.8, 0.4, 0.0, /* orange */
0.8, 0.001, 0.001, /* red */
0.001, 0.8, 0.001, /* green */
0.05, 0.001, 0.4, /* blue */
0.9, 0.9, 0.9, /* white */
0.9, 0.9, 0.001, /* yellow */
0.8, 0.4, 0.001, /* orange */
};

int i, j;
@@ -317,24 +357,101 @@ int main(int argc, char *argv[])
{
double p[3];
/* FIXME: Imlib fucks up the RGB order. */
p[2] = srcdata[4 * (j * w + i)];
p[1] = srcdata[4 * (j * w + i) + 1];
p[0] = srcdata[4 * (j * w + i) + 2];
p[B] = srcdata[4 * (j * w + i)];
p[G] = srcdata[4 * (j * w + i) + 1];
p[R] = srcdata[4 * (j * w + i) + 2];
double gray = BRIGHT(p);

double dx = (p[0] - gray * y[0]) * u[0]
+ (p[1] - gray * y[1]) * u[1]
+ (p[2] - gray * y[2]) * u[2];
double dy = (p[0] - gray * y[0]) * v[0]
+ (p[1] - gray * y[1]) * v[1]
+ (p[2] - gray * y[2]) * v[2];

dstdata[4 * (j * w + i)] = myhull->bary[(int)(gray * STEPS)][2]
+ dx * u[2] * .3 + dy * v[2] * .3;
dstdata[4 * (j * w + i) + 1] = myhull->bary[(int)(gray * STEPS)][1]
+ dx * u[1] * .3 + dy * v[1] * .3;
dstdata[4 * (j * w + i) + 2] = myhull->bary[(int)(gray * STEPS)][0]
+ dx * u[0] * .3 + dy * v[0] * .3;
/* Convert to 2-D */
double xp = (p[R] - gray) * u[R]
+ (p[G] - gray) * u[G]
+ (p[B] - gray) * u[B];
double yp = (p[R] - gray) * v[R]
+ (p[G] - gray) * v[G]
+ (p[B] - gray) * v[B];
printf("PIX(%g,%g,%g): gray %g + P(%g,%g)\n", p[R], p[G], p[B], gray, xp, yp);

int pt = (int)(gray * STEPS + 0.5);

/* 1. find the excentricity in RGB space. There is an easier
* way to do this, which is to find the intersection of our
* line with the RGB cube itself, but we'd lose the possibility
* of having an original colour space other than RGB. */

/* First, find the relevant triangle. */
int n, count = rgbhull->hullsize[pt];
double angle = atan2(yp, xp);
printf(" slice %i angle %g\n", pt, angle);
for(n = 0; n < count; n++)
{
double a1 = rgbhull->pts[pt][n][A];
double a2 = rgbhull->pts[pt][(n + 1) % count][A];
if(a1 > a2)
{
if(angle >= a1)
a2 += 2 * M_PI;
else
a1 -= 2 * M_PI;
}
if(angle >= a1 && angle <= a2)
break;
}

/* Now compute the distance to the triangle's edge. If the edge
* intersection is M, then t is such as P = t.M (can be zero) */
double xa = rgbhull->pts[pt][n % count][X];
double ya = rgbhull->pts[pt][n % count][Y];
double xb = rgbhull->pts[pt][(n + 1) % count][X];
double yb = rgbhull->pts[pt][(n + 1) % count][Y];
double t = (xp * (yb - ya) - yp * (xb - xa)) / (xa * yb - xb * ya);
printf(" best is %g %g (%g) -- %g %g (%g)\n", xa, ya, rgbhull->pts[pt][n % count][A], xb, yb, rgbhull->pts[pt][(n + 1) % count][A]);

/* 2. apply the excentricity in reduced space. */

count = myhull->hullsize[pt];
for(n = 0; n < count; n++)
{
double a1 = myhull->pts[pt][n][A];
double a2 = myhull->pts[pt][(n + 1) % count][A];
if(a1 > a2)
{
if(angle >= a1)
a2 += 2 * M_PI;
else
a1 -= 2 * M_PI;
}
if(angle >= a1 && angle <= a2)
break;
}

/* If the edge intersection is M', s is such as P = s.M'. We
* want P' = t.M' = t.P/s */
xa = myhull->pts[pt][n % count][X];
ya = myhull->pts[pt][n % count][Y];
xb = myhull->pts[pt][(n + 1) % count][X];
yb = myhull->pts[pt][(n + 1) % count][Y];
double s = (xp * (yb - ya) - yp * (xb - xa)) / (xa * yb - xb * ya);
printf(" apply to %g %g (%g) -- %g %g (%g)\n", xa, ya, myhull->pts[pt][n % count][A], xb, yb, myhull->pts[pt][(n + 1) % count][A]);

if(s > 0)
{
xp *= t / s;
yp *= t / s;
}
printf(" t = %g s = %g -> (%g,%g)\n", t, s, xp, yp);

dstdata[4 * (j * w + i)] = myhull->bary[(int)(gray * STEPS)][B]
+ xp * u[B] + yp * v[B];
dstdata[4 * (j * w + i) + 1] = myhull->bary[(int)(gray * STEPS)][G]
+ xp * u[G] + yp * v[G];
dstdata[4 * (j * w + i) + 2] = myhull->bary[(int)(gray * STEPS)][R]
+ xp * u[R] + yp * v[R];
if(dstdata[4 * (j * w + i)] < 0.) dstdata[4 * (j * w + i)] = 0.;
if(dstdata[4 * (j * w + i)] > 1.) dstdata[4 * (j * w + i)] = 1.;
if(dstdata[4 * (j * w + i) + 1] < 0.) dstdata[4 * (j * w + i) + 1] = 0.;
if(dstdata[4 * (j * w + i) + 1] > 1.) dstdata[4 * (j * w + i) + 1] = 1.;
if(dstdata[4 * (j * w + i) + 2] < 0.) dstdata[4 * (j * w + i) + 2] = 0.;
if(dstdata[4 * (j * w + i) + 2] > 1.) dstdata[4 * (j * w + i) + 2] = 1.;
}

pipi_save(dst, argv[2]);


Write Preview
Loading…
Cancel
Save