From 9e9fb5864435d7376381ea142d7f8a2d40255f1d Mon Sep 17 00:00:00 2001 From: John Cupitt Date: Thu, 30 Oct 2008 13:25:26 +0000 Subject: [PATCH] oops, forgot yafrnohalo --- libsrc/mosaicing/gegl-sampler-yafr-nohalo.c | 1 + libsrc/mosaicing/yafrnohalo.c | 1087 +++++++++++++++++++ 2 files changed, 1088 insertions(+) create mode 100644 libsrc/mosaicing/yafrnohalo.c diff --git a/libsrc/mosaicing/gegl-sampler-yafr-nohalo.c b/libsrc/mosaicing/gegl-sampler-yafr-nohalo.c index 13b93fd6..c120ef10 100644 --- a/libsrc/mosaicing/gegl-sampler-yafr-nohalo.c +++ b/libsrc/mosaicing/gegl-sampler-yafr-nohalo.c @@ -1019,6 +1019,7 @@ gegl_sampler_yafr_get ( GeglSampler *self, left_width_times_down_hight, rite_width_times_down_hight, sampler_bptr++); + newval[1] = bilinear_yafr (c_horizo_left___up, c_horizo_left_down, c_horizo_rite___up, diff --git a/libsrc/mosaicing/yafrnohalo.c b/libsrc/mosaicing/yafrnohalo.c new file mode 100644 index 00000000..b16b5c1a --- /dev/null +++ b/libsrc/mosaicing/yafrnohalo.c @@ -0,0 +1,1087 @@ +/* yafrnohalo ... yafr-nohalo as a vips interpolate class + */ + +/* + + This file is part of VIPS. + + VIPS is free software; you can redistribute it and/or modify + it under the terms of the GNU Lesser General Public License as published by + the Free Software Foundation; either version 2 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU Lesser General Public License for more details. + + You should have received a copy of the GNU Lesser General Public License + along with this program; if not, write to the Free Software + Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA + + */ + +/* + + These files are distributed with VIPS - http://www.vips.ecs.soton.ac.uk + + */ + +/* +#define DEBUG + */ + +#ifdef HAVE_CONFIG_H +#include +#endif /*HAVE_CONFIG_H*/ +#include + +#include +#include + +#include +#include + +#ifdef WITH_DMALLOC +#include +#endif /*WITH_DMALLOC*/ + +/* "fast" floor() ... on my laptop, anyway. + */ +#define FLOOR( V ) ((V) >= 0 ? (int)(V) : (int)((V) - 1)) + +static VipsInterpolateClass *vips_interpolate_yafrnohalo_parent_class = NULL; + +/* Copy-paste of gegl-sampler-yafr-nohalo.c starts + */ + +/* + * 2008 (c) Nicolas Robidoux (developer of Yet Another Fast + * Resampler). + * + * Acknowledgement: N. Robidoux's research on YAFR funded in part by + * an NSERC (National Science and Engineering Research Council of + * Canada) Discovery Grant. + */ + +#ifndef restrict +#ifdef __restrict +#define restrict __restrict +#else +#ifdef __restrict__ +#define restrict __restrict__ +#else +#define restrict +#endif +#endif +#endif + +#ifndef unlikely +#ifdef __builtin_expect +#define unlikely(x) __builtin_expect((x),0) +#else +#define unlikely(x) (x) +#endif +#endif + +/* + * YAFR = Yet Another Fast Resampler + * + * Yet Another Fast Resampler is a nonlinear resampler which consists + * of a linear scheme (in this version, bilinear) plus a nonlinear + * sharpening correction the purpose of which is the straightening of + * diagonal interfaces between flat colour areas. + * + * Key properties: + * + * YAFR (nohalo) is a parameterized method: + * + * The key parameter is called "sharpening." When sharpening = 0, YAFR + * (nohalo) is standard bilinear interpolation, method which is + * coconvex (and comonotone). More generally, YAFR (nohalo) is a + * coconvex (and comonotone) method when sharpening <= .5, When + * sharpening <= 1, YAFR (nohalo) is comonotone (but not coconvex), + * and consequently does not add haloing where it is not present. For + * the default value of sharpening, namely 4/3, YAFR (nohalo) does add + * some haloing, so little that it is generally not noticeable; + * however, for this sharpening value, YAFR (nohalo) correctly fixes + * key values near a sharp diagonal interface between flat color + * areas. In some situations (for example, photographs of text), + * values of sharpening as large as 2 are recommended. (Negative + * values of sharpening gives "artistic" results, as do very large + * ones). + * + * YAFR (nohalo) is interpolatory: + * + * If asked for the value at the center of an input pixel, it will + * return the corresponding value, unchanged. + * + * YAFR (nohalo) preserves local averages: + * + * The average of the reconstructed intensity surface over any region + * is the same as the average of the piecewise constant surface with + * values over pixel areas equal to the input pixel values (the + * "nearest neighbour" surface), except for a small amount of blur at + * the boundary of the region. More precicely: YAFR (nohalo) is a box + * filtered exact area method. + * + * Main weaknesses of YAFR (nohalo): + * + * YAFR (nohalo) improves on the standard bilinear method only for + * images with at least a little bit of smoothness. + * + * Like bilinear, YAFR (nohalo) often suffers from noticeable Mach + * banding. When the images are reasonably smooth, however, the Mach + * banding produced by YAFR for moderate values of the sharpening + * parameter is much less noticeable than the one for bilinear. + */ + +/* Pointers to write to / read from, how much to add to move right a pixel, + * how much to add to move down a line. + */ + +static inline void +bilinear_yafrnohalo (float* restrict out, const float* restrict in, + const int channels, + const int pixels_per_buffer_row, + const float user_sharpening, + + const gfloat c_horizo_left___up, + const gfloat c_horizo_left_down, + const gfloat c_horizo_rite___up, + const gfloat c_horizo_rite_down, + const gfloat c_vertic_left___up, + const gfloat c_vertic_rite___up, + const gfloat c_vertic_left_down, + const gfloat c_vertic_rite_down, + const gfloat c_settin_left___up, + const gfloat c_settin_rite___up, + const gfloat c_settin_left_down, + const gfloat c_settin_rite_down, + const gfloat c_rising_left___up, + const gfloat c_rising_rite___up, + const gfloat c_rising_left_down, + const gfloat c_rising_rite_down, + const gfloat left_width_times___up_hight, + const gfloat rite_width_times___up_hight, + const gfloat left_width_times_down_hight, + const gfloat rite_width_times_down_hight) +{ + + /* + * "sharpening" is a continuous method parameter which is + * proportional to the amount of "diagonal straightening" which the + * nonlinear correction part of the method may add to the underlying + * linear scheme. You may also think of it as a sharpening + * parameter: higher values correspond to more sharpening, and + * negative values lead to strange looking effects. + * + * The default value is sharpening = 4/3 when the scheme being + * "straightened" is bilinear---as is the case here. This value + * fixes key pixel values near the diagonal boundary between two + * monochrome regions (the diagonal boundary pixel values being set + * to the halfway colour). + * + * If resampling seems to add unwanted texture artifacts, push + * sharpening toward 0. It is not generally not recommended to set + * sharpening to a value much larger than 2. + * + * In order to simplify interfacing with users, the parameter which + * should be set by the user is normalized so that user_sharpening = + * 1 when sharpening is equal to the recommended value. Consistently + * with the above discussion, values of user_sharpening between 0 + * and about 1.5 give good results. + */ + const gfloat sharpening = user_sharpening * ( 4.f / 3.f ); + + /* + * The input pixel values are described by the following stencil. + * Spanish abbreviations are used to label positions from top to + * bottom, English ones to label positions from left to right,: + * + * (ix-1,iy-1) (ix,iy-1) (ix+1,iy-1) (ix+2,iy-1) + * =uno_one =uno_two =uno_thr = uno_fou + * + * (ix-1,iy) (ix,iy) (ix+1,iy) (ix+2,iy) + * =dos_one =dos_two =dos_thr = dos_fou + * + * (ix-1,iy+1) (ix,iy+1) (ix+1,iy+1) (ix+2,iy+1) + * =tre_one =tre_two =tre_thr = tre_fou + * + * (ix-1,iy+2) (ix,iy+2) (ix+1,iy+2) (ix+2,iy+2) + * =qua_one =qua_two =qua_thr = qua_fou + */ + + /* + * Load the useful pixel values for the channel under + * consideration. The in pointer is assumed + * to point to uno_one when bilinear_yafrnohalo is entered. + */ + const float uno_one = in[ 0 ]; + const float uno_two = in[ channels ]; + const float uno_thr = in[ 2 * channels ]; + const float uno_fou = in[ 3 * channels ]; + + const float dos_one = in[ pixels_per_buffer_row * channels ]; + const float dos_two = in[ ( 1 + pixels_per_buffer_row ) * channels ]; + const float dos_thr = in[ ( 2 + pixels_per_buffer_row ) * channels ]; + const float dos_fou = in[ ( 3 + pixels_per_buffer_row ) * channels ]; + + const float tre_one = in[ 2 * pixels_per_buffer_row * channels ]; + const float tre_two = in[ ( 1 + 2 * pixels_per_buffer_row ) * channels ]; + const float tre_thr = in[ ( 2 + 2 * pixels_per_buffer_row ) * channels ]; + const float tre_fou = in[ ( 3 + 2 * pixels_per_buffer_row ) * channels ]; + + const float qua_one = in[ 3 * pixels_per_buffer_row * channels ]; + const float qua_two = in[ ( 1 + 3 * pixels_per_buffer_row ) * channels ]; + const float qua_thr = in[ ( 2 + 3 * pixels_per_buffer_row ) * channels ]; + const float qua_fou = in[ ( 3 + 3 * pixels_per_buffer_row ) * channels ]; + + /* + * Bilinear (piecewise constant histopolant pieces) baseline + * contribution: + */ + const gfloat bilinear = + left_width_times___up_hight * dos_two + + + rite_width_times___up_hight * dos_thr + + + left_width_times_down_hight * tre_two + + + rite_width_times_down_hight * tre_thr; + + /* + * Computation of the YAFR correction: + * + * Basically, if two consecutive pixel value differences have the + * same sign, the smallest one (in absolute value) is taken to be + * the corresponding slope. If they don't have the same sign, the + * corresponding slope is set to 0. + * + * For each of the four overlapped pixels, four slopes are thus + * computed: horizontal, vertical, "rising" (45 degree angle) and + * "setting" (-45 degree angle). + */ + + /* + * Beginning of the computation of the "up" horizontal differences + * (left to right): + */ + const gfloat prem___up_horizo = dos_two - dos_one; + const gfloat deux___up_horizo = dos_thr - dos_two; + const gfloat troi___up_horizo = dos_fou - dos_thr; + /* + * "down" horizontal differences: + */ + const gfloat prem_down_horizo = tre_two - tre_one; + const gfloat deux_down_horizo = tre_thr - tre_two; + const gfloat troi_down_horizo = tre_fou - tre_thr; + /* + * "left" vertical differences (top to bottom): + */ + const gfloat prem_left_vertic = dos_two - uno_two; + const gfloat deux_left_vertic = tre_two - dos_two; + const gfloat troi_left_vertic = qua_two - tre_two; + /* + * "right" vertical differences: + */ + const gfloat prem_rite_vertic = dos_thr - uno_thr; + const gfloat deux_rite_vertic = tre_thr - dos_thr; + const gfloat troi_rite_vertic = qua_thr - tre_thr; + /* + * "left rising" diagonal differences (bottom left to upper right): + */ + const gfloat prem_left_rising = dos_two - tre_one; + const gfloat deux_left_rising = uno_thr - dos_two; + /* + * "middle rising" diagonal differences: + */ + const gfloat prem_midl_rising = tre_two - qua_one; + const gfloat deux_midl_rising = dos_thr - tre_two; + const gfloat troi_midl_rising = uno_fou - dos_thr; + /* + * "right rising" diagonal differences: + */ + const gfloat prem_rite_rising = tre_thr - qua_two; + const gfloat deux_rite_rising = dos_fou - tre_thr; + /* + * "left setting" diagonal differences (top left to bottom right): + */ + const gfloat prem_left_settin = tre_two - dos_one; + const gfloat deux_left_settin = qua_thr - tre_two; + /* + * "middle setting" diagonal differences: + */ + const gfloat prem_midl_settin = dos_two - uno_one; + const gfloat deux_midl_settin = tre_thr - dos_two; + const gfloat troi_midl_settin = qua_fou - tre_thr; + /* + * "right" setting diagonal differences: + */ + const gfloat prem_rite_settin = dos_thr - uno_two; + const gfloat deux_rite_settin = tre_fou - dos_thr; + + /* + * Back to "up": + */ + const gfloat prem___up_horizo_squared = prem___up_horizo * prem___up_horizo; + const gfloat deux___up_horizo_squared = deux___up_horizo * deux___up_horizo; + const gfloat troi___up_horizo_squared = troi___up_horizo * troi___up_horizo; + /* + * Back to "down": + */ + const gfloat prem_down_horizo_squared = prem_down_horizo * prem_down_horizo; + const gfloat deux_down_horizo_squared = deux_down_horizo * deux_down_horizo; + const gfloat troi_down_horizo_squared = troi_down_horizo * troi_down_horizo; + /* + * Back to "left": + */ + const gfloat prem_left_vertic_squared = prem_left_vertic * prem_left_vertic; + const gfloat deux_left_vertic_squared = deux_left_vertic * deux_left_vertic; + const gfloat troi_left_vertic_squared = troi_left_vertic * troi_left_vertic; + /* + * Back to "right": + */ + const gfloat prem_rite_vertic_squared = prem_rite_vertic * prem_rite_vertic; + const gfloat deux_rite_vertic_squared = deux_rite_vertic * deux_rite_vertic; + const gfloat troi_rite_vertic_squared = troi_rite_vertic * troi_rite_vertic; + /* + * Back to "left rising": + */ + const gfloat prem_left_rising_squared = prem_left_rising * prem_left_rising; + const gfloat deux_left_rising_squared = deux_left_rising * deux_left_rising; + /* + * Back to "middle rising": + */ + const gfloat prem_midl_rising_squared = prem_midl_rising * prem_midl_rising; + const gfloat deux_midl_rising_squared = deux_midl_rising * deux_midl_rising; + const gfloat troi_midl_rising_squared = troi_midl_rising * troi_midl_rising; + /* + * Back to "right rising": + */ + const gfloat prem_rite_rising_squared = prem_rite_rising * prem_rite_rising; + const gfloat deux_rite_rising_squared = deux_rite_rising * deux_rite_rising; + /* + * Back to "left setting": + */ + const gfloat prem_left_settin_squared = prem_left_settin * prem_left_settin; + const gfloat deux_left_settin_squared = deux_left_settin * deux_left_settin; + /* + * Back to "middle setting": + */ + const gfloat prem_midl_settin_squared = prem_midl_settin * prem_midl_settin; + const gfloat deux_midl_settin_squared = deux_midl_settin * deux_midl_settin; + const gfloat troi_midl_settin_squared = troi_midl_settin * troi_midl_settin; + /* + * Back to "right setting": + */ + const gfloat prem_rite_settin_squared = prem_rite_settin * prem_rite_settin; + const gfloat deux_rite_settin_squared = deux_rite_settin * deux_rite_settin; + + /* + * "up": + */ + const gfloat prem___up_horizo_times_deux___up_horizo = + prem___up_horizo * deux___up_horizo; + const gfloat deux___up_horizo_times_troi___up_horizo = + deux___up_horizo * troi___up_horizo; + /* + * "down": + */ + const gfloat prem_down_horizo_times_deux_down_horizo = + prem_down_horizo * deux_down_horizo; + const gfloat deux_down_horizo_times_troi_down_horizo = + deux_down_horizo * troi_down_horizo; + /* + * "left": + */ + const gfloat prem_left_vertic_times_deux_left_vertic = + prem_left_vertic * deux_left_vertic; + const gfloat deux_left_vertic_times_troi_left_vertic = + deux_left_vertic * troi_left_vertic; + /* + * "right": + */ + const gfloat prem_rite_vertic_times_deux_rite_vertic = + prem_rite_vertic * deux_rite_vertic; + const gfloat deux_rite_vertic_times_troi_rite_vertic = + deux_rite_vertic * troi_rite_vertic; + /* + * "left rising": + */ + const gfloat prem_left_rising_times_deux_left_rising = + prem_left_rising * deux_left_rising; + /* + * "middle rising": + */ + const gfloat prem_midl_rising_times_deux_midl_rising = + prem_midl_rising * deux_midl_rising; + const gfloat deux_midl_rising_times_troi_midl_rising = + deux_midl_rising * troi_midl_rising; + /* + * "right rising": + */ + const gfloat prem_rite_rising_times_deux_rite_rising = + prem_rite_rising * deux_rite_rising; + /* + * "left setting": + */ + const gfloat prem_left_settin_times_deux_left_settin = + prem_left_settin * deux_left_settin; + /* + * "middle setting": + */ + const gfloat prem_midl_settin_times_deux_midl_settin = + prem_midl_settin * deux_midl_settin; + const gfloat deux_midl_settin_times_troi_midl_settin = + deux_midl_settin * troi_midl_settin; + /* + * "right setting": + */ + const gfloat prem_rite_settin_times_deux_rite_settin = + prem_rite_settin * deux_rite_settin; + + /* + * Branching parts of the computation of the YAFR correction (could + * be unbranched using arithmetic branching and C99 math intrinsics, + * although the compiler may be smart enough to remove the branching + * on its own): + */ + /* + * "up": + */ + const gfloat prem___up_horizo_vs_deux___up_horizo = + prem___up_horizo_squared < deux___up_horizo_squared + ? prem___up_horizo + : deux___up_horizo; + const gfloat deux___up_horizo_vs_troi___up_horizo = + deux___up_horizo_squared < troi___up_horizo_squared + ? deux___up_horizo + : troi___up_horizo; + /* + * "down": + */ + const gfloat prem_down_horizo_vs_deux_down_horizo = + prem_down_horizo_squared < deux_down_horizo_squared + ? prem_down_horizo + : deux_down_horizo; + const gfloat deux_down_horizo_vs_troi_down_horizo = + deux_down_horizo_squared < troi_down_horizo_squared + ? deux_down_horizo + : troi_down_horizo; + /* + * "left": + */ + const gfloat prem_left_vertic_vs_deux_left_vertic = + prem_left_vertic_squared < deux_left_vertic_squared + ? prem_left_vertic + : deux_left_vertic; + const gfloat deux_left_vertic_vs_troi_left_vertic = + deux_left_vertic_squared < troi_left_vertic_squared + ? deux_left_vertic + : troi_left_vertic; + /* + * "right": + */ + const gfloat prem_rite_vertic_vs_deux_rite_vertic = + prem_rite_vertic_squared < deux_rite_vertic_squared + ? prem_rite_vertic + : deux_rite_vertic; + const gfloat deux_rite_vertic_vs_troi_rite_vertic = + deux_rite_vertic_squared < troi_rite_vertic_squared + ? deux_rite_vertic + : troi_rite_vertic; + /* + * "left rising": + */ + const gfloat prem_left_rising_vs_deux_left_rising = + prem_left_rising_squared < deux_left_rising_squared + ? prem_left_rising + : deux_left_rising; + /* + * "middle rising": + */ + const gfloat prem_midl_rising_vs_deux_midl_rising = + prem_midl_rising_squared < deux_midl_rising_squared + ? prem_midl_rising + : deux_midl_rising; + const gfloat deux_midl_rising_vs_troi_midl_rising = + deux_midl_rising_squared < troi_midl_rising_squared + ? deux_midl_rising + : troi_midl_rising; + /* + * "right rising": + */ + const gfloat prem_rite_rising_vs_deux_rite_rising = + prem_rite_rising_squared < deux_rite_rising_squared + ? prem_rite_rising + : deux_rite_rising; + /* + * "left setting": + */ + const gfloat prem_left_settin_vs_deux_left_settin = + prem_left_settin_squared < deux_left_settin_squared + ? prem_left_settin + : deux_left_settin; + /* + * "middle setting": + */ + const gfloat prem_midl_settin_vs_deux_midl_settin = + prem_midl_settin_squared < deux_midl_settin_squared + ? prem_midl_settin + : deux_midl_settin; + const gfloat deux_midl_settin_vs_troi_midl_settin = + deux_midl_settin_squared < troi_midl_settin_squared + ? deux_midl_settin + : troi_midl_settin; + /* + * "right setting": + */ + const gfloat prem_rite_settin_vs_deux_rite_settin = + prem_rite_settin_squared < deux_rite_settin_squared + ? prem_rite_settin + : deux_rite_settin; + + /* + * Computation of the YAFR slopes. + */ + /* + * "up": + */ + const gfloat slope_horizo_left___up = + prem___up_horizo_times_deux___up_horizo < 0.f + ? 0.f + : prem___up_horizo_vs_deux___up_horizo; + const gfloat slope_horizo_rite___up = + deux___up_horizo_times_troi___up_horizo < 0.f + ? 0.f + : deux___up_horizo_vs_troi___up_horizo; + /* + * "down": + */ + const gfloat slope_horizo_left_down = + prem_down_horizo_times_deux_down_horizo < 0.f + ? 0.f + : prem_down_horizo_vs_deux_down_horizo; + const gfloat slope_horizo_rite_down = + deux_down_horizo_times_troi_down_horizo < 0.f + ? 0.f + : deux_down_horizo_vs_troi_down_horizo; + /* + * "left": + */ + const gfloat slope_vertic_left___up = + prem_left_vertic_times_deux_left_vertic < 0.f + ? 0.f + : prem_left_vertic_vs_deux_left_vertic; + const gfloat slope_vertic_left_down = + deux_left_vertic_times_troi_left_vertic < 0.f + ? 0.f + : deux_left_vertic_vs_troi_left_vertic; + /* + * "down": + */ + const gfloat slope_vertic_rite___up = + prem_rite_vertic_times_deux_rite_vertic < 0.f + ? 0.f + : prem_rite_vertic_vs_deux_rite_vertic; + const gfloat slope_vertic_rite_down = + deux_rite_vertic_times_troi_rite_vertic < 0.f + ? 0.f + : deux_rite_vertic_vs_troi_rite_vertic; + /* + * "left rising": + */ + const gfloat slope_rising_left___up = + prem_left_rising_times_deux_left_rising < 0.f + ? 0.f + : prem_left_rising_vs_deux_left_rising; + /* + * "middle rising": + */ + const gfloat slope_rising_left_down = + prem_midl_rising_times_deux_midl_rising < 0.f + ? 0.f + : prem_midl_rising_vs_deux_midl_rising; + const gfloat slope_rising_rite___up = + deux_midl_rising_times_troi_midl_rising < 0.f + ? 0.f + : deux_midl_rising_vs_troi_midl_rising; + /* + * "right rising": + */ + const gfloat slope_rising_rite_down = + prem_rite_rising_times_deux_rite_rising < 0.f + ? 0.f + : prem_rite_rising_vs_deux_rite_rising; + /* + * "left setting": + */ + const gfloat slope_settin_left_down = + prem_left_settin_times_deux_left_settin < 0.f + ? 0.f + : prem_left_settin_vs_deux_left_settin; + /* + * "middle setting": + */ + const gfloat slope_settin_left___up = + prem_midl_settin_times_deux_midl_settin < 0.f + ? 0.f + : prem_midl_settin_vs_deux_midl_settin; + const gfloat slope_settin_rite_down = + deux_midl_settin_times_troi_midl_settin < 0.f + ? 0.f + : deux_midl_settin_vs_troi_midl_settin; + /* + * "right setting": + */ + const gfloat slope_settin_rite___up = + prem_rite_settin_times_deux_rite_settin < 0.f + ? 0.f + : prem_rite_settin_vs_deux_rite_settin; + + /* + * Assemble the unweighted YAFR correction: + */ + const gfloat unweighted_yafr_correction = + c_horizo_left___up * slope_horizo_left___up + + + c_horizo_left_down * slope_horizo_left_down + + + c_horizo_rite___up * slope_horizo_rite___up + + + c_horizo_rite_down * slope_horizo_rite_down + + + c_vertic_left___up * slope_vertic_left___up + + + c_vertic_rite___up * slope_vertic_rite___up + + + c_vertic_left_down * slope_vertic_left_down + + + c_vertic_rite_down * slope_vertic_rite_down + + + c_settin_left___up * slope_settin_left___up + + + c_settin_rite___up * slope_settin_rite___up + + + c_settin_left_down * slope_settin_left_down + + + c_settin_rite_down * slope_settin_rite_down + + + c_rising_left___up * slope_rising_left___up + + + c_rising_rite___up * slope_rising_rite___up + + + c_rising_left_down * slope_rising_left_down + + + c_rising_rite_down * slope_rising_rite_down; + + /* + * Add the bilinear baseline and the weighted YAFR correction: + */ + const gfloat newval = sharpening * unweighted_yafr_correction + bilinear; + + *out = newval; +} + +static void +vips_interpolate_yafrnohalo_interpolate( VipsInterpolate *interpolate, + REGION *out, REGION *in, + int out_x, int out_y, double x, double y ) +{ + VipsInterpolateYafrnohalo *yafrnohalo = + VIPS_INTERPOLATE_YAFRNOHALO( interpolate ); + + /* + * Note: The computation is structured to foster software + * pipelining. + */ + + /* + * x is understood to increase from left to right, y, from top to + * bottom. Consequently, ix and iy are the indices of the pixel + * located at or to the left, and at or above. the sampling point. + * + * floor is used to make sure that the transition through 0 is + * smooth. If it is known that negative x and y will never be used, + * cast (which truncates) could be used instead. + */ + const gint ix = FLOOR (x); + const gint iy = FLOOR (y); + + /* + * Each (channel's) output pixel value is obtained by combining four + * "pieces," each piece corresponding to the set of points which are + * closest to the four pixels closest to the (x,y) position, pixel + * positions which have coordinates and labels as follows: + * + * (ix,iy) (ix+1,iy) + * =left__up =rite__up + * + * <- (x,y) is somewhere in the convex hull + * + * (ix,iy+1) (ix+1,iy+1) + * =left_dow =rite_dow + */ + /* + * rite_width is the width of the overlaps of the unit averaging box + * (which is centered at the position where an interpolated value is + * desired), with the closest unit pixel areas to the right. + * + * left_width is the width of the overlaps of the unit averaging box + * (which is centered at the position where an interpolated value is + * desired), with the closest unit pixel areas to the left. + */ + const gfloat rite_width = x - ix; + const gfloat down_hight = y - iy; + const gfloat left_width = 1.f - rite_width; + const gfloat up___hight = 1.f - down_hight; + + const gfloat rite_width_minus_half = rite_width - .5f; + const gfloat down_hight_minus_half = down_hight - .5f; + const gfloat left_width_minus_half = left_width - .5f; + const gfloat up___hight_minus_half = up___hight - .5f; + + const gfloat rite_rite_overlap = + rite_width_minus_half < 0.f ? 0.f : rite_width_minus_half; + const gfloat down_down_overlap = + down_hight_minus_half < 0.f ? 0.f : down_hight_minus_half; + /* + * Note that if rite_width_minus_half is (non)negative, + * left_width_minus_half is (non)positive. Consequently, the + * following two branching tests could be merged with the above two. + * The duplications are kept in the hope that the compiler translate + * them into "max" machine instructions: + */ + const gfloat left_left_overlap = + left_width_minus_half < 0.f ? 0.f : left_width_minus_half; + const gfloat up_____up_overlap = + up___hight_minus_half < 0.f ? 0.f : up___hight_minus_half; + /* + * The computation of quantities useful for the YAFR correction + * resumes after the computation of the needed bilinear baseline + * weights. + */ + + /* + * Bilinear (the plain vanilla standard bilinear) weights: + */ + const gfloat left_width_times___up_hight = left_width * up___hight; + const gfloat rite_width_times___up_hight = rite_width * up___hight; + const gfloat left_width_times_down_hight = left_width * down_hight; + const gfloat rite_width_times_down_hight = rite_width * down_hight; + + /* + * Back to quantities useful for the YAFR correction: + */ + /* + * Remaining overlap lengths: + */ + const gfloat rite_left_overlap = left_width - left_left_overlap; + const gfloat left_rite_overlap = rite_width - rite_rite_overlap; + const gfloat down___up_overlap = up___hight - up_____up_overlap; + const gfloat up___down_overlap = down_hight - down_down_overlap; + + const gfloat one_minus_left_left_overlap = 1.f - left_left_overlap; + const gfloat one_minus_rite_rite_overlap = 1.f - rite_rite_overlap; + const gfloat one_minus_up_____up_overlap = 1.f - up_____up_overlap; + const gfloat one_minus_down_down_overlap = 1.f - down_down_overlap; + const gfloat one_minus_rite_left_overlap = 1.f - rite_left_overlap; + const gfloat one_minus_left_rite_overlap = 1.f - left_rite_overlap; + const gfloat one_minus_down___up_overlap = 1.f - down___up_overlap; + const gfloat one_minus_up___down_overlap = 1.f - up___down_overlap; + + /* + * Recyclable products: + */ + const gfloat left_left_overlap_times_left_left_overlap = + left_left_overlap * left_left_overlap; + const gfloat rite_rite_overlap_times_rite_rite_overlap = + rite_rite_overlap * rite_rite_overlap; + const gfloat up_____up_overlap_times_up_____up_overlap = + up_____up_overlap * up_____up_overlap; + const gfloat down_down_overlap_times_down_down_overlap = + down_down_overlap * down_down_overlap; + + const gfloat rite_left_overlap_times_one_minus_rite_left_overlap = + rite_left_overlap * one_minus_rite_left_overlap; + const gfloat left_rite_overlap_times_one_minus_left_rite_overlap = + left_rite_overlap * one_minus_left_rite_overlap; + const gfloat down___up_overlap_times_one_minus_down___up_overlap = + down___up_overlap * one_minus_down___up_overlap; + const gfloat up___down_overlap_times_one_minus_up___down_overlap = + up___down_overlap * one_minus_up___down_overlap; + + /* + * "Cardinal" contributions of the various YAFR slopes. Each of the + * (up to) four contributing pixels contributes four slopes: one + * horizontal slope, one vertical slope, and two diagonal slopes + * (one "rising", and one "setting"). Consequently, sixteen + * "cardinal" basis contributions need to be computed. + */ + /* + * "Cardinal" contributions of the horizontal and vertical slopes: + */ + /* + * Common factors: + */ + const gfloat horizo___up = + up_____up_overlap /* height of the overlap */ + * + one_minus_up_____up_overlap /* bilinear coefficient */ + + /* second interpolated value */ + down___up_overlap /* height of the overlap */ + * + down___up_overlap; /* bilinear coefficient */ + const gfloat horizo_down = + down_down_overlap /* height of the overlap */ + * + one_minus_down_down_overlap /* bilinear coefficient */ + + /* second interpolated value */ + up___down_overlap /* height of the overlap */ + * + up___down_overlap; /* bilinear coefficient */ + const gfloat horizo_left = + rite_left_overlap_times_one_minus_rite_left_overlap + - + left_left_overlap_times_left_left_overlap; + const gfloat horizo_rite = + rite_rite_overlap_times_rite_rite_overlap + - + left_rite_overlap_times_one_minus_left_rite_overlap; + const gfloat vertic_left = + left_left_overlap + * + one_minus_left_left_overlap + + + rite_left_overlap + * + rite_left_overlap; + const gfloat vertic_rite = + rite_rite_overlap + * + one_minus_rite_rite_overlap + + + left_rite_overlap + * + left_rite_overlap; + const gfloat vertic___up = + down___up_overlap_times_one_minus_down___up_overlap + - + up_____up_overlap_times_up_____up_overlap; + const gfloat vertic_down = + down_down_overlap_times_down_down_overlap + - + up___down_overlap_times_one_minus_up___down_overlap; + /* + * "Cardinal" contribution of the left horizontal slopes: + */ + const gfloat c_horizo_left___up = horizo_left * horizo___up; + const gfloat c_horizo_left_down = horizo_left * horizo_down; + /* + * "Cardinal" contribution of the right horizontal slopes: + */ + const gfloat c_horizo_rite___up = horizo_rite * horizo___up; + const gfloat c_horizo_rite_down = horizo_rite * horizo_down; + /* + * "Cardinal" contribution of the up vertical slopes: + */ + const gfloat c_vertic_left___up = vertic___up * vertic_left; + const gfloat c_vertic_rite___up = vertic___up * vertic_rite; + /* + * "Cardinal" contribution of the down vertical slopes: + */ + const gfloat c_vertic_left_down = vertic_down * vertic_left; + const gfloat c_vertic_rite_down = vertic_down * vertic_rite; + + /* + * "Cardinal" contributions of the diagonal slopes: + */ + const gfloat c_settin_left___up = + rite_left_overlap_times_one_minus_rite_left_overlap + * + down___up_overlap_times_one_minus_down___up_overlap + - + left_left_overlap_times_left_left_overlap + * + up_____up_overlap_times_up_____up_overlap; + const gfloat c_settin_rite___up = + rite_rite_overlap_times_rite_rite_overlap + * + down___up_overlap_times_one_minus_down___up_overlap + - + left_rite_overlap_times_one_minus_left_rite_overlap + * + up_____up_overlap_times_up_____up_overlap; + const gfloat c_settin_left_down = + rite_left_overlap_times_one_minus_rite_left_overlap + * + down_down_overlap_times_down_down_overlap + - + left_left_overlap_times_left_left_overlap + * + up___down_overlap_times_one_minus_up___down_overlap; + const gfloat c_settin_rite_down = + rite_rite_overlap_times_rite_rite_overlap + * + down_down_overlap_times_down_down_overlap + - + left_rite_overlap_times_one_minus_left_rite_overlap + * + up___down_overlap_times_one_minus_up___down_overlap; + const gfloat c_rising_left___up = + rite_left_overlap_times_one_minus_rite_left_overlap + * + up_____up_overlap_times_up_____up_overlap + - + left_left_overlap_times_left_left_overlap + * + down___up_overlap_times_one_minus_down___up_overlap; + const gfloat c_rising_rite___up = + rite_rite_overlap_times_rite_rite_overlap + * + up_____up_overlap_times_up_____up_overlap + - + left_rite_overlap_times_one_minus_left_rite_overlap + * + down___up_overlap_times_one_minus_down___up_overlap; + const gfloat c_rising_left_down = + rite_left_overlap_times_one_minus_rite_left_overlap + * + up___down_overlap_times_one_minus_up___down_overlap + - + left_left_overlap_times_left_left_overlap + * + down_down_overlap_times_down_down_overlap; + const gfloat c_rising_rite_down = + rite_rite_overlap_times_rite_rite_overlap + * + up___down_overlap_times_one_minus_up___down_overlap + - + left_rite_overlap_times_one_minus_left_rite_overlap + * + down_down_overlap_times_down_down_overlap; + + /* + * Set the tile pointer to the first relevant value. Since the + * pointer initially points to dos_two, we need to rewind it one + * tile row, then go back one additional pixel. + */ + const PEL *p = (PEL *) IM_REGION_ADDR( in, ix - 1, iy - 1 ); + + /* Pel size and line size. + */ + const int channels = in->im->Bands; + const int pixels_per_buffer_row = + IM_REGION_LSKIP( in ) / (sizeof( float ) * channels); + + /* Where we write the result. + */ + PEL *q = (PEL *) IM_REGION_ADDR( out, out_x, out_y ); + int z; + + for( z = 0; z < channels; z++ ) + bilinear_yafrnohalo ((float *) q + z, (float *) p + z, + channels, pixels_per_buffer_row, + yafrnohalo->sharpening, + c_horizo_left___up, + c_horizo_left_down, + c_horizo_rite___up, + c_horizo_rite_down, + c_vertic_left___up, + c_vertic_rite___up, + c_vertic_left_down, + c_vertic_rite_down, + c_settin_left___up, + c_settin_rite___up, + c_settin_left_down, + c_settin_rite_down, + c_rising_left___up, + c_rising_rite___up, + c_rising_left_down, + c_rising_rite_down, + left_width_times___up_hight, + rite_width_times___up_hight, + left_width_times_down_hight, + rite_width_times_down_hight); +} + +static void +vips_interpolate_yafrnohalo_class_init( VipsInterpolateYafrnohaloClass *class ) +{ + VipsInterpolateClass *interpolate_class = + VIPS_INTERPOLATE_CLASS( class ); + + vips_interpolate_yafrnohalo_parent_class = + g_type_class_peek_parent( class ); + + interpolate_class->interpolate = + vips_interpolate_yafrnohalo_interpolate; + interpolate_class->window_size = 4; +} + +static void +vips_interpolate_yafrnohalo_init( VipsInterpolateYafrnohalo *yafrnohalo ) +{ +#ifdef DEBUG + printf( "vips_interpolate_yafrnohalo_init: " ); + vips_object_print( VIPS_OBJECT( yafrnohalo ) ); +#endif /*DEBUG*/ + + yafrnohalo->sharpening = 1.0; +} + +GType +vips_interpolate_yafrnohalo_get_type( void ) +{ + static GType type = 0; + + if( !type ) { + static const GTypeInfo info = { + sizeof( VipsInterpolateYafrnohaloClass ), + NULL, /* base_init */ + NULL, /* base_finalize */ + (GClassInitFunc) vips_interpolate_yafrnohalo_class_init, + NULL, /* class_finalize */ + NULL, /* class_data */ + sizeof( VipsInterpolateYafrnohalo ), + 32, /* n_preallocs */ + (GInstanceInitFunc) vips_interpolate_yafrnohalo_init, + }; + + type = g_type_register_static( VIPS_TYPE_INTERPOLATE, + "VipsInterpolateYafrnohalo", &info, 0 ); + } + + return( type ); +} + +VipsInterpolate * +vips_interpolate_yafrnohalo_new( void ) +{ + return( VIPS_INTERPOLATE( g_object_new( + VIPS_TYPE_INTERPOLATE_YAFRNOHALO, NULL ) ) ); +} + +void +vips_interpolate_yafrnohalo_set_sharpening( + VipsInterpolateYafrnohalo *yafrnohalo, + double sharpening ) +{ + yafrnohalo->sharpening = sharpening; +} + +/* Convenience: return a static yafrnohalo you don't need to free. + */ +VipsInterpolate * +vips_interpolate_yafrnohalo_static( void ) +{ + static VipsInterpolate *interpolate = NULL; + + if( !interpolate ) + interpolate = vips_interpolate_yafrnohalo_new(); + + return( interpolate ); +} + +