Fixed some bad math in pixblock-transform

diff --git a/src/display/pixblock-scaler.cpp b/src/display/pixblock-scaler.cpp
index cefbfe4fbe5db4774d3f86b087fa6606267b1f57..84274d0ff30158d9493ee929e418ee5b57181947 100644 (file)
--- a/src/display/pixblock-scaler.cpp
+++ b/src/display/pixblock-scaler.cpp
@@ -6,7 +6,7 @@
  * Author:
  *   Niko Kiirala <niko@kiirala.com>
  *
- * Copyright (C) 2006 Niko Kiirala
+ * Copyright (C) 2006,2009 Niko Kiirala
  *
  * Released under GNU GPL, read the file 'COPYING' for more information
  */

diff --git a/src/display/pixblock-transform.cpp b/src/display/pixblock-transform.cpp
index d0ba898063fe0be69cdb2055987c7c80dfe8f3f4..73b467d5a8f4da558f0af1bd25d9373b74e15839 100644 (file)
--- a/src/display/pixblock-transform.cpp
+++ b/src/display/pixblock-transform.cpp
@@ -6,7 +6,7 @@
  * Author:
  *   Niko Kiirala <niko@kiirala.com>
  *
- * Copyright (C) 2006 Niko Kiirala
+ * Copyright (C) 2006,2009 Niko Kiirala
  *
  * Released under GNU GPL, read the file 'COPYING' for more information
  */
@@ -28,6 +28,9 @@ using std::floor;
 namespace NR {
 
 struct RGBA {
+    double r, g, b, a;
+};
+struct RGBAi {
     int r, g, b, a;
 };
 
@@ -72,14 +75,14 @@ void transform_nearest(NRPixBlock *to, NRPixBlock *from, Geom::Matrix const &tra
     // Loop through every pixel of destination image, a line at a time
     for (int to_y = 0 ; to_y < to_height ; to_y++) {
         for (int to_x = 0 ; to_x < to_width ; to_x++) {
-            RGBA result = {0,0,0,0};
+            RGBAi result = {0,0,0,0};
 
-            int from_x = (int)round(itrans[0] * (to_x + to->area.x0)
-                                    + itrans[2] * (to_y + to->area.y0)
+            int from_x = (int)floor(itrans[0] * (to_x + 0.5 + to->area.x0)
+                                    + itrans[2] * (to_y + 0.5 + to->area.y0)
                                     + itrans[4]);
             from_x -= from->area.x0;
-            int from_y = (int)round(itrans[1] * (to_x + to->area.x0)
-                                    + itrans[3] * (to_y + to->area.y0)
+            int from_y = (int)floor(itrans[1] * (to_x + 0.5 + to->area.x0)
+                                    + itrans[3] * (to_y + 0.5 + to->area.y0)
                                     + itrans[5]);
             from_y -= from->area.y0;
 
@@ -106,53 +109,28 @@ void transform_nearest(NRPixBlock *to, NRPixBlock *from, Geom::Matrix const &tra
 }
 
 /** Calculates cubically interpolated value of the four given pixel values.
- * The pixel values should be from four vertically adjacent pixels.
- * If we are calculating a pixel, whose y-coordinate in source image is
- * i, these pixel values a, b, c and d should come from lines
- * floor(i) - 1, floor(i), floor(i) + 1, floor(i) + 2, respectively.
- * Parameter len should be set to i.
- * Returns the interpolated value in fixed point format with 8 bit
- * decimal part. (24.8 assuming 32-bit int)
+ * The pixel values should be from four adjacent pixels in source image or
+ * four adjacent interpolated values. len should be the x- or y-coordinate
+ * (depending on interpolation direction) of the center of the target pixel
+ * in source image coordinates.
  */
 __attribute__ ((const))
-inline static int sampley(unsigned const char a, unsigned const char b,
-                   unsigned const char c, unsigned const char d,
-                   const double len)
+inline static double sample(double const a, double const b,
+                            double const c, double const d,
+                            double const len)
 {
-    double lenf = len - floor(len);
-    int sum = 0;
-    sum += (int)((((-1.0 / 3.0) * lenf + 4.0 / 5.0) * lenf - 7.0 / 15.0)
-                 * lenf * 256 * a);
-    sum += (int)((((lenf - 9.0 / 5.0) * lenf - 1.0 / 5.0) * lenf + 1.0)
-                 * 256 * b);
-    sum += (int)(((((1 - lenf) - 9.0 / 5.0) * (1 - lenf) - 1.0 / 5.0)
-                  * (1 - lenf) + 1.0) * 256 * c);
-    sum += (int)((((-1.0 / 3.0) * (1 - lenf) + 4.0 / 5.0) * (1 - lenf)
-                  - 7.0 / 15.0) * (1 - lenf) * 256 * d);
-    return sum;
-}
+    double lena = 1.5 + (len - round(len));
+    double lenb = 0.5 + (len - round(len));
+    double lenc = 0.5 - (len - round(len));
+    double lend = 1.5 - (len - round(len));
+    double const f = -0.5; // corresponds to cubic Hermite spline
+    double sum = 0;
+    sum += ((((f * lena) - 5.0 * f) * lena + 8.0 * f) * lena - 4 * f) * a;
+    sum += (((f + 2.0) * lenb - (f + 3.0)) * lenb * lenb + 1.0) * b;
+    sum += (((f + 2.0) * lenc - (f + 3.0)) * lenc * lenc + 1.0) * c;
+    sum += ((((f * lend) - 5.0 * f) * lend + 8.0 * f) * lend - 4 * f) * d;
 
-/** Calculates cubically interpolated value of the four given pixel values.
- * The pixel values should be interpolated values from sampley, from four
- * horizontally adjacent vertical lines. The parameters a, b, c and d
- * should be in fixed point format with 8-bit decimal part.
- * If we are calculating a pixel, whose x-coordinate in source image is
- * i, these vertical  lines from where a, b, c and d are calculated, should be
- * floor(i) - 1, floor(i), floor(i) + 1, floor(i) + 2, respectively.
- * Parameter len should be set to i.
- * Returns the interpolated value in 8-bit format, ready to be written
- * to output buffer.
- */
-inline static int samplex(const int a, const int b, const int c, const int d, const double len) {
-    double lenf = len - floor(len);
-    int sum = 0;
-    sum += (int)(a * (((-1.0 / 3.0) * lenf + 4.0 / 5.0) * lenf - 7.0 / 15.0) * lenf);
-    sum += (int)(b * (((lenf - 9.0 / 5.0) * lenf - 1.0 / 5.0) * lenf + 1.0));
-    sum += (int)(c * ((((1 - lenf) - 9.0 / 5.0) * (1 - lenf) - 1.0 / 5.0) * (1 - lenf) + 1.0));
-    sum += (int)(d * (((-1.0 / 3.0) * (1 - lenf) + 4.0 / 5.0) * (1 - lenf) - 7.0 / 15.0) * (1 - lenf));
-    //if (sum < 0) sum = 0;
-    //if (sum > 255 * 256) sum = 255 * 256;
-    return sum / 256;
+    return sum;
 }
 
 void transform_bicubic(NRPixBlock *to, NRPixBlock *from, Geom::Matrix const &trans)
@@ -182,11 +160,11 @@ void transform_bicubic(NRPixBlock *to, NRPixBlock *from, Geom::Matrix const &tra
     // Loop through every pixel of destination image, a line at a time
     for (int to_y = 0 ; to_y < to_height ; to_y++) {
         for (int to_x = 0 ; to_x < to_width ; to_x++) {
-            double from_x = itrans[0] * (to_x + to->area.x0)
-                + itrans[2] * (to_y + to->area.y0)
+            double from_x = itrans[0] * (to_x + 0.5 + to->area.x0)
+                + itrans[2] * (to_y + 0.5 + to->area.y0)
                 + itrans[4] - from->area.x0;
-            double from_y = itrans[1] * (to_x + to->area.x0)
-                + itrans[3] * (to_y + to->area.y0)
+            double from_y = itrans[1] * (to_x + 0.5 + to->area.x0)
+                + itrans[3] * (to_y + 0.5 + to->area.y0)
                 + itrans[5] - from->area.y0;
 
             if (from_x < 0 || from_x >= from_width ||
@@ -198,9 +176,10 @@ void transform_bicubic(NRPixBlock *to, NRPixBlock *from, Geom::Matrix const &tra
 
             int from_line[4];
             for (int i = 0 ; i < 4 ; i++) {
-                if ((int)floor(from_y) + i - 1 >= 0) {
-                    if ((int)floor(from_y) + i - 1 < from_height) {
-                        from_line[i] = ((int)floor(from_y) + i - 1) * from->rs;
+                int fy_line = (int)round(from_y) + i - 2;
+                if (fy_line >= 0) {
+                    if (fy_line < from_height) {
+                        from_line[i] = fy_line * from->rs;
                     } else {
                         from_line[i] = (from_height - 1) * from->rs;
                     }
@@ -210,7 +189,7 @@ void transform_bicubic(NRPixBlock *to, NRPixBlock *from, Geom::Matrix const &tra
             }
 
             for (int i = 0 ; i < 4 ; i++) {
-                int k = (int)floor(from_x) + i - 1;
+                int k = (int)round(from_x) + i - 2;
                 if (k < 0) k = 0;
                 if (k >= from_width) k = from_width - 1;
                 k *= 4;
@@ -218,36 +197,36 @@ void transform_bicubic(NRPixBlock *to, NRPixBlock *from, Geom::Matrix const &tra
                 _check_index(from, from_line[1] + k, __LINE__);
                 _check_index(from, from_line[2] + k, __LINE__);
                 _check_index(from, from_line[3] + k, __LINE__);
-                line[i].r = sampley(NR_PIXBLOCK_PX(from)[from_line[0] + k],
-                                    NR_PIXBLOCK_PX(from)[from_line[1] + k],
-                                    NR_PIXBLOCK_PX(from)[from_line[2] + k],
-                                    NR_PIXBLOCK_PX(from)[from_line[3] + k],
-                                    from_y);
-                line[i].g = sampley(NR_PIXBLOCK_PX(from)[from_line[0] + k + 1],
-                                    NR_PIXBLOCK_PX(from)[from_line[1] + k + 1],
-                                    NR_PIXBLOCK_PX(from)[from_line[2] + k + 1],
-                                    NR_PIXBLOCK_PX(from)[from_line[3] + k + 1],
-                                    from_y);
-                line[i].b = sampley(NR_PIXBLOCK_PX(from)[from_line[0] + k + 2],
-                                    NR_PIXBLOCK_PX(from)[from_line[1] + k + 2],
-                                    NR_PIXBLOCK_PX(from)[from_line[2] + k + 2],
-                                    NR_PIXBLOCK_PX(from)[from_line[3] + k + 2],
-                                    from_y);
-                line[i].a = sampley(NR_PIXBLOCK_PX(from)[from_line[0] + k + 3],
-                                    NR_PIXBLOCK_PX(from)[from_line[1] + k + 3],
-                                    NR_PIXBLOCK_PX(from)[from_line[2] + k + 3],
-                                    NR_PIXBLOCK_PX(from)[from_line[3] + k + 3],
-                                    from_y);
+                line[i].r = sample(NR_PIXBLOCK_PX(from)[from_line[0] + k],
+                                   NR_PIXBLOCK_PX(from)[from_line[1] + k],
+                                   NR_PIXBLOCK_PX(from)[from_line[2] + k],
+                                   NR_PIXBLOCK_PX(from)[from_line[3] + k],
+                                   from_y);
+                line[i].g = sample(NR_PIXBLOCK_PX(from)[from_line[0] + k + 1],
+                                   NR_PIXBLOCK_PX(from)[from_line[1] + k + 1],
+                                   NR_PIXBLOCK_PX(from)[from_line[2] + k + 1],
+                                   NR_PIXBLOCK_PX(from)[from_line[3] + k + 1],
+                                   from_y);
+                line[i].b = sample(NR_PIXBLOCK_PX(from)[from_line[0] + k + 2],
+                                   NR_PIXBLOCK_PX(from)[from_line[1] + k + 2],
+                                   NR_PIXBLOCK_PX(from)[from_line[2] + k + 2],
+                                   NR_PIXBLOCK_PX(from)[from_line[3] + k + 2],
+                                   from_y);
+                line[i].a = sample(NR_PIXBLOCK_PX(from)[from_line[0] + k + 3],
+                                   NR_PIXBLOCK_PX(from)[from_line[1] + k + 3],
+                                   NR_PIXBLOCK_PX(from)[from_line[2] + k + 3],
+                                   NR_PIXBLOCK_PX(from)[from_line[3] + k + 3],
+                                   from_y);
             }
             RGBA result;
-            result.r = samplex(line[0].r, line[1].r, line[2].r, line[3].r,
-                               from_x);
-            result.g = samplex(line[0].g, line[1].g, line[2].g, line[3].g,
-                               from_x);
-            result.b = samplex(line[0].b, line[1].b, line[2].b, line[3].b,
-                               from_x);
-            result.a = samplex(line[0].a, line[1].a, line[2].a, line[3].a,
-                               from_x);
+            result.r = round(sample(line[0].r, line[1].r, line[2].r, line[3].r,
+                                    from_x));
+            result.g = round(sample(line[0].g, line[1].g, line[2].g, line[3].g,
+                                    from_x));
+            result.b = round(sample(line[0].b, line[1].b, line[2].b, line[3].b,
+                                    from_x));
+            result.a = round(sample(line[0].a, line[1].a, line[2].a, line[3].a,
+                                    from_x));
 
             using Inkscape::Filters::clamp;
             using Inkscape::Filters::clamp_alpha;
@@ -255,20 +234,24 @@ void transform_bicubic(NRPixBlock *to, NRPixBlock *from, Geom::Matrix const &tra
             if (to->mode == NR_PIXBLOCK_MODE_R8G8B8A8P) {
                 /* Make sure, none of the RGB channels exceeds 100% intensity
                  * in premultiplied output */
-                result.a = clamp(result.a);
+                int const alpha = clamp((int)result.a);
                 NR_PIXBLOCK_PX(to)[to_y * to->rs + to_x * 4] = 
-                    clamp_alpha(result.r, result.a);
+                    clamp_alpha((int)result.r, alpha);
                 NR_PIXBLOCK_PX(to)[to_y * to->rs + to_x * 4 + 1] = 
-                    clamp_alpha(result.g, result.a);
+                    clamp_alpha((int)result.g, alpha);
                 NR_PIXBLOCK_PX(to)[to_y * to->rs + to_x * 4 + 2] = 
-                    clamp_alpha(result.b, result.a);
-                NR_PIXBLOCK_PX(to)[to_y * to->rs + to_x * 4 + 3] = result.a;
+                    clamp_alpha((int)result.b, alpha);
+                NR_PIXBLOCK_PX(to)[to_y * to->rs + to_x * 4 + 3] = alpha;
             } else {
                 /* Clamp the output to unsigned char range */
-                NR_PIXBLOCK_PX(to)[to_y * to->rs + to_x * 4] = clamp(result.r);
-                NR_PIXBLOCK_PX(to)[to_y * to->rs + to_x * 4 + 1] = clamp(result.g);
-                NR_PIXBLOCK_PX(to)[to_y * to->rs + to_x * 4 + 2] = clamp(result.b);
-                NR_PIXBLOCK_PX(to)[to_y * to->rs + to_x * 4 + 3] = clamp(result.a);
+                NR_PIXBLOCK_PX(to)[to_y * to->rs + to_x * 4]
+                    = clamp((int)result.r);
+                NR_PIXBLOCK_PX(to)[to_y * to->rs + to_x * 4 + 1]
+                    = clamp((int)result.g);
+                NR_PIXBLOCK_PX(to)[to_y * to->rs + to_x * 4 + 2]
+                    = clamp((int)result.b);
+                NR_PIXBLOCK_PX(to)[to_y * to->rs + to_x * 4 + 3]
+                    = clamp((int)result.a);
             }
         }
     }