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

   1 /*
   2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   3 %                                                                             %
   4 %                                                                             %
   5 %                                                                             %
   6 %        SSSSS  TTTTT   AAA   TTTTT  IIIII  SSSSS  TTTTT  IIIII   CCCC        %
   7 %        SS       T    A   A    T      I    SS       T      I    C            %
   8 %         SSS     T    AAAAA    T      I     SSS     T      I    C            %
   9 %           SS    T    A   A    T      I       SS    T      I    C            %
  10 %        SSSSS    T    A   A    T    IIIII  SSSSS    T    IIIII   CCCC        %
  11 %                                                                             %
  12 %                                                                             %
  13 %                          MagickCore Image Methods                           %
  14 %                                                                             %
  15 %                              Software Design                                %
  16 %                                John Cristy                                  %
  17 %                                 July 1992                                   %
  18 %                                                                             %
  19 %                                                                             %
  20 %  Copyright 1999-2009 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 */
  39 \f
  40 /*
  41   Include declarations.
  42 */
  43 #include "magick/studio.h"
  44 #include "magick/property.h"
  45 #include "magick/animate.h"
  46 #include "magick/blob.h"
  47 #include "magick/blob-private.h"
  48 #include "magick/cache.h"
  49 #include "magick/cache-private.h"
  50 #include "magick/cache-view.h"
  51 #include "magick/client.h"
  52 #include "magick/color.h"
  53 #include "magick/color-private.h"
  54 #include "magick/colorspace.h"
  55 #include "magick/colorspace-private.h"
  56 #include "magick/composite.h"
  57 #include "magick/composite-private.h"
  58 #include "magick/compress.h"
  59 #include "magick/constitute.h"
  60 #include "magick/deprecate.h"
  61 #include "magick/display.h"
  62 #include "magick/draw.h"
  63 #include "magick/enhance.h"
  64 #include "magick/exception.h"
  65 #include "magick/exception-private.h"
  66 #include "magick/gem.h"
  67 #include "magick/geometry.h"
  68 #include "magick/list.h"
  69 #include "magick/image-private.h"
  70 #include "magick/magic.h"
  71 #include "magick/magick.h"
  72 #include "magick/memory_.h"
  73 #include "magick/module.h"
  74 #include "magick/monitor.h"
  75 #include "magick/option.h"
  76 #include "magick/paint.h"
  77 #include "magick/pixel-private.h"
  78 #include "magick/profile.h"
  79 #include "magick/quantize.h"
  80 #include "magick/random_.h"
  81 #include "magick/segment.h"
  82 #include "magick/semaphore.h"
  83 #include "magick/signature-private.h"
  84 #include "magick/statistic.h"
  85 #include "magick/string_.h"
  86 #include "magick/thread-private.h"
  87 #include "magick/timer.h"
  88 #include "magick/utility.h"
  89 #include "magick/version.h"
  90 \f
  91 /*
  92 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  93 %                                                                             %
  94 %                                                                             %
  95 %                                                                             %
  96 +   G e t I m a g e B o u n d i n g B o x                                     %
  97 %                                                                             %
  98 %                                                                             %
  99 %                                                                             %
 100 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 101 %
 102 %  GetImageBoundingBox() returns the bounding box of an image canvas.
 103 %
 104 %  The format of the GetImageBoundingBox method is:
 105 %
 106 %      RectangleInfo GetImageBoundingBox(const Image *image,
 107 %        ExceptionInfo *exception)
 108 %
 109 %  A description of each parameter follows:
 110 %
 111 %    o bounds: Method GetImageBoundingBox returns the bounding box of an
 112 %      image canvas.
 113 %
 114 %    o image: the image.
 115 %
 116 %    o exception: return any errors or warnings in this structure.
 117 %
 118 */
 119 MagickExport RectangleInfo GetImageBoundingBox(const Image *image,
 120   ExceptionInfo *exception)
 121 {
 122   long
 123     y;
 124
 125   MagickBooleanType
 126     status;
 127
 128   MagickPixelPacket
 129     target[3],
 130     zero;
 131
 132   RectangleInfo
 133     bounds;
 134
 135   register const PixelPacket
 136     *p;
 137
 138   CacheView
 139     *image_view;
 140
 141   assert(image != (Image *) NULL);
 142   assert(image->signature == MagickSignature);
 143   if (image->debug != MagickFalse)
 144     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
 145   bounds.width=0;
 146   bounds.height=0;
 147   bounds.x=(long) image->columns;
 148   bounds.y=(long) image->rows;
 149   GetMagickPixelPacket(image,&target[0]);
 150   image_view=AcquireCacheView(image);
 151   p=GetCacheViewVirtualPixels(image_view,0,0,1,1,exception);
 152   if (p == (const PixelPacket *) NULL)
 153     {
 154       image_view=DestroyCacheView(image_view);
 155       return(bounds);
 156     }
 157   SetMagickPixelPacket(image,p,GetCacheViewAuthenticIndexQueue(image_view),
 158     &target[0]);
 159   GetMagickPixelPacket(image,&target[1]);
 160   p=GetCacheViewVirtualPixels(image_view,(long) image->columns-1,0,1,1,
 161     exception);
 162   SetMagickPixelPacket(image,p,GetCacheViewAuthenticIndexQueue(image_view),
 163     &target[1]);
 164   GetMagickPixelPacket(image,&target[2]);
 165   p=GetCacheViewVirtualPixels(image_view,0,(long) image->rows-1,1,1,exception);
 166   SetMagickPixelPacket(image,p,GetCacheViewAuthenticIndexQueue(image_view),
 167     &target[2]);
 168   status=MagickTrue;
 169   GetMagickPixelPacket(image,&zero);
 170 #if defined(MAGICKCORE_OPENMP_SUPPORT)
 171   #pragma omp parallel for schedule(dynamic,4) shared(status)
 172 #endif
 173   for (y=0; y < (long) image->rows; y++)
 174   {
 175     MagickPixelPacket
 176       pixel;
 177
 178     RectangleInfo
 179       bounding_box;
 180
 181     register const IndexPacket
 182       *__restrict indexes;
 183
 184     register const PixelPacket
 185       *__restrict p;
 186
 187     register long
 188       x;
 189
 190     if (status == MagickFalse)
 191       continue;
 192 #if defined(HAVE_OPENMP)
 193 #  pragma omp critical (MagickCore_GetImageBoundingBox)
 194 #endif
 195     bounding_box=bounds;
 196     p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
 197     if (p == (const PixelPacket *) NULL)
 198       {
 199         status=MagickFalse;
 200         continue;
 201       }
 202     indexes=GetCacheViewVirtualIndexQueue(image_view);
 203     pixel=zero;
 204     for (x=0; x < (long) image->columns; x++)
 205     {
 206       SetMagickPixelPacket(image,p,indexes+x,&pixel);
 207       if ((x < bounding_box.x) &&
 208           (IsMagickColorSimilar(&pixel,&target[0]) == MagickFalse))
 209         bounding_box.x=x;
 210       if ((x > (long) bounding_box.width) &&
 211           (IsMagickColorSimilar(&pixel,&target[1]) == MagickFalse))
 212         bounding_box.width=(unsigned long) x;
 213       if ((y < bounding_box.y) &&
 214           (IsMagickColorSimilar(&pixel,&target[0]) == MagickFalse))
 215         bounding_box.y=y;
 216       if ((y > (long) bounding_box.height) &&
 217           (IsMagickColorSimilar(&pixel,&target[2]) == MagickFalse))
 218         bounding_box.height=(unsigned long) y;
 219       p++;
 220     }
 221 #if defined(HAVE_OPENMP)
 222 #  pragma omp critical (MagickCore_GetImageBoundingBox)
 223 #endif
 224     {
 225       if (bounding_box.x < bounds.x)
 226         bounds.x=bounding_box.x;
 227       if (bounding_box.y < bounds.y)
 228         bounds.y=bounding_box.y;
 229       if (bounding_box.width > bounds.width)
 230         bounds.width=bounding_box.width;
 231       if (bounding_box.height > bounds.height)
 232         bounds.height=bounding_box.height;
 233     }
 234   }
 235   image_view=DestroyCacheView(image_view);
 236   if ((bounds.width == 0) || (bounds.height == 0))
 237     (void) ThrowMagickException(exception,GetMagickModule(),OptionWarning,
 238       "GeometryDoesNotContainImage","`%s'",image->filename);
 239   else
 240     {
 241       bounds.width-=(bounds.x-1);
 242       bounds.height-=(bounds.y-1);
 243     }
 244   return(bounds);
 245 }
 246 \f
 247 /*
 248 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 249 %                                                                             %
 250 %                                                                             %
 251 %                                                                             %
 252 %   G e t I m a g e C h a n n e l D e p t h                                   %
 253 %                                                                             %
 254 %                                                                             %
 255 %                                                                             %
 256 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 257 %
 258 %  GetImageChannelDepth() returns the depth of a particular image channel.
 259 %
 260 %  The format of the GetImageChannelDepth method is:
 261 %
 262 %      unsigned long GetImageDepth(const Image *image,ExceptionInfo *exception)
 263 %      unsigned long GetImageChannelDepth(const Image *image,
 264 %        const ChannelType channel,ExceptionInfo *exception)
 265 %
 266 %  A description of each parameter follows:
 267 %
 268 %    o image: the image.
 269 %
 270 %    o channel: the channel.
 271 %
 272 %    o exception: return any errors or warnings in this structure.
 273 %
 274 */
 275
 276 MagickExport unsigned long GetImageDepth(const Image *image,
 277   ExceptionInfo *exception)
 278 {
 279   return(GetImageChannelDepth(image,AllChannels,exception));
 280 }
 281
 282 MagickExport unsigned long GetImageChannelDepth(const Image *image,
 283   const ChannelType channel,ExceptionInfo *exception)
 284 {
 285   long
 286     y;
 287
 288   MagickBooleanType
 289     status;
 290
 291   register long
 292     id;
 293
 294   unsigned long
 295     *current_depth,
 296     depth,
 297     number_threads;
 298
 299   CacheView
 300     *image_view;
 301
 302   /*
 303     Compute image depth.
 304   */
 305   assert(image != (Image *) NULL);
 306   assert(image->signature == MagickSignature);
 307   if (image->debug != MagickFalse)
 308     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
 309   number_threads=GetOpenMPMaximumThreads();
 310   current_depth=(unsigned long *) AcquireQuantumMemory(number_threads,
 311     sizeof(*current_depth));
 312   if (current_depth == (unsigned long *) NULL)
 313     ThrowFatalException(ResourceLimitFatalError,"MemoryAllocationFailed");
 314   status=MagickTrue;
 315   for (id=0; id < (long) number_threads; id++)
 316     current_depth[id]=1;
 317   if ((image->storage_class == PseudoClass) && (image->matte == MagickFalse))
 318     {
 319       register const PixelPacket
 320         *__restrict p;
 321
 322       register long
 323         i;
 324
 325       p=image->colormap;
 326 #if defined(MAGICKCORE_OPENMP_SUPPORT)
 327   #pragma omp parallel for schedule(dynamic,4) shared(status)
 328 #endif
 329       for (i=0; i < (long) image->colors; i++)
 330       {
 331         if (status == MagickFalse)
 332           continue;
 333         id=GetOpenMPThreadId();
 334         while (current_depth[id] < MAGICKCORE_QUANTUM_DEPTH)
 335         {
 336           MagickStatusType
 337             status;
 338
 339           QuantumAny
 340             range;
 341
 342           status=0;
 343           range=GetQuantumRange(current_depth[id]);
 344           if ((channel & RedChannel) != 0)
 345             status|=p->red != ScaleAnyToQuantum(ScaleQuantumToAny(p->red,
 346               range),range);
 347           if ((channel & GreenChannel) != 0)
 348             status|=p->green != ScaleAnyToQuantum(ScaleQuantumToAny(p->green,
 349               range),range);
 350           if ((channel & BlueChannel) != 0)
 351             status|=p->blue != ScaleAnyToQuantum(ScaleQuantumToAny(p->blue,
 352               range),range);
 353           if (status == 0)
 354             break;
 355           current_depth[id]++;
 356         }
 357         p++;
 358       }
 359       depth=current_depth[0];
 360       for (id=1; id < (long) number_threads; id++)
 361         if (depth < current_depth[id])
 362           depth=current_depth[id];
 363       current_depth=(unsigned long *) RelinquishMagickMemory(current_depth);
 364       return(depth);
 365     }
 366   image_view=AcquireCacheView(image);
 367 #if defined(MAGICKCORE_OPENMP_SUPPORT)
 368   #pragma omp parallel for schedule(dynamic,4) shared(status)
 369 #endif
 370   for (y=0; y < (long) image->rows; y++)
 371   {
 372     register const IndexPacket
 373       *__restrict indexes;
 374
 375     register const PixelPacket
 376       *__restrict p;
 377
 378     register long
 379       id,
 380       x;
 381
 382     if (status == MagickFalse)
 383       continue;
 384     id=GetOpenMPThreadId();
 385     p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
 386     if (p == (const PixelPacket *) NULL)
 387       continue;
 388     indexes=GetCacheViewVirtualIndexQueue(image_view);
 389     for (x=0; x < (long) image->columns; x++)
 390     {
 391       while (current_depth[id] < MAGICKCORE_QUANTUM_DEPTH)
 392       {
 393         MagickStatusType
 394           status;
 395
 396         QuantumAny
 397           range;
 398
 399         status=0;
 400         range=GetQuantumRange(current_depth[id]);
 401         if ((channel & RedChannel) != 0)
 402           status|=p->red != ScaleAnyToQuantum(ScaleQuantumToAny(p->red,range),
 403             range);
 404         if ((channel & GreenChannel) != 0)
 405           status|=p->green != ScaleAnyToQuantum(ScaleQuantumToAny(p->green,
 406             range),range);
 407         if ((channel & BlueChannel) != 0)
 408           status|=p->blue != ScaleAnyToQuantum(ScaleQuantumToAny(p->blue,range),
 409             range);
 410         if (((channel & OpacityChannel) != 0) && (image->matte != MagickFalse))
 411           status|=p->opacity != ScaleAnyToQuantum(ScaleQuantumToAny(p->opacity,
 412             range),range);
 413         if (((channel & IndexChannel) != 0) &&
 414             (image->colorspace == CMYKColorspace))
 415           status|=indexes[x] != ScaleAnyToQuantum(ScaleQuantumToAny(indexes[x],
 416             range),range);
 417         if (status == 0)
 418           break;
 419         current_depth[id]++;
 420       }
 421       p++;
 422     }
 423     if (current_depth[id] == MAGICKCORE_QUANTUM_DEPTH)
 424       status=MagickFalse;
 425   }
 426   image_view=DestroyCacheView(image_view);
 427   depth=current_depth[0];
 428   for (id=1; id < (long) number_threads; id++)
 429     if (depth < current_depth[id])
 430       depth=current_depth[id];
 431   current_depth=(unsigned long *) RelinquishMagickMemory(current_depth);
 432   return(depth);
 433 }
 434 \f
 435 /*
 436 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 437 %                                                                             %
 438 %                                                                             %
 439 %                                                                             %
 440 +   G e t I m a g e C h a n n e l E x t r e m a                               %
 441 %                                                                             %
 442 %                                                                             %
 443 %                                                                             %
 444 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 445 %
 446 %  GetImageChannelExtrema() returns the extrema of one or more image channels.
 447 %
 448 %  The format of the GetImageChannelExtrema method is:
 449 %
 450 %      MagickBooleanType GetImageChannelExtrema(const Image *image,
 451 %        const ChannelType channel,unsigned long *minima,unsigned long *maxima,
 452 %        ExceptionInfo *exception)
 453 %
 454 %  A description of each parameter follows:
 455 %
 456 %    o image: the image.
 457 %
 458 %    o channel: the channel.
 459 %
 460 %    o minima: the minimum value in the channel.
 461 %
 462 %    o maxima: the maximum value in the channel.
 463 %
 464 %    o exception: return any errors or warnings in this structure.
 465 %
 466 */
 467
 468 MagickExport MagickBooleanType GetImageExtrema(const Image *image,
 469   unsigned long *minima,unsigned long *maxima,ExceptionInfo *exception)
 470 {
 471   return(GetImageChannelExtrema(image,AllChannels,minima,maxima,exception));
 472 }
 473
 474 MagickExport MagickBooleanType GetImageChannelExtrema(const Image *image,
 475   const ChannelType channel,unsigned long *minima,unsigned long *maxima,
 476   ExceptionInfo *exception)
 477 {
 478   double
 479     max,
 480     min;
 481
 482   MagickBooleanType
 483     status;
 484
 485   assert(image != (Image *) NULL);
 486   assert(image->signature == MagickSignature);
 487   if (image->debug != MagickFalse)
 488     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
 489   status=GetImageChannelRange(image,channel,&min,&max,exception);
 490   *minima=(unsigned long) (min+0.5);
 491   *maxima=(unsigned long) (max+0.5);
 492   return(status);
 493 }
 494 \f
 495 /*
 496 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 497 %                                                                             %
 498 %                                                                             %
 499 %                                                                             %
 500 %   G e t I m a g e C h a n n e l M e a n                                     %
 501 %                                                                             %
 502 %                                                                             %
 503 %                                                                             %
 504 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 505 %
 506 %  GetImageChannelMean() returns the mean and standard deviation of one or more
 507 %  image channels.
 508 %
 509 %  The format of the GetImageChannelMean method is:
 510 %
 511 %      MagickBooleanType GetImageChannelMean(const Image *image,
 512 %        const ChannelType channel,double *mean,double *standard_deviation,
 513 %        ExceptionInfo *exception)
 514 %
 515 %  A description of each parameter follows:
 516 %
 517 %    o image: the image.
 518 %
 519 %    o channel: the channel.
 520 %
 521 %    o mean: the average value in the channel.
 522 %
 523 %    o standard_deviation: the standard deviation of the channel.
 524 %
 525 %    o exception: return any errors or warnings in this structure.
 526 %
 527 */
 528
 529 MagickExport MagickBooleanType GetImageMean(const Image *image,double *mean,
 530   double *standard_deviation,ExceptionInfo *exception)
 531 {
 532   MagickBooleanType
 533     status;
 534
 535   status=GetImageChannelMean(image,AllChannels,mean,standard_deviation,
 536     exception);
 537   return(status);
 538 }
 539
 540 MagickExport MagickBooleanType GetImageChannelMean(const Image *image,
 541   const ChannelType channel,double *mean,double *standard_deviation,
 542   ExceptionInfo *exception)
 543 {
 544   double
 545     area;
 546
 547   long
 548     y;
 549
 550   assert(image != (Image *) NULL);
 551   assert(image->signature == MagickSignature);
 552   if (image->debug != MagickFalse)
 553     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
 554   *mean=0.0;
 555   *standard_deviation=0.0;
 556   area=0.0;
 557   for (y=0; y < (long) image->rows; y++)
 558   {
 559     register const IndexPacket
 560       *__restrict indexes;
 561
 562     register const PixelPacket
 563       *__restrict p;
 564
 565     register long
 566       x;
 567
 568     p=GetVirtualPixels(image,0,y,image->columns,1,exception);
 569     if (p == (const PixelPacket *) NULL)
 570       break;
 571     indexes=GetVirtualIndexQueue(image);
 572     for (x=0; x < (long) image->columns; x++)
 573     {
 574       if ((channel & RedChannel) != 0)
 575         {
 576           *mean+=p->red;
 577           *standard_deviation+=(double) p->red*p->red;
 578           area++;
 579         }
 580       if ((channel & GreenChannel) != 0)
 581         {
 582           *mean+=p->green;
 583           *standard_deviation+=(double) p->green*p->green;
 584           area++;
 585         }
 586       if ((channel & BlueChannel) != 0)
 587         {
 588           *mean+=p->blue;
 589           *standard_deviation+=(double) p->blue*p->blue;
 590           area++;
 591         }
 592       if ((channel & OpacityChannel) != 0)
 593         {
 594           *mean+=p->opacity;
 595           *standard_deviation+=(double) p->opacity*p->opacity;
 596           area++;
 597         }
 598       if (((channel & IndexChannel) != 0) &&
 599           (image->colorspace == CMYKColorspace))
 600         {
 601           *mean+=indexes[x];
 602           *standard_deviation+=(double) indexes[x]*indexes[x];
 603           area++;
 604         }
 605       p++;
 606     }
 607   }
 608   if (y < (long) image->rows)
 609     return(MagickFalse);
 610   if (area != 0)
 611     {
 612       *mean/=area;
 613       *standard_deviation/=area;
 614     }
 615   *standard_deviation=sqrt(*standard_deviation-(*mean*(*mean)));
 616   return(y == (long) image->rows ? MagickTrue : MagickFalse);
 617 }
 618 \f
 619 /*
 620 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 621 %                                                                             %
 622 %                                                                             %
 623 %                                                                             %
 624 %   G e t I m a g e C h a n n e l K u r t o s i s                             %
 625 %                                                                             %
 626 %                                                                             %
 627 %                                                                             %
 628 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 629 %
 630 %  GetImageChannelKurtosis() returns the kurtosis and skewness of one or more
 631 %  image channels.
 632 %
 633 %  The format of the GetImageChannelKurtosis method is:
 634 %
 635 %      MagickBooleanType GetImageChannelKurtosis(const Image *image,
 636 %        const ChannelType channel,double *kurtosis,double *skewness,
 637 %        ExceptionInfo *exception)
 638 %
 639 %  A description of each parameter follows:
 640 %
 641 %    o image: the image.
 642 %
 643 %    o channel: the channel.
 644 %
 645 %    o kurtosis: the kurtosis of the channel.
 646 %
 647 %    o skewness: the skewness of the channel.
 648 %
 649 %    o exception: return any errors or warnings in this structure.
 650 %
 651 */
 652
 653 MagickExport MagickBooleanType GetImageKurtosis(const Image *image,
 654   double *kurtosis,double *skewness,ExceptionInfo *exception)
 655 {
 656   MagickBooleanType
 657     status;
 658
 659   status=GetImageChannelKurtosis(image,AllChannels,kurtosis,skewness,
 660     exception);
 661   return(status);
 662 }
 663
 664 MagickExport MagickBooleanType GetImageChannelKurtosis(const Image *image,
 665   const ChannelType channel,double *kurtosis,double *skewness,
 666   ExceptionInfo *exception)
 667 {
 668   double
 669     area,
 670     mean,
 671     standard_deviation,
 672     sum_squares,
 673     sum_cubes,
 674     sum_fourth_power;
 675
 676   long
 677     y;
 678
 679   assert(image != (Image *) NULL);
 680   assert(image->signature == MagickSignature);
 681   if (image->debug != MagickFalse)
 682     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
 683   *kurtosis=0.0;
 684   *skewness=0.0;
 685   area=0.0;
 686   mean=0.0;
 687   standard_deviation=0.0;
 688   sum_squares=0.0;
 689   sum_cubes=0.0;
 690   sum_fourth_power=0.0;
 691   for (y=0; y < (long) image->rows; y++)
 692   {
 693     register const IndexPacket
 694       *__restrict indexes;
 695
 696     register const PixelPacket
 697       *__restrict p;
 698
 699     register long
 700       x;
 701
 702     p=GetVirtualPixels(image,0,y,image->columns,1,exception);
 703     if (p == (const PixelPacket *) NULL)
 704       break;
 705     indexes=GetVirtualIndexQueue(image);
 706     for (x=0; x < (long) image->columns; x++)
 707     {
 708       if ((channel & RedChannel) != 0)
 709         {
 710           mean+=p->red;
 711           sum_squares+=(double) p->red*p->red;
 712           sum_cubes+=(double) p->red*p->red*p->red;
 713           sum_fourth_power+=(double) p->red*p->red*p->red*p->red;
 714           area++;
 715         }
 716       if ((channel & GreenChannel) != 0)
 717         {
 718           mean+=p->green;
 719           sum_squares+=(double) p->green*p->green;
 720           sum_cubes+=(double) p->green*p->green*p->green;
 721           sum_fourth_power+=(double) p->green*p->green*p->green*p->green;
 722           area++;
 723         }
 724       if ((channel & BlueChannel) != 0)
 725         {
 726           mean+=p->blue;
 727           sum_squares+=(double) p->blue*p->blue;
 728           sum_cubes+=(double) p->blue*p->blue*p->blue;
 729           sum_fourth_power+=(double) p->blue*p->blue*p->blue*p->blue;
 730           area++;
 731         }
 732       if ((channel & OpacityChannel) != 0)
 733         {
 734           mean+=p->opacity;
 735           sum_squares+=(double) p->opacity*p->opacity;
 736           sum_cubes+=(double) p->opacity*p->opacity*p->opacity;
 737           sum_fourth_power+=(double) p->opacity*p->opacity*p->opacity*
 738             p->opacity;
 739           area++;
 740         }
 741       if (((channel & IndexChannel) != 0) &&
 742           (image->colorspace == CMYKColorspace))
 743         {
 744           mean+=indexes[x];
 745           sum_squares+=(double) indexes[x]*indexes[x];
 746           sum_cubes+=(double) indexes[x]*indexes[x]*indexes[x];
 747           sum_fourth_power+=(double) indexes[x]*indexes[x]*indexes[x]*
 748             indexes[x];
 749           area++;
 750         }
 751       p++;
 752     }
 753   }
 754   if (y < (long) image->rows)
 755     return(MagickFalse);
 756   if (area != 0.0)
 757     {
 758       mean/=area;
 759       sum_squares/=area;
 760       sum_cubes/=area;
 761       sum_fourth_power/=area;
 762     }
 763   standard_deviation=sqrt(sum_squares-(mean*mean));
 764   if (standard_deviation != 0.0)
 765     {
 766       *kurtosis=sum_fourth_power-4.0*mean*sum_cubes+6.0*mean*mean*sum_squares-
 767         3.0*mean*mean*mean*mean;
 768       *kurtosis/=standard_deviation*standard_deviation*standard_deviation*
 769         standard_deviation;
 770       *kurtosis-=3.0;
 771       *skewness=sum_cubes-3.0*mean*sum_squares+2.0*mean*mean*mean;
 772       *skewness/=standard_deviation*standard_deviation*standard_deviation;
 773     }
 774   return(y == (long) image->rows ? MagickTrue : MagickFalse);
 775 }
 776
 777 /*
 778 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 779 %                                                                             %
 780 %                                                                             %
 781 %                                                                             %
 782 %   G e t I m a g e C h a n n e l R a n g e                                   %
 783 %                                                                             %
 784 %                                                                             %
 785 %                                                                             %
 786 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 787 %
 788 %  GetImageChannelRange() returns the range of one or more image channels.
 789 %
 790 %  The format of the GetImageChannelRange method is:
 791 %
 792 %      MagickBooleanType GetImageChannelRange(const Image *image,
 793 %        const ChannelType channel,double *minima,double *maxima,
 794 %        ExceptionInfo *exception)
 795 %
 796 %  A description of each parameter follows:
 797 %
 798 %    o image: the image.
 799 %
 800 %    o channel: the channel.
 801 %
 802 %    o minima: the minimum value in the channel.
 803 %
 804 %    o maxima: the maximum value in the channel.
 805 %
 806 %    o exception: return any errors or warnings in this structure.
 807 %
 808 */
 809
 810 MagickExport MagickBooleanType GetImageRange(const Image *image,
 811   double *minima,double *maxima,ExceptionInfo *exception)
 812 {
 813   return(GetImageChannelRange(image,AllChannels,minima,maxima,exception));
 814 }
 815
 816 MagickExport MagickBooleanType GetImageChannelRange(const Image *image,
 817   const ChannelType channel,double *minima,double *maxima,
 818   ExceptionInfo *exception)
 819 {
 820   long
 821     y;
 822
 823   MagickPixelPacket
 824     pixel;
 825
 826   assert(image != (Image *) NULL);
 827   assert(image->signature == MagickSignature);
 828   if (image->debug != MagickFalse)
 829     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
 830   *maxima=(-1.0E-37);
 831   *minima=1.0E+37;
 832   GetMagickPixelPacket(image,&pixel);
 833   for (y=0; y < (long) image->rows; y++)
 834   {
 835     register const IndexPacket
 836       *__restrict indexes;
 837
 838     register const PixelPacket
 839       *__restrict p;
 840
 841     register long
 842       x;
 843
 844     p=GetVirtualPixels(image,0,y,image->columns,1,exception);
 845     if (p == (const PixelPacket *) NULL)
 846       break;
 847     indexes=GetVirtualIndexQueue(image);
 848     for (x=0; x < (long) image->columns; x++)
 849     {
 850       SetMagickPixelPacket(image,p,indexes+x,&pixel);
 851       if ((channel & RedChannel) != 0)
 852         {
 853           if (pixel.red < *minima)
 854             *minima=(double) pixel.red;
 855           if (pixel.red > *maxima)
 856             *maxima=(double) pixel.red;
 857         }
 858       if ((channel & GreenChannel) != 0)
 859         {
 860           if (pixel.green < *minima)
 861             *minima=(double) pixel.green;
 862           if (pixel.green > *maxima)
 863             *maxima=(double) pixel.green;
 864         }
 865       if ((channel & BlueChannel) != 0)
 866         {
 867           if (pixel.blue < *minima)
 868             *minima=(double) pixel.blue;
 869           if (pixel.blue > *maxima)
 870             *maxima=(double) pixel.blue;
 871         }
 872       if ((channel & OpacityChannel) != 0)
 873         {
 874           if (pixel.opacity < *minima)
 875             *minima=(double) pixel.opacity;
 876           if (pixel.opacity > *maxima)
 877             *maxima=(double) pixel.opacity;
 878         }
 879       if (((channel & IndexChannel) != 0) &&
 880           (image->colorspace == CMYKColorspace))
 881         {
 882           if ((double) indexes[x] < *minima)
 883             *minima=(double) indexes[x];
 884           if ((double) indexes[x] > *maxima)
 885             *maxima=(double) indexes[x];
 886         }
 887       p++;
 888     }
 889   }
 890   return(y == (long) image->rows ? MagickTrue : MagickFalse);
 891 }
 892 \f
 893 /*
 894 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 895 %                                                                             %
 896 %                                                                             %
 897 %                                                                             %
 898 %   G e t I m a g e C h a n n e l S t a t i s t i c s                         %
 899 %                                                                             %
 900 %                                                                             %
 901 %                                                                             %
 902 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 903 %
 904 %  GetImageChannelStatistics() returns statistics for each channel in the
 905 %  image.  The statistics include the channel depth, its minima, maxima, mean,
 906 %  standard deviation, kurtosis and skewness.  You can access the red channel
 907 %  mean, for example, like this:
 908 %
 909 %      channel_statistics=GetImageChannelStatistics(image,excepton);
 910 %      red_mean=channel_statistics[RedChannel].mean;
 911 %
 912 %  Use MagickRelinquishMemory() to free the statistics buffer.
 913 %
 914 %  The format of the GetImageChannelStatistics method is:
 915 %
 916 %      ChannelStatistics *GetImageChannelStatistics(const Image *image,
 917 %        ExceptionInfo *exception)
 918 %
 919 %  A description of each parameter follows:
 920 %
 921 %    o image: the image.
 922 %
 923 %    o exception: return any errors or warnings in this structure.
 924 %
 925 */
 926
 927 static inline double MagickMax(const double x,const double y)
 928 {
 929   if (x > y)
 930     return(x);
 931   return(y);
 932 }
 933
 934 static inline double MagickMin(const double x,const double y)
 935 {
 936   if (x < y)
 937     return(x);
 938   return(y);
 939 }
 940
 941 MagickExport ChannelStatistics *GetImageChannelStatistics(const Image *image,
 942   ExceptionInfo *exception)
 943 {
 944   ChannelStatistics
 945     *channel_statistics;
 946
 947   double
 948     area,
 949     sum_squares,
 950     sum_cubes;
 951
 952   long
 953     y;
 954
 955   MagickStatusType
 956     status;
 957
 958   QuantumAny
 959     range;
 960
 961   register long
 962     i;
 963
 964   size_t
 965     length;
 966
 967   unsigned long
 968     channels,
 969     depth;
 970
 971   assert(image != (Image *) NULL);
 972   assert(image->signature == MagickSignature);
 973   if (image->debug != MagickFalse)
 974     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
 975   length=AllChannels+1UL;
 976   channel_statistics=(ChannelStatistics *) AcquireQuantumMemory(length,
 977     sizeof(*channel_statistics));
 978   if (channel_statistics == (ChannelStatistics *) NULL)
 979     ThrowFatalException(ResourceLimitFatalError,"MemoryAllocationFailed");
 980   (void) ResetMagickMemory(channel_statistics,0,length*
 981     sizeof(*channel_statistics));
 982   for (i=0; i <= AllChannels; i++)
 983   {
 984     channel_statistics[i].depth=1;
 985     channel_statistics[i].maxima=(-1.0E-37);
 986     channel_statistics[i].minima=1.0E+37;
 987     channel_statistics[i].mean=0.0;
 988     channel_statistics[i].standard_deviation=0.0;
 989     channel_statistics[i].kurtosis=0.0;
 990     channel_statistics[i].skewness=0.0;
 991   }
 992   for (y=0; y < (long) image->rows; y++)
 993   {
 994     register const IndexPacket
 995       *__restrict indexes;
 996
 997     register const PixelPacket
 998       *__restrict p;
 999
1000     register long
1001       x;
1002
1003     p=GetVirtualPixels(image,0,y,image->columns,1,exception);
1004     if (p == (const PixelPacket *) NULL)
1005       break;
1006     indexes=GetVirtualIndexQueue(image);
1007     for (x=0; x < (long) image->columns; )
1008     {
1009       if (channel_statistics[RedChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1010         {
1011           depth=channel_statistics[RedChannel].depth;
1012           range=GetQuantumRange(depth);
1013           status=p->red != ScaleAnyToQuantum(ScaleQuantumToAny(p->red,range),
1014             range) ? MagickTrue : MagickFalse;
1015           if (status != MagickFalse)
1016             {
1017               channel_statistics[RedChannel].depth++;
1018               continue;
1019             }
1020         }
1021       if (channel_statistics[GreenChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1022         {
1023           depth=channel_statistics[GreenChannel].depth;
1024           range=GetQuantumRange(depth);
1025           status=p->green != ScaleAnyToQuantum(ScaleQuantumToAny(p->green,
1026             range),range) ? MagickTrue : MagickFalse;
1027           if (status != MagickFalse)
1028             {
1029               channel_statistics[GreenChannel].depth++;
1030               continue;
1031             }
1032         }
1033       if (channel_statistics[BlueChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1034         {
1035           depth=channel_statistics[BlueChannel].depth;
1036           range=GetQuantumRange(depth);
1037           status=p->blue != ScaleAnyToQuantum(ScaleQuantumToAny(p->blue,
1038             range),range) ? MagickTrue : MagickFalse;
1039           if (status != MagickFalse)
1040             {
1041               channel_statistics[BlueChannel].depth++;
1042               continue;
1043             }
1044         }
1045       if (image->matte != MagickFalse)
1046         {
1047           if (channel_statistics[OpacityChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1048             {
1049               depth=channel_statistics[OpacityChannel].depth;
1050               range=GetQuantumRange(depth);
1051               status=p->opacity != ScaleAnyToQuantum(ScaleQuantumToAny(
1052                 p->opacity,range),range) ? MagickTrue : MagickFalse;
1053               if (status != MagickFalse)
1054                 {
1055                   channel_statistics[OpacityChannel].depth++;
1056                   continue;
1057                 }
1058             }
1059           }
1060       if (image->colorspace == CMYKColorspace)
1061         {
1062           if (channel_statistics[BlackChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1063             {
1064               depth=channel_statistics[BlackChannel].depth;
1065               range=GetQuantumRange(depth);
1066               status=indexes[x] != ScaleAnyToQuantum(ScaleQuantumToAny(
1067                 indexes[x],range),range) ? MagickTrue : MagickFalse;
1068               if (status != MagickFalse)
1069                 {
1070                   channel_statistics[BlackChannel].depth++;
1071                   continue;
1072                 }
1073             }
1074         }
1075       if ((double) p->red < channel_statistics[RedChannel].minima)
1076         channel_statistics[RedChannel].minima=(double) p->red;
1077       if ((double) p->red > channel_statistics[RedChannel].maxima)
1078         channel_statistics[RedChannel].maxima=(double) p->red;
1079       channel_statistics[RedChannel].mean+=p->red;
1080       channel_statistics[RedChannel].standard_deviation+=(double) p->red*p->red;
1081       channel_statistics[RedChannel].kurtosis+=(double) p->red*p->red*
1082         p->red*p->red;
1083       channel_statistics[RedChannel].skewness+=(double) p->red*p->red*p->red;
1084       if ((double) p->green < channel_statistics[GreenChannel].minima)
1085         channel_statistics[GreenChannel].minima=(double) p->green;
1086       if ((double) p->green > channel_statistics[GreenChannel].maxima)
1087         channel_statistics[GreenChannel].maxima=(double) p->green;
1088       channel_statistics[GreenChannel].mean+=p->green;
1089       channel_statistics[GreenChannel].standard_deviation+=(double) p->green*
1090         p->green;
1091       channel_statistics[GreenChannel].kurtosis+=(double) p->green*p->green*
1092         p->green*p->green;
1093       channel_statistics[GreenChannel].skewness+=(double) p->green*p->green*
1094         p->green;
1095       if ((double) p->blue < channel_statistics[BlueChannel].minima)
1096         channel_statistics[BlueChannel].minima=(double) p->blue;
1097       if ((double) p->blue > channel_statistics[BlueChannel].maxima)
1098         channel_statistics[BlueChannel].maxima=(double) p->blue;
1099       channel_statistics[BlueChannel].mean+=p->blue;
1100       channel_statistics[BlueChannel].standard_deviation+=(double) p->blue*
1101         p->blue;
1102       channel_statistics[BlueChannel].kurtosis+=(double) p->blue*p->blue*
1103         p->blue*p->blue;
1104       channel_statistics[BlueChannel].skewness+=(double) p->blue*p->blue*
1105         p->blue;
1106       if (image->matte != MagickFalse)
1107         {
1108           if ((double) p->opacity < channel_statistics[OpacityChannel].minima)
1109             channel_statistics[OpacityChannel].minima=(double) p->opacity;
1110           if ((double) p->opacity > channel_statistics[OpacityChannel].maxima)
1111             channel_statistics[OpacityChannel].maxima=(double) p->opacity;
1112           channel_statistics[OpacityChannel].mean+=p->opacity;
1113           channel_statistics[OpacityChannel].standard_deviation+=(double)
1114             p->opacity*p->opacity;
1115           channel_statistics[OpacityChannel].kurtosis+=(double) p->opacity*
1116             p->opacity*p->opacity*p->opacity;
1117           channel_statistics[OpacityChannel].skewness+=(double) p->opacity*
1118             p->opacity*p->opacity;
1119         }
1120       if (image->colorspace == CMYKColorspace)
1121         {
1122           if ((double) indexes[x] < channel_statistics[BlackChannel].minima)
1123             channel_statistics[BlackChannel].minima=(double) indexes[x];
1124           if ((double) indexes[x] > channel_statistics[BlackChannel].maxima)
1125             channel_statistics[BlackChannel].maxima=(double) indexes[x];
1126           channel_statistics[BlackChannel].mean+=indexes[x];
1127           channel_statistics[BlackChannel].standard_deviation+=(double)
1128             indexes[x]*indexes[x];
1129           channel_statistics[BlackChannel].kurtosis+=(double) indexes[x]*
1130             indexes[x]*indexes[x]*indexes[x];
1131           channel_statistics[BlackChannel].skewness+=(double) indexes[x]*
1132             indexes[x]*indexes[x];
1133         }
1134       x++;
1135       p++;
1136     }
1137   }
1138   area=(double) image->columns*image->rows;
1139   for (i=0; i < AllChannels; i++)
1140   {
1141     channel_statistics[i].mean/=area;
1142     channel_statistics[i].standard_deviation/=area;
1143     channel_statistics[i].kurtosis/=area;
1144     channel_statistics[i].skewness/=area;
1145   }
1146   for (i=0; i < AllChannels; i++)
1147   {
1148     channel_statistics[AllChannels].depth=(unsigned long) MagickMax((double)
1149       channel_statistics[AllChannels].depth,(double)
1150       channel_statistics[i].depth);
1151     channel_statistics[AllChannels].minima=MagickMin(
1152       channel_statistics[AllChannels].minima,channel_statistics[i].minima);
1153     channel_statistics[AllChannels].maxima=MagickMax(
1154       channel_statistics[AllChannels].maxima,channel_statistics[i].maxima);
1155     channel_statistics[AllChannels].mean+=channel_statistics[i].mean;
1156     channel_statistics[AllChannels].standard_deviation+=
1157       channel_statistics[i].standard_deviation;
1158     channel_statistics[AllChannels].kurtosis+=channel_statistics[i].kurtosis;
1159     channel_statistics[AllChannels].skewness+=channel_statistics[i].skewness;
1160   }
1161   channels=4;
1162   if (image->colorspace == CMYKColorspace)
1163     channels++;
1164   channel_statistics[AllChannels].mean/=channels;
1165   channel_statistics[AllChannels].standard_deviation/=channels;
1166   channel_statistics[AllChannels].kurtosis/=channels;
1167   channel_statistics[AllChannels].skewness/=channels;
1168   for (i=0; i <= AllChannels; i++)
1169   {
1170     sum_squares=0.0;
1171     sum_squares=channel_statistics[i].standard_deviation;
1172     sum_cubes=0.0;
1173     sum_cubes=channel_statistics[i].skewness;
1174     channel_statistics[i].standard_deviation=sqrt(
1175       channel_statistics[i].standard_deviation-
1176        (channel_statistics[i].mean*channel_statistics[i].mean));
1177     if (channel_statistics[i].standard_deviation == 0.0)
1178       {
1179         channel_statistics[i].kurtosis=0.0;
1180         channel_statistics[i].skewness=0.0;
1181       }
1182     else
1183       {
1184         channel_statistics[i].skewness=(channel_statistics[i].skewness-
1185           3.0*channel_statistics[i].mean*sum_squares+
1186           2.0*channel_statistics[i].mean*channel_statistics[i].mean*
1187           channel_statistics[i].mean)/
1188           (channel_statistics[i].standard_deviation*
1189            channel_statistics[i].standard_deviation*
1190            channel_statistics[i].standard_deviation);
1191         channel_statistics[i].kurtosis=(channel_statistics[i].kurtosis-
1192           4.0*channel_statistics[i].mean*sum_cubes+
1193           6.0*channel_statistics[i].mean*channel_statistics[i].mean*sum_squares-
1194           3.0*channel_statistics[i].mean*channel_statistics[i].mean*
1195           1.0*channel_statistics[i].mean*channel_statistics[i].mean)/
1196           (channel_statistics[i].standard_deviation*
1197            channel_statistics[i].standard_deviation*
1198            channel_statistics[i].standard_deviation*
1199            channel_statistics[i].standard_deviation)-3.0;
1200       }
1201   }
1202   return(channel_statistics);
1203 }
1204 \f
1205 /*
1206 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1207 %                                                                             %
1208 %                                                                             %
1209 %                                                                             %
1210 %   G e t I m a g e Q u a n t u m D e p t h                                   %
1211 %                                                                             %
1212 %                                                                             %
1213 %                                                                             %
1214 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1215 %
1216 %  GetImageQuantumDepth() returns the depth of the image rounded to a legal
1217 %  quantum depth: 8, 16, or 32.
1218 %
1219 %  The format of the GetImageQuantumDepth method is:
1220 %
1221 %      unsigned long GetImageQuantumDepth(const Image *image,
1222 %        const MagickBooleanType constrain)
1223 %
1224 %  A description of each parameter follows:
1225 %
1226 %    o image: the image.
1227 %
1228 %    o constrain: A value other than MagickFalse, constrains the depth to
1229 %      a maximum of MAGICKCORE_QUANTUM_DEPTH.
1230 %
1231 */
1232 MagickExport unsigned long GetImageQuantumDepth(const Image *image,
1233   const MagickBooleanType constrain)
1234 {
1235   unsigned long
1236     depth;
1237
1238   depth=image->depth;
1239   if (depth <= 8)
1240     depth=8;
1241   else
1242     if (depth <= 16)
1243       depth=16;
1244     else
1245       if (depth <= 32)
1246         depth=32;
1247       else
1248         if (depth <= 64)
1249           depth=64;
1250   if (constrain != MagickFalse)
1251     depth=(unsigned long) MagickMin((double) depth,(double)
1252       MAGICKCORE_QUANTUM_DEPTH);
1253   return(depth);
1254 }
1255 \f
1256 /*
1257 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1258 %                                                                             %
1259 %                                                                             %
1260 %                                                                             %
1261 %   S e t I m a g e C h a n n e l D e p t h                                   %
1262 %                                                                             %
1263 %                                                                             %
1264 %                                                                             %
1265 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1266 %
1267 %  SetImageChannelDepth() sets the depth of the image.
1268 %
1269 %  The format of the SetImageChannelDepth method is:
1270 %
1271 %      MagickBooleanType SetImageDepth(Image *image,const unsigned long depth)
1272 %      MagickBooleanType SetImageChannelDepth(Image *image,
1273 %        const ChannelType channel,const unsigned long depth)
1274 %
1275 %  A description of each parameter follows:
1276 %
1277 %    o image: the image.
1278 %
1279 %    o channel: the channel.
1280 %
1281 %    o depth: the image depth.
1282 %
1283 */
1284
1285 MagickExport MagickBooleanType SetImageDepth(Image *image,
1286   const unsigned long depth)
1287 {
1288   return(SetImageChannelDepth(image,AllChannels,depth));
1289 }
1290
1291 MagickExport MagickBooleanType SetImageChannelDepth(Image *image,
1292   const ChannelType channel,const unsigned long depth)
1293 {
1294   ExceptionInfo
1295     *exception;
1296
1297   long
1298     y;
1299
1300   MagickBooleanType
1301     status;
1302
1303   QuantumAny
1304     range;
1305
1306   CacheView
1307     *image_view;
1308
1309   assert(image != (Image *) NULL);
1310   if (image->debug != MagickFalse)
1311     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"...");
1312   assert(image->signature == MagickSignature);
1313   if (GetImageDepth(image,&image->exception) <= (unsigned long)
1314       MagickMin((double) depth,(double) MAGICKCORE_QUANTUM_DEPTH))
1315     {
1316       image->depth=depth;
1317       return(MagickTrue);
1318     }
1319   /*
1320     Scale pixels to desired depth.
1321   */
1322   status=MagickTrue;
1323   range=GetQuantumRange(depth);
1324   exception=(&image->exception);
1325   image_view=AcquireCacheView(image);
1326 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1327   #pragma omp parallel for schedule(dynamic,4) shared(status)
1328 #endif
1329   for (y=0; y < (long) image->rows; y++)
1330   {
1331     register IndexPacket
1332       *__restrict indexes;
1333
1334     register long
1335       x;
1336
1337     register PixelPacket
1338       *__restrict q;
1339
1340     if (status == MagickFalse)
1341       continue;
1342     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,
1343       exception);
1344     if (q == (PixelPacket *) NULL)
1345       {
1346         status=MagickFalse;
1347         continue;
1348       }
1349     indexes=GetCacheViewAuthenticIndexQueue(image_view);
1350     for (x=0; x < (long) image->columns; x++)
1351     {
1352       if ((channel & RedChannel) != 0)
1353         q->red=ScaleAnyToQuantum(ScaleQuantumToAny(q->red,range),range);
1354       if ((channel & GreenChannel) != 0)
1355         q->green=ScaleAnyToQuantum(ScaleQuantumToAny(q->green,range),range);
1356       if ((channel & BlueChannel) != 0)
1357         q->blue=ScaleAnyToQuantum(ScaleQuantumToAny(q->blue,range),range);
1358       if (((channel & OpacityChannel) != 0) && (image->matte != MagickFalse))
1359         q->opacity=ScaleAnyToQuantum(ScaleQuantumToAny(q->opacity,range),range);
1360       if (((channel & IndexChannel) != 0) &&
1361           (image->colorspace == CMYKColorspace))
1362         indexes[x]=ScaleAnyToQuantum(ScaleQuantumToAny(indexes[x],range),range);
1363       q++;
1364     }
1365     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
1366       {
1367         status=MagickFalse;
1368         continue;
1369       }
1370   }
1371   image_view=DestroyCacheView(image_view);
1372   if (image->storage_class == PseudoClass)
1373     {
1374       QuantumAny
1375         range;
1376
1377       register long
1378         i;
1379
1380       register PixelPacket
1381         *__restrict p;
1382
1383       p=image->colormap;
1384       range=GetQuantumRange(depth);
1385 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1386   #pragma omp parallel for schedule(dynamic,4) shared(status)
1387 #endif
1388       for (i=0; i < (long) image->colors; i++)
1389       {
1390         if ((channel & RedChannel) != 0)
1391           p->red=ScaleAnyToQuantum(ScaleQuantumToAny(p->red,range),range);
1392         if ((channel & GreenChannel) != 0)
1393           p->green=ScaleAnyToQuantum(ScaleQuantumToAny(p->green,range),range);
1394         if ((channel & BlueChannel) != 0)
1395           p->blue=ScaleAnyToQuantum(ScaleQuantumToAny(p->blue,range),range);
1396         if ((channel & OpacityChannel) != 0)
1397           p->opacity=ScaleAnyToQuantum(ScaleQuantumToAny(p->opacity,range),
1398             range);
1399         p++;
1400       }
1401     }
1402   image->depth=depth;
1403   return(status);
1404 }