2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
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 %
13 % MagickCore Image Statistical Methods %
20 % Copyright 1999-2010 ImageMagick Studio LLC, a non-profit organization %
21 % dedicated to making software imaging solutions freely available. %
23 % You may not use this file except in compliance with the License. You may %
24 % obtain a copy of the License at %
26 % http://www.imagemagick.org/script/license.php %
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. %
34 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
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/monitor-private.h"
76 #include "magick/option.h"
77 #include "magick/paint.h"
78 #include "magick/pixel-private.h"
79 #include "magick/profile.h"
80 #include "magick/quantize.h"
81 #include "magick/random_.h"
82 #include "magick/random-private.h"
83 #include "magick/segment.h"
84 #include "magick/semaphore.h"
85 #include "magick/signature-private.h"
86 #include "magick/statistic.h"
87 #include "magick/string_.h"
88 #include "magick/thread-private.h"
89 #include "magick/timer.h"
90 #include "magick/utility.h"
91 #include "magick/version.h"
94 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
98 % E v a l u a t e I m a g e %
102 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
104 % EvaluateImage() applies a value to the image with an arithmetic, relational,
105 % or logical operator to an image. Use these operations to lighten or darken
106 % an image, to increase or decrease contrast in an image, or to produce the
107 % "negative" of an image.
109 % The format of the EvaluateImageChannel method is:
111 % MagickBooleanType EvaluateImage(Image *image,
112 % const MagickEvaluateOperator op,const double value,
113 % ExceptionInfo *exception)
114 % MagickBooleanType EvaluateImages(Image *images,
115 % const MagickEvaluateOperator op,const double value,
116 % ExceptionInfo *exception)
117 % MagickBooleanType EvaluateImageChannel(Image *image,
118 % const ChannelType channel,const MagickEvaluateOperator op,
119 % const double value,ExceptionInfo *exception)
121 % A description of each parameter follows:
123 % o image: the image.
125 % o channel: the channel.
127 % o op: A channel op.
129 % o value: A value value.
131 % o exception: return any errors or warnings in this structure.
135 static MagickPixelPacket **DestroyPixelThreadSet(MagickPixelPacket **pixels)
140 assert(pixels != (MagickPixelPacket **) NULL);
141 for (i=0; i < (ssize_t) GetOpenMPMaximumThreads(); i++)
142 if (pixels[i] != (MagickPixelPacket *) NULL)
143 pixels[i]=(MagickPixelPacket *) RelinquishMagickMemory(pixels[i]);
144 pixels=(MagickPixelPacket **) RelinquishMagickMemory(pixels);
148 static MagickPixelPacket **AcquirePixelThreadSet(const Image *image)
160 number_threads=GetOpenMPMaximumThreads();
161 pixels=(MagickPixelPacket **) AcquireQuantumMemory(number_threads,
163 if (pixels == (MagickPixelPacket **) NULL)
164 return((MagickPixelPacket **) NULL);
165 (void) ResetMagickMemory(pixels,0,number_threads*sizeof(*pixels));
166 for (i=0; i < (ssize_t) number_threads; i++)
168 pixels[i]=(MagickPixelPacket *) AcquireQuantumMemory(image->columns,
170 if (pixels[i] == (MagickPixelPacket *) NULL)
171 return(DestroyPixelThreadSet(pixels));
172 for (j=0; j < (ssize_t) image->columns; j++)
173 GetMagickPixelPacket(image,&pixels[i][j]);
178 static inline double MagickMax(const double x,const double y)
185 static inline double MagickMin(const double x,const double y)
192 static MagickRealType ApplyEvaluateOperator(RandomInfo *random_info,
193 Quantum pixel,const MagickEvaluateOperator op,const MagickRealType value)
201 case UndefinedEvaluateOperator:
203 case AbsEvaluateOperator:
205 result=(MagickRealType) fabs((double) (pixel+value));
208 case AddEvaluateOperator:
210 result=(MagickRealType) (pixel+value);
213 case AddModulusEvaluateOperator:
216 This returns a 'floored modulus' of the addition which is a
217 positive result. It differs from % or fmod() which returns a
218 'truncated modulus' result, where floor() is replaced by trunc()
219 and could return a negative result (which is clipped).
222 result-=(QuantumRange+1.0)*floor((double) result/(QuantumRange+1.0));
225 case AndEvaluateOperator:
227 result=(MagickRealType) ((size_t) pixel & (size_t) (value+0.5));
230 case CosineEvaluateOperator:
232 result=(MagickRealType) (QuantumRange*(0.5*cos((double) (2.0*MagickPI*
233 QuantumScale*pixel*value))+0.5));
236 case DivideEvaluateOperator:
238 result=pixel/(value == 0.0 ? 1.0 : value);
241 case ExponentialEvaluateOperator:
243 result=(MagickRealType) (QuantumRange*exp((double) (value*QuantumScale*
247 case GaussianNoiseEvaluateOperator:
249 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
250 GaussianNoise,value);
253 case ImpulseNoiseEvaluateOperator:
255 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
259 case LaplacianNoiseEvaluateOperator:
261 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
262 LaplacianNoise,value);
265 case LeftShiftEvaluateOperator:
267 result=(MagickRealType) ((size_t) pixel << (size_t) (value+0.5));
270 case LogEvaluateOperator:
272 result=(MagickRealType) (QuantumRange*log((double) (QuantumScale*value*
273 pixel+1.0))/log((double) (value+1.0)));
276 case MaxEvaluateOperator:
278 result=(MagickRealType) MagickMax((double) pixel,value);
281 case MeanEvaluateOperator:
283 result=(MagickRealType) (pixel+value);
286 case MinEvaluateOperator:
288 result=(MagickRealType) MagickMin((double) pixel,value);
291 case MultiplicativeNoiseEvaluateOperator:
293 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
294 MultiplicativeGaussianNoise,value);
297 case MultiplyEvaluateOperator:
299 result=(MagickRealType) (value*pixel);
302 case OrEvaluateOperator:
304 result=(MagickRealType) ((size_t) pixel | (size_t) (value+0.5));
307 case PoissonNoiseEvaluateOperator:
309 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
313 case PowEvaluateOperator:
315 result=(MagickRealType) (QuantumRange*pow((double) (QuantumScale*pixel),
319 case RightShiftEvaluateOperator:
321 result=(MagickRealType) ((size_t) pixel >> (size_t) (value+0.5));
324 case SetEvaluateOperator:
329 case SineEvaluateOperator:
331 result=(MagickRealType) (QuantumRange*(0.5*sin((double) (2.0*MagickPI*
332 QuantumScale*pixel*value))+0.5));
335 case SubtractEvaluateOperator:
337 result=(MagickRealType) (pixel-value);
340 case ThresholdEvaluateOperator:
342 result=(MagickRealType) (((MagickRealType) pixel <= value) ? 0 :
346 case ThresholdBlackEvaluateOperator:
348 result=(MagickRealType) (((MagickRealType) pixel <= value) ? 0 : pixel);
351 case ThresholdWhiteEvaluateOperator:
353 result=(MagickRealType) (((MagickRealType) pixel > value) ? QuantumRange :
357 case UniformNoiseEvaluateOperator:
359 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
363 case XorEvaluateOperator:
365 result=(MagickRealType) ((size_t) pixel ^ (size_t) (value+0.5));
372 MagickExport MagickBooleanType EvaluateImage(Image *image,
373 const MagickEvaluateOperator op,const double value,ExceptionInfo *exception)
378 status=EvaluateImageChannel(image,AllChannels,op,value,exception);
382 MagickExport Image *EvaluateImages(const Image *images,
383 const MagickEvaluateOperator op,ExceptionInfo *exception)
385 #define EvaluateImageTag "Evaluate/Image"
403 **restrict evaluate_pixels,
407 **restrict random_info;
416 Ensure the image are the same size.
418 assert(images != (Image *) NULL);
419 assert(images->signature == MagickSignature);
420 if (images->debug != MagickFalse)
421 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",images->filename);
422 assert(exception != (ExceptionInfo *) NULL);
423 assert(exception->signature == MagickSignature);
424 for (next=images; next != (Image *) NULL; next=GetNextImageInList(next))
425 if ((next->columns != images->columns) || (next->rows != images->rows))
427 (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
428 "ImageWidthsOrHeightsDiffer","`%s'",images->filename);
429 return((Image *) NULL);
432 Initialize evaluate next attributes.
434 evaluate_image=CloneImage(images,images->columns,images->rows,MagickTrue,
436 if (evaluate_image == (Image *) NULL)
437 return((Image *) NULL);
438 if (SetImageStorageClass(evaluate_image,DirectClass) == MagickFalse)
440 InheritException(exception,&evaluate_image->exception);
441 evaluate_image=DestroyImage(evaluate_image);
442 return((Image *) NULL);
444 evaluate_pixels=AcquirePixelThreadSet(images);
445 if (evaluate_pixels == (MagickPixelPacket **) NULL)
447 evaluate_image=DestroyImage(evaluate_image);
448 (void) ThrowMagickException(exception,GetMagickModule(),
449 ResourceLimitError,"MemoryAllocationFailed","`%s'",images->filename);
450 return((Image *) NULL);
453 Evaluate image pixels.
457 GetMagickPixelPacket(images,&zero);
458 random_info=AcquireRandomInfoThreadSet();
459 number_images=GetImageListLength(images);
460 evaluate_view=AcquireCacheView(evaluate_image);
461 #if defined(MAGICKCORE_OPENMP_SUPPORT)
462 #pragma omp parallel for schedule(dynamic) shared(progress,status)
464 for (y=0; y < (ssize_t) evaluate_image->rows; y++)
473 id = GetOpenMPThreadId();
479 *restrict evaluate_indexes;
485 register MagickPixelPacket
491 if (status == MagickFalse)
493 q=QueueCacheViewAuthenticPixels(evaluate_view,0,y,evaluate_image->columns,1,
495 if (q == (PixelPacket *) NULL)
500 evaluate_indexes=GetCacheViewAuthenticIndexQueue(evaluate_view);
502 evaluate_pixel=evaluate_pixels[id];
503 for (x=0; x < (ssize_t) evaluate_image->columns; x++)
504 evaluate_pixel[x]=zero;
506 for (i=0; i < (ssize_t) number_images; i++)
508 register const IndexPacket
511 register const PixelPacket
514 image_view=AcquireCacheView(next);
515 p=GetCacheViewVirtualPixels(image_view,0,y,next->columns,1,exception);
516 if (p == (const PixelPacket *) NULL)
518 image_view=DestroyCacheView(image_view);
521 indexes=GetCacheViewVirtualIndexQueue(image_view);
522 for (x=0; x < (ssize_t) next->columns; x++)
524 evaluate_pixel[x].red=ApplyEvaluateOperator(random_info[id],p->red,
525 i == 0 ? AddEvaluateOperator : op,evaluate_pixel[x].red);
526 evaluate_pixel[x].green=ApplyEvaluateOperator(random_info[id],p->green,
527 i == 0 ? AddEvaluateOperator : op,evaluate_pixel[x].green);
528 evaluate_pixel[x].blue=ApplyEvaluateOperator(random_info[id],p->blue,
529 i == 0 ? AddEvaluateOperator : op,evaluate_pixel[x].blue);
530 evaluate_pixel[x].opacity=ApplyEvaluateOperator(random_info[id],
531 p->opacity,i == 0 ? AddEvaluateOperator : op,
532 evaluate_pixel[x].opacity);
533 if (evaluate_image->colorspace == CMYKColorspace)
534 evaluate_pixel[x].index=ApplyEvaluateOperator(random_info[id],
535 indexes[x],i == 0 ? AddEvaluateOperator : op,
536 evaluate_pixel[x].index);
539 image_view=DestroyCacheView(image_view);
540 next=GetNextImageInList(next);
542 if (op == MeanEvaluateOperator)
543 for (x=0; x < (ssize_t) evaluate_image->columns; x++)
545 evaluate_pixel[x].red/=number_images;
546 evaluate_pixel[x].green/=number_images;
547 evaluate_pixel[x].blue/=number_images;
548 evaluate_pixel[x].opacity/=number_images;
549 evaluate_pixel[x].index/=number_images;
551 for (x=0; x < (ssize_t) evaluate_image->columns; x++)
553 q->red=ClampToQuantum(evaluate_pixel[x].red);
554 q->green=ClampToQuantum(evaluate_pixel[x].green);
555 q->blue=ClampToQuantum(evaluate_pixel[x].blue);
556 if (evaluate_image->matte == MagickFalse)
557 q->opacity=ClampToQuantum(evaluate_pixel[x].opacity);
559 q->opacity=ClampToQuantum(QuantumRange-evaluate_pixel[x].opacity);
560 if (evaluate_image->colorspace == CMYKColorspace)
561 evaluate_indexes[x]=ClampToQuantum(evaluate_pixel[x].index);
564 if (SyncCacheViewAuthenticPixels(evaluate_view,exception) == MagickFalse)
566 if (images->progress_monitor != (MagickProgressMonitor) NULL)
571 #if defined(MAGICKCORE_OPENMP_SUPPORT)
572 #pragma omp critical (MagickCore_EvaluateImages)
574 proceed=SetImageProgress(images,EvaluateImageTag,progress++,
575 evaluate_image->rows);
576 if (proceed == MagickFalse)
580 evaluate_view=DestroyCacheView(evaluate_view);
581 evaluate_pixels=DestroyPixelThreadSet(evaluate_pixels);
582 random_info=DestroyRandomInfoThreadSet(random_info);
583 if (status == MagickFalse)
584 evaluate_image=DestroyImage(evaluate_image);
585 return(evaluate_image);
588 MagickExport MagickBooleanType EvaluateImageChannel(Image *image,
589 const ChannelType channel,const MagickEvaluateOperator op,const double value,
590 ExceptionInfo *exception)
602 **restrict random_info;
607 assert(image != (Image *) NULL);
608 assert(image->signature == MagickSignature);
609 if (image->debug != MagickFalse)
610 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
611 assert(exception != (ExceptionInfo *) NULL);
612 assert(exception->signature == MagickSignature);
613 if (SetImageStorageClass(image,DirectClass) == MagickFalse)
615 InheritException(exception,&image->exception);
620 random_info=AcquireRandomInfoThreadSet();
621 image_view=AcquireCacheView(image);
622 #if defined(MAGICKCORE_OPENMP_SUPPORT)
623 #pragma omp parallel for schedule(dynamic,4) shared(progress,status)
625 for (y=0; y < (ssize_t) image->rows; y++)
628 id = GetOpenMPThreadId();
639 if (status == MagickFalse)
641 q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
642 if (q == (PixelPacket *) NULL)
647 indexes=GetCacheViewAuthenticIndexQueue(image_view);
648 for (x=0; x < (ssize_t) image->columns; x++)
650 if ((channel & RedChannel) != 0)
651 q->red=ClampToQuantum(ApplyEvaluateOperator(random_info[id],q->red,op,
653 if ((channel & GreenChannel) != 0)
654 q->green=ClampToQuantum(ApplyEvaluateOperator(random_info[id],q->green,
656 if ((channel & BlueChannel) != 0)
657 q->blue=ClampToQuantum(ApplyEvaluateOperator(random_info[id],q->blue,op,
659 if ((channel & OpacityChannel) != 0)
661 if (image->matte == MagickFalse)
662 q->opacity=ClampToQuantum(ApplyEvaluateOperator(random_info[id],
663 q->opacity,op,value));
665 q->opacity=ClampToQuantum(QuantumRange-ApplyEvaluateOperator(
666 random_info[id],(Quantum) GetAlphaPixelComponent(q),op,value));
668 if (((channel & IndexChannel) != 0) && (indexes != (IndexPacket *) NULL))
669 indexes[x]=(IndexPacket) ClampToQuantum(ApplyEvaluateOperator(
670 random_info[id],indexes[x],op,value));
673 if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
675 if (image->progress_monitor != (MagickProgressMonitor) NULL)
680 #if defined(MAGICKCORE_OPENMP_SUPPORT)
681 #pragma omp critical (MagickCore_EvaluateImageChannel)
683 proceed=SetImageProgress(image,EvaluateImageTag,progress++,image->rows);
684 if (proceed == MagickFalse)
688 image_view=DestroyCacheView(image_view);
689 random_info=DestroyRandomInfoThreadSet(random_info);
694 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
698 % F u n c t i o n I m a g e %
702 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
704 % FunctionImage() applies a value to the image with an arithmetic, relational,
705 % or logical operator to an image. Use these operations to lighten or darken
706 % an image, to increase or decrease contrast in an image, or to produce the
707 % "negative" of an image.
709 % The format of the FunctionImageChannel method is:
711 % MagickBooleanType FunctionImage(Image *image,
712 % const MagickFunction function,const ssize_t number_parameters,
713 % const double *parameters,ExceptionInfo *exception)
714 % MagickBooleanType FunctionImageChannel(Image *image,
715 % const ChannelType channel,const MagickFunction function,
716 % const ssize_t number_parameters,const double *argument,
717 % ExceptionInfo *exception)
719 % A description of each parameter follows:
721 % o image: the image.
723 % o channel: the channel.
725 % o function: A channel function.
727 % o parameters: one or more parameters.
729 % o exception: return any errors or warnings in this structure.
733 static Quantum ApplyFunction(Quantum pixel,const MagickFunction function,
734 const size_t number_parameters,const double *parameters,
735 ExceptionInfo *exception)
747 case PolynomialFunction:
751 * Parameters: polynomial constants, highest to lowest order
752 * For example: c0*x^3 + c1*x^2 + c2*x + c3
755 for (i=0; i < (ssize_t) number_parameters; i++)
756 result = result*QuantumScale*pixel + parameters[i];
757 result *= QuantumRange;
760 case SinusoidFunction:
763 * Parameters: Freq, Phase, Ampl, bias
765 double freq,phase,ampl,bias;
766 freq = ( number_parameters >= 1 ) ? parameters[0] : 1.0;
767 phase = ( number_parameters >= 2 ) ? parameters[1] : 0.0;
768 ampl = ( number_parameters >= 3 ) ? parameters[2] : 0.5;
769 bias = ( number_parameters >= 4 ) ? parameters[3] : 0.5;
770 result=(MagickRealType) (QuantumRange*(ampl*sin((double) (2.0*MagickPI*
771 (freq*QuantumScale*pixel + phase/360.0) )) + bias ) );
776 /* Arcsin Function (peged at range limits for invalid results)
777 * Parameters: Width, Center, Range, Bias
779 double width,range,center,bias;
780 width = ( number_parameters >= 1 ) ? parameters[0] : 1.0;
781 center = ( number_parameters >= 2 ) ? parameters[1] : 0.5;
782 range = ( number_parameters >= 3 ) ? parameters[2] : 1.0;
783 bias = ( number_parameters >= 4 ) ? parameters[3] : 0.5;
784 result = 2.0/width*(QuantumScale*pixel - center);
785 if ( result <= -1.0 )
786 result = bias - range/2.0;
787 else if ( result >= 1.0 )
788 result = bias + range/2.0;
790 result=range/MagickPI*asin((double)result) + bias;
791 result *= QuantumRange;
797 * Parameters: Slope, Center, Range, Bias
799 double slope,range,center,bias;
800 slope = ( number_parameters >= 1 ) ? parameters[0] : 1.0;
801 center = ( number_parameters >= 2 ) ? parameters[1] : 0.5;
802 range = ( number_parameters >= 3 ) ? parameters[2] : 1.0;
803 bias = ( number_parameters >= 4 ) ? parameters[3] : 0.5;
804 result = MagickPI*slope*(QuantumScale*pixel - center);
805 result=(MagickRealType) (QuantumRange*(range/MagickPI*atan((double)
809 case UndefinedFunction:
812 return(ClampToQuantum(result));
815 MagickExport MagickBooleanType FunctionImage(Image *image,
816 const MagickFunction function,const size_t number_parameters,
817 const double *parameters,ExceptionInfo *exception)
822 status=FunctionImageChannel(image,AllChannels,function,number_parameters,
823 parameters,exception);
827 MagickExport MagickBooleanType FunctionImageChannel(Image *image,
828 const ChannelType channel,const MagickFunction function,
829 const size_t number_parameters,const double *parameters,
830 ExceptionInfo *exception)
832 #define FunctionImageTag "Function/Image "
846 assert(image != (Image *) NULL);
847 assert(image->signature == MagickSignature);
848 if (image->debug != MagickFalse)
849 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
850 assert(exception != (ExceptionInfo *) NULL);
851 assert(exception->signature == MagickSignature);
852 if (SetImageStorageClass(image,DirectClass) == MagickFalse)
854 InheritException(exception,&image->exception);
859 image_view=AcquireCacheView(image);
860 #if defined(MAGICKCORE_OPENMP_SUPPORT)
861 #pragma omp parallel for schedule(dynamic,4) shared(progress,status)
863 for (y=0; y < (ssize_t) image->rows; y++)
874 if (status == MagickFalse)
876 q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
877 if (q == (PixelPacket *) NULL)
882 indexes=GetCacheViewAuthenticIndexQueue(image_view);
883 for (x=0; x < (ssize_t) image->columns; x++)
885 if ((channel & RedChannel) != 0)
886 q->red=ApplyFunction(q->red,function,number_parameters,parameters,
888 if ((channel & GreenChannel) != 0)
889 q->green=ApplyFunction(q->green,function,number_parameters,parameters,
891 if ((channel & BlueChannel) != 0)
892 q->blue=ApplyFunction(q->blue,function,number_parameters,parameters,
894 if ((channel & OpacityChannel) != 0)
896 if (image->matte == MagickFalse)
897 q->opacity=ApplyFunction(q->opacity,function,number_parameters,
898 parameters,exception);
900 q->opacity=(Quantum) QuantumRange-ApplyFunction((Quantum)
901 GetAlphaPixelComponent(q),function,number_parameters,parameters,
904 if (((channel & IndexChannel) != 0) && (indexes != (IndexPacket *) NULL))
905 indexes[x]=(IndexPacket) ApplyFunction(indexes[x],function,
906 number_parameters,parameters,exception);
909 if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
911 if (image->progress_monitor != (MagickProgressMonitor) NULL)
916 #if defined(MAGICKCORE_OPENMP_SUPPORT)
917 #pragma omp critical (MagickCore_FunctionImageChannel)
919 proceed=SetImageProgress(image,FunctionImageTag,progress++,image->rows);
920 if (proceed == MagickFalse)
924 image_view=DestroyCacheView(image_view);
929 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
933 + G e t I m a g e C h a n n e l E x t r e m a %
937 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
939 % GetImageChannelExtrema() returns the extrema of one or more image channels.
941 % The format of the GetImageChannelExtrema method is:
943 % MagickBooleanType GetImageChannelExtrema(const Image *image,
944 % const ChannelType channel,size_t *minima,size_t *maxima,
945 % ExceptionInfo *exception)
947 % A description of each parameter follows:
949 % o image: the image.
951 % o channel: the channel.
953 % o minima: the minimum value in the channel.
955 % o maxima: the maximum value in the channel.
957 % o exception: return any errors or warnings in this structure.
961 MagickExport MagickBooleanType GetImageExtrema(const Image *image,
962 size_t *minima,size_t *maxima,ExceptionInfo *exception)
964 return(GetImageChannelExtrema(image,AllChannels,minima,maxima,exception));
967 MagickExport MagickBooleanType GetImageChannelExtrema(const Image *image,
968 const ChannelType channel,size_t *minima,size_t *maxima,
969 ExceptionInfo *exception)
978 assert(image != (Image *) NULL);
979 assert(image->signature == MagickSignature);
980 if (image->debug != MagickFalse)
981 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
982 status=GetImageChannelRange(image,channel,&min,&max,exception);
983 *minima=(size_t) ceil(min-0.5);
984 *maxima=(size_t) floor(max+0.5);
989 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
993 % G e t I m a g e C h a n n e l M e a n %
997 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
999 % GetImageChannelMean() returns the mean and standard deviation of one or more
1002 % The format of the GetImageChannelMean method is:
1004 % MagickBooleanType GetImageChannelMean(const Image *image,
1005 % const ChannelType channel,double *mean,double *standard_deviation,
1006 % ExceptionInfo *exception)
1008 % A description of each parameter follows:
1010 % o image: the image.
1012 % o channel: the channel.
1014 % o mean: the average value in the channel.
1016 % o standard_deviation: the standard deviation of the channel.
1018 % o exception: return any errors or warnings in this structure.
1022 MagickExport MagickBooleanType GetImageMean(const Image *image,double *mean,
1023 double *standard_deviation,ExceptionInfo *exception)
1028 status=GetImageChannelMean(image,AllChannels,mean,standard_deviation,
1033 MagickExport MagickBooleanType GetImageChannelMean(const Image *image,
1034 const ChannelType channel,double *mean,double *standard_deviation,
1035 ExceptionInfo *exception)
1038 *channel_statistics;
1043 assert(image != (Image *) NULL);
1044 assert(image->signature == MagickSignature);
1045 if (image->debug != MagickFalse)
1046 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1047 channel_statistics=GetImageChannelStatistics(image,exception);
1048 if (channel_statistics == (ChannelStatistics *) NULL)
1049 return(MagickFalse);
1051 channel_statistics[AllChannels].mean=0.0;
1052 channel_statistics[AllChannels].standard_deviation=0.0;
1053 if ((channel & RedChannel) != 0)
1055 channel_statistics[AllChannels].mean+=
1056 channel_statistics[RedChannel].mean;
1057 channel_statistics[AllChannels].standard_deviation+=
1058 channel_statistics[RedChannel].variance-
1059 channel_statistics[RedChannel].mean*
1060 channel_statistics[RedChannel].mean;
1063 if ((channel & GreenChannel) != 0)
1065 channel_statistics[AllChannels].mean+=
1066 channel_statistics[GreenChannel].mean;
1067 channel_statistics[AllChannels].standard_deviation+=
1068 channel_statistics[GreenChannel].variance-
1069 channel_statistics[GreenChannel].mean*
1070 channel_statistics[GreenChannel].mean;
1073 if ((channel & BlueChannel) != 0)
1075 channel_statistics[AllChannels].mean+=
1076 channel_statistics[BlueChannel].mean;
1077 channel_statistics[AllChannels].standard_deviation+=
1078 channel_statistics[BlueChannel].variance-
1079 channel_statistics[BlueChannel].mean*
1080 channel_statistics[BlueChannel].mean;
1083 if (((channel & OpacityChannel) != 0) &&
1084 (image->matte != MagickFalse))
1086 channel_statistics[AllChannels].mean+=
1087 channel_statistics[OpacityChannel].mean;
1088 channel_statistics[AllChannels].standard_deviation+=
1089 channel_statistics[OpacityChannel].variance-
1090 channel_statistics[OpacityChannel].mean*
1091 channel_statistics[OpacityChannel].mean;
1094 if (((channel & IndexChannel) != 0) &&
1095 (image->colorspace == CMYKColorspace))
1097 channel_statistics[AllChannels].mean+=
1098 channel_statistics[BlackChannel].mean;
1099 channel_statistics[AllChannels].standard_deviation+=
1100 channel_statistics[BlackChannel].variance-
1101 channel_statistics[BlackChannel].mean*
1102 channel_statistics[BlackChannel].mean;
1105 channel_statistics[AllChannels].mean/=channels;
1106 channel_statistics[AllChannels].standard_deviation=
1107 sqrt(channel_statistics[AllChannels].standard_deviation/channels);
1108 *mean=channel_statistics[AllChannels].mean;
1109 *standard_deviation=channel_statistics[AllChannels].standard_deviation;
1110 channel_statistics=(ChannelStatistics *) RelinquishMagickMemory(
1111 channel_statistics);
1116 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1120 % G e t I m a g e C h a n n e l K u r t o s i s %
1124 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1126 % GetImageChannelKurtosis() returns the kurtosis and skewness of one or more
1129 % The format of the GetImageChannelKurtosis method is:
1131 % MagickBooleanType GetImageChannelKurtosis(const Image *image,
1132 % const ChannelType channel,double *kurtosis,double *skewness,
1133 % ExceptionInfo *exception)
1135 % A description of each parameter follows:
1137 % o image: the image.
1139 % o channel: the channel.
1141 % o kurtosis: the kurtosis of the channel.
1143 % o skewness: the skewness of the channel.
1145 % o exception: return any errors or warnings in this structure.
1149 MagickExport MagickBooleanType GetImageKurtosis(const Image *image,
1150 double *kurtosis,double *skewness,ExceptionInfo *exception)
1155 status=GetImageChannelKurtosis(image,AllChannels,kurtosis,skewness,
1160 MagickExport MagickBooleanType GetImageChannelKurtosis(const Image *image,
1161 const ChannelType channel,double *kurtosis,double *skewness,
1162 ExceptionInfo *exception)
1175 assert(image != (Image *) NULL);
1176 assert(image->signature == MagickSignature);
1177 if (image->debug != MagickFalse)
1178 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1183 standard_deviation=0.0;
1186 sum_fourth_power=0.0;
1187 for (y=0; y < (ssize_t) image->rows; y++)
1189 register const IndexPacket
1192 register const PixelPacket
1198 p=GetVirtualPixels(image,0,y,image->columns,1,exception);
1199 if (p == (const PixelPacket *) NULL)
1201 indexes=GetVirtualIndexQueue(image);
1202 for (x=0; x < (ssize_t) image->columns; x++)
1204 if ((channel & RedChannel) != 0)
1206 mean+=GetRedPixelComponent(p);
1207 sum_squares+=(double) p->red*GetRedPixelComponent(p);
1208 sum_cubes+=(double) p->red*p->red*GetRedPixelComponent(p);
1209 sum_fourth_power+=(double) p->red*p->red*p->red*
1210 GetRedPixelComponent(p);
1213 if ((channel & GreenChannel) != 0)
1215 mean+=GetGreenPixelComponent(p);
1216 sum_squares+=(double) p->green*GetGreenPixelComponent(p);
1217 sum_cubes+=(double) p->green*p->green*GetGreenPixelComponent(p);
1218 sum_fourth_power+=(double) p->green*p->green*p->green*
1219 GetGreenPixelComponent(p);
1222 if ((channel & BlueChannel) != 0)
1224 mean+=GetBluePixelComponent(p);
1225 sum_squares+=(double) p->blue*GetBluePixelComponent(p);
1226 sum_cubes+=(double) p->blue*p->blue*GetBluePixelComponent(p);
1227 sum_fourth_power+=(double) p->blue*p->blue*p->blue*
1228 GetBluePixelComponent(p);
1231 if ((channel & OpacityChannel) != 0)
1233 mean+=GetOpacityPixelComponent(p);
1234 sum_squares+=(double) p->opacity*GetOpacityPixelComponent(p);
1235 sum_cubes+=(double) p->opacity*p->opacity*GetOpacityPixelComponent(p);
1236 sum_fourth_power+=(double) p->opacity*p->opacity*p->opacity*
1237 GetOpacityPixelComponent(p);
1240 if (((channel & IndexChannel) != 0) &&
1241 (image->colorspace == CMYKColorspace))
1244 sum_squares+=(double) indexes[x]*indexes[x];
1245 sum_cubes+=(double) indexes[x]*indexes[x]*indexes[x];
1246 sum_fourth_power+=(double) indexes[x]*indexes[x]*indexes[x]*
1253 if (y < (ssize_t) image->rows)
1254 return(MagickFalse);
1260 sum_fourth_power/=area;
1262 standard_deviation=sqrt(sum_squares-(mean*mean));
1263 if (standard_deviation != 0.0)
1265 *kurtosis=sum_fourth_power-4.0*mean*sum_cubes+6.0*mean*mean*sum_squares-
1266 3.0*mean*mean*mean*mean;
1267 *kurtosis/=standard_deviation*standard_deviation*standard_deviation*
1270 *skewness=sum_cubes-3.0*mean*sum_squares+2.0*mean*mean*mean;
1271 *skewness/=standard_deviation*standard_deviation*standard_deviation;
1273 return(y == (ssize_t) image->rows ? MagickTrue : MagickFalse);
1277 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1281 % G e t I m a g e C h a n n e l R a n g e %
1285 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1287 % GetImageChannelRange() returns the range of one or more image channels.
1289 % The format of the GetImageChannelRange method is:
1291 % MagickBooleanType GetImageChannelRange(const Image *image,
1292 % const ChannelType channel,double *minima,double *maxima,
1293 % ExceptionInfo *exception)
1295 % A description of each parameter follows:
1297 % o image: the image.
1299 % o channel: the channel.
1301 % o minima: the minimum value in the channel.
1303 % o maxima: the maximum value in the channel.
1305 % o exception: return any errors or warnings in this structure.
1309 MagickExport MagickBooleanType GetImageRange(const Image *image,
1310 double *minima,double *maxima,ExceptionInfo *exception)
1312 return(GetImageChannelRange(image,AllChannels,minima,maxima,exception));
1315 MagickExport MagickBooleanType GetImageChannelRange(const Image *image,
1316 const ChannelType channel,double *minima,double *maxima,
1317 ExceptionInfo *exception)
1325 assert(image != (Image *) NULL);
1326 assert(image->signature == MagickSignature);
1327 if (image->debug != MagickFalse)
1328 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1331 GetMagickPixelPacket(image,&pixel);
1332 for (y=0; y < (ssize_t) image->rows; y++)
1334 register const IndexPacket
1337 register const PixelPacket
1343 p=GetVirtualPixels(image,0,y,image->columns,1,exception);
1344 if (p == (const PixelPacket *) NULL)
1346 indexes=GetVirtualIndexQueue(image);
1347 for (x=0; x < (ssize_t) image->columns; x++)
1349 SetMagickPixelPacket(image,p,indexes+x,&pixel);
1350 if ((channel & RedChannel) != 0)
1352 if (pixel.red < *minima)
1353 *minima=(double) pixel.red;
1354 if (pixel.red > *maxima)
1355 *maxima=(double) pixel.red;
1357 if ((channel & GreenChannel) != 0)
1359 if (pixel.green < *minima)
1360 *minima=(double) pixel.green;
1361 if (pixel.green > *maxima)
1362 *maxima=(double) pixel.green;
1364 if ((channel & BlueChannel) != 0)
1366 if (pixel.blue < *minima)
1367 *minima=(double) pixel.blue;
1368 if (pixel.blue > *maxima)
1369 *maxima=(double) pixel.blue;
1371 if ((channel & OpacityChannel) != 0)
1373 if (pixel.opacity < *minima)
1374 *minima=(double) pixel.opacity;
1375 if (pixel.opacity > *maxima)
1376 *maxima=(double) pixel.opacity;
1378 if (((channel & IndexChannel) != 0) &&
1379 (image->colorspace == CMYKColorspace))
1381 if ((double) indexes[x] < *minima)
1382 *minima=(double) indexes[x];
1383 if ((double) indexes[x] > *maxima)
1384 *maxima=(double) indexes[x];
1389 return(y == (ssize_t) image->rows ? MagickTrue : MagickFalse);
1393 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1397 % 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 %
1401 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1403 % GetImageChannelStatistics() returns statistics for each channel in the
1404 % image. The statistics include the channel depth, its minima, maxima, mean,
1405 % standard deviation, kurtosis and skewness. You can access the red channel
1406 % mean, for example, like this:
1408 % channel_statistics=GetImageChannelStatistics(image,excepton);
1409 % red_mean=channel_statistics[RedChannel].mean;
1411 % Use MagickRelinquishMemory() to free the statistics buffer.
1413 % The format of the GetImageChannelStatistics method is:
1415 % ChannelStatistics *GetImageChannelStatistics(const Image *image,
1416 % ExceptionInfo *exception)
1418 % A description of each parameter follows:
1420 % o image: the image.
1422 % o exception: return any errors or warnings in this structure.
1425 MagickExport ChannelStatistics *GetImageChannelStatistics(const Image *image,
1426 ExceptionInfo *exception)
1429 *channel_statistics;
1453 assert(image != (Image *) NULL);
1454 assert(image->signature == MagickSignature);
1455 if (image->debug != MagickFalse)
1456 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1457 length=AllChannels+1UL;
1458 channel_statistics=(ChannelStatistics *) AcquireQuantumMemory(length,
1459 sizeof(*channel_statistics));
1460 if (channel_statistics == (ChannelStatistics *) NULL)
1461 ThrowFatalException(ResourceLimitFatalError,"MemoryAllocationFailed");
1462 (void) ResetMagickMemory(channel_statistics,0,length*
1463 sizeof(*channel_statistics));
1464 for (i=0; i <= AllChannels; i++)
1466 channel_statistics[i].depth=1;
1467 channel_statistics[i].maxima=(-1.0E-37);
1468 channel_statistics[i].minima=1.0E+37;
1470 for (y=0; y < (ssize_t) image->rows; y++)
1472 register const IndexPacket
1475 register const PixelPacket
1481 p=GetVirtualPixels(image,0,y,image->columns,1,exception);
1482 if (p == (const PixelPacket *) NULL)
1484 indexes=GetVirtualIndexQueue(image);
1485 for (x=0; x < (ssize_t) image->columns; )
1487 if (channel_statistics[RedChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1489 depth=channel_statistics[RedChannel].depth;
1490 range=GetQuantumRange(depth);
1491 status=p->red != ScaleAnyToQuantum(ScaleQuantumToAny(p->red,range),
1492 range) ? MagickTrue : MagickFalse;
1493 if (status != MagickFalse)
1495 channel_statistics[RedChannel].depth++;
1499 if (channel_statistics[GreenChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1501 depth=channel_statistics[GreenChannel].depth;
1502 range=GetQuantumRange(depth);
1503 status=p->green != ScaleAnyToQuantum(ScaleQuantumToAny(p->green,
1504 range),range) ? MagickTrue : MagickFalse;
1505 if (status != MagickFalse)
1507 channel_statistics[GreenChannel].depth++;
1511 if (channel_statistics[BlueChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1513 depth=channel_statistics[BlueChannel].depth;
1514 range=GetQuantumRange(depth);
1515 status=p->blue != ScaleAnyToQuantum(ScaleQuantumToAny(p->blue,
1516 range),range) ? MagickTrue : MagickFalse;
1517 if (status != MagickFalse)
1519 channel_statistics[BlueChannel].depth++;
1523 if (image->matte != MagickFalse)
1525 if (channel_statistics[OpacityChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1527 depth=channel_statistics[OpacityChannel].depth;
1528 range=GetQuantumRange(depth);
1529 status=p->opacity != ScaleAnyToQuantum(ScaleQuantumToAny(
1530 p->opacity,range),range) ? MagickTrue : MagickFalse;
1531 if (status != MagickFalse)
1533 channel_statistics[OpacityChannel].depth++;
1538 if (image->colorspace == CMYKColorspace)
1540 if (channel_statistics[BlackChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1542 depth=channel_statistics[BlackChannel].depth;
1543 range=GetQuantumRange(depth);
1544 status=indexes[x] != ScaleAnyToQuantum(ScaleQuantumToAny(
1545 indexes[x],range),range) ? MagickTrue : MagickFalse;
1546 if (status != MagickFalse)
1548 channel_statistics[BlackChannel].depth++;
1553 if ((double) p->red < channel_statistics[RedChannel].minima)
1554 channel_statistics[RedChannel].minima=(double) GetRedPixelComponent(p);
1555 if ((double) p->red > channel_statistics[RedChannel].maxima)
1556 channel_statistics[RedChannel].maxima=(double) GetRedPixelComponent(p);
1557 channel_statistics[RedChannel].sum+=GetRedPixelComponent(p);
1558 channel_statistics[RedChannel].sum_squared+=(double) p->red*
1559 GetRedPixelComponent(p);
1560 channel_statistics[RedChannel].sum_cubed+=(double) p->red*p->red*
1561 GetRedPixelComponent(p);
1562 channel_statistics[RedChannel].sum_fourth_power+=(double) p->red*p->red*
1563 p->red*GetRedPixelComponent(p);
1564 if ((double) p->green < channel_statistics[GreenChannel].minima)
1565 channel_statistics[GreenChannel].minima=(double)
1566 GetGreenPixelComponent(p);
1567 if ((double) p->green > channel_statistics[GreenChannel].maxima)
1568 channel_statistics[GreenChannel].maxima=(double)
1569 GetGreenPixelComponent(p);
1570 channel_statistics[GreenChannel].sum+=GetGreenPixelComponent(p);
1571 channel_statistics[GreenChannel].sum_squared+=(double) p->green*
1572 GetGreenPixelComponent(p);
1573 channel_statistics[GreenChannel].sum_cubed+=(double) p->green*p->green*
1574 GetGreenPixelComponent(p);
1575 channel_statistics[GreenChannel].sum_fourth_power+=(double) p->green*
1576 p->green*p->green*GetGreenPixelComponent(p);
1577 if ((double) p->blue < channel_statistics[BlueChannel].minima)
1578 channel_statistics[BlueChannel].minima=(double)
1579 GetBluePixelComponent(p);
1580 if ((double) p->blue > channel_statistics[BlueChannel].maxima)
1581 channel_statistics[BlueChannel].maxima=(double)
1582 GetBluePixelComponent(p);
1583 channel_statistics[BlueChannel].sum+=GetBluePixelComponent(p);
1584 channel_statistics[BlueChannel].sum_squared+=(double) p->blue*
1585 GetBluePixelComponent(p);
1586 channel_statistics[BlueChannel].sum_cubed+=(double) p->blue*p->blue*
1587 GetBluePixelComponent(p);
1588 channel_statistics[BlueChannel].sum_fourth_power+=(double) p->blue*
1589 p->blue*p->blue*GetBluePixelComponent(p);
1590 if (image->matte != MagickFalse)
1592 if ((double) p->opacity < channel_statistics[OpacityChannel].minima)
1593 channel_statistics[OpacityChannel].minima=(double)
1594 GetOpacityPixelComponent(p);
1595 if ((double) p->opacity > channel_statistics[OpacityChannel].maxima)
1596 channel_statistics[OpacityChannel].maxima=(double)
1597 GetOpacityPixelComponent(p);
1598 channel_statistics[OpacityChannel].sum+=GetOpacityPixelComponent(p);
1599 channel_statistics[OpacityChannel].sum_squared+=(double)
1600 p->opacity*GetOpacityPixelComponent(p);
1601 channel_statistics[OpacityChannel].sum_cubed+=(double) p->opacity*
1602 p->opacity*GetOpacityPixelComponent(p);
1603 channel_statistics[OpacityChannel].sum_fourth_power+=(double)
1604 p->opacity*p->opacity*p->opacity*GetOpacityPixelComponent(p);
1606 if (image->colorspace == CMYKColorspace)
1608 if ((double) indexes[x] < channel_statistics[BlackChannel].minima)
1609 channel_statistics[BlackChannel].minima=(double) indexes[x];
1610 if ((double) indexes[x] > channel_statistics[BlackChannel].maxima)
1611 channel_statistics[BlackChannel].maxima=(double) indexes[x];
1612 channel_statistics[BlackChannel].sum+=indexes[x];
1613 channel_statistics[BlackChannel].sum_squared+=(double)
1614 indexes[x]*indexes[x];
1615 channel_statistics[BlackChannel].sum_cubed+=(double) indexes[x]*
1616 indexes[x]*indexes[x];
1617 channel_statistics[BlackChannel].sum_fourth_power+=(double)
1618 indexes[x]*indexes[x]*indexes[x]*indexes[x];
1624 area=(double) image->columns*image->rows;
1625 for (i=0; i < AllChannels; i++)
1627 channel_statistics[i].sum/=area;
1628 channel_statistics[i].sum_squared/=area;
1629 channel_statistics[i].sum_cubed/=area;
1630 channel_statistics[i].sum_fourth_power/=area;
1631 channel_statistics[i].mean=channel_statistics[i].sum;
1632 channel_statistics[i].variance=channel_statistics[i].sum_squared;
1633 channel_statistics[i].standard_deviation=sqrt(
1634 channel_statistics[i].variance-(channel_statistics[i].mean*
1635 channel_statistics[i].mean));
1637 for (i=0; i < AllChannels; i++)
1639 channel_statistics[AllChannels].depth=(size_t) MagickMax((double)
1640 channel_statistics[AllChannels].depth,(double)
1641 channel_statistics[i].depth);
1642 channel_statistics[AllChannels].minima=MagickMin(
1643 channel_statistics[AllChannels].minima,channel_statistics[i].minima);
1644 channel_statistics[AllChannels].maxima=MagickMax(
1645 channel_statistics[AllChannels].maxima,channel_statistics[i].maxima);
1646 channel_statistics[AllChannels].sum+=channel_statistics[i].sum;
1647 channel_statistics[AllChannels].sum_squared+=
1648 channel_statistics[i].sum_squared;
1649 channel_statistics[AllChannels].sum_cubed+=channel_statistics[i].sum_cubed;
1650 channel_statistics[AllChannels].sum_fourth_power+=
1651 channel_statistics[i].sum_fourth_power;
1652 channel_statistics[AllChannels].mean+=channel_statistics[i].mean;
1653 channel_statistics[AllChannels].variance+=channel_statistics[i].variance-
1654 channel_statistics[i].mean*channel_statistics[i].mean;
1655 channel_statistics[AllChannels].standard_deviation+=
1656 channel_statistics[i].variance-channel_statistics[i].mean*
1657 channel_statistics[i].mean;
1660 if (image->matte != MagickFalse)
1662 if (image->colorspace == CMYKColorspace)
1664 channel_statistics[AllChannels].sum/=channels;
1665 channel_statistics[AllChannels].sum_squared/=channels;
1666 channel_statistics[AllChannels].sum_cubed/=channels;
1667 channel_statistics[AllChannels].sum_fourth_power/=channels;
1668 channel_statistics[AllChannels].mean/=channels;
1669 channel_statistics[AllChannels].variance/=channels;
1670 channel_statistics[AllChannels].standard_deviation=
1671 sqrt(channel_statistics[AllChannels].standard_deviation/channels);
1672 channel_statistics[AllChannels].kurtosis/=channels;
1673 channel_statistics[AllChannels].skewness/=channels;
1674 for (i=0; i <= AllChannels; i++)
1676 if (channel_statistics[i].standard_deviation == 0.0)
1678 channel_statistics[i].skewness=(channel_statistics[i].sum_cubed-
1679 3.0*channel_statistics[i].mean*channel_statistics[i].sum_squared+
1680 2.0*channel_statistics[i].mean*channel_statistics[i].mean*
1681 channel_statistics[i].mean)/(channel_statistics[i].standard_deviation*
1682 channel_statistics[i].standard_deviation*
1683 channel_statistics[i].standard_deviation);
1684 channel_statistics[i].kurtosis=(channel_statistics[i].sum_fourth_power-
1685 4.0*channel_statistics[i].mean*channel_statistics[i].sum_cubed+
1686 6.0*channel_statistics[i].mean*channel_statistics[i].mean*
1687 channel_statistics[i].sum_squared-3.0*channel_statistics[i].mean*
1688 channel_statistics[i].mean*1.0*channel_statistics[i].mean*
1689 channel_statistics[i].mean)/(channel_statistics[i].standard_deviation*
1690 channel_statistics[i].standard_deviation*
1691 channel_statistics[i].standard_deviation*
1692 channel_statistics[i].standard_deviation)-3.0;
1694 return(channel_statistics);