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