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