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