granicus.if.org Git - imagemagick/blob - magick/resize.c

   1 /*
   2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   3 %                                                                             %
   4 %                                                                             %
   5 %                                                                             %
   6 %                 RRRR   EEEEE  SSSSS  IIIII  ZZZZZ  EEEEE                    %
   7 %                 R   R  E      SS       I       ZZ  E                        %
   8 %                 RRRR   EEE     SSS     I     ZZZ   EEE                      %
   9 %                 R R    E         SS    I    ZZ     E                        %
  10 %                 R  R   EEEEE  SSSSS  IIIII  ZZZZZ  EEEEE                    %
  11 %                                                                             %
  12 %                                                                             %
  13 %                      MagickCore Image Resize Methods                        %
  14 %                                                                             %
  15 %                              Software Design                                %
  16 %                                John Cristy                                  %
  17 %                                 July 1992                                   %
  18 %                                                                             %
  19 %                                                                             %
  20 %  Copyright 1999-2010 ImageMagick Studio LLC, a non-profit organization      %
  21 %  dedicated to making software imaging solutions freely available.           %
  22 %                                                                             %
  23 %  You may not use this file except in compliance with the License.  You may  %
  24 %  obtain a copy of the License at                                            %
  25 %                                                                             %
  26 %    http://www.imagemagick.org/script/license.php                            %
  27 %                                                                             %
  28 %  Unless required by applicable law or agreed to in writing, software        %
  29 %  distributed under the License is distributed on an "AS IS" BASIS,          %
  30 %  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.   %
  31 %  See the License for the specific language governing permissions and        %
  32 %  limitations under the License.                                             %
  33 %                                                                             %
  34 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  35 %
  36 %
  37 */
  38 \f
  39 /*
  40   Include declarations.
  41 */
  42 #include "magick/studio.h"
  43 #include "magick/artifact.h"
  44 #include "magick/blob.h"
  45 #include "magick/cache.h"
  46 #include "magick/cache-view.h"
  47 #include "magick/color.h"
  48 #include "magick/color-private.h"
  49 #include "magick/draw.h"
  50 #include "magick/exception.h"
  51 #include "magick/exception-private.h"
  52 #include "magick/gem.h"
  53 #include "magick/image.h"
  54 #include "magick/image-private.h"
  55 #include "magick/list.h"
  56 #include "magick/memory_.h"
  57 #include "magick/magick.h"
  58 #include "magick/pixel-private.h"
  59 #include "magick/property.h"
  60 #include "magick/monitor.h"
  61 #include "magick/monitor-private.h"
  62 #include "magick/pixel.h"
  63 #include "magick/option.h"
  64 #include "magick/resample.h"
  65 #include "magick/resize.h"
  66 #include "magick/resize-private.h"
  67 #include "magick/string_.h"
  68 #include "magick/string-private.h"
  69 #include "magick/thread-private.h"
  70 #include "magick/utility.h"
  71 #include "magick/version.h"
  72 #if defined(MAGICKCORE_LQR_DELEGATE)
  73 #include <lqr.h>
  74 #endif
  75 \f
  76 /*
  77   Typedef declarations.
  78 */
  79 struct _ResizeFilter
  80 {
  81   MagickRealType
  82     (*filter)(const MagickRealType,const ResizeFilter *),
  83     (*window)(const MagickRealType,const ResizeFilter *),
  84     support,        /* filter region of support - the filter support limit */
  85     window_support, /* window support, usally equal to support (expert only) */
  86     scale,          /* dimension scaling to fit window support (usally 1.0) */
  87     blur,           /* x-scale (blur-sharpen) */
  88     cubic[8];       /* cubic coefficents for smooth Cubic filters */
  89
  90   size_t
  91     signature;
  92 };
  93 \f
  94 /*
  95   Forward declaractions.
  96 */
  97 static MagickRealType
  98   I0(MagickRealType x),
  99   BesselOrderOne(MagickRealType),
 100   Sinc(const MagickRealType, const ResizeFilter *),
 101   SincFast(const MagickRealType, const ResizeFilter *);
 102 \f
 103 /*
 104 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 105 %                                                                             %
 106 %                                                                             %
 107 %                                                                             %
 108 +   F i l t e r F u n c t i o n s                                             %
 109 %                                                                             %
 110 %                                                                             %
 111 %                                                                             %
 112 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 113 %
 114 %  These are the various filter and windowing functions that are provided.
 115 %
 116 %  They are internal to this module only.  See AcquireResizeFilterInfo() for
 117 %  details of the access to these functions, via the GetResizeFilterSupport()
 118 %  and GetResizeFilterWeight() API interface.
 119 %
 120 %  The individual filter functions have this format...
 121 %
 122 %     static MagickRealtype *FilterName(const MagickRealType x,
 123 %        const MagickRealType support)
 124 %
 125 %  A description of each parameter follows:
 126 %
 127 %    o x: the distance from the sampling point generally in the range of 0 to
 128 %      support.  The GetResizeFilterWeight() ensures this a positive value.
 129 %
 130 %    o resize_filter: current filter information.  This allows function to
 131 %      access support, and possibly other pre-calculated information defining
 132 %      the functions.
 133 %
 134 */
 135
 136 #define MagickPIL ((MagickRealType) 3.14159265358979323846264338327950288420L)
 137
 138 static MagickRealType Jinc(const MagickRealType x,
 139   const ResizeFilter *magick_unused(resize_filter))
 140 {
 141   /*
 142     See Pratt "Digital Image Processing" p.97 for Jinc/Bessel functions.
 143     http://mathworld.wolfram.com/JincFunction.html and page 11 of
 144     http://www.ph.ed.ac.uk/%7ewjh/teaching/mo/slides/lens/lens.pdf
 145
 146     The original "zoom" program by Paul Heckbert called this "Bessel"
 147     But really its is more accuritally named "Jinc".
 148   */
 149   if (x == 0.0)
 150     return(0.25*MagickPIL);
 151   return(BesselOrderOne(MagickPIL*x)/(x+x));
 152 }
 153
 154 static MagickRealType Blackman(const MagickRealType x,
 155   const ResizeFilter *magick_unused(resize_filter))
 156 {
 157   /*
 158     Blackman: 2nd order cosine windowing function:
 159       0.42 + 0.5 cos(pi x) + 0.08 cos(2pi x)
 160     Refactored by Chantal Racette and Nicolas Robidoux to one trig
 161     call and five flops.
 162   */
 163   const MagickRealType cospix = cos((double) (MagickPIL*x));
 164   return(0.34+cospix*(0.5+cospix*0.16));
 165 }
 166
 167 static MagickRealType Bohman(const MagickRealType x,
 168   const ResizeFilter *magick_unused(resize_filter))
 169 {
 170   /*
 171     Bohman: 2rd Order cosine windowing function:
 172       (1-x) cos(pi x) + sin(pi x) / pi.
 173     Refactored by Nicolas Robidoux to one trig call, one sqrt call,
 174     and 7 flops, taking advantage of the fact that the support of
 175     Bohman is 1 (so that we know that sin(pi x) >= 0).
 176   */
 177   const double cospix = cos((double) (MagickPIL*x));
 178   const double sinpix = sqrt(1.0-cospix*cospix);
 179   return((1.0-x)*cospix+(1.0/MagickPIL)*sinpix);
 180 }
 181
 182 static MagickRealType Box(const MagickRealType magick_unused(x),
 183   const ResizeFilter *magick_unused(resize_filter))
 184 {
 185   /*
 186     A Box filter is a equal weighting function (all weights equal).
 187     DO NOT LIMIT results by support or resize point sampling will work
 188     as it requests points beyond its normal 0.0 support size.
 189   */
 190   return(1.0);
 191 }
 192
 193 static MagickRealType CubicBC(const MagickRealType x,
 194   const ResizeFilter *resize_filter)
 195 {
 196   /*
 197     Cubic Filters using B,C determined values:
 198        Mitchell-Netravali  B=1/3 C=1/3   Qualitively ideal Cubic Filter
 199        Catmull-Rom         B= 0  C=1/2   Cublic Interpolation Function
 200        Cubic B-Spline      B= 1  C= 0    Spline Approximation of Gaussian
 201        Hermite             B= 0  C= 0    Quadratic Spline (support = 1)
 202
 203     See paper by Mitchell and Netravali, Reconstruction Filters in Computer
 204     Graphics Computer Graphics, Volume 22, Number 4, August 1988
 205     http://www.cs.utexas.edu/users/fussell/courses/cs384g/lectures/mitchell/
 206     Mitchell.pdf.
 207
 208     Coefficents are determined from B,C values:
 209        P0 = (  6 - 2*B       )/6
 210        P1 =         0
 211        P2 = (-18 +12*B + 6*C )/6
 212        P3 = ( 12 - 9*B - 6*C )/6
 213        Q0 = (      8*B +24*C )/6
 214        Q1 = (    -12*B -48*C )/6
 215        Q2 = (      6*B +30*C )/6
 216        Q3 = (    - 1*B - 6*C )/6
 217
 218     which are used to define the filter:
 219
 220        P0 + P1*x + P2*x^2 + P3*x^3      0 <= x < 1
 221        Q0 + Q1*x + Q2*x^2 + Q3*x^3      1 <= x <= 2
 222
 223     which ensures function is continuous in value and derivative (slope).
 224   */
 225   if (x < 1.0)
 226     return(resize_filter->cubic[0]+x*(resize_filter->cubic[1]+x*
 227       (resize_filter->cubic[2]+x*resize_filter->cubic[3])));
 228   if (x < 2.0)
 229     return(resize_filter->cubic[4]+x*(resize_filter->cubic[5]+x*
 230       (resize_filter->cubic[6]+x*resize_filter->cubic[7])));
 231   return(0.0);
 232 }
 233
 234 static MagickRealType Gaussian(const MagickRealType x,
 235   const ResizeFilter *magick_unused(resize_filter))
 236 {
 237   /*
 238     1D Gaussian with sigma=1/2:
 239       exp(-2 x^2)/sqrt(pi/2))
 240   */
 241   /*const MagickRealType alpha = (MagickRealType) (2.0/MagickSQ2PI);*/
 242   return(exp((double) (-2.0*x*x)));
 243 }
 244
 245 static MagickRealType Hanning(const MagickRealType x,
 246   const ResizeFilter *magick_unused(resize_filter))
 247 {
 248   /*
 249     Cosine window function:
 250       .5+.5cos(pi x).
 251   */
 252   const MagickRealType cospix = cos((double) (MagickPIL*x));
 253   return(0.5+0.5*cospix);
 254 }
 255
 256 static MagickRealType Hamming(const MagickRealType x,
 257   const ResizeFilter *magick_unused(resize_filter))
 258 {
 259   /*
 260     Offset cosine window function:
 261      .54 + .46 cos(pi x).
 262   */
 263   const MagickRealType cospix = cos((double) (MagickPIL*x));
 264   return(0.54+0.46*cospix);
 265 }
 266
 267 static MagickRealType Kaiser(const MagickRealType x,
 268   const ResizeFilter *magick_unused(resize_filter))
 269 {
 270 #define Alpha  6.5
 271 #define I0A  (1.0/I0(Alpha))
 272
 273   /*
 274     Kaiser Windowing Function (bessel windowing): Alpha is a free
 275     value from 5 to 8 (currently hardcoded to 6.5).
 276     Future: make alpha the IOA pre-calculation, an 'expert' setting.
 277   */
 278   return(I0A*I0(Alpha*sqrt((double) (1.0-x*x))));
 279 }
 280
 281 static MagickRealType Lagrange(const MagickRealType x,
 282   const ResizeFilter *resize_filter)
 283 {
 284   MagickRealType
 285     value;
 286
 287   register ssize_t
 288     i;
 289
 290   ssize_t
 291     n,
 292     order;
 293
 294   /*
 295     Lagrange piecewise polynomial fit of sinc: N is the 'order' of the
 296     lagrange function and depends on the overall support window size
 297     of the filter. That is: for a support of 2, it gives a lagrange-4
 298     (piecewise cubic function).
 299
 300     "n" identifies the piece of the piecewise polynomial.
 301
 302     See Survey: Interpolation Methods, IEEE Transactions on Medical
 303     Imaging, Vol 18, No 11, November 1999, p1049-1075, -- Equation 27
 304     on p1064.
 305   */
 306   if (x > resize_filter->support)
 307     return(0.0);
 308   order=(ssize_t) (2.0*resize_filter->window_support);  /* number of pieces */
 309   /*n=(ssize_t)((1.0*order)/2.0+x);   --  which piece does x belong to */
 310   n = (ssize_t)(resize_filter->window_support + x);
 311   value=1.0f;
 312   for (i=0; i < order; i++)
 313     if (i != n)
 314       value*=(n-i-x)/(n-i);
 315   return(value);
 316 }
 317
 318 static MagickRealType Quadratic(const MagickRealType x,
 319   const ResizeFilter *magick_unused(resize_filter))
 320 {
 321   /*
 322     2rd order (quadratic) B-Spline approximation of Gaussian.
 323   */
 324   if (x < 0.5)
 325     return(0.75-x*x);
 326   if (x < 1.5)
 327     return(0.5*(x-1.5)*(x-1.5));
 328   return(0.0);
 329 }
 330
 331 static MagickRealType Sinc(const MagickRealType x,
 332   const ResizeFilter *magick_unused(resize_filter))
 333 {
 334   /*
 335     Scaled sinc(x) function using a trig call:
 336       sinc(x) == sin(pi x)/(pi x).
 337   */
 338   if (x != 0.0)
 339   {
 340     const MagickRealType pix = (MagickRealType) (MagickPIL*x);
 341     return(sin((double) pix)/pix);
 342   }
 343   return((MagickRealType) 1.0);
 344 }
 345
 346 static MagickRealType SincFast(const MagickRealType x,
 347   const ResizeFilter *magick_unused(resize_filter))
 348 {
 349   /*
 350     Approximations of the sinc function sin(pi x)/(pi x) over the
 351     interval [-4,4] constructed by Nicolas Robidoux and Chantal
 352     Racette with funding from the Natural Sciences and Engineering
 353     Research Council of Canada.
 354
 355     Although the approximations are polynomials (for low order of
 356     approximation) and quotients of polynomials (for higher order of
 357     approximation) and consequently are similar in form to Taylor
 358     polynomials/Pade approximants, the approximations are computed
 359     with a completely different technique.
 360
 361     Summary: These approximations are "the best" in terms of bang
 362     (accuracy) for the buck (flops). More specifically: Among the
 363     polynomial quotients that can be computed using a fixed number of
 364     flops (with a given "+ - * / budget"), the chosen polynomial
 365     quotient is the one closest to the approximated function with
 366     respect to maximum absolute relative error over the given
 367     interval.
 368
 369     The Remez algorithm, as implemented in the boost library's minimax
 370     package, is the key to the construction:
 371     http://www.boost.org/doc/libs/1_36_0/libs/math/doc/...
 372     ...sf_and_dist/html/math_toolkit/backgrounders/remez.html
 373   */
 374   /*
 375     If outside of the interval of approximation, use the standard trig
 376     formula.
 377   */
 378   if (x > 4.0)
 379     {
 380       const MagickRealType pix = (MagickRealType) (MagickPIL*x);
 381       return(sin((double) pix)/pix);
 382     }
 383   {
 384     /*
 385       The approximations only depend on x^2 (sinc is an even
 386       function).
 387     */
 388     const MagickRealType xx = x*x;
 389 #if MAGICKCORE_QUANTUM_DEPTH <= 8
 390     /*
 391       Maximum absolute relative error 6.3e-6 < 1/2^17.
 392     */
 393     const MagickRealType c0 = 0.173610016489197553621906385078711564924e-2L;
 394     const MagickRealType c1 = -0.384186115075660162081071290162149315834e-3L;
 395     const MagickRealType c2 = 0.393684603287860108352720146121813443561e-4L;
 396     const MagickRealType c3 = -0.248947210682259168029030370205389323899e-5L;
 397     const MagickRealType c4 = 0.107791837839662283066379987646635416692e-6L;
 398     const MagickRealType c5 = -0.324874073895735800961260474028013982211e-8L;
 399     const MagickRealType c6 = 0.628155216606695311524920882748052490116e-10L;
 400     const MagickRealType c7 = -0.586110644039348333520104379959307242711e-12L;
 401     const MagickRealType p =
 402       c0+xx*(c1+xx*(c2+xx*(c3+xx*(c4+xx*(c5+xx*(c6+xx*c7))))));
 403     return((xx-1.0)*(xx-4.0)*(xx-9.0)*(xx-16.0)*p);
 404 #elif MAGICKCORE_QUANTUM_DEPTH <= 16
 405     /*
 406       Max. abs. rel. error 2.2e-8 < 1/2^25.
 407     */
 408     const MagickRealType c0 = 0.173611107357320220183368594093166520811e-2L;
 409     const MagickRealType c1 = -0.384240921114946632192116762889211361285e-3L;
 410     const MagickRealType c2 = 0.394201182359318128221229891724947048771e-4L;
 411     const MagickRealType c3 = -0.250963301609117217660068889165550534856e-5L;
 412     const MagickRealType c4 = 0.111902032818095784414237782071368805120e-6L;
 413     const MagickRealType c5 = -0.372895101408779549368465614321137048875e-8L;
 414     const MagickRealType c6 = 0.957694196677572570319816780188718518330e-10L;
 415     const MagickRealType c7 = -0.187208577776590710853865174371617338991e-11L;
 416     const MagickRealType c8 = 0.253524321426864752676094495396308636823e-13L;
 417     const MagickRealType c9 = -0.177084805010701112639035485248501049364e-15L;
 418     const MagickRealType p =
 419       c0+xx*(c1+xx*(c2+xx*(c3+xx*(c4+xx*(c5+xx*(c6+xx*(c7+xx*(c8+xx*c9))))))));
 420     return((xx-1.0)*(xx-4.0)*(xx-9.0)*(xx-16.0)*p);
 421 #else
 422     /*
 423       Max. abs. rel. error 1.2e-12 < 1/2^39.
 424     */
 425     const MagickRealType c0 = 0.173611111110910715186413700076827593074e-2L;
 426     const MagickRealType c1 = -0.289105544717893415815859968653611245425e-3L;
 427     const MagickRealType c2 = 0.206952161241815727624413291940849294025e-4L;
 428     const MagickRealType c3 = -0.834446180169727178193268528095341741698e-6L;
 429     const MagickRealType c4 = 0.207010104171026718629622453275917944941e-7L;
 430     const MagickRealType c5 = -0.319724784938507108101517564300855542655e-9L;
 431     const MagickRealType c6 = 0.288101675249103266147006509214934493930e-11L;
 432     const MagickRealType c7 = -0.118218971804934245819960233886876537953e-13L;
 433     const MagickRealType p =
 434       c0+xx*(c1+xx*(c2+xx*(c3+xx*(c4+xx*(c5+xx*(c6+xx*c7))))));
 435     const MagickRealType d0 = 1.0L;
 436     const MagickRealType d1 = 0.547981619622284827495856984100563583948e-1L;
 437     const MagickRealType d2 = 0.134226268835357312626304688047086921806e-2L;
 438     const MagickRealType d3 = 0.178994697503371051002463656833597608689e-4L;
 439     const MagickRealType d4 = 0.114633394140438168641246022557689759090e-6L;
 440     const MagickRealType q = d0+xx*(d1+xx*(d2+xx*(d3+xx*d4)));
 441     return((xx-1.0)*(xx-4.0)*(xx-9.0)*(xx-16.0)/q*p);
 442 #endif
 443   }
 444 }
 445
 446 static MagickRealType Triangle(const MagickRealType x,
 447   const ResizeFilter *magick_unused(resize_filter))
 448 {
 449   /*
 450     1st order (linear) B-Spline, bilinear interpolation, Tent 1D
 451     filter, or a Bartlett 2D Cone filter.
 452   */
 453   if (x < 1.0)
 454     return(1.0-x);
 455   return(0.0);
 456 }
 457
 458 static MagickRealType Welsh(const MagickRealType x,
 459   const ResizeFilter *magick_unused(resize_filter))
 460 {
 461   /*
 462     Welsh parabolic windowing filter.
 463   */
 464   if (x < 1.0)
 465     return(1.0-x*x);
 466   return(0.0);
 467 }
 468 \f
 469 /*
 470 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 471 %                                                                             %
 472 %                                                                             %
 473 %                                                                             %
 474 +   A c q u i r e R e s i z e F i l t e r                                     %
 475 %                                                                             %
 476 %                                                                             %
 477 %                                                                             %
 478 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 479 %
 480 %  AcquireResizeFilter() allocates the ResizeFilter structure.  Choose
 481 %  from these filters:
 482 %
 483 %  FIR (Finite impulse Response) Filters
 484 %      Box         Triangle   Quadratic
 485 %      Cubic       Hermite    Catrom
 486 %      Mitchell
 487 %
 488 %  IIR (Infinite impulse Response) Filters
 489 %      Gaussian     Sinc        Jinc (Bessel)
 490 %
 491 %  Windowed Sinc/Jinc Filters
 492 %      Blackman     Hanning     Hamming
 493 %      Kaiser       Lanczos
 494 %
 495 %  Special purpose Filters
 496 %      SincFast  Lanczos2D
 497 %
 498 %  The users "-filter" selection is used to lookup the default 'expert'
 499 %  settings for that filter from a internal table.  However any provided
 500 %  'expert' settings (see below) may override this selection.
 501 %
 502 %  FIR filters are used as is, and are limited to that filters support
 503 %  window (unless over-ridden).  'Gaussian' while classed as an IIR
 504 %  filter, is also simply clipped by its support size (currently 1.5
 505 %  ro approximatally 3*sigma as recommended by many references)
 506 %
 507 %  The selection is typically either a windowed Sinc, or interpolated
 508 %  filter, for use by functions such as ResizeImage().  However if a
 509 %  'cylindrical' filter flag is requested, any default Sinc weighting
 510 %  and windowing functions will be promoted to cylindrical Jinc form of
 511 %  function.
 512 %
 513 %  Directly requesting 'Sinc' or 'Jinc' will force the use of that
 514 %  filter function without any windowing. This is not recommended,
 515 %  except by image processing experts or in expert options.  Selecting a
 516 %  window filtering version of these functions is better.
 517 %
 518 %  Lanczos is a special case of a Sinc-windowed Sinc, (or Jinc-Jinc for
 519 %  the cylindrical case) but defaulting to 3-lobe support, rather that
 520 %  the default 4 lobe support of the other windowed sinc/jinc filters.
 521 %
 522 %  Two forms of the 'Sinc' function are available: Sinc and SincFast.
 523 %  Sinc is computed using the traditional sin(pi*x)/(pi*x); it is
 524 %  selected if the user specifically specifies the use of a Sinc
 525 %  filter. SincFast uses highly accurate (and fast) polynomial (low Q)
 526 %  and rational (high Q) approximations, and will be used by default in
 527 %  most cases.
 528 %
 529 %  The Lanczos2D filter is just a normal 2-lobed Lanczos filter. But if
 530 %  selected as a cylindrical (radial) filter, the 2-lobed Jinc windowed
 531 %  Jinc will be modified by a blur factor to negate the effect of windowing
 532 %  the Jinc function, and greatly reduce resulting blur.
 533 %
 534 %  Special 'expert' options can be used to override any and all filter
 535 %  settings. This is not advised unless you have expert knowledge of
 536 %  the use of resampling filtered techniques.  Check on the results of
 537 %  your selections using the "filter:verbose" setting to make sure you
 538 %  get the exact filter that you are tring to achieve.
 539 %
 540 %    "filter:filter"    Select the main function associated with
 541 %        this filter name, as the weighting function of the filter.
 542 %        This can be used to set a windowing function as a weighting
 543 %        function, for special purposes, such as graphing.
 544 %
 545 %        If a "filter:window" operation has not been provided, then a 'Box'
 546 %        windowing function will be set to denote that no windowing function
 547 %        is being used.
 548 %
 549 %    "filter:window"   Select this windowing function for the filter.
 550 %        While any filter could be used as a windowing function, using the
 551 %        'first lobe' of that filter over the whole support window, using a
 552 %        non-windowing function is not advisible. If no weighting filter
 553 %        function is specifed a 'SincFast' filter will be used.
 554 %
 555 %    "filter:lobes"    Number of lobes to use for the Sinc/Jinc filter.
 556 %        This a simpler method of setting filter support size that will
 557 %        correctly handle the Sinc/Jinc switch for an operators filtering
 558 %        requirements.  Only integers should be given.
 559 %
 560 %    "filter:support"  Set the support size for filtering to the size given
 561 %        This not recommended for Sinc/Jinc windowed filters (lobes should
 562 %        be used instead).  This will override any 'filter:lobes' option.
 563 %
 564 %    "filter:win-support"  Scale windowing function to this size instead.
 565 %        This causes the windowing (or self-windowing Lagrange filter) to act
 566 %        is if the support window it much much larger than what is actually
 567 %        supplied to the calling operator.  The filter however is still
 568 %        clipped to the real support size given, by the support range suppiled
 569 %        to the caller.  If unset this will equal the normal filter support
 570 %        size.
 571 %
 572 %    "filter:blur"     Scale the filter and support window by this amount.
 573 %        A value >1 will generally result in a more burred image with
 574 %        more ringing effects, while a value <1 will sharpen the
 575 %        resulting image with more aliasing and Morie effects.
 576 %
 577 %    "filter:b"
 578 %    "filter:c"    Override the preset B,C values for a Cubic type of filter
 579 %         If only one of these are given it is assumes to be a 'Keys'
 580 %         type of filter such that B+2C=1, where Keys 'alpha' value = C
 581 %
 582 %    "filter:verbose"   Output the exact results of the filter selections
 583 %         made, as well as plotting data for graphing the resulting filter
 584 %         over support range (blur adjusted).
 585 %
 586 %  Set a true un-windowed Sinc filter with 10 lobes (very slow)
 587 %     -define filter:filter=Sinc
 588 %     -define filter:lobes=8
 589 %
 590 %  For example force an 8 lobe Lanczos (Sinc or Jinc) filter...
 591 %     -filter Lanczos
 592 %     -define filter:lobes=8
 593 %
 594 %  The format of the AcquireResizeFilter method is:
 595 %
 596 %      ResizeFilter *AcquireResizeFilter(const Image *image,
 597 %        const FilterTypes filter_type, const MagickBooleanType radial,
 598 %        ExceptionInfo *exception)
 599 %
 600 %  A description of each parameter follows:
 601 %
 602 %    o image: the image.
 603 %
 604 %    o filter: the filter type, defining a preset filter, window and
 605 %      support.  The artifact settings listed above will override
 606 %      those selections.
 607 %
 608 %    o blur: blur the filter by this amount, use 1.0 if unknown.  Image
 609 %      artifact "filter:blur" will override this API call usage, including
 610 %      any internal change (such as for cylindrical usage).
 611 %
 612 %    o radial: use a 1D orthogonal filter (Sinc) or 2D cylindrical
 613 %      (radial) filter (Jinc)
 614 %
 615 %    o exception: return any errors or warnings in this structure.
 616 %
 617 */
 618 MagickExport ResizeFilter *AcquireResizeFilter(const Image *image,
 619   const FilterTypes filter,const MagickRealType blur,
 620   const MagickBooleanType cylindrical,ExceptionInfo *exception)
 621 {
 622   const char
 623     *artifact;
 624
 625   FilterTypes
 626     filter_type,
 627     window_type;
 628
 629   MagickRealType
 630     B,
 631     C;
 632
 633   register ResizeFilter
 634     *resize_filter;
 635
 636   ssize_t
 637     option;
 638
 639   /*
 640     Table Mapping given Filter, into Weighting and Windowing functions.
 641     A 'Box' windowing function means its a simble non-windowed filter.
 642     An 'SincFast' filter function could be upgraded to a 'Jinc' filter
 643     if a "cylindrical", unless a 'Sinc' or 'SincFast' filter was
 644     specifically requested.
 645
 646     WARNING: The order of this tabel must match the order of the
 647     FilterTypes enumeration specified in "resample.h", or the filter
 648     names will not match the filter being setup.
 649
 650     You can check filter setups with the "filter:verbose" setting.
 651   */
 652   static struct
 653   {
 654     FilterTypes
 655       filter,
 656       window;
 657   } const mapping[SentinelFilter] =
 658   {
 659     { UndefinedFilter, BoxFilter },      /* Undefined (default to Box)       */
 660     { PointFilter,     BoxFilter },      /* SPECIAL: Nearest neighbour       */
 661     { BoxFilter,       BoxFilter },      /* Box averaging filter             */
 662     { TriangleFilter,  BoxFilter },      /* Linear interpolation filter      */
 663     { HermiteFilter,   BoxFilter },      /* Hermite interpolation filter     */
 664     { SincFastFilter,  HanningFilter },  /* Hanning -- cosine-sinc           */
 665     { SincFastFilter,  HammingFilter },  /* Hamming --      ''    variation  */
 666     { SincFastFilter,  BlackmanFilter }, /* Blackman -- 2*cosine-sinc        */
 667     { GaussianFilter,  BoxFilter },      /* Gaussian blur filter             */
 668     { QuadraticFilter, BoxFilter },      /* Quadratic Gaussian approximation */
 669     { CubicFilter,     BoxFilter },      /* Cubic B-Spline                   */
 670     { CatromFilter,    BoxFilter },      /* Cubic interpolator               */
 671     { MitchellFilter,  BoxFilter },      /* 'Ideal' cubic filter             */
 672     { LanczosFilter,   SincFastFilter }, /* SPECIAL: 3-lobed sinc-sinc       */
 673     { JincFilter,      BoxFilter },      /* Raw 3-lobed Jinc function        */
 674     { SincFilter,      BoxFilter },      /* Raw 4-lobed Sinc function        */
 675     { SincFastFilter,  KaiserFilter },   /* Kaiser -- square root-sinc       */
 676     { SincFastFilter,  WelshFilter },    /* Welsh -- parabolic-sinc          */
 677     { SincFastFilter,  CubicFilter },    /* Parzen -- cubic-sinc             */
 678     { LagrangeFilter,  BoxFilter },      /* Lagrange self-windowing filter   */
 679     { SincFastFilter,  BohmanFilter },   /* Bohman -- 2*cosine-sinc          */
 680     { SincFastFilter,  TriangleFilter }, /* Bartlett -- triangle-sinc        */
 681     { SincFastFilter,  BoxFilter },      /* Raw fast sinc ("Pade"-type)      */
 682     { Lanczos2DFilter, JincFilter },     /* SPECIAL: 2-lobed jinc-jinc       */
 683   };
 684   /*
 685     Table mapping the filter/window from the above table to an actual
 686     function.  The default support size for that filter as a weighting
 687     function, the range to scale with to use that function as a sinc
 688     windowing function, (typ 1.0).
 689
 690     Note that the filter_type -> function is 1 to 1 except for Sinc(),
 691     SincFast(), and CubicBC() functions, which may have multiple
 692     filter to function associations.
 693
 694     See "filter:verbose" handling below for the function -> filter
 695     mapping.
 696   */
 697   static struct
 698   {
 699     MagickRealType
 700       (*function)(const MagickRealType, const ResizeFilter*),
 701       lobes,    /* default lobes/support size of the weighting filter */
 702       scale,    /* windowing function range, for scaling windowing function */
 703       B,C;      /* Cubic Filter factors for a CubicBC function, else ignored */
 704   } const filters[SentinelFilter] =
 705   {
 706     { Box,       0.5, 0.5,     0.0, 0.0 }, /* Undefined (default to Box)  */
 707     { Box,       0.0, 0.5,     0.0, 0.0 }, /* Point (special handling)    */
 708     { Box,       0.5, 0.5,     0.0, 0.0 }, /* Box                         */
 709     { Triangle,  1.0, 1.0,     0.0, 0.0 }, /* Triangle                    */
 710     { CubicBC,   1.0, 1.0,     0.0, 0.0 }, /* Hermite (cubic  B=C=0)      */
 711     { Hanning,   1.0, 1.0,     0.0, 0.0 }, /* Hanning, cosine window      */
 712     { Hamming,   1.0, 1.0,     0.0, 0.0 }, /* Hamming, '' variation       */
 713     { Blackman,  1.0, 1.0,     0.0, 0.0 }, /* Blackman, 2*cosine window   */
 714     { Gaussian,  1.5, 1.5,     0.0, 0.0 }, /* Gaussian                    */
 715     { Quadratic, 1.5, 1.5,     0.0, 0.0 }, /* Quadratic gaussian          */
 716     { CubicBC,   2.0, 2.0,     1.0, 0.0 }, /* Cubic B-Spline (B=1,C=0)    */
 717     { CubicBC,   2.0, 1.0,     0.0, 0.5 }, /* Catmull-Rom    (B=0,C=1/2)  */
 718     { CubicBC,   2.0, 1.0, 1./3., 1./3. }, /* Mitchell      (B=C=1/3)    */
 719     { SincFast,  3.0, 1.0,     0.0, 0.0 }, /* Lanczos, 3-lobed sinc-sinc  */
 720     { Jinc,      3.0, 1.21967, 0.0, 0.0 }, /* Raw 3-lobed Jinc          */
 721     { Sinc,      4.0, 1.0,     0.0, 0.0 }, /* Raw 4-lobed Sinc            */
 722     { Kaiser,    1.0, 1.0,     0.0, 0.0 }, /* Kaiser (square root window) */
 723     { Welsh,     1.0, 1.0,     0.0, 0.0 }, /* Welsh (parabolic window)    */
 724     { CubicBC,   2.0, 2.0,     1.0, 0.0 }, /* Parzen (B-Spline window)    */
 725     { Lagrange,  2.0, 1.0,     0.0, 0.0 }, /* Lagrange sinc approximation */
 726     { Bohman,    1.0, 1.0,     0.0, 0.0 }, /* Bohman, 2*Cosine window     */
 727     { Triangle,  1.0, 1.0,     0.0, 0.0 }, /* Bartlett (triangle window)  */
 728     { SincFast,  4.0, 1.0,     0.0, 0.0 }, /* Raw fast sinc ("Pade"-type) */
 729     { Jinc,      2.0, 1.0,     0.0, 0.0 }, /* Lanczos2D, adjusted below   */
 730   };
 731   /*
 732     The known zero crossings of the Jinc() or more accuritally the Jinc(x*PI)
 733     function being used as a filter. It is used by the "filter:lobes" for of
 734     support selection, so users do not have to deal with the highly irrational
 735     sizes of the 'lobes' of the Jinc filter.
 736
 737     Values were sourced from
 738       http://cose.math.bas.bg/webMathematica/webComputing/BesselZeros.jsp
 739     Using Jv-function with v=1, then divided by PI.
 740   */
 741   static MagickRealType
 742     jinc_zeros[16] =
 743     {
 744       1.21966989126651,
 745       2.23313059438153,
 746       3.23831548416624,
 747       4.24106286379607,
 748       5.24276437687019,
 749       6.24392168986449,
 750       7.24475986871996,
 751       8.24539491395205,
 752       9.24589268494948,
 753       10.2462933487549,
 754       11.2466227948779,
 755       12.2468984611381,
 756       13.2471325221811,
 757       14.2473337358069,
 758       15.2475085630373,
 759       16.247661874701
 760    };
 761
 762   /*
 763     Allocate resize filter.
 764   */
 765   assert(image != (const Image *) NULL);
 766   assert(image->signature == MagickSignature);
 767   if (image->debug != MagickFalse)
 768     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
 769   assert(UndefinedFilter < filter && filter < SentinelFilter);
 770   assert(exception != (ExceptionInfo *) NULL);
 771   assert(exception->signature == MagickSignature);
 772   resize_filter=(ResizeFilter *) AcquireMagickMemory(sizeof(*resize_filter));
 773   if (resize_filter == (ResizeFilter *) NULL)
 774     ThrowFatalException(ResourceLimitFatalError,"MemoryAllocationFailed");
 775   /*
 776     Defaults for the requested filter.
 777   */
 778   filter_type=mapping[filter].filter;
 779   window_type=mapping[filter].window;
 780   /* Cylindrical Filters should use Jinc instead of Sinc */
 781   if (cylindrical != MagickFalse)
 782     switch (filter_type)
 783     {
 784       case SincFilter:
 785         /* Promote 1D Sinc Filter to a 2D Jinc filter. */
 786         if ( filter != SincFilter )
 787           filter_type=JincFilter;
 788         break;
 789       case SincFastFilter:
 790         /* Ditto for SincFast variant */
 791         if ( filter != SincFastFilter )
 792           filter_type=JincFilter;
 793         break;
 794
 795       case LanczosFilter:
 796         /* Promote Lanczos from a Sinc-Sinc to a Jinc-Jinc */
 797         filter_type=JincFilter;
 798         window_type=JincFilter;
 799         break;
 800       default:
 801         break;
 802     }
 803   else
 804     switch (filter_type)
 805     {
 806       case Lanczos2DFilter:
 807         /* depromote to a 2-lobe Sinc-Sinc for orthoginal use */
 808         window_type=SincFastFilter;
 809         break;
 810       default:
 811         break;
 812     }
 813
 814   artifact=GetImageArtifact(image,"filter:filter");
 815   if (artifact != (const char *) NULL)
 816     {
 817       option=ParseMagickOption(MagickFilterOptions,MagickFalse,artifact);
 818       if ((UndefinedFilter < option) && (option < SentinelFilter))
 819         { /* Raw filter request - no window function. */
 820           filter_type=(FilterTypes) option;
 821           window_type=BoxFilter;
 822         }
 823       if (option == LanczosFilter)
 824         { /* Lanczos is not a real filter but a self windowing Sinc/Jinc. */
 825           filter_type=cylindrical != MagickFalse ? JincFilter : Lanczos2DFilter;
 826           window_type=cylindrical != MagickFalse ? JincFilter : SincFastFilter;
 827         }
 828       /* Filter override with a specific window function. */
 829       artifact=GetImageArtifact(image,"filter:window");
 830       if (artifact != (const char *) NULL)
 831         {
 832           option=ParseMagickOption(MagickFilterOptions,MagickFalse,artifact);
 833           if ((UndefinedFilter < option) && (option < SentinelFilter))
 834             {
 835               if (option != LanczosFilter)
 836                 window_type=(FilterTypes) option;
 837               else
 838                 window_type=cylindrical != MagickFalse ? JincFilter :
 839                   SincFastFilter;
 840             }
 841         }
 842     }
 843   else
 844     {
 845       /* Window specified, but no filter function?  Assume Sinc/Jinc. */
 846       artifact=GetImageArtifact(image,"filter:window");
 847       if (artifact != (const char *) NULL)
 848         {
 849           option=ParseMagickOption(MagickFilterOptions,MagickFalse,
 850             artifact);
 851           if ((UndefinedFilter < option) && (option < SentinelFilter))
 852             {
 853               filter_type=cylindrical != MagickFalse ?
 854                          JincFilter : SincFastFilter;
 855               window_type=(FilterTypes) option;
 856             }
 857         }
 858     }
 859   /* Assign the real functions to use for the filters selected. */
 860   resize_filter->filter=filters[filter_type].function;
 861   resize_filter->support=filters[filter_type].lobes;
 862   resize_filter->window=filters[window_type].function;
 863   resize_filter->scale=filters[window_type].scale;
 864   resize_filter->signature=MagickSignature;
 865
 866   /* Filter blur -- scaling both filter and support window.  */
 867   resize_filter->blur=blur;
 868   artifact=GetImageArtifact(image,"filter:blur");
 869   if (artifact != (const char *) NULL)
 870     resize_filter->blur=StringToDouble(artifact);
 871   if (resize_filter->blur < MagickEpsilon)
 872     resize_filter->blur=(MagickRealType) MagickEpsilon;
 873
 874   if (cylindrical != MagickFalse)
 875     switch (filter_type)
 876     {
 877       case PointFilter:
 878       case BoxFilter:
 879         /* Support for Cylindrical Box should be sqrt(2)/2 */
 880         resize_filter->support=(MagickRealType) MagickSQ1_2;
 881         break;
 882       case GaussianFilter:
 883         /* Cylindrical Gaussian should have a sigma of sqrt(2)/2
 884          * and not the default sigma of 1/2 - so use blur to enlarge
 885          * and adjust support so actual practical support = 2.0 by default
 886          */
 887         resize_filter->blur *= MagickSQ2;
 888         resize_filter->support = (MagickRealType) MagickSQ2; /* which times blur => 2.0 */
 889         break;
 890       case Lanczos2DFilter:
 891         /* Special 2 lobed cylindrical Jinc-Jinc filter,
 892          * with a special blur adjustment to remove the blurring
 893          * effect of the windowing of the Jinc function (in the 2
 894          * lobed case only).  To be used as the default filter for EWA
 895          * Resampling and Distorts.
 896          *
 897          * Derivation: Set the scaling s=1/blur of the Lanczos2D
 898          * filter function so that
 899          *   Lanczos2D(s)=-2*Lanczos2D(s*sqrt(2))-Lanczos2D(s*2)
 900          * which implies that with a no-op a single vertical line is
 901          * not amplified (although it has negative ripples), and also
 902          * so that the high frequency checkerboard mode, as well as
 903          * the high frequency vertical, and horizontal, stripe modes
 904          * are slightly damped.
 905          *
 906          * (The value which preserves the high frequency vertical, and
 907          * horizontal, stripe modes (and slightly dampens the
 908          * checkerboard mode) satisfies
 909          *   Lanczos2D(s)=-2*Lanczos2D(s*sqrt(2))
 910          * which gives 0.9549921738 instead of 0.958033808.)
 911          *
 912          * Note from Nicolas: It is still not totally clear what
 913          * scaling of the Lanczos2D kernel is optimal for Clamped-EWA
 914          * upsampling.
 915          */
 916         resize_filter->blur *= (MagickRealType) 0.958033808;
 917       default:
 918         break;
 919     }
 920   else
 921     switch (filter_type)
 922     {
 923       case Lanczos2DFilter:
 924         /* depromote to a 2-lobe Sinc-Sinc for orthoginal use */
 925         resize_filter->filter=SincFast;
 926         break;
 927       default:
 928         break;
 929     }
 930
 931   /* Filter support overrides. */
 932   artifact=GetImageArtifact(image,"filter:lobes");
 933   if (artifact != (const char *) NULL)
 934     {
 935       ssize_t
 936         lobes;
 937
 938       lobes=(ssize_t) StringToLong(artifact);
 939       if (lobes < 1)
 940         lobes=1;
 941       resize_filter->support=(MagickRealType) lobes;
 942     }
 943   /* convert Jinc lobes to a real support value */
 944   if (resize_filter->filter == Jinc)
 945     {
 946       if (resize_filter->support > 16)
 947         resize_filter->support=jinc_zeros[15];  /* largest entry in table */
 948       else
 949         resize_filter->support = jinc_zeros[((long)resize_filter->support)-1];
 950     }
 951   /* expert override of the support setting */
 952   artifact=GetImageArtifact(image,"filter:support");
 953   if (artifact != (const char *) NULL)
 954     resize_filter->support=fabs(StringToDouble(artifact));
 955   /*
 956     Scale windowing function separatally to the support 'clipping'
 957     window that calling operator is planning to actually use. (Expert
 958     override)
 959   */
 960   resize_filter->window_support=resize_filter->support; /* default */
 961   artifact=GetImageArtifact(image,"filter:win-support");
 962   if (artifact != (const char *) NULL)
 963     resize_filter->window_support=fabs(StringToDouble(artifact));
 964   /*
 965     Adjust window function scaling to the windowing support for
 966     weighting function.  This avoids a division on every filter call.
 967   */
 968   resize_filter->scale /= resize_filter->window_support;
 969   /*
 970     Set Cubic Spline B,C values, calculate Cubic coefficients.
 971   */
 972   B=0.0;
 973   C=0.0;
 974   if ((filters[filter_type].function == CubicBC) ||
 975       (filters[window_type].function == CubicBC))
 976     {
 977       B=filters[filter_type].B;
 978       C=filters[filter_type].C;
 979       if (filters[window_type].function == CubicBC)
 980         {
 981           B=filters[window_type].B;
 982           C=filters[window_type].C;
 983         }
 984       artifact=GetImageArtifact(image,"filter:b");
 985       if (artifact != (const char *) NULL)
 986         {
 987           B=StringToDouble(artifact);
 988           C=(1.0-B)/2.0; /* Calculate C as if it is a Keys cubic filter. */
 989           artifact=GetImageArtifact(image,"filter:c");
 990           if (artifact != (const char *) NULL)
 991             C=StringToDouble(artifact);
 992         }
 993       else
 994         {
 995           artifact=GetImageArtifact(image,"filter:c");
 996           if (artifact != (const char *) NULL)
 997             {
 998               C=StringToDouble(artifact);
 999               B=1.0-2.0*C;  /* Calculate B as if it is a Keys cubic filter. */
1000             }
1001         }
1002     /*
1003       Convert B,C values into Cubic Coefficents.  See CubicBC().
1004     */
1005     resize_filter->cubic[0]=(6.0-2.0*B)/6.0;
1006     resize_filter->cubic[1]=0.0;
1007     resize_filter->cubic[2]=(-18.0+12.0*B+6.0*C)/6.0;
1008     resize_filter->cubic[3]=(12.0-9.0*B-6.0*C)/6.0;
1009     resize_filter->cubic[4]=(8.0*B+24.0*C)/6.0;
1010     resize_filter->cubic[5]=(-12.0*B-48.0*C)/6.0;
1011     resize_filter->cubic[6]=(6.0*B+30.0*C)/6.0;
1012     resize_filter->cubic[7]=(-B-6.0*C)/6.0;
1013   }
1014   /*
1015     Expert Option Request for verbose details of the resulting filter.
1016   */
1017 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1018   if( GetOpenMPThreadId() == 0 ) {
1019 #endif
1020     artifact=GetImageArtifact(image,"filter:verbose");
1021     if (artifact != (const char *) NULL)
1022       {
1023         double
1024           support,
1025           x;
1026
1027         /*
1028           Set the weighting function properly when the weighting
1029           function may not exactly match the filter of the same name.
1030           EG: a Point filter really uses a Box weighting function
1031           with a different support than is typically used.
1032
1033         */
1034         if (resize_filter->filter == Box)       filter_type=BoxFilter;
1035         if (resize_filter->filter == Sinc)      filter_type=SincFilter;
1036         if (resize_filter->filter == SincFast)  filter_type=SincFastFilter;
1037         if (resize_filter->filter == Jinc)      filter_type=JincFilter;
1038         if (resize_filter->filter == CubicBC)   filter_type=CubicFilter;
1039         /*
1040           Report Filter Details.
1041         */
1042         support=GetResizeFilterSupport(resize_filter); /* support range */
1043         (void) fprintf(stdout,"# Resize Filter (for graphing)\n#\n");
1044         (void) fprintf(stdout,"# filter = %s\n",MagickOptionToMnemonic(
1045            MagickFilterOptions,filter_type));
1046         (void) fprintf(stdout,"# window = %s\n",MagickOptionToMnemonic(
1047            MagickFilterOptions, window_type));
1048         (void) fprintf(stdout,"# support = %.*g\n",GetMagickPrecision(),
1049            (double) resize_filter->support);
1050         (void) fprintf(stdout,"# win-support = %.*g\n",GetMagickPrecision(),
1051            (double) resize_filter->window_support);
1052         (void) fprintf(stdout,"# blur = %.*g\n",GetMagickPrecision(),
1053            (double) resize_filter->blur);
1054        (void) fprintf(stdout,"# blurred_support = %.*g\n",GetMagickPrecision(),
1055            (double) support);
1056         (void) fprintf(stdout,"# B,C = %.*g,%.*g\n",GetMagickPrecision(),
1057            (double) B,GetMagickPrecision(),(double) C);
1058         (void) fprintf(stdout,"\n");
1059         /*
1060           Output values of resulting filter graph -- for graphing
1061           filter result.
1062         */
1063         for (x=0.0; x <= support; x+=0.01f)
1064           (void) fprintf(stdout,"%5.2lf\t%.*g\n",x,GetMagickPrecision(),
1065             (double) GetResizeFilterWeight(resize_filter,x));
1066         /* A final value so gnuplot can graph the 'stop' properly. */
1067         (void) fprintf(stdout,"%5.2lf\t%.*g\n",support,GetMagickPrecision(),
1068           0.0);
1069       }
1070 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1071     }
1072 #endif
1073   return(resize_filter);
1074 }
1075 \f
1076 /*
1077 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1078 %                                                                             %
1079 %                                                                             %
1080 %                                                                             %
1081 %   A d a p t i v e R e s i z e I m a g e                                     %
1082 %                                                                             %
1083 %                                                                             %
1084 %                                                                             %
1085 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1086 %
1087 %  AdaptiveResizeImage() adaptively resize image with pixel resampling.
1088 %
1089 %  The format of the AdaptiveResizeImage method is:
1090 %
1091 %      Image *AdaptiveResizeImage(const Image *image,const size_t columns,
1092 %        const size_t rows,ExceptionInfo *exception)
1093 %
1094 %  A description of each parameter follows:
1095 %
1096 %    o image: the image.
1097 %
1098 %    o columns: the number of columns in the resized image.
1099 %
1100 %    o rows: the number of rows in the resized image.
1101 %
1102 %    o exception: return any errors or warnings in this structure.
1103 %
1104 */
1105 MagickExport Image *AdaptiveResizeImage(const Image *image,
1106   const size_t columns,const size_t rows,ExceptionInfo *exception)
1107 {
1108 #define AdaptiveResizeImageTag  "Resize/Image"
1109
1110   CacheView
1111     *resize_view;
1112
1113   Image
1114     *resize_image;
1115
1116   MagickBooleanType
1117     proceed;
1118
1119   MagickPixelPacket
1120     pixel;
1121
1122   PointInfo
1123     offset;
1124
1125   ResampleFilter
1126     *resample_filter;
1127
1128   ssize_t
1129     y;
1130
1131   /*
1132     Adaptively resize image.
1133   */
1134   assert(image != (const Image *) NULL);
1135   assert(image->signature == MagickSignature);
1136   if (image->debug != MagickFalse)
1137     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1138   assert(exception != (ExceptionInfo *) NULL);
1139   assert(exception->signature == MagickSignature);
1140   if ((columns == 0) || (rows == 0))
1141     return((Image *) NULL);
1142   if ((columns == image->columns) && (rows == image->rows))
1143     return(CloneImage(image,0,0,MagickTrue,exception));
1144   resize_image=CloneImage(image,columns,rows,MagickTrue,exception);
1145   if (resize_image == (Image *) NULL)
1146     return((Image *) NULL);
1147   if (SetImageStorageClass(resize_image,DirectClass) == MagickFalse)
1148     {
1149       InheritException(exception,&resize_image->exception);
1150       resize_image=DestroyImage(resize_image);
1151       return((Image *) NULL);
1152     }
1153   GetMagickPixelPacket(image,&pixel);
1154   resample_filter=AcquireResampleFilter(image,exception);
1155   (void) SetResampleFilter(resample_filter,PointFilter,1.0);
1156   (void) SetResampleFilterInterpolateMethod(resample_filter,
1157     MeshInterpolatePixel);
1158   resize_view=AcquireCacheView(resize_image);
1159   for (y=0; y < (ssize_t) resize_image->rows; y++)
1160   {
1161     register IndexPacket
1162       *restrict resize_indexes;
1163
1164     register ssize_t
1165       x;
1166
1167     register PixelPacket
1168       *restrict q;
1169
1170     q=QueueCacheViewAuthenticPixels(resize_view,0,y,resize_image->columns,1,
1171       exception);
1172     if (q == (PixelPacket *) NULL)
1173       break;
1174     resize_indexes=GetCacheViewAuthenticIndexQueue(resize_view);
1175     offset.y=((MagickRealType) y*image->rows/resize_image->rows);
1176     for (x=0; x < (ssize_t) resize_image->columns; x++)
1177     {
1178       offset.x=((MagickRealType) x*image->columns/resize_image->columns);
1179       (void) ResamplePixelColor(resample_filter,offset.x-0.5,offset.y-0.5,
1180         &pixel);
1181       SetPixelPacket(resize_image,&pixel,q,resize_indexes+x);
1182       q++;
1183     }
1184     if (SyncCacheViewAuthenticPixels(resize_view,exception) == MagickFalse)
1185       break;
1186     proceed=SetImageProgress(image,AdaptiveResizeImageTag,(MagickOffsetType) y,
1187       image->rows);
1188     if (proceed == MagickFalse)
1189       break;
1190   }
1191   resample_filter=DestroyResampleFilter(resample_filter);
1192   resize_view=DestroyCacheView(resize_view);
1193   return(resize_image);
1194 }
1195 \f
1196 /*
1197 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1198 %                                                                             %
1199 %                                                                             %
1200 %                                                                             %
1201 +   B e s s e l O r d e r O n e                                               %
1202 %                                                                             %
1203 %                                                                             %
1204 %                                                                             %
1205 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1206 %
1207 %  BesselOrderOne() computes the Bessel function of x of the first kind of
1208 %  order 0.  This is used to create the Jinc() filter function below.
1209 %
1210 %    Reduce x to |x| since j1(x)= -j1(-x), and for x in (0,8]
1211 %
1212 %       j1(x) = x*j1(x);
1213 %
1214 %    For x in (8,inf)
1215 %
1216 %       j1(x) = sqrt(2/(pi*x))*(p1(x)*cos(x1)-q1(x)*sin(x1))
1217 %
1218 %    where x1 = x-3*pi/4. Compute sin(x1) and cos(x1) as follow:
1219 %
1220 %       cos(x1) =  cos(x)cos(3pi/4)+sin(x)sin(3pi/4)
1221 %               =  1/sqrt(2) * (sin(x) - cos(x))
1222 %       sin(x1) =  sin(x)cos(3pi/4)-cos(x)sin(3pi/4)
1223 %               = -1/sqrt(2) * (sin(x) + cos(x))
1224 %
1225 %  The format of the BesselOrderOne method is:
1226 %
1227 %      MagickRealType BesselOrderOne(MagickRealType x)
1228 %
1229 %  A description of each parameter follows:
1230 %
1231 %    o x: MagickRealType value.
1232 %
1233 */
1234
1235 #undef I0
1236 static MagickRealType I0(MagickRealType x)
1237 {
1238   MagickRealType
1239     sum,
1240     t,
1241     y;
1242
1243   register ssize_t
1244     i;
1245
1246   /*
1247     Zeroth order Bessel function of the first kind.
1248   */
1249   sum=1.0;
1250   y=x*x/4.0;
1251   t=y;
1252   for (i=2; t > MagickEpsilon; i++)
1253   {
1254     sum+=t;
1255     t*=y/((MagickRealType) i*i);
1256   }
1257   return(sum);
1258 }
1259
1260 #undef J1
1261 static MagickRealType J1(MagickRealType x)
1262 {
1263   MagickRealType
1264     p,
1265     q;
1266
1267   register ssize_t
1268     i;
1269
1270   static const double
1271     Pone[] =
1272     {
1273        0.581199354001606143928050809e+21,
1274       -0.6672106568924916298020941484e+20,
1275        0.2316433580634002297931815435e+19,
1276       -0.3588817569910106050743641413e+17,
1277        0.2908795263834775409737601689e+15,
1278       -0.1322983480332126453125473247e+13,
1279        0.3413234182301700539091292655e+10,
1280       -0.4695753530642995859767162166e+7,
1281        0.270112271089232341485679099e+4
1282     },
1283     Qone[] =
1284     {
1285       0.11623987080032122878585294e+22,
1286       0.1185770712190320999837113348e+20,
1287       0.6092061398917521746105196863e+17,
1288       0.2081661221307607351240184229e+15,
1289       0.5243710262167649715406728642e+12,
1290       0.1013863514358673989967045588e+10,
1291       0.1501793594998585505921097578e+7,
1292       0.1606931573481487801970916749e+4,
1293       0.1e+1
1294     };
1295
1296   p=Pone[8];
1297   q=Qone[8];
1298   for (i=7; i >= 0; i--)
1299   {
1300     p=p*x*x+Pone[i];
1301     q=q*x*x+Qone[i];
1302   }
1303   return(p/q);
1304 }
1305
1306 #undef P1
1307 static MagickRealType P1(MagickRealType x)
1308 {
1309   MagickRealType
1310     p,
1311     q;
1312
1313   register ssize_t
1314     i;
1315
1316   static const double
1317     Pone[] =
1318     {
1319       0.352246649133679798341724373e+5,
1320       0.62758845247161281269005675e+5,
1321       0.313539631109159574238669888e+5,
1322       0.49854832060594338434500455e+4,
1323       0.2111529182853962382105718e+3,
1324       0.12571716929145341558495e+1
1325     },
1326     Qone[] =
1327     {
1328       0.352246649133679798068390431e+5,
1329       0.626943469593560511888833731e+5,
1330       0.312404063819041039923015703e+5,
1331       0.4930396490181088979386097e+4,
1332       0.2030775189134759322293574e+3,
1333       0.1e+1
1334     };
1335
1336   p=Pone[5];
1337   q=Qone[5];
1338   for (i=4; i >= 0; i--)
1339   {
1340     p=p*(8.0/x)*(8.0/x)+Pone[i];
1341     q=q*(8.0/x)*(8.0/x)+Qone[i];
1342   }
1343   return(p/q);
1344 }
1345
1346 #undef Q1
1347 static MagickRealType Q1(MagickRealType x)
1348 {
1349   MagickRealType
1350     p,
1351     q;
1352
1353   register ssize_t
1354     i;
1355
1356   static const double
1357     Pone[] =
1358     {
1359       0.3511751914303552822533318e+3,
1360       0.7210391804904475039280863e+3,
1361       0.4259873011654442389886993e+3,
1362       0.831898957673850827325226e+2,
1363       0.45681716295512267064405e+1,
1364       0.3532840052740123642735e-1
1365     },
1366     Qone[] =
1367     {
1368       0.74917374171809127714519505e+4,
1369       0.154141773392650970499848051e+5,
1370       0.91522317015169922705904727e+4,
1371       0.18111867005523513506724158e+4,
1372       0.1038187585462133728776636e+3,
1373       0.1e+1
1374     };
1375
1376   p=Pone[5];
1377   q=Qone[5];
1378   for (i=4; i >= 0; i--)
1379   {
1380     p=p*(8.0/x)*(8.0/x)+Pone[i];
1381     q=q*(8.0/x)*(8.0/x)+Qone[i];
1382   }
1383   return(p/q);
1384 }
1385
1386 static MagickRealType BesselOrderOne(MagickRealType x)
1387 {
1388   MagickRealType
1389     p,
1390     q;
1391
1392   if (x == 0.0)
1393     return(0.0);
1394   p=x;
1395   if (x < 0.0)
1396     x=(-x);
1397   if (x < 8.0)
1398     return(p*J1(x));
1399   q=sqrt((double) (2.0/(MagickPI*x)))*(P1(x)*(1.0/sqrt(2.0)*(sin((double) x)-
1400     cos((double) x)))-8.0/x*Q1(x)*(-1.0/sqrt(2.0)*(sin((double) x)+
1401     cos((double) x))));
1402   if (p < 0.0)
1403     q=(-q);
1404   return(q);
1405 }
1406 \f
1407 /*
1408 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1409 %                                                                             %
1410 %                                                                             %
1411 %                                                                             %
1412 +   D e s t r o y R e s i z e F i l t e r                                     %
1413 %                                                                             %
1414 %                                                                             %
1415 %                                                                             %
1416 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1417 %
1418 %  DestroyResizeFilter() destroy the resize filter.
1419 %
1420 %  The format of the DestroyResizeFilter method is:
1421 %
1422 %      ResizeFilter *DestroyResizeFilter(ResizeFilter *resize_filter)
1423 %
1424 %  A description of each parameter follows:
1425 %
1426 %    o resize_filter: the resize filter.
1427 %
1428 */
1429 MagickExport ResizeFilter *DestroyResizeFilter(ResizeFilter *resize_filter)
1430 {
1431   assert(resize_filter != (ResizeFilter *) NULL);
1432   assert(resize_filter->signature == MagickSignature);
1433   resize_filter->signature=(~MagickSignature);
1434   resize_filter=(ResizeFilter *) RelinquishMagickMemory(resize_filter);
1435   return(resize_filter);
1436 }
1437 \f
1438 /*
1439 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1440 %                                                                             %
1441 %                                                                             %
1442 %                                                                             %
1443 +   G e t R e s i z e F i l t e r S u p p o r t                               %
1444 %                                                                             %
1445 %                                                                             %
1446 %                                                                             %
1447 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1448 %
1449 %  GetResizeFilterSupport() return the current support window size for this
1450 %  filter.  Note that this may have been enlarged by filter:blur factor.
1451 %
1452 %  The format of the GetResizeFilterSupport method is:
1453 %
1454 %      MagickRealType GetResizeFilterSupport(const ResizeFilter *resize_filter)
1455 %
1456 %  A description of each parameter follows:
1457 %
1458 %    o filter: Image filter to use.
1459 %
1460 */
1461 MagickExport MagickRealType GetResizeFilterSupport(
1462   const ResizeFilter *resize_filter)
1463 {
1464   assert(resize_filter != (ResizeFilter *) NULL);
1465   assert(resize_filter->signature == MagickSignature);
1466   return(resize_filter->support*resize_filter->blur);
1467 }
1468 \f
1469 /*
1470 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1471 %                                                                             %
1472 %                                                                             %
1473 %                                                                             %
1474 +   G e t R e s i z e F i l t e r W e i g h t                                 %
1475 %                                                                             %
1476 %                                                                             %
1477 %                                                                             %
1478 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1479 %
1480 %  GetResizeFilterWeight evaluates the specified resize filter at the point x
1481 %  which usally lies between zero and the filters current 'support' and
1482 %  returns the weight of the filter function at that point.
1483 %
1484 %  The format of the GetResizeFilterWeight method is:
1485 %
1486 %      MagickRealType GetResizeFilterWeight(const ResizeFilter *resize_filter,
1487 %        const MagickRealType x)
1488 %
1489 %  A description of each parameter follows:
1490 %
1491 %    o filter: the filter type.
1492 %
1493 %    o x: the point.
1494 %
1495 */
1496 MagickExport MagickRealType GetResizeFilterWeight(
1497   const ResizeFilter *resize_filter,const MagickRealType x)
1498 {
1499   MagickRealType
1500     scale,
1501     x_blur;
1502
1503   /*
1504     Windowing function - scale the weighting filter by this amount.
1505   */
1506   assert(resize_filter != (ResizeFilter *) NULL);
1507   assert(resize_filter->signature == MagickSignature);
1508   x_blur=fabs((double) x)/resize_filter->blur;  /* X offset with blur scaling */
1509   if ((resize_filter->window_support < MagickEpsilon) ||
1510       (resize_filter->window == Box))
1511     scale=1.0;  /* Point or Box Filter -- avoid division by zero */
1512   else
1513     {
1514       scale=resize_filter->scale;
1515       scale=resize_filter->window(x_blur*scale,resize_filter);
1516     }
1517   return(scale*resize_filter->filter(x_blur,resize_filter));
1518 }
1519 \f
1520 /*
1521 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1522 %                                                                             %
1523 %                                                                             %
1524 %                                                                             %
1525 %   M a g n i f y I m a g e                                                   %
1526 %                                                                             %
1527 %                                                                             %
1528 %                                                                             %
1529 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1530 %
1531 %  MagnifyImage() is a convenience method that scales an image proportionally
1532 %  to twice its size.
1533 %
1534 %  The format of the MagnifyImage method is:
1535 %
1536 %      Image *MagnifyImage(const Image *image,ExceptionInfo *exception)
1537 %
1538 %  A description of each parameter follows:
1539 %
1540 %    o image: the image.
1541 %
1542 %    o exception: return any errors or warnings in this structure.
1543 %
1544 */
1545 MagickExport Image *MagnifyImage(const Image *image,ExceptionInfo *exception)
1546 {
1547   Image
1548     *magnify_image;
1549
1550   assert(image != (Image *) NULL);
1551   assert(image->signature == MagickSignature);
1552   if (image->debug != MagickFalse)
1553     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1554   assert(exception != (ExceptionInfo *) NULL);
1555   assert(exception->signature == MagickSignature);
1556   magnify_image=ResizeImage(image,2*image->columns,2*image->rows,CubicFilter,
1557     1.0,exception);
1558   return(magnify_image);
1559 }
1560 \f
1561 /*
1562 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1563 %                                                                             %
1564 %                                                                             %
1565 %                                                                             %
1566 %   M i n i f y I m a g e                                                     %
1567 %                                                                             %
1568 %                                                                             %
1569 %                                                                             %
1570 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1571 %
1572 %  MinifyImage() is a convenience method that scales an image proportionally
1573 %  to half its size.
1574 %
1575 %  The format of the MinifyImage method is:
1576 %
1577 %      Image *MinifyImage(const Image *image,ExceptionInfo *exception)
1578 %
1579 %  A description of each parameter follows:
1580 %
1581 %    o image: the image.
1582 %
1583 %    o exception: return any errors or warnings in this structure.
1584 %
1585 */
1586 MagickExport Image *MinifyImage(const Image *image,ExceptionInfo *exception)
1587 {
1588   Image
1589     *minify_image;
1590
1591   assert(image != (Image *) NULL);
1592   assert(image->signature == MagickSignature);
1593   if (image->debug != MagickFalse)
1594     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1595   assert(exception != (ExceptionInfo *) NULL);
1596   assert(exception->signature == MagickSignature);
1597   minify_image=ResizeImage(image,image->columns/2,image->rows/2,CubicFilter,
1598     1.0,exception);
1599   return(minify_image);
1600 }
1601 \f
1602 /*
1603 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1604 %                                                                             %
1605 %                                                                             %
1606 %                                                                             %
1607 %   R e s a m p l e I m a g e                                                 %
1608 %                                                                             %
1609 %                                                                             %
1610 %                                                                             %
1611 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1612 %
1613 %  ResampleImage() resize image in terms of its pixel size, so that when
1614 %  displayed at the given resolution it will be the same size in terms of
1615 %  real world units as the original image at the original resolution.
1616 %
1617 %  The format of the ResampleImage method is:
1618 %
1619 %      Image *ResampleImage(Image *image,const double x_resolution,
1620 %        const double y_resolution,const FilterTypes filter,const double blur,
1621 %        ExceptionInfo *exception)
1622 %
1623 %  A description of each parameter follows:
1624 %
1625 %    o image: the image to be resized to fit the given resolution.
1626 %
1627 %    o x_resolution: the new image x resolution.
1628 %
1629 %    o y_resolution: the new image y resolution.
1630 %
1631 %    o filter: Image filter to use.
1632 %
1633 %    o blur: the blur factor where > 1 is blurry, < 1 is sharp.
1634 %
1635 */
1636 MagickExport Image *ResampleImage(const Image *image,const double x_resolution,
1637   const double y_resolution,const FilterTypes filter,const double blur,
1638   ExceptionInfo *exception)
1639 {
1640 #define ResampleImageTag  "Resample/Image"
1641
1642   Image
1643     *resample_image;
1644
1645   size_t
1646     height,
1647     width;
1648
1649   /*
1650     Initialize sampled image attributes.
1651   */
1652   assert(image != (const Image *) NULL);
1653   assert(image->signature == MagickSignature);
1654   if (image->debug != MagickFalse)
1655     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1656   assert(exception != (ExceptionInfo *) NULL);
1657   assert(exception->signature == MagickSignature);
1658   width=(size_t) (x_resolution*image->columns/(image->x_resolution == 0.0 ?
1659     72.0 : image->x_resolution)+0.5);
1660   height=(size_t) (y_resolution*image->rows/(image->y_resolution == 0.0 ?
1661     72.0 : image->y_resolution)+0.5);
1662   resample_image=ResizeImage(image,width,height,filter,blur,exception);
1663   if (resample_image != (Image *) NULL)
1664     {
1665       resample_image->x_resolution=x_resolution;
1666       resample_image->y_resolution=y_resolution;
1667     }
1668   return(resample_image);
1669 }
1670 #if defined(MAGICKCORE_LQR_DELEGATE)
1671 \f
1672 /*
1673 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1674 %                                                                             %
1675 %                                                                             %
1676 %                                                                             %
1677 %   L i q u i d R e s c a l e I m a g e                                       %
1678 %                                                                             %
1679 %                                                                             %
1680 %                                                                             %
1681 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1682 %
1683 %  LiquidRescaleImage() rescales image with seam carving.
1684 %
1685 %  The format of the LiquidRescaleImage method is:
1686 %
1687 %      Image *LiquidRescaleImage(const Image *image,
1688 %        const size_t columns,const size_t rows,
1689 %        const double delta_x,const double rigidity,ExceptionInfo *exception)
1690 %
1691 %  A description of each parameter follows:
1692 %
1693 %    o image: the image.
1694 %
1695 %    o columns: the number of columns in the rescaled image.
1696 %
1697 %    o rows: the number of rows in the rescaled image.
1698 %
1699 %    o delta_x: maximum seam transversal step (0 means straight seams).
1700 %
1701 %    o rigidity: introduce a bias for non-straight seams (typically 0).
1702 %
1703 %    o exception: return any errors or warnings in this structure.
1704 %
1705 */
1706 MagickExport Image *LiquidRescaleImage(const Image *image,const size_t columns,
1707   const size_t rows,const double delta_x,const double rigidity,
1708   ExceptionInfo *exception)
1709 {
1710 #define LiquidRescaleImageTag  "Rescale/Image"
1711
1712   CacheView
1713     *rescale_view;
1714
1715   const char
1716     *map;
1717
1718   guchar
1719     *packet;
1720
1721   Image
1722     *rescale_image;
1723
1724   int
1725     x,
1726     y;
1727
1728   LqrCarver
1729     *carver;
1730
1731   LqrRetVal
1732     lqr_status;
1733
1734   MagickBooleanType
1735     status;
1736
1737   MagickPixelPacket
1738     pixel;
1739
1740   unsigned char
1741     *pixels;
1742
1743   /*
1744     Liquid rescale image.
1745   */
1746   assert(image != (const Image *) NULL);
1747   assert(image->signature == MagickSignature);
1748   if (image->debug != MagickFalse)
1749     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1750   assert(exception != (ExceptionInfo *) NULL);
1751   assert(exception->signature == MagickSignature);
1752   if ((columns == 0) || (rows == 0))
1753     return((Image *) NULL);
1754   if ((columns == image->columns) && (rows == image->rows))
1755     return(CloneImage(image,0,0,MagickTrue,exception));
1756   if ((columns <= 2) || (rows <= 2))
1757     return(ResizeImage(image,columns,rows,image->filter,image->blur,exception));
1758   if ((columns >= (2*image->columns)) || (rows >= (2*image->rows)))
1759     {
1760       Image
1761         *resize_image;
1762
1763       size_t
1764         height,
1765         width;
1766
1767       /*
1768         Honor liquid resize size limitations.
1769       */
1770       for (width=image->columns; columns >= (2*width-1); width*=2);
1771       for (height=image->rows; rows >= (2*height-1); height*=2);
1772       resize_image=ResizeImage(image,width,height,image->filter,image->blur,
1773         exception);
1774       if (resize_image == (Image *) NULL)
1775         return((Image *) NULL);
1776       rescale_image=LiquidRescaleImage(resize_image,columns,rows,delta_x,
1777         rigidity,exception);
1778       resize_image=DestroyImage(resize_image);
1779       return(rescale_image);
1780     }
1781   map="RGB";
1782   if (image->matte == MagickFalse)
1783     map="RGBA";
1784   if (image->colorspace == CMYKColorspace)
1785     {
1786       map="CMYK";
1787       if (image->matte == MagickFalse)
1788         map="CMYKA";
1789     }
1790   pixels=(unsigned char *) AcquireQuantumMemory(image->columns,image->rows*
1791     strlen(map)*sizeof(*pixels));
1792   if (pixels == (unsigned char *) NULL)
1793     return((Image *) NULL);
1794   status=ExportImagePixels(image,0,0,image->columns,image->rows,map,CharPixel,
1795     pixels,exception);
1796   if (status == MagickFalse)
1797     {
1798       pixels=(unsigned char *) RelinquishMagickMemory(pixels);
1799       ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
1800     }
1801   carver=lqr_carver_new(pixels,image->columns,image->rows,strlen(map));
1802   if (carver == (LqrCarver *) NULL)
1803     {
1804       pixels=(unsigned char *) RelinquishMagickMemory(pixels);
1805       ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
1806     }
1807   lqr_status=lqr_carver_init(carver,(int) delta_x,rigidity);
1808   lqr_status=lqr_carver_resize(carver,columns,rows);
1809   rescale_image=CloneImage(image,lqr_carver_get_width(carver),
1810     lqr_carver_get_height(carver),MagickTrue,exception);
1811   if (rescale_image == (Image *) NULL)
1812     {
1813       pixels=(unsigned char *) RelinquishMagickMemory(pixels);
1814       return((Image *) NULL);
1815     }
1816   if (SetImageStorageClass(rescale_image,DirectClass) == MagickFalse)
1817     {
1818       InheritException(exception,&rescale_image->exception);
1819       rescale_image=DestroyImage(rescale_image);
1820       return((Image *) NULL);
1821     }
1822   GetMagickPixelPacket(rescale_image,&pixel);
1823   (void) lqr_carver_scan_reset(carver);
1824   rescale_view=AcquireCacheView(rescale_image);
1825   while (lqr_carver_scan(carver,&x,&y,&packet) != 0)
1826   {
1827     register IndexPacket
1828       *restrict rescale_indexes;
1829
1830     register PixelPacket
1831       *restrict q;
1832
1833     q=QueueCacheViewAuthenticPixels(rescale_view,x,y,1,1,exception);
1834     if (q == (PixelPacket *) NULL)
1835       break;
1836     rescale_indexes=GetCacheViewAuthenticIndexQueue(rescale_view);
1837     pixel.red=QuantumRange*(packet[0]/255.0);
1838     pixel.green=QuantumRange*(packet[1]/255.0);
1839     pixel.blue=QuantumRange*(packet[2]/255.0);
1840     if (image->colorspace != CMYKColorspace)
1841       {
1842         if (image->matte == MagickFalse)
1843           pixel.opacity=QuantumRange*(packet[3]/255.0);
1844       }
1845     else
1846       {
1847         pixel.index=QuantumRange*(packet[3]/255.0);
1848         if (image->matte == MagickFalse)
1849           pixel.opacity=QuantumRange*(packet[4]/255.0);
1850       }
1851     SetPixelPacket(rescale_image,&pixel,q,rescale_indexes);
1852     if (SyncCacheViewAuthenticPixels(rescale_view,exception) == MagickFalse)
1853       break;
1854   }
1855   rescale_view=DestroyCacheView(rescale_view);
1856   /*
1857     Relinquish resources.
1858   */
1859   lqr_carver_destroy(carver);
1860   return(rescale_image);
1861 }
1862 #else
1863 MagickExport Image *LiquidRescaleImage(const Image *image,
1864   const size_t magick_unused(columns),const size_t magick_unused(rows),
1865   const double magick_unused(delta_x),const double magick_unused(rigidity),
1866   ExceptionInfo *exception)
1867 {
1868   assert(image != (const Image *) NULL);
1869   assert(image->signature == MagickSignature);
1870   if (image->debug != MagickFalse)
1871     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1872   assert(exception != (ExceptionInfo *) NULL);
1873   assert(exception->signature == MagickSignature);
1874   (void) ThrowMagickException(exception,GetMagickModule(),MissingDelegateError,
1875     "DelegateLibrarySupportNotBuiltIn","`%s' (LQR)",image->filename);
1876   return((Image *) NULL);
1877 }
1878 #endif
1879 \f
1880 /*
1881 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1882 %                                                                             %
1883 %                                                                             %
1884 %                                                                             %
1885 %   R e s i z e I m a g e                                                     %
1886 %                                                                             %
1887 %                                                                             %
1888 %                                                                             %
1889 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1890 %
1891 %  ResizeImage() scales an image to the desired dimensions, using the given
1892 %  filter (see AcquireFilterInfo()).
1893 %
1894 %  If an undefined filter is given the filter defaults to Mitchell for a
1895 %  colormapped image, a image with a matte channel, or if the image is
1896 %  enlarged.  Otherwise the filter defaults to a Lanczos.
1897 %
1898 %  ResizeImage() was inspired by Paul Heckbert's "zoom" program.
1899 %
1900 %  The format of the ResizeImage method is:
1901 %
1902 %      Image *ResizeImage(Image *image,const size_t columns,
1903 %        const size_t rows,const FilterTypes filter,const double blur,
1904 %        ExceptionInfo *exception)
1905 %
1906 %  A description of each parameter follows:
1907 %
1908 %    o image: the image.
1909 %
1910 %    o columns: the number of columns in the scaled image.
1911 %
1912 %    o rows: the number of rows in the scaled image.
1913 %
1914 %    o filter: Image filter to use.
1915 %
1916 %    o blur: the blur factor where > 1 is blurry, < 1 is sharp.  Typically set
1917 %      this to 1.0.
1918 %
1919 %    o exception: return any errors or warnings in this structure.
1920 %
1921 */
1922
1923 typedef struct _ContributionInfo
1924 {
1925   MagickRealType
1926     weight;
1927
1928   ssize_t
1929     pixel;
1930 } ContributionInfo;
1931
1932 static ContributionInfo **DestroyContributionThreadSet(
1933   ContributionInfo **contribution)
1934 {
1935   register ssize_t
1936     i;
1937
1938   assert(contribution != (ContributionInfo **) NULL);
1939   for (i=0; i < (ssize_t) GetOpenMPMaximumThreads(); i++)
1940     if (contribution[i] != (ContributionInfo *) NULL)
1941       contribution[i]=(ContributionInfo *) RelinquishMagickMemory(
1942         contribution[i]);
1943   contribution=(ContributionInfo **) RelinquishMagickMemory(contribution);
1944   return(contribution);
1945 }
1946
1947 static ContributionInfo **AcquireContributionThreadSet(const size_t count)
1948 {
1949   register ssize_t
1950     i;
1951
1952   ContributionInfo
1953     **contribution;
1954
1955   size_t
1956     number_threads;
1957
1958   number_threads=GetOpenMPMaximumThreads();
1959   contribution=(ContributionInfo **) AcquireQuantumMemory(number_threads,
1960     sizeof(*contribution));
1961   if (contribution == (ContributionInfo **) NULL)
1962     return((ContributionInfo **) NULL);
1963   (void) ResetMagickMemory(contribution,0,number_threads*sizeof(*contribution));
1964   for (i=0; i < (ssize_t) number_threads; i++)
1965   {
1966     contribution[i]=(ContributionInfo *) AcquireQuantumMemory(count,
1967       sizeof(**contribution));
1968     if (contribution[i] == (ContributionInfo *) NULL)
1969       return(DestroyContributionThreadSet(contribution));
1970   }
1971   return(contribution);
1972 }
1973
1974 static inline double MagickMax(const double x,const double y)
1975 {
1976   if (x > y)
1977     return(x);
1978   return(y);
1979 }
1980
1981 static inline double MagickMin(const double x,const double y)
1982 {
1983   if (x < y)
1984     return(x);
1985   return(y);
1986 }
1987
1988 static MagickBooleanType HorizontalFilter(const ResizeFilter *resize_filter,
1989   const Image *image,Image *resize_image,const MagickRealType x_factor,
1990   const MagickSizeType span,MagickOffsetType *offset,ExceptionInfo *exception)
1991 {
1992 #define ResizeImageTag  "Resize/Image"
1993
1994   CacheView
1995     *image_view,
1996     *resize_view;
1997
1998   ClassType
1999     storage_class;
2000
2001   ContributionInfo
2002     **restrict contributions;
2003
2004   MagickBooleanType
2005     status;
2006
2007   MagickPixelPacket
2008     zero;
2009
2010   MagickRealType
2011     scale,
2012     support;
2013
2014   ssize_t
2015     x;
2016
2017   /*
2018     Apply filter to resize horizontally from image to resize image.
2019   */
2020   scale=MagickMax(1.0/x_factor+MagickEpsilon,1.0);
2021   support=scale*GetResizeFilterSupport(resize_filter);
2022   storage_class=support > 0.5 ? DirectClass : image->storage_class;
2023   if (SetImageStorageClass(resize_image,storage_class) == MagickFalse)
2024     {
2025       InheritException(exception,&resize_image->exception);
2026       return(MagickFalse);
2027     }
2028   if (support < 0.5)
2029     {
2030       /*
2031         Support too small even for nearest neighbour: Reduce to point
2032         sampling.
2033       */
2034       support=(MagickRealType) 0.5;
2035       scale=1.0;
2036     }
2037   contributions=AcquireContributionThreadSet((size_t) (2.0*support+3.0));
2038   if (contributions == (ContributionInfo **) NULL)
2039     {
2040       (void) ThrowMagickException(exception,GetMagickModule(),
2041         ResourceLimitError,"MemoryAllocationFailed","`%s'",image->filename);
2042       return(MagickFalse);
2043     }
2044   status=MagickTrue;
2045   scale=1.0/scale;
2046   (void) ResetMagickMemory(&zero,0,sizeof(zero));
2047   image_view=AcquireCacheView(image);
2048   resize_view=AcquireCacheView(resize_image);
2049 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2050   #pragma omp parallel for shared(status)
2051 #endif
2052   for (x=0; x < (ssize_t) resize_image->columns; x++)
2053   {
2054     MagickRealType
2055       center,
2056       density;
2057
2058     register const IndexPacket
2059       *restrict indexes;
2060
2061     register const PixelPacket
2062       *restrict p;
2063
2064     register ContributionInfo
2065       *restrict contribution;
2066
2067     register IndexPacket
2068       *restrict resize_indexes;
2069
2070     register PixelPacket
2071       *restrict q;
2072
2073     register ssize_t
2074       y;
2075
2076     ssize_t
2077       n,
2078       start,
2079       stop;
2080
2081     if (status == MagickFalse)
2082       continue;
2083     center=(MagickRealType) (x+0.5)/x_factor;
2084     start=(ssize_t) MagickMax(center-support+0.5,0.0);
2085     stop=(ssize_t) MagickMin(center+support+0.5,(double) image->columns);
2086     density=0.0;
2087     contribution=contributions[GetOpenMPThreadId()];
2088     for (n=0; n < (stop-start); n++)
2089     {
2090       contribution[n].pixel=start+n;
2091       contribution[n].weight=GetResizeFilterWeight(resize_filter,scale*
2092         ((MagickRealType) (start+n)-center+0.5));
2093       density+=contribution[n].weight;
2094     }
2095     if ((density != 0.0) && (density != 1.0))
2096       {
2097         register ssize_t
2098           i;
2099
2100         /*
2101           Normalize.
2102         */
2103         density=1.0/density;
2104         for (i=0; i < n; i++)
2105           contribution[i].weight*=density;
2106       }
2107     p=GetCacheViewVirtualPixels(image_view,contribution[0].pixel,0,(size_t)
2108       (contribution[n-1].pixel-contribution[0].pixel+1),image->rows,exception);
2109     q=QueueCacheViewAuthenticPixels(resize_view,x,0,1,resize_image->rows,
2110       exception);
2111     if ((p == (const PixelPacket *) NULL) || (q == (PixelPacket *) NULL))
2112       {
2113         status=MagickFalse;
2114         continue;
2115       }
2116     indexes=GetCacheViewVirtualIndexQueue(image_view);
2117     resize_indexes=GetCacheViewAuthenticIndexQueue(resize_view);
2118     for (y=0; y < (ssize_t) resize_image->rows; y++)
2119     {
2120       MagickPixelPacket
2121         pixel;
2122
2123       MagickRealType
2124         alpha;
2125
2126       register ssize_t
2127         i;
2128
2129       ssize_t
2130         j;
2131
2132       pixel=zero;
2133       if (image->matte == MagickFalse)
2134         {
2135           for (i=0; i < n; i++)
2136           {
2137             j=y*(contribution[n-1].pixel-contribution[0].pixel+1)+
2138               (contribution[i].pixel-contribution[0].pixel);
2139             alpha=contribution[i].weight;
2140             pixel.red+=alpha*(p+j)->red;
2141             pixel.green+=alpha*(p+j)->green;
2142             pixel.blue+=alpha*(p+j)->blue;
2143             pixel.opacity+=alpha*(p+j)->opacity;
2144           }
2145           SetRedPixelComponent(q,ClampRedPixelComponent(&pixel));
2146           SetGreenPixelComponent(q,ClampGreenPixelComponent(&pixel));
2147           SetBluePixelComponent(q,ClampBluePixelComponent(&pixel));
2148           SetOpacityPixelComponent(q,ClampOpacityPixelComponent(&pixel));
2149           if ((image->colorspace == CMYKColorspace) &&
2150               (resize_image->colorspace == CMYKColorspace))
2151             {
2152               for (i=0; i < n; i++)
2153               {
2154                 j=y*(contribution[n-1].pixel-contribution[0].pixel+1)+
2155                   (contribution[i].pixel-contribution[0].pixel);
2156                 alpha=contribution[i].weight;
2157                 pixel.index+=alpha*indexes[j];
2158               }
2159               resize_indexes[y]=(IndexPacket) ClampToQuantum(pixel.index);
2160             }
2161         }
2162       else
2163         {
2164           MagickRealType
2165             gamma;
2166
2167           gamma=0.0;
2168           for (i=0; i < n; i++)
2169           {
2170             j=y*(contribution[n-1].pixel-contribution[0].pixel+1)+
2171               (contribution[i].pixel-contribution[0].pixel);
2172             alpha=contribution[i].weight*QuantumScale*
2173               GetAlphaPixelComponent(p+j);
2174             pixel.red+=alpha*(p+j)->red;
2175             pixel.green+=alpha*(p+j)->green;
2176             pixel.blue+=alpha*(p+j)->blue;
2177             pixel.opacity+=contribution[i].weight*(p+j)->opacity;
2178             gamma+=alpha;
2179           }
2180           gamma=1.0/(fabs((double) gamma) <= MagickEpsilon ? 1.0 : gamma);
2181           q->red=ClampToQuantum(gamma*GetRedPixelComponent(&pixel));
2182           q->green=ClampToQuantum(gamma*GetGreenPixelComponent(&pixel));
2183           q->blue=ClampToQuantum(gamma*GetBluePixelComponent(&pixel));
2184           SetOpacityPixelComponent(q,ClampOpacityPixelComponent(&pixel));
2185           if ((image->colorspace == CMYKColorspace) &&
2186               (resize_image->colorspace == CMYKColorspace))
2187             {
2188               for (i=0; i < n; i++)
2189               {
2190                 j=y*(contribution[n-1].pixel-contribution[0].pixel+1)+
2191                   (contribution[i].pixel-contribution[0].pixel);
2192                 alpha=contribution[i].weight*QuantumScale*
2193                   GetAlphaPixelComponent(p+j);
2194                 pixel.index+=alpha*indexes[j];
2195               }
2196               resize_indexes[y]=(IndexPacket) ClampToQuantum(gamma*
2197                GetIndexPixelComponent(&pixel));
2198             }
2199         }
2200       if ((resize_image->storage_class == PseudoClass) &&
2201           (image->storage_class == PseudoClass))
2202         {
2203           i=(ssize_t) (MagickMin(MagickMax(center,(double) start),(double) stop-
2204             1.0)+0.5);
2205           j=y*(contribution[n-1].pixel-contribution[0].pixel+1)+
2206             (contribution[i-start].pixel-contribution[0].pixel);
2207           resize_indexes[y]=indexes[j];
2208         }
2209       q++;
2210     }
2211     if (SyncCacheViewAuthenticPixels(resize_view,exception) == MagickFalse)
2212       status=MagickFalse;
2213     if (image->progress_monitor != (MagickProgressMonitor) NULL)
2214       {
2215         MagickBooleanType
2216           proceed;
2217
2218 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2219   #pragma omp critical (MagickCore_HorizontalFilter)
2220 #endif
2221         proceed=SetImageProgress(image,ResizeImageTag,(*offset)++,span);
2222         if (proceed == MagickFalse)
2223           status=MagickFalse;
2224       }
2225   }
2226   resize_view=DestroyCacheView(resize_view);
2227   image_view=DestroyCacheView(image_view);
2228   contributions=DestroyContributionThreadSet(contributions);
2229   return(status);
2230 }
2231
2232 static MagickBooleanType VerticalFilter(const ResizeFilter *resize_filter,
2233   const Image *image,Image *resize_image,const MagickRealType y_factor,
2234   const MagickSizeType span,MagickOffsetType *offset,ExceptionInfo *exception)
2235 {
2236   CacheView
2237     *image_view,
2238     *resize_view;
2239
2240   ClassType
2241     storage_class;
2242
2243   ContributionInfo
2244     **restrict contributions;
2245
2246   MagickBooleanType
2247     status;
2248
2249   MagickPixelPacket
2250     zero;
2251
2252   MagickRealType
2253     scale,
2254     support;
2255
2256   ssize_t
2257     y;
2258
2259   /*
2260     Apply filter to resize vertically from image to resize image.
2261   */
2262   scale=MagickMax(1.0/y_factor+MagickEpsilon,1.0);
2263   support=scale*GetResizeFilterSupport(resize_filter);
2264   storage_class=support > 0.5 ? DirectClass : image->storage_class;
2265   if (SetImageStorageClass(resize_image,storage_class) == MagickFalse)
2266     {
2267       InheritException(exception,&resize_image->exception);
2268       return(MagickFalse);
2269     }
2270   if (support < 0.5)
2271     {
2272       /*
2273         Support too small even for nearest neighbour: Reduce to point
2274         sampling.
2275       */
2276       support=(MagickRealType) 0.5;
2277       scale=1.0;
2278     }
2279   contributions=AcquireContributionThreadSet((size_t) (2.0*support+3.0));
2280   if (contributions == (ContributionInfo **) NULL)
2281     {
2282       (void) ThrowMagickException(exception,GetMagickModule(),
2283         ResourceLimitError,"MemoryAllocationFailed","`%s'",image->filename);
2284       return(MagickFalse);
2285     }
2286   status=MagickTrue;
2287   scale=1.0/scale;
2288   (void) ResetMagickMemory(&zero,0,sizeof(zero));
2289   image_view=AcquireCacheView(image);
2290   resize_view=AcquireCacheView(resize_image);
2291 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2292   #pragma omp parallel for shared(status)
2293 #endif
2294   for (y=0; y < (ssize_t) resize_image->rows; y++)
2295   {
2296     MagickRealType
2297       center,
2298       density;
2299
2300     register const IndexPacket
2301       *restrict indexes;
2302
2303     register const PixelPacket
2304       *restrict p;
2305
2306     register ContributionInfo
2307       *restrict contribution;
2308
2309     register IndexPacket
2310       *restrict resize_indexes;
2311
2312     register PixelPacket
2313       *restrict q;
2314
2315     register ssize_t
2316       x;
2317
2318     ssize_t
2319       n,
2320       start,
2321       stop;
2322
2323     if (status == MagickFalse)
2324       continue;
2325     center=(MagickRealType) (y+0.5)/y_factor;
2326     start=(ssize_t) MagickMax(center-support+0.5,0.0);
2327     stop=(ssize_t) MagickMin(center+support+0.5,(double) image->rows);
2328     density=0.0;
2329     contribution=contributions[GetOpenMPThreadId()];
2330     for (n=0; n < (stop-start); n++)
2331     {
2332       contribution[n].pixel=start+n;
2333       contribution[n].weight=GetResizeFilterWeight(resize_filter,scale*
2334         ((MagickRealType) (start+n)-center+0.5));
2335       density+=contribution[n].weight;
2336     }
2337     if ((density != 0.0) && (density != 1.0))
2338       {
2339         register ssize_t
2340           i;
2341
2342         /*
2343           Normalize.
2344         */
2345         density=1.0/density;
2346         for (i=0; i < n; i++)
2347           contribution[i].weight*=density;
2348       }
2349     p=GetCacheViewVirtualPixels(image_view,0,contribution[0].pixel,
2350       image->columns,(size_t) (contribution[n-1].pixel-contribution[0].pixel+1),
2351       exception);
2352     q=QueueCacheViewAuthenticPixels(resize_view,0,y,resize_image->columns,1,
2353       exception);
2354     if ((p == (const PixelPacket *) NULL) || (q == (PixelPacket *) NULL))
2355       {
2356         status=MagickFalse;
2357         continue;
2358       }
2359     indexes=GetCacheViewVirtualIndexQueue(image_view);
2360     resize_indexes=GetCacheViewAuthenticIndexQueue(resize_view);
2361     for (x=0; x < (ssize_t) resize_image->columns; x++)
2362     {
2363       MagickPixelPacket
2364         pixel;
2365
2366       MagickRealType
2367         alpha;
2368
2369       register ssize_t
2370         i;
2371
2372       ssize_t
2373         j;
2374
2375       pixel=zero;
2376       if (image->matte == MagickFalse)
2377         {
2378           for (i=0; i < n; i++)
2379           {
2380             j=(ssize_t) ((contribution[i].pixel-contribution[0].pixel)*
2381               image->columns+x);
2382             alpha=contribution[i].weight;
2383             pixel.red+=alpha*(p+j)->red;
2384             pixel.green+=alpha*(p+j)->green;
2385             pixel.blue+=alpha*(p+j)->blue;
2386             pixel.opacity+=alpha*(p+j)->opacity;
2387           }
2388           SetRedPixelComponent(q,ClampRedPixelComponent(&pixel));
2389           SetGreenPixelComponent(q,ClampGreenPixelComponent(&pixel));
2390           SetBluePixelComponent(q,ClampBluePixelComponent(&pixel));
2391           SetOpacityPixelComponent(q,ClampOpacityPixelComponent(&pixel));
2392           if ((image->colorspace == CMYKColorspace) &&
2393               (resize_image->colorspace == CMYKColorspace))
2394             {
2395               for (i=0; i < n; i++)
2396               {
2397                 j=(ssize_t) ((contribution[i].pixel-contribution[0].pixel)*
2398                   image->columns+x);
2399                 alpha=contribution[i].weight;
2400                 pixel.index+=alpha*indexes[j];
2401               }
2402               resize_indexes[x]=(IndexPacket) ClampToQuantum(pixel.index);
2403             }
2404         }
2405       else
2406         {
2407           MagickRealType
2408             gamma;
2409
2410           gamma=0.0;
2411           for (i=0; i < n; i++)
2412           {
2413             j=(ssize_t) ((contribution[i].pixel-contribution[0].pixel)*
2414               image->columns+x);
2415             alpha=contribution[i].weight*QuantumScale*
2416               GetAlphaPixelComponent(p+j);
2417             pixel.red+=alpha*(p+j)->red;
2418             pixel.green+=alpha*(p+j)->green;
2419             pixel.blue+=alpha*(p+j)->blue;
2420             pixel.opacity+=contribution[i].weight*(p+j)->opacity;
2421             gamma+=alpha;
2422           }
2423           gamma=1.0/(fabs((double) gamma) <= MagickEpsilon ? 1.0 : gamma);
2424           q->red=ClampToQuantum(gamma*GetRedPixelComponent(&pixel));
2425           q->green=ClampToQuantum(gamma*GetGreenPixelComponent(&pixel));
2426           q->blue=ClampToQuantum(gamma*GetBluePixelComponent(&pixel));
2427           SetOpacityPixelComponent(q,ClampOpacityPixelComponent(&pixel));
2428           if ((image->colorspace == CMYKColorspace) &&
2429               (resize_image->colorspace == CMYKColorspace))
2430             {
2431               for (i=0; i < n; i++)
2432               {
2433                 j=(ssize_t) ((contribution[i].pixel-contribution[0].pixel)*
2434                   image->columns+x);
2435                 alpha=contribution[i].weight*QuantumScale*
2436                   GetAlphaPixelComponent(p+j);
2437                 pixel.index+=alpha*indexes[j];
2438               }
2439               resize_indexes[x]=(IndexPacket) ClampToQuantum(gamma*
2440                 GetIndexPixelComponent(&pixel));
2441             }
2442         }
2443       if ((resize_image->storage_class == PseudoClass) &&
2444           (image->storage_class == PseudoClass))
2445         {
2446           i=(ssize_t) (MagickMin(MagickMax(center,(double) start),(double) stop-
2447             1.0)+0.5);
2448           j=(ssize_t) ((contribution[i-start].pixel-contribution[0].pixel)*
2449             image->columns+x);
2450           resize_indexes[x]=indexes[j];
2451         }
2452       q++;
2453     }
2454     if (SyncCacheViewAuthenticPixels(resize_view,exception) == MagickFalse)
2455       status=MagickFalse;
2456     if (image->progress_monitor != (MagickProgressMonitor) NULL)
2457       {
2458         MagickBooleanType
2459           proceed;
2460
2461 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2462   #pragma omp critical (MagickCore_VerticalFilter)
2463 #endif
2464         proceed=SetImageProgress(image,ResizeImageTag,(*offset)++,span);
2465         if (proceed == MagickFalse)
2466           status=MagickFalse;
2467       }
2468   }
2469   resize_view=DestroyCacheView(resize_view);
2470   image_view=DestroyCacheView(image_view);
2471   contributions=DestroyContributionThreadSet(contributions);
2472   return(status);
2473 }
2474
2475 MagickExport Image *ResizeImage(const Image *image,const size_t columns,
2476   const size_t rows,const FilterTypes filter,const double blur,
2477   ExceptionInfo *exception)
2478 {
2479 #define WorkLoadFactor  0.265
2480
2481   FilterTypes
2482     filter_type;
2483
2484   Image
2485     *filter_image,
2486     *resize_image;
2487
2488   MagickOffsetType
2489     offset;
2490
2491   MagickRealType
2492     x_factor,
2493     y_factor;
2494
2495   MagickSizeType
2496     span;
2497
2498   MagickStatusType
2499     status;
2500
2501   ResizeFilter
2502     *resize_filter;
2503
2504   /*
2505     Acquire resize image.
2506   */
2507   assert(image != (Image *) NULL);
2508   assert(image->signature == MagickSignature);
2509   if (image->debug != MagickFalse)
2510     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2511   assert(exception != (ExceptionInfo *) NULL);
2512   assert(exception->signature == MagickSignature);
2513   if ((columns == 0) || (rows == 0))
2514     ThrowImageException(ImageError,"NegativeOrZeroImageSize");
2515   if ((columns == image->columns) && (rows == image->rows) &&
2516       (filter == UndefinedFilter) && (blur == 1.0))
2517     return(CloneImage(image,0,0,MagickTrue,exception));
2518   resize_image=CloneImage(image,columns,rows,MagickTrue,exception);
2519   if (resize_image == (Image *) NULL)
2520     return(resize_image);
2521   /*
2522     Acquire resize filter.
2523   */
2524   x_factor=(MagickRealType) columns/(MagickRealType) image->columns;
2525   y_factor=(MagickRealType) rows/(MagickRealType) image->rows;
2526   if ((x_factor*y_factor) > WorkLoadFactor)
2527     filter_image=CloneImage(image,columns,image->rows,MagickTrue,exception);
2528   else
2529     filter_image=CloneImage(image,image->columns,rows,MagickTrue,exception);
2530   if (filter_image == (Image *) NULL)
2531     return(DestroyImage(resize_image));
2532   filter_type=LanczosFilter;
2533   if (filter != UndefinedFilter)
2534     filter_type=filter;
2535   else
2536     if ((x_factor == 1.0) && (y_factor == 1.0))
2537       filter_type=PointFilter;
2538     else
2539       if ((image->storage_class == PseudoClass) ||
2540           (image->matte != MagickFalse) || ((x_factor*y_factor) > 1.0))
2541         filter_type=MitchellFilter;
2542   resize_filter=AcquireResizeFilter(image,filter_type,blur,MagickFalse,
2543     exception);
2544   /*
2545     Resize image.
2546   */
2547   offset=0;
2548   if ((x_factor*y_factor) > WorkLoadFactor)
2549     {
2550       span=(MagickSizeType) (filter_image->columns+rows);
2551       status=HorizontalFilter(resize_filter,image,filter_image,x_factor,span,
2552         &offset,exception);
2553       status&=VerticalFilter(resize_filter,filter_image,resize_image,y_factor,
2554         span,&offset,exception);
2555     }
2556   else
2557     {
2558       span=(MagickSizeType) (filter_image->rows+columns);
2559       status=VerticalFilter(resize_filter,image,filter_image,y_factor,span,
2560         &offset,exception);
2561       status&=HorizontalFilter(resize_filter,filter_image,resize_image,x_factor,
2562         span,&offset,exception);
2563     }
2564   /*
2565     Free resources.
2566   */
2567   filter_image=DestroyImage(filter_image);
2568   resize_filter=DestroyResizeFilter(resize_filter);
2569   if ((status == MagickFalse) || (resize_image == (Image *) NULL))
2570     return((Image *) NULL);
2571   resize_image->type=image->type;
2572   return(resize_image);
2573 }
2574 \f
2575 /*
2576 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2577 %                                                                             %
2578 %                                                                             %
2579 %                                                                             %
2580 %   S a m p l e I m a g e                                                     %
2581 %                                                                             %
2582 %                                                                             %
2583 %                                                                             %
2584 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2585 %
2586 %  SampleImage() scales an image to the desired dimensions with pixel
2587 %  sampling.  Unlike other scaling methods, this method does not introduce
2588 %  any additional color into the scaled image.
2589 %
2590 %  The format of the SampleImage method is:
2591 %
2592 %      Image *SampleImage(const Image *image,const size_t columns,
2593 %        const size_t rows,ExceptionInfo *exception)
2594 %
2595 %  A description of each parameter follows:
2596 %
2597 %    o image: the image.
2598 %
2599 %    o columns: the number of columns in the sampled image.
2600 %
2601 %    o rows: the number of rows in the sampled image.
2602 %
2603 %    o exception: return any errors or warnings in this structure.
2604 %
2605 */
2606 MagickExport Image *SampleImage(const Image *image,const size_t columns,
2607   const size_t rows,ExceptionInfo *exception)
2608 {
2609 #define SampleImageTag  "Sample/Image"
2610
2611   CacheView
2612     *image_view,
2613     *sample_view;
2614
2615   Image
2616     *sample_image;
2617
2618   MagickBooleanType
2619     status;
2620
2621   MagickOffsetType
2622     progress;
2623
2624   register ssize_t
2625     x;
2626
2627   ssize_t
2628     *x_offset,
2629     y;
2630
2631   /*
2632     Initialize sampled image attributes.
2633   */
2634   assert(image != (const Image *) NULL);
2635   assert(image->signature == MagickSignature);
2636   if (image->debug != MagickFalse)
2637     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2638   assert(exception != (ExceptionInfo *) NULL);
2639   assert(exception->signature == MagickSignature);
2640   if ((columns == 0) || (rows == 0))
2641     ThrowImageException(ImageError,"NegativeOrZeroImageSize");
2642   if ((columns == image->columns) && (rows == image->rows))
2643     return(CloneImage(image,0,0,MagickTrue,exception));
2644   sample_image=CloneImage(image,columns,rows,MagickTrue,exception);
2645   if (sample_image == (Image *) NULL)
2646     return((Image *) NULL);
2647   /*
2648     Allocate scan line buffer and column offset buffers.
2649   */
2650   x_offset=(ssize_t *) AcquireQuantumMemory((size_t) sample_image->columns,
2651     sizeof(*x_offset));
2652   if (x_offset == (ssize_t *) NULL)
2653     {
2654       sample_image=DestroyImage(sample_image);
2655       ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
2656     }
2657   for (x=0; x < (ssize_t) sample_image->columns; x++)
2658     x_offset[x]=(ssize_t) (((MagickRealType) x+0.5)*image->columns/
2659       sample_image->columns);
2660   /*
2661     Sample each row.
2662   */
2663   status=MagickTrue;
2664   progress=0;
2665   image_view=AcquireCacheView(image);
2666   sample_view=AcquireCacheView(sample_image);
2667 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2668   #pragma omp parallel for schedule(dynamic,4) shared(progress,status)
2669 #endif
2670   for (y=0; y < (ssize_t) sample_image->rows; y++)
2671   {
2672     register const IndexPacket
2673       *restrict indexes;
2674
2675     register const PixelPacket
2676       *restrict p;
2677
2678     register IndexPacket
2679       *restrict sample_indexes;
2680
2681     register PixelPacket
2682       *restrict q;
2683
2684     register ssize_t
2685       x;
2686
2687     ssize_t
2688       y_offset;
2689
2690     if (status == MagickFalse)
2691       continue;
2692     y_offset=(ssize_t) (((MagickRealType) y+0.5)*image->rows/
2693       sample_image->rows);
2694     p=GetCacheViewVirtualPixels(image_view,0,y_offset,image->columns,1,
2695       exception);
2696     q=QueueCacheViewAuthenticPixels(sample_view,0,y,sample_image->columns,1,
2697       exception);
2698     if ((p == (const PixelPacket *) NULL) || (q == (PixelPacket *) NULL))
2699       {
2700         status=MagickFalse;
2701         continue;
2702       }
2703     indexes=GetCacheViewAuthenticIndexQueue(image_view);
2704     sample_indexes=GetCacheViewAuthenticIndexQueue(sample_view);
2705     /*
2706       Sample each column.
2707     */
2708     for (x=0; x < (ssize_t) sample_image->columns; x++)
2709       *q++=p[x_offset[x]];
2710     if ((image->storage_class == PseudoClass) ||
2711         (image->colorspace == CMYKColorspace))
2712       for (x=0; x < (ssize_t) sample_image->columns; x++)
2713         sample_indexes[x]=indexes[x_offset[x]];
2714     if (SyncCacheViewAuthenticPixels(sample_view,exception) == MagickFalse)
2715       status=MagickFalse;
2716     if (image->progress_monitor != (MagickProgressMonitor) NULL)
2717       {
2718         MagickBooleanType
2719           proceed;
2720
2721 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2722         #pragma omp critical (MagickCore_SampleImage)
2723 #endif
2724         proceed=SetImageProgress(image,SampleImageTag,progress++,image->rows);
2725         if (proceed == MagickFalse)
2726           status=MagickFalse;
2727       }
2728   }
2729   image_view=DestroyCacheView(image_view);
2730   sample_view=DestroyCacheView(sample_view);
2731   x_offset=(ssize_t *) RelinquishMagickMemory(x_offset);
2732   sample_image->type=image->type;
2733   return(sample_image);
2734 }
2735 \f
2736 /*
2737 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2738 %                                                                             %
2739 %                                                                             %
2740 %                                                                             %
2741 %   S c a l e I m a g e                                                       %
2742 %                                                                             %
2743 %                                                                             %
2744 %                                                                             %
2745 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2746 %
2747 %  ScaleImage() changes the size of an image to the given dimensions.
2748 %
2749 %  The format of the ScaleImage method is:
2750 %
2751 %      Image *ScaleImage(const Image *image,const size_t columns,
2752 %        const size_t rows,ExceptionInfo *exception)
2753 %
2754 %  A description of each parameter follows:
2755 %
2756 %    o image: the image.
2757 %
2758 %    o columns: the number of columns in the scaled image.
2759 %
2760 %    o rows: the number of rows in the scaled image.
2761 %
2762 %    o exception: return any errors or warnings in this structure.
2763 %
2764 */
2765 MagickExport Image *ScaleImage(const Image *image,const size_t columns,
2766   const size_t rows,ExceptionInfo *exception)
2767 {
2768 #define ScaleImageTag  "Scale/Image"
2769
2770   CacheView
2771     *image_view,
2772     *scale_view;
2773
2774   Image
2775     *scale_image;
2776
2777   MagickBooleanType
2778     next_column,
2779     next_row,
2780     proceed;
2781
2782   MagickPixelPacket
2783     pixel,
2784     *scale_scanline,
2785     *scanline,
2786     *x_vector,
2787     *y_vector,
2788     zero;
2789
2790   PointInfo
2791     scale,
2792     span;
2793
2794   register ssize_t
2795     i;
2796
2797   ssize_t
2798     number_rows,
2799     y;
2800
2801   /*
2802     Initialize scaled image attributes.
2803   */
2804   assert(image != (const Image *) NULL);
2805   assert(image->signature == MagickSignature);
2806   if (image->debug != MagickFalse)
2807     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2808   assert(exception != (ExceptionInfo *) NULL);
2809   assert(exception->signature == MagickSignature);
2810   if ((columns == 0) || (rows == 0))
2811     return((Image *) NULL);
2812   if ((columns == image->columns) && (rows == image->rows))
2813     return(CloneImage(image,0,0,MagickTrue,exception));
2814   scale_image=CloneImage(image,columns,rows,MagickTrue,exception);
2815   if (scale_image == (Image *) NULL)
2816     return((Image *) NULL);
2817   if (SetImageStorageClass(scale_image,DirectClass) == MagickFalse)
2818     {
2819       InheritException(exception,&scale_image->exception);
2820       scale_image=DestroyImage(scale_image);
2821       return((Image *) NULL);
2822     }
2823   /*
2824     Allocate memory.
2825   */
2826   x_vector=(MagickPixelPacket *) AcquireQuantumMemory((size_t) image->columns,
2827     sizeof(*x_vector));
2828   scanline=x_vector;
2829   if (image->rows != scale_image->rows)
2830     scanline=(MagickPixelPacket *) AcquireQuantumMemory((size_t) image->columns,
2831       sizeof(*scanline));
2832   scale_scanline=(MagickPixelPacket *) AcquireQuantumMemory((size_t)
2833     scale_image->columns,sizeof(*scale_scanline));
2834   y_vector=(MagickPixelPacket *) AcquireQuantumMemory((size_t) image->columns,
2835     sizeof(*y_vector));
2836   if ((scanline == (MagickPixelPacket *) NULL) ||
2837       (scale_scanline == (MagickPixelPacket *) NULL) ||
2838       (x_vector == (MagickPixelPacket *) NULL) ||
2839       (y_vector == (MagickPixelPacket *) NULL))
2840     {
2841       scale_image=DestroyImage(scale_image);
2842       ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
2843     }
2844   /*
2845     Scale image.
2846   */
2847   number_rows=0;
2848   next_row=MagickTrue;
2849   span.y=1.0;
2850   scale.y=(double) scale_image->rows/(double) image->rows;
2851   (void) ResetMagickMemory(y_vector,0,(size_t) image->columns*
2852     sizeof(*y_vector));
2853   GetMagickPixelPacket(image,&pixel);
2854   (void) ResetMagickMemory(&zero,0,sizeof(zero));
2855   i=0;
2856   image_view=AcquireCacheView(image);
2857   scale_view=AcquireCacheView(scale_image);
2858   for (y=0; y < (ssize_t) scale_image->rows; y++)
2859   {
2860     register const IndexPacket
2861       *restrict indexes;
2862
2863     register const PixelPacket
2864       *restrict p;
2865
2866     register IndexPacket
2867       *restrict scale_indexes;
2868
2869     register MagickPixelPacket
2870       *restrict s,
2871       *restrict t;
2872
2873     register PixelPacket
2874       *restrict q;
2875
2876     register ssize_t
2877       x;
2878
2879     q=QueueCacheViewAuthenticPixels(scale_view,0,y,scale_image->columns,1,
2880       exception);
2881     if (q == (PixelPacket *) NULL)
2882       break;
2883     scale_indexes=GetAuthenticIndexQueue(scale_image);
2884     if (scale_image->rows == image->rows)
2885       {
2886         /*
2887           Read a new scanline.
2888         */
2889         p=GetCacheViewVirtualPixels(image_view,0,i++,image->columns,1,
2890           exception);
2891         if (p == (const PixelPacket *) NULL)
2892           break;
2893         indexes=GetCacheViewVirtualIndexQueue(image_view);
2894         for (x=0; x < (ssize_t) image->columns; x++)
2895         {
2896           x_vector[x].red=(MagickRealType) GetRedPixelComponent(p);
2897           x_vector[x].green=(MagickRealType) GetGreenPixelComponent(p);
2898           x_vector[x].blue=(MagickRealType) GetBluePixelComponent(p);
2899           if (image->matte != MagickFalse)
2900             x_vector[x].opacity=(MagickRealType) GetOpacityPixelComponent(p);
2901           if (indexes != (IndexPacket *) NULL)
2902             x_vector[x].index=(MagickRealType) indexes[x];
2903           p++;
2904         }
2905       }
2906     else
2907       {
2908         /*
2909           Scale Y direction.
2910         */
2911         while (scale.y < span.y)
2912         {
2913           if ((next_row != MagickFalse) &&
2914               (number_rows < (ssize_t) image->rows))
2915             {
2916               /*
2917                 Read a new scanline.
2918               */
2919               p=GetCacheViewVirtualPixels(image_view,0,i++,image->columns,1,
2920                 exception);
2921               if (p == (const PixelPacket *) NULL)
2922                 break;
2923               indexes=GetCacheViewVirtualIndexQueue(image_view);
2924               for (x=0; x < (ssize_t) image->columns; x++)
2925               {
2926                 x_vector[x].red=(MagickRealType) GetRedPixelComponent(p);
2927                 x_vector[x].green=(MagickRealType) GetGreenPixelComponent(p);
2928                 x_vector[x].blue=(MagickRealType) GetBluePixelComponent(p);
2929                 if (image->matte != MagickFalse)
2930                   x_vector[x].opacity=(MagickRealType)
2931                     GetOpacityPixelComponent(p);
2932                 if (indexes != (IndexPacket *) NULL)
2933                   x_vector[x].index=(MagickRealType) indexes[x];
2934                 p++;
2935               }
2936               number_rows++;
2937             }
2938           for (x=0; x < (ssize_t) image->columns; x++)
2939           {
2940             y_vector[x].red+=scale.y*x_vector[x].red;
2941             y_vector[x].green+=scale.y*x_vector[x].green;
2942             y_vector[x].blue+=scale.y*x_vector[x].blue;
2943             if (scale_image->matte != MagickFalse)
2944               y_vector[x].opacity+=scale.y*x_vector[x].opacity;
2945             if (scale_indexes != (IndexPacket *) NULL)
2946               y_vector[x].index+=scale.y*x_vector[x].index;
2947           }
2948           span.y-=scale.y;
2949           scale.y=(double) scale_image->rows/(double) image->rows;
2950           next_row=MagickTrue;
2951         }
2952         if ((next_row != MagickFalse) && (number_rows < (ssize_t) image->rows))
2953           {
2954             /*
2955               Read a new scanline.
2956             */
2957             p=GetCacheViewVirtualPixels(image_view,0,i++,image->columns,1,
2958               exception);
2959             if (p == (const PixelPacket *) NULL)
2960               break;
2961             indexes=GetCacheViewVirtualIndexQueue(image_view);
2962             for (x=0; x < (ssize_t) image->columns; x++)
2963             {
2964               x_vector[x].red=(MagickRealType) GetRedPixelComponent(p);
2965               x_vector[x].green=(MagickRealType) GetGreenPixelComponent(p);
2966               x_vector[x].blue=(MagickRealType) GetBluePixelComponent(p);
2967               if (image->matte != MagickFalse)
2968                 x_vector[x].opacity=(MagickRealType)
2969                   GetOpacityPixelComponent(p);
2970               if (indexes != (IndexPacket *) NULL)
2971                 x_vector[x].index=(MagickRealType) indexes[x];
2972               p++;
2973             }
2974             number_rows++;
2975             next_row=MagickFalse;
2976           }
2977         s=scanline;
2978         for (x=0; x < (ssize_t) image->columns; x++)
2979         {
2980           pixel.red=y_vector[x].red+span.y*x_vector[x].red;
2981           pixel.green=y_vector[x].green+span.y*x_vector[x].green;
2982           pixel.blue=y_vector[x].blue+span.y*x_vector[x].blue;
2983           if (image->matte != MagickFalse)
2984             pixel.opacity=y_vector[x].opacity+span.y*x_vector[x].opacity;
2985           if (scale_indexes != (IndexPacket *) NULL)
2986             pixel.index=y_vector[x].index+span.y*x_vector[x].index;
2987           s->red=pixel.red;
2988           s->green=pixel.green;
2989           s->blue=pixel.blue;
2990           if (scale_image->matte != MagickFalse)
2991             s->opacity=pixel.opacity;
2992           if (scale_indexes != (IndexPacket *) NULL)
2993             s->index=pixel.index;
2994           s++;
2995           y_vector[x]=zero;
2996         }
2997         scale.y-=span.y;
2998         if (scale.y <= 0)
2999           {
3000             scale.y=(double) scale_image->rows/(double) image->rows;
3001             next_row=MagickTrue;
3002           }
3003         span.y=1.0;
3004       }
3005     if (scale_image->columns == image->columns)
3006       {
3007         /*
3008           Transfer scanline to scaled image.
3009         */
3010         s=scanline;
3011         for (x=0; x < (ssize_t) scale_image->columns; x++)
3012         {
3013           q->red=ClampToQuantum(s->red);
3014           q->green=ClampToQuantum(s->green);
3015           q->blue=ClampToQuantum(s->blue);
3016           if (scale_image->matte != MagickFalse)
3017             q->opacity=ClampToQuantum(s->opacity);
3018           if (scale_indexes != (IndexPacket *) NULL)
3019             scale_indexes[x]=(IndexPacket) ClampToQuantum(s->index);
3020           q++;
3021           s++;
3022         }
3023       }
3024     else
3025       {
3026         /*
3027           Scale X direction.
3028         */
3029         pixel=zero;
3030         next_column=MagickFalse;
3031         span.x=1.0;
3032         s=scanline;
3033         t=scale_scanline;
3034         for (x=0; x < (ssize_t) image->columns; x++)
3035         {
3036           scale.x=(double) scale_image->columns/(double) image->columns;
3037           while (scale.x >= span.x)
3038           {
3039             if (next_column != MagickFalse)
3040               {
3041                 pixel=zero;
3042                 t++;
3043               }
3044             pixel.red+=span.x*s->red;
3045             pixel.green+=span.x*s->green;
3046             pixel.blue+=span.x*s->blue;
3047             if (image->matte != MagickFalse)
3048               pixel.opacity+=span.x*s->opacity;
3049             if (scale_indexes != (IndexPacket *) NULL)
3050               pixel.index+=span.x*s->index;
3051             t->red=pixel.red;
3052             t->green=pixel.green;
3053             t->blue=pixel.blue;
3054             if (scale_image->matte != MagickFalse)
3055               t->opacity=pixel.opacity;
3056             if (scale_indexes != (IndexPacket *) NULL)
3057               t->index=pixel.index;
3058             scale.x-=span.x;
3059             span.x=1.0;
3060             next_column=MagickTrue;
3061           }
3062         if (scale.x > 0)
3063           {
3064             if (next_column != MagickFalse)
3065               {
3066                 pixel=zero;
3067                 next_column=MagickFalse;
3068                 t++;
3069               }
3070             pixel.red+=scale.x*s->red;
3071             pixel.green+=scale.x*s->green;
3072             pixel.blue+=scale.x*s->blue;
3073             if (scale_image->matte != MagickFalse)
3074               pixel.opacity+=scale.x*s->opacity;
3075             if (scale_indexes != (IndexPacket *) NULL)
3076               pixel.index+=scale.x*s->index;
3077             span.x-=scale.x;
3078           }
3079         s++;
3080       }
3081       if (span.x > 0)
3082         {
3083           s--;
3084           pixel.red+=span.x*s->red;
3085           pixel.green+=span.x*s->green;
3086           pixel.blue+=span.x*s->blue;
3087           if (scale_image->matte != MagickFalse)
3088             pixel.opacity+=span.x*s->opacity;
3089           if (scale_indexes != (IndexPacket *) NULL)
3090             pixel.index+=span.x*s->index;
3091         }
3092       if ((next_column == MagickFalse) &&
3093           ((ssize_t) (t-scale_scanline) < (ssize_t) scale_image->columns))
3094         {
3095           t->red=pixel.red;
3096           t->green=pixel.green;
3097           t->blue=pixel.blue;
3098           if (scale_image->matte != MagickFalse)
3099             t->opacity=pixel.opacity;
3100           if (scale_indexes != (IndexPacket *) NULL)
3101             t->index=pixel.index;
3102         }
3103       /*
3104         Transfer scanline to scaled image.
3105       */
3106       t=scale_scanline;
3107       for (x=0; x < (ssize_t) scale_image->columns; x++)
3108       {
3109         q->red=ClampToQuantum(t->red);
3110         q->green=ClampToQuantum(t->green);
3111         q->blue=ClampToQuantum(t->blue);
3112         if (scale_image->matte != MagickFalse)
3113           q->opacity=ClampToQuantum(t->opacity);
3114         if (scale_indexes != (IndexPacket *) NULL)
3115           scale_indexes[x]=(IndexPacket) ClampToQuantum(t->index);
3116         t++;
3117         q++;
3118       }
3119     }
3120     if (SyncCacheViewAuthenticPixels(scale_view,exception) == MagickFalse)
3121       break;
3122     proceed=SetImageProgress(image,ScaleImageTag,(MagickOffsetType) y,
3123       image->rows);
3124     if (proceed == MagickFalse)
3125       break;
3126   }
3127   scale_view=DestroyCacheView(scale_view);
3128   image_view=DestroyCacheView(image_view);
3129   /*
3130     Free allocated memory.
3131   */
3132   y_vector=(MagickPixelPacket *) RelinquishMagickMemory(y_vector);
3133   scale_scanline=(MagickPixelPacket *) RelinquishMagickMemory(scale_scanline);
3134   if (scale_image->rows != image->rows)
3135     scanline=(MagickPixelPacket *) RelinquishMagickMemory(scanline);
3136   x_vector=(MagickPixelPacket *) RelinquishMagickMemory(x_vector);
3137   scale_image->type=image->type;
3138   return(scale_image);
3139 }
3140 \f
3141 /*
3142 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3143 %                                                                             %
3144 %                                                                             %
3145 %                                                                             %
3146 +   S e t R e s i z e F i l t e r S u p p o r t                               %
3147 %                                                                             %
3148 %                                                                             %
3149 %                                                                             %
3150 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3151 %
3152 %  SetResizeFilterSupport() specifies which IR filter to use to window
3153 %
3154 %  The format of the SetResizeFilterSupport method is:
3155 %
3156 %      void SetResizeFilterSupport(ResizeFilter *resize_filter,
3157 %        const MagickRealType support)
3158 %
3159 %  A description of each parameter follows:
3160 %
3161 %    o resize_filter: the resize filter.
3162 %
3163 %    o support: the filter spport radius.
3164 %
3165 */
3166 MagickExport void SetResizeFilterSupport(ResizeFilter *resize_filter,
3167   const MagickRealType support)
3168 {
3169   assert(resize_filter != (ResizeFilter *) NULL);
3170   assert(resize_filter->signature == MagickSignature);
3171   resize_filter->support=support;
3172 }
3173 \f
3174 /*
3175 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3176 %                                                                             %
3177 %                                                                             %
3178 %                                                                             %
3179 %   T h u m b n a i l I m a g e                                               %
3180 %                                                                             %
3181 %                                                                             %
3182 %                                                                             %
3183 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3184 %
3185 %  ThumbnailImage() changes the size of an image to the given dimensions and
3186 %  removes any associated profiles.  The goal is to produce small low cost
3187 %  thumbnail images suited for display on the Web.
3188 %
3189 %  The format of the ThumbnailImage method is:
3190 %
3191 %      Image *ThumbnailImage(const Image *image,const size_t columns,
3192 %        const size_t rows,ExceptionInfo *exception)
3193 %
3194 %  A description of each parameter follows:
3195 %
3196 %    o image: the image.
3197 %
3198 %    o columns: the number of columns in the scaled image.
3199 %
3200 %    o rows: the number of rows in the scaled image.
3201 %
3202 %    o exception: return any errors or warnings in this structure.
3203 %
3204 */
3205 MagickExport Image *ThumbnailImage(const Image *image,const size_t columns,
3206   const size_t rows,ExceptionInfo *exception)
3207 {
3208 #define SampleFactor  5
3209
3210   char
3211     value[MaxTextExtent];
3212
3213   const char
3214     *name;
3215
3216   Image
3217     *thumbnail_image;
3218
3219   MagickRealType
3220     x_factor,
3221     y_factor;
3222
3223   size_t
3224     version;
3225
3226   struct stat
3227     attributes;
3228
3229   assert(image != (Image *) NULL);
3230   assert(image->signature == MagickSignature);
3231   if (image->debug != MagickFalse)
3232     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
3233   assert(exception != (ExceptionInfo *) NULL);
3234   assert(exception->signature == MagickSignature);
3235   x_factor=(MagickRealType) columns/(MagickRealType) image->columns;
3236   y_factor=(MagickRealType) rows/(MagickRealType) image->rows;
3237   if ((x_factor*y_factor) > 0.1)
3238     thumbnail_image=ResizeImage(image,columns,rows,image->filter,image->blur,
3239       exception);
3240   else
3241     if (((SampleFactor*columns) < 128) || ((SampleFactor*rows) < 128))
3242       thumbnail_image=ResizeImage(image,columns,rows,image->filter,
3243         image->blur,exception);
3244     else
3245       {
3246         Image
3247           *sample_image;
3248
3249         sample_image=SampleImage(image,SampleFactor*columns,SampleFactor*rows,
3250           exception);
3251         if (sample_image == (Image *) NULL)
3252           return((Image *) NULL);
3253         thumbnail_image=ResizeImage(sample_image,columns,rows,image->filter,
3254           image->blur,exception);
3255         sample_image=DestroyImage(sample_image);
3256       }
3257   if (thumbnail_image == (Image *) NULL)
3258     return(thumbnail_image);
3259   (void) ParseAbsoluteGeometry("0x0+0+0",&thumbnail_image->page);
3260   if (thumbnail_image->matte == MagickFalse)
3261     (void) SetImageAlphaChannel(thumbnail_image,OpaqueAlphaChannel);
3262   thumbnail_image->depth=8;
3263   thumbnail_image->interlace=NoInterlace;
3264   /*
3265     Strip all profiles except color profiles.
3266   */
3267   ResetImageProfileIterator(thumbnail_image);
3268   for (name=GetNextImageProfile(thumbnail_image); name != (const char *) NULL; )
3269   {
3270     if ((LocaleCompare(name,"icc") != 0) && (LocaleCompare(name,"icm") != 0))
3271      {
3272        (void) DeleteImageProfile(thumbnail_image,name);
3273        ResetImageProfileIterator(thumbnail_image);
3274      }
3275     name=GetNextImageProfile(thumbnail_image);
3276   }
3277   (void) DeleteImageProperty(thumbnail_image,"comment");
3278   (void) CopyMagickString(value,image->magick_filename,MaxTextExtent);
3279   if (strstr(image->magick_filename,"//") == (char *) NULL)
3280     (void) FormatMagickString(value,MaxTextExtent,"file://%s",
3281       image->magick_filename);
3282   (void) SetImageProperty(thumbnail_image,"Thumb::URI",value);
3283   (void) CopyMagickString(value,image->magick_filename,MaxTextExtent);
3284   if (GetPathAttributes(image->filename,&attributes) != MagickFalse)
3285     {
3286       (void) FormatMagickString(value,MaxTextExtent,"%.20g",(double)
3287         attributes.st_mtime);
3288       (void) SetImageProperty(thumbnail_image,"Thumb::MTime",value);
3289     }
3290   (void) FormatMagickString(value,MaxTextExtent,"%.20g",(double)
3291     attributes.st_mtime);
3292   (void) FormatMagickSize(GetBlobSize(image),MagickFalse,value);
3293   (void) ConcatenateMagickString(value,"B",MaxTextExtent);
3294   (void) SetImageProperty(thumbnail_image,"Thumb::Size",value);
3295   (void) FormatMagickString(value,MaxTextExtent,"image/%s",image->magick);
3296   LocaleLower(value);
3297   (void) SetImageProperty(thumbnail_image,"Thumb::Mimetype",value);
3298   (void) SetImageProperty(thumbnail_image,"software",
3299     GetMagickVersion(&version));
3300   (void) FormatMagickString(value,MaxTextExtent,"%.20g",(double)
3301     image->magick_columns);
3302   (void) SetImageProperty(thumbnail_image,"Thumb::Image::Width",value);
3303   (void) FormatMagickString(value,MaxTextExtent,"%.20g",(double)
3304     image->magick_rows);
3305   (void) SetImageProperty(thumbnail_image,"Thumb::Image::height",value);
3306   (void) FormatMagickString(value,MaxTextExtent,"%.20g",(double)
3307     GetImageListLength(image));
3308   (void) SetImageProperty(thumbnail_image,"Thumb::Document::Pages",value);
3309   return(thumbnail_image);
3310 }