Shepards distortion
Posted: 2012-12-10T08:47:32-07:00
IM has two methods of creating a Shepards distortion. They produce similar but different results, and I believe one is wrong. (Perhaps the bug is in my understanding.) The methods are:
(1) With "-distort Shepards", eg:
(2) With "-sparse-color Shepards" to create a displacement file, followed by a "composite -displace". Eg:
I have written software that creates a displacement file, emulating IM's "-sparse-color Shepards", "Inverse" and "Voronoi". It gives identical results to the IM "-sparse-color" methods. Thus, I believe IM's method (1) is incorrect.
The guts of my software is:
In my software, the variable DistPower is the power the distance to each control point is to be raised. A value of 1 is equivalent to IM's "-sparse-color Inverse", while 2 is "-sparse-color Shepards". For my experiments, which is about matching photographs, larger powers give better results. I would like this feature to be added to ImageMagick, something like "-define ShepardsDistancePower:3.5". Ideally I'd like both "-distort Shepards" and "-sparse-color Shepards" to respect this, but only after the bug (if it is a bug) is fixed.
EDIT: should have said: IM 6.7.9 on Windows 7.
Perhaps I should also say that my simple examples above have only 2 control points, but my real-world images have 35 or more.
(1) With "-distort Shepards", eg:
Code: Select all
"%IMG%convert" rose: -depth 16 rose.tiff
"%IMG%identify" rose.tiff
rem Check that it really is 16 bits.
rem Method (1).
"%IMG%convert" rose16.tiff ^
-distort Shepards "10,10 20,10 60,10 60,5" ^
b_d_shep.tiff
Code: Select all
rem Method (2).
"%IMG%convert" rose16.tiff ^
-sparse-color Shepards "10 10 #7FF680000000 60 10 #800080050000" ^
ac.tiff
"%IMG%composite" ac.tiff rose16.tiff -displace 32768x32768 b_shep.tiff
The guts of my software is:
Code: Select all
static void WrPpm (void)
{
FILE * f;
if (strcmp(PpmFile, "-") == 0) f = stdout;
else f = fopen (PpmFile, "wt");
if (!f) FatalError ("Can't open %s", PpmFile);
fprintf (f, "P3\n");
fprintf (f, "%i %i\n", PpmWidth, PpmHeight);
fprintf (f, "%i\n", 65535);
int const mid_val = 32768;
int x, y, i;
float mult;
float max_dist = (PpmWidth + PpmHeight) * 2;
float ClosestDist;
for (y = 0; y < PpmHeight; y++) {
for (x = 0; x < PpmWidth; x++) {
float shiftX = 0.0, shiftY = 0.0;
float sigMult = 0.0;
ClosestDist = max_dist;
int ndxClosest;
for (i = 0; i < nValues; i++) {
ValueT * pv = &Values[i];
if (pv->include) {
float dx = (pv->x0 - x);
float dy = (pv->y0 - y);
float distSq = dx * dx + dy * dy;
if (distSq==0) distSq = 1;
if (DoVoronoi) {
float dist = sqrt(distSq);
if (ClosestDist > dist) {
ClosestDist = dist;
ndxClosest = i;
}
} else {
// For speed, treat DistPower 1 and 2 as special cases.
if (DistPower == 1) mult = 1.0/sqrt(distSq);
else if (DistPower == 2) mult = 1.0/distSq;
else mult = 1.0 / pow (distSq, DistPower/2);
shiftX += pv->dx * mult;
shiftY += pv->dy * mult;
sigMult += mult;
}
}
}
int dx, dy;
if (DoVoronoi) {
dx = Values[ndxClosest].dx + mid_val;
dy = Values[ndxClosest].dy + mid_val;
} else {
dx = (shiftX / sigMult + 0.5) + mid_val;
dy = (shiftY / sigMult + 0.5) + mid_val;
}
fprintf (f, "%i %i 0\n", dx, dy);
}
}
if (strcmp(PpmFile, "-") != 0) fclose (f);
if (gProcesses) printf ("Created %s\n", PpmFile);
}
EDIT: should have said: IM 6.7.9 on Windows 7.
Perhaps I should also say that my simple examples above have only 2 control points, but my real-world images have 35 or more.