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 **) RelinquishAlignedMemory(pixels);
148 static MagickPixelPacket **AcquirePixelThreadSet(const Image *image)
160 number_threads=GetOpenMPMaximumThreads();
161 pixels=(MagickPixelPacket **) AcquireAlignedMemory(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 AddEvaluateOperator:
205 result=(MagickRealType) (pixel+value);
208 case AddModulusEvaluateOperator:
211 This returns a 'floored modulus' of the addition which is a
212 positive result. It differs from % or fmod() which returns a
213 'truncated modulus' result, where floor() is replaced by trunc()
214 and could return a negative result (which is clipped).
217 result-=(QuantumRange+1.0)*floor((double) result/(QuantumRange+1.0));
220 case AndEvaluateOperator:
222 result=(MagickRealType) ((size_t) pixel & (size_t)
226 case CosineEvaluateOperator:
228 result=(MagickRealType) (QuantumRange*(0.5*cos((double) (2.0*MagickPI*
229 QuantumScale*pixel*value))+0.5));
232 case DivideEvaluateOperator:
234 result=pixel/(value == 0.0 ? 1.0 : value);
237 case GaussianNoiseEvaluateOperator:
239 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
240 GaussianNoise,value);
243 case ImpulseNoiseEvaluateOperator:
245 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
249 case LaplacianNoiseEvaluateOperator:
251 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
252 LaplacianNoise,value);
255 case LeftShiftEvaluateOperator:
257 result=(MagickRealType) ((size_t) pixel << (size_t)
261 case LogEvaluateOperator:
263 result=(MagickRealType) (QuantumRange*log((double) (QuantumScale*value*
264 pixel+1.0))/log((double) (value+1.0)));
267 case MaxEvaluateOperator:
269 result=(MagickRealType) MagickMax((double) pixel,value);
272 case MeanEvaluateOperator:
274 result=(MagickRealType) (pixel+value);
277 case MinEvaluateOperator:
279 result=(MagickRealType) MagickMin((double) pixel,value);
282 case MultiplicativeNoiseEvaluateOperator:
284 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
285 MultiplicativeGaussianNoise,value);
288 case MultiplyEvaluateOperator:
290 result=(MagickRealType) (value*pixel);
293 case OrEvaluateOperator:
295 result=(MagickRealType) ((size_t) pixel | (size_t)
299 case PoissonNoiseEvaluateOperator:
301 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
305 case PowEvaluateOperator:
307 result=(MagickRealType) (QuantumRange*pow((double) (QuantumScale*pixel),
311 case RightShiftEvaluateOperator:
313 result=(MagickRealType) ((size_t) pixel >> (size_t)
317 case SetEvaluateOperator:
322 case SineEvaluateOperator:
324 result=(MagickRealType) (QuantumRange*(0.5*sin((double) (2.0*MagickPI*
325 QuantumScale*pixel*value))+0.5));
328 case SubtractEvaluateOperator:
330 result=(MagickRealType) (pixel-value);
333 case ThresholdEvaluateOperator:
335 result=(MagickRealType) (((MagickRealType) pixel <= value) ? 0 :
339 case ThresholdBlackEvaluateOperator:
341 result=(MagickRealType) (((MagickRealType) pixel <= value) ? 0 : pixel);
344 case ThresholdWhiteEvaluateOperator:
346 result=(MagickRealType) (((MagickRealType) pixel > value) ? QuantumRange :
350 case UniformNoiseEvaluateOperator:
352 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
356 case XorEvaluateOperator:
358 result=(MagickRealType) ((size_t) pixel ^ (size_t)
366 MagickExport MagickBooleanType EvaluateImage(Image *image,
367 const MagickEvaluateOperator op,const double value,ExceptionInfo *exception)
372 status=EvaluateImageChannel(image,AllChannels,op,value,exception);
376 MagickExport Image *EvaluateImages(const Image *images,
377 const MagickEvaluateOperator op,ExceptionInfo *exception)
379 #define EvaluateImageTag "Evaluate/Image"
397 **restrict evaluate_pixels,
401 **restrict random_info;
410 Ensure the image are the same size.
412 assert(images != (Image *) NULL);
413 assert(images->signature == MagickSignature);
414 if (images->debug != MagickFalse)
415 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",images->filename);
416 assert(exception != (ExceptionInfo *) NULL);
417 assert(exception->signature == MagickSignature);
418 for (next=images; next != (Image *) NULL; next=GetNextImageInList(next))
419 if ((next->columns != images->columns) || (next->rows != images->rows))
421 (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
422 "ImageWidthsOrHeightsDiffer","`%s'",images->filename);
423 return((Image *) NULL);
426 Initialize evaluate next attributes.
428 evaluate_image=CloneImage(images,images->columns,images->rows,MagickTrue,
430 if (evaluate_image == (Image *) NULL)
431 return((Image *) NULL);
432 if (SetImageStorageClass(evaluate_image,DirectClass) == MagickFalse)
434 InheritException(exception,&evaluate_image->exception);
435 evaluate_image=DestroyImage(evaluate_image);
436 return((Image *) NULL);
438 evaluate_pixels=AcquirePixelThreadSet(images);
439 if (evaluate_pixels == (MagickPixelPacket **) NULL)
441 evaluate_image=DestroyImage(evaluate_image);
442 (void) ThrowMagickException(exception,GetMagickModule(),
443 ResourceLimitError,"MemoryAllocationFailed","`%s'",images->filename);
444 return((Image *) NULL);
447 Evaluate image pixels.
451 GetMagickPixelPacket(images,&zero);
452 random_info=AcquireRandomInfoThreadSet();
453 number_images=GetImageListLength(images);
454 evaluate_view=AcquireCacheView(evaluate_image);
455 #if defined(MAGICKCORE_OPENMP_SUPPORT)
456 #pragma omp parallel for schedule(dynamic) shared(progress,status)
458 for (y=0; y < (ssize_t) evaluate_image->rows; y++)
473 *restrict evaluate_indexes;
479 register MagickPixelPacket
485 if (status == MagickFalse)
487 q=QueueCacheViewAuthenticPixels(evaluate_view,0,y,evaluate_image->columns,1,
489 if (q == (PixelPacket *) NULL)
494 evaluate_indexes=GetCacheViewAuthenticIndexQueue(evaluate_view);
496 id=GetOpenMPThreadId();
497 evaluate_pixel=evaluate_pixels[id];
498 for (x=0; x < (ssize_t) evaluate_image->columns; x++)
499 evaluate_pixel[x]=zero;
501 for (i=0; i < (ssize_t) number_images; i++)
503 register const IndexPacket
506 register const PixelPacket
509 image_view=AcquireCacheView(next);
510 p=GetCacheViewVirtualPixels(image_view,0,y,next->columns,1,exception);
511 if (p == (const PixelPacket *) NULL)
513 image_view=DestroyCacheView(image_view);
516 indexes=GetCacheViewVirtualIndexQueue(image_view);
517 for (x=0; x < (ssize_t) next->columns; x++)
519 evaluate_pixel[x].red=ApplyEvaluateOperator(random_info[id],p->red,
520 i == 0 ? AddEvaluateOperator : op,evaluate_pixel[x].red);
521 evaluate_pixel[x].green=ApplyEvaluateOperator(random_info[id],p->green,
522 i == 0 ? AddEvaluateOperator : op,evaluate_pixel[x].green);
523 evaluate_pixel[x].blue=ApplyEvaluateOperator(random_info[id],p->blue,
524 i == 0 ? AddEvaluateOperator : op,evaluate_pixel[x].blue);
525 evaluate_pixel[x].opacity=ApplyEvaluateOperator(random_info[id],
526 p->opacity,i == 0 ? AddEvaluateOperator : op,
527 evaluate_pixel[x].opacity);
528 if (evaluate_image->colorspace == CMYKColorspace)
529 evaluate_pixel[x].index=ApplyEvaluateOperator(random_info[id],
530 indexes[x],i == 0 ? AddEvaluateOperator : op,
531 evaluate_pixel[x].index);
534 image_view=DestroyCacheView(image_view);
535 next=GetNextImageInList(next);
537 if (op == MeanEvaluateOperator)
538 for (x=0; x < (ssize_t) evaluate_image->columns; x++)
540 evaluate_pixel[x].red/=number_images;
541 evaluate_pixel[x].green/=number_images;
542 evaluate_pixel[x].blue/=number_images;
543 evaluate_pixel[x].opacity/=number_images;
544 evaluate_pixel[x].index/=number_images;
546 for (x=0; x < (ssize_t) evaluate_image->columns; x++)
548 q->red=ClampToQuantum(evaluate_pixel[x].red);
549 q->green=ClampToQuantum(evaluate_pixel[x].green);
550 q->blue=ClampToQuantum(evaluate_pixel[x].blue);
551 if (evaluate_image->matte == MagickFalse)
552 q->opacity=ClampToQuantum(evaluate_pixel[x].opacity);
554 q->opacity=ClampToQuantum(QuantumRange-evaluate_pixel[x].opacity);
555 if (evaluate_image->colorspace == CMYKColorspace)
556 evaluate_indexes[x]=ClampToQuantum(evaluate_pixel[x].index);
559 if (SyncCacheViewAuthenticPixels(evaluate_view,exception) == MagickFalse)
561 if (images->progress_monitor != (MagickProgressMonitor) NULL)
566 #if defined(MAGICKCORE_OPENMP_SUPPORT)
567 #pragma omp critical (MagickCore_EvaluateImages)
569 proceed=SetImageProgress(images,EvaluateImageTag,progress++,
570 evaluate_image->rows);
571 if (proceed == MagickFalse)
575 evaluate_view=DestroyCacheView(evaluate_view);
576 evaluate_pixels=DestroyPixelThreadSet(evaluate_pixels);
577 random_info=DestroyRandomInfoThreadSet(random_info);
578 if (status == MagickFalse)
579 evaluate_image=DestroyImage(evaluate_image);
580 return(evaluate_image);
583 MagickExport MagickBooleanType EvaluateImageChannel(Image *image,
584 const ChannelType channel,const MagickEvaluateOperator op,const double value,
585 ExceptionInfo *exception)
597 **restrict random_info;
602 assert(image != (Image *) NULL);
603 assert(image->signature == MagickSignature);
604 if (image->debug != MagickFalse)
605 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
606 assert(exception != (ExceptionInfo *) NULL);
607 assert(exception->signature == MagickSignature);
608 if (SetImageStorageClass(image,DirectClass) == MagickFalse)
610 InheritException(exception,&image->exception);
615 random_info=AcquireRandomInfoThreadSet();
616 image_view=AcquireCacheView(image);
617 #if defined(MAGICKCORE_OPENMP_SUPPORT)
618 #pragma omp parallel for schedule(dynamic,4) shared(progress,status)
620 for (y=0; y < (ssize_t) image->rows; y++)
634 if (status == MagickFalse)
636 q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
637 if (q == (PixelPacket *) NULL)
642 indexes=GetCacheViewAuthenticIndexQueue(image_view);
643 id=GetOpenMPThreadId();
644 for (x=0; x < (ssize_t) image->columns; x++)
646 if ((channel & RedChannel) != 0)
647 q->red=ClampToQuantum(ApplyEvaluateOperator(random_info[id],q->red,op,
649 if ((channel & GreenChannel) != 0)
650 q->green=ClampToQuantum(ApplyEvaluateOperator(random_info[id],q->green,
652 if ((channel & BlueChannel) != 0)
653 q->blue=ClampToQuantum(ApplyEvaluateOperator(random_info[id],q->blue,op,
655 if ((channel & OpacityChannel) != 0)
657 if (image->matte == MagickFalse)
658 q->opacity=ClampToQuantum(ApplyEvaluateOperator(random_info[id],
659 q->opacity,op,value));
661 q->opacity=ClampToQuantum(QuantumRange-ApplyEvaluateOperator(
662 random_info[id],(Quantum) GetAlphaPixelComponent(q),op,value));
664 if (((channel & IndexChannel) != 0) && (indexes != (IndexPacket *) NULL))
665 indexes[x]=(IndexPacket) ClampToQuantum(ApplyEvaluateOperator(
666 random_info[id],indexes[x],op,value));
669 if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
671 if (image->progress_monitor != (MagickProgressMonitor) NULL)
676 #if defined(MAGICKCORE_OPENMP_SUPPORT)
677 #pragma omp critical (MagickCore_EvaluateImageChannel)
679 proceed=SetImageProgress(image,EvaluateImageTag,progress++,image->rows);
680 if (proceed == MagickFalse)
684 image_view=DestroyCacheView(image_view);
685 random_info=DestroyRandomInfoThreadSet(random_info);
690 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
694 % F u n c t i o n I m a g e %
698 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
700 % FunctionImage() applies a value to the image with an arithmetic, relational,
701 % or logical operator to an image. Use these operations to lighten or darken
702 % an image, to increase or decrease contrast in an image, or to produce the
703 % "negative" of an image.
705 % The format of the FunctionImageChannel method is:
707 % MagickBooleanType FunctionImage(Image *image,
708 % const MagickFunction function,const ssize_t number_parameters,
709 % const double *parameters,ExceptionInfo *exception)
710 % MagickBooleanType FunctionImageChannel(Image *image,
711 % const ChannelType channel,const MagickFunction function,
712 % const ssize_t number_parameters,const double *argument,
713 % ExceptionInfo *exception)
715 % A description of each parameter follows:
717 % o image: the image.
719 % o channel: the channel.
721 % o function: A channel function.
723 % o parameters: one or more parameters.
725 % o exception: return any errors or warnings in this structure.
729 static Quantum ApplyFunction(Quantum pixel,const MagickFunction function,
730 const size_t number_parameters,const double *parameters,
731 ExceptionInfo *exception)
743 case PolynomialFunction:
747 * Parameters: polynomial constants, highest to lowest order
748 * For example: c0*x^3 + c1*x^2 + c2*x + c3
751 for (i=0; i < (ssize_t) number_parameters; i++)
752 result = result*QuantumScale*pixel + parameters[i];
753 result *= QuantumRange;
756 case SinusoidFunction:
759 * Parameters: Freq, Phase, Ampl, bias
761 double freq,phase,ampl,bias;
762 freq = ( number_parameters >= 1 ) ? parameters[0] : 1.0;
763 phase = ( number_parameters >= 2 ) ? parameters[1] : 0.0;
764 ampl = ( number_parameters >= 3 ) ? parameters[2] : 0.5;
765 bias = ( number_parameters >= 4 ) ? parameters[3] : 0.5;
766 result=(MagickRealType) (QuantumRange*(ampl*sin((double) (2.0*MagickPI*
767 (freq*QuantumScale*pixel + phase/360.0) )) + bias ) );
772 /* Arcsin Function (peged at range limits for invalid results)
773 * Parameters: Width, Center, Range, Bias
775 double width,range,center,bias;
776 width = ( number_parameters >= 1 ) ? parameters[0] : 1.0;
777 center = ( number_parameters >= 2 ) ? parameters[1] : 0.5;
778 range = ( number_parameters >= 3 ) ? parameters[2] : 1.0;
779 bias = ( number_parameters >= 4 ) ? parameters[3] : 0.5;
780 result = 2.0/width*(QuantumScale*pixel - center);
781 if ( result <= -1.0 )
782 result = bias - range/2.0;
783 else if ( result >= 1.0 )
784 result = bias + range/2.0;
786 result=range/MagickPI*asin((double)result) + bias;
787 result *= QuantumRange;
793 * Parameters: Slope, Center, Range, Bias
795 double slope,range,center,bias;
796 slope = ( number_parameters >= 1 ) ? parameters[0] : 1.0;
797 center = ( number_parameters >= 2 ) ? parameters[1] : 0.5;
798 range = ( number_parameters >= 3 ) ? parameters[2] : 1.0;
799 bias = ( number_parameters >= 4 ) ? parameters[3] : 0.5;
800 result = MagickPI*slope*(QuantumScale*pixel - center);
801 result=(MagickRealType) (QuantumRange*(range/MagickPI*atan((double)
805 case UndefinedFunction:
808 return(ClampToQuantum(result));
811 MagickExport MagickBooleanType FunctionImage(Image *image,
812 const MagickFunction function,const size_t number_parameters,
813 const double *parameters,ExceptionInfo *exception)
818 status=FunctionImageChannel(image,AllChannels,function,number_parameters,
819 parameters,exception);
823 MagickExport MagickBooleanType FunctionImageChannel(Image *image,
824 const ChannelType channel,const MagickFunction function,
825 const size_t number_parameters,const double *parameters,
826 ExceptionInfo *exception)
828 #define FunctionImageTag "Function/Image "
842 assert(image != (Image *) NULL);
843 assert(image->signature == MagickSignature);
844 if (image->debug != MagickFalse)
845 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
846 assert(exception != (ExceptionInfo *) NULL);
847 assert(exception->signature == MagickSignature);
848 if (SetImageStorageClass(image,DirectClass) == MagickFalse)
850 InheritException(exception,&image->exception);
855 image_view=AcquireCacheView(image);
856 #if defined(MAGICKCORE_OPENMP_SUPPORT)
857 #pragma omp parallel for schedule(dynamic,4) shared(progress,status)
859 for (y=0; y < (ssize_t) image->rows; y++)
870 if (status == MagickFalse)
872 q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
873 if (q == (PixelPacket *) NULL)
878 indexes=GetCacheViewAuthenticIndexQueue(image_view);
879 for (x=0; x < (ssize_t) image->columns; x++)
881 if ((channel & RedChannel) != 0)
882 q->red=ApplyFunction(q->red,function,number_parameters,parameters,
884 if ((channel & GreenChannel) != 0)
885 q->green=ApplyFunction(q->green,function,number_parameters,parameters,
887 if ((channel & BlueChannel) != 0)
888 q->blue=ApplyFunction(q->blue,function,number_parameters,parameters,
890 if ((channel & OpacityChannel) != 0)
892 if (image->matte == MagickFalse)
893 q->opacity=ApplyFunction(q->opacity,function,number_parameters,
894 parameters,exception);
896 q->opacity=(Quantum) QuantumRange-ApplyFunction((Quantum)
897 GetAlphaPixelComponent(q),function,number_parameters,parameters,
900 if (((channel & IndexChannel) != 0) && (indexes != (IndexPacket *) NULL))
901 indexes[x]=(IndexPacket) ApplyFunction(indexes[x],function,
902 number_parameters,parameters,exception);
905 if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
907 if (image->progress_monitor != (MagickProgressMonitor) NULL)
912 #if defined(MAGICKCORE_OPENMP_SUPPORT)
913 #pragma omp critical (MagickCore_FunctionImageChannel)
915 proceed=SetImageProgress(image,FunctionImageTag,progress++,image->rows);
916 if (proceed == MagickFalse)
920 image_view=DestroyCacheView(image_view);
925 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
929 + G e t I m a g e C h a n n e l E x t r e m a %
933 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
935 % GetImageChannelExtrema() returns the extrema of one or more image channels.
937 % The format of the GetImageChannelExtrema method is:
939 % MagickBooleanType GetImageChannelExtrema(const Image *image,
940 % const ChannelType channel,size_t *minima,size_t *maxima,
941 % ExceptionInfo *exception)
943 % A description of each parameter follows:
945 % o image: the image.
947 % o channel: the channel.
949 % o minima: the minimum value in the channel.
951 % o maxima: the maximum value in the channel.
953 % o exception: return any errors or warnings in this structure.
957 MagickExport MagickBooleanType GetImageExtrema(const Image *image,
958 size_t *minima,size_t *maxima,ExceptionInfo *exception)
960 return(GetImageChannelExtrema(image,AllChannels,minima,maxima,exception));
963 MagickExport MagickBooleanType GetImageChannelExtrema(const Image *image,
964 const ChannelType channel,size_t *minima,size_t *maxima,
965 ExceptionInfo *exception)
974 assert(image != (Image *) NULL);
975 assert(image->signature == MagickSignature);
976 if (image->debug != MagickFalse)
977 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
978 status=GetImageChannelRange(image,channel,&min,&max,exception);
979 *minima=(size_t) ceil(min-0.5);
980 *maxima=(size_t) floor(max+0.5);
985 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
989 % G e t I m a g e C h a n n e l M e a n %
993 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
995 % GetImageChannelMean() returns the mean and standard deviation of one or more
998 % The format of the GetImageChannelMean method is:
1000 % MagickBooleanType GetImageChannelMean(const Image *image,
1001 % const ChannelType channel,double *mean,double *standard_deviation,
1002 % ExceptionInfo *exception)
1004 % A description of each parameter follows:
1006 % o image: the image.
1008 % o channel: the channel.
1010 % o mean: the average value in the channel.
1012 % o standard_deviation: the standard deviation of the channel.
1014 % o exception: return any errors or warnings in this structure.
1018 MagickExport MagickBooleanType GetImageMean(const Image *image,double *mean,
1019 double *standard_deviation,ExceptionInfo *exception)
1024 status=GetImageChannelMean(image,AllChannels,mean,standard_deviation,
1029 MagickExport MagickBooleanType GetImageChannelMean(const Image *image,
1030 const ChannelType channel,double *mean,double *standard_deviation,
1031 ExceptionInfo *exception)
1034 *channel_statistics;
1039 assert(image != (Image *) NULL);
1040 assert(image->signature == MagickSignature);
1041 if (image->debug != MagickFalse)
1042 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1043 channel_statistics=GetImageChannelStatistics(image,exception);
1044 if (channel_statistics == (ChannelStatistics *) NULL)
1045 return(MagickFalse);
1047 channel_statistics[AllChannels].mean=0.0;
1048 channel_statistics[AllChannels].standard_deviation=0.0;
1049 if ((channel & RedChannel) != 0)
1051 channel_statistics[AllChannels].mean+=
1052 channel_statistics[RedChannel].mean;
1053 channel_statistics[AllChannels].standard_deviation+=
1054 channel_statistics[RedChannel].variance-
1055 channel_statistics[RedChannel].mean*
1056 channel_statistics[RedChannel].mean;
1059 if ((channel & GreenChannel) != 0)
1061 channel_statistics[AllChannels].mean+=
1062 channel_statistics[GreenChannel].mean;
1063 channel_statistics[AllChannels].standard_deviation+=
1064 channel_statistics[GreenChannel].variance-
1065 channel_statistics[GreenChannel].mean*
1066 channel_statistics[GreenChannel].mean;
1069 if ((channel & BlueChannel) != 0)
1071 channel_statistics[AllChannels].mean+=
1072 channel_statistics[BlueChannel].mean;
1073 channel_statistics[AllChannels].standard_deviation+=
1074 channel_statistics[BlueChannel].variance-
1075 channel_statistics[BlueChannel].mean*
1076 channel_statistics[BlueChannel].mean;
1079 if (((channel & OpacityChannel) != 0) &&
1080 (image->matte != MagickFalse))
1082 channel_statistics[AllChannels].mean+=
1083 channel_statistics[OpacityChannel].mean;
1084 channel_statistics[AllChannels].standard_deviation+=
1085 channel_statistics[OpacityChannel].variance-
1086 channel_statistics[OpacityChannel].mean*
1087 channel_statistics[OpacityChannel].mean;
1090 if (((channel & IndexChannel) != 0) &&
1091 (image->colorspace == CMYKColorspace))
1093 channel_statistics[AllChannels].mean+=
1094 channel_statistics[BlackChannel].mean;
1095 channel_statistics[AllChannels].standard_deviation+=
1096 channel_statistics[BlackChannel].variance-
1097 channel_statistics[BlackChannel].mean*
1098 channel_statistics[BlackChannel].mean;
1101 channel_statistics[AllChannels].mean/=channels;
1102 channel_statistics[AllChannels].standard_deviation=
1103 sqrt(channel_statistics[AllChannels].standard_deviation/channels);
1104 *mean=channel_statistics[AllChannels].mean;
1105 *standard_deviation=channel_statistics[AllChannels].standard_deviation;
1106 channel_statistics=(ChannelStatistics *) RelinquishMagickMemory(
1107 channel_statistics);
1112 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1116 % G e t I m a g e C h a n n e l K u r t o s i s %
1120 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1122 % GetImageChannelKurtosis() returns the kurtosis and skewness of one or more
1125 % The format of the GetImageChannelKurtosis method is:
1127 % MagickBooleanType GetImageChannelKurtosis(const Image *image,
1128 % const ChannelType channel,double *kurtosis,double *skewness,
1129 % ExceptionInfo *exception)
1131 % A description of each parameter follows:
1133 % o image: the image.
1135 % o channel: the channel.
1137 % o kurtosis: the kurtosis of the channel.
1139 % o skewness: the skewness of the channel.
1141 % o exception: return any errors or warnings in this structure.
1145 MagickExport MagickBooleanType GetImageKurtosis(const Image *image,
1146 double *kurtosis,double *skewness,ExceptionInfo *exception)
1151 status=GetImageChannelKurtosis(image,AllChannels,kurtosis,skewness,
1156 MagickExport MagickBooleanType GetImageChannelKurtosis(const Image *image,
1157 const ChannelType channel,double *kurtosis,double *skewness,
1158 ExceptionInfo *exception)
1171 assert(image != (Image *) NULL);
1172 assert(image->signature == MagickSignature);
1173 if (image->debug != MagickFalse)
1174 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1179 standard_deviation=0.0;
1182 sum_fourth_power=0.0;
1183 for (y=0; y < (ssize_t) image->rows; y++)
1185 register const IndexPacket
1188 register const PixelPacket
1194 p=GetVirtualPixels(image,0,y,image->columns,1,exception);
1195 if (p == (const PixelPacket *) NULL)
1197 indexes=GetVirtualIndexQueue(image);
1198 for (x=0; x < (ssize_t) image->columns; x++)
1200 if ((channel & RedChannel) != 0)
1202 mean+=GetRedPixelComponent(p);
1203 sum_squares+=(double) p->red*GetRedPixelComponent(p);
1204 sum_cubes+=(double) p->red*p->red*GetRedPixelComponent(p);
1205 sum_fourth_power+=(double) p->red*p->red*p->red*
1206 GetRedPixelComponent(p);
1209 if ((channel & GreenChannel) != 0)
1211 mean+=GetGreenPixelComponent(p);
1212 sum_squares+=(double) p->green*GetGreenPixelComponent(p);
1213 sum_cubes+=(double) p->green*p->green*GetGreenPixelComponent(p);
1214 sum_fourth_power+=(double) p->green*p->green*p->green*
1215 GetGreenPixelComponent(p);
1218 if ((channel & BlueChannel) != 0)
1220 mean+=GetBluePixelComponent(p);
1221 sum_squares+=(double) p->blue*GetBluePixelComponent(p);
1222 sum_cubes+=(double) p->blue*p->blue*GetBluePixelComponent(p);
1223 sum_fourth_power+=(double) p->blue*p->blue*p->blue*
1224 GetBluePixelComponent(p);
1227 if ((channel & OpacityChannel) != 0)
1229 mean+=GetOpacityPixelComponent(p);
1230 sum_squares+=(double) p->opacity*GetOpacityPixelComponent(p);
1231 sum_cubes+=(double) p->opacity*p->opacity*GetOpacityPixelComponent(p);
1232 sum_fourth_power+=(double) p->opacity*p->opacity*p->opacity*
1233 GetOpacityPixelComponent(p);
1236 if (((channel & IndexChannel) != 0) &&
1237 (image->colorspace == CMYKColorspace))
1240 sum_squares+=(double) indexes[x]*indexes[x];
1241 sum_cubes+=(double) indexes[x]*indexes[x]*indexes[x];
1242 sum_fourth_power+=(double) indexes[x]*indexes[x]*indexes[x]*
1249 if (y < (ssize_t) image->rows)
1250 return(MagickFalse);
1256 sum_fourth_power/=area;
1258 standard_deviation=sqrt(sum_squares-(mean*mean));
1259 if (standard_deviation != 0.0)
1261 *kurtosis=sum_fourth_power-4.0*mean*sum_cubes+6.0*mean*mean*sum_squares-
1262 3.0*mean*mean*mean*mean;
1263 *kurtosis/=standard_deviation*standard_deviation*standard_deviation*
1266 *skewness=sum_cubes-3.0*mean*sum_squares+2.0*mean*mean*mean;
1267 *skewness/=standard_deviation*standard_deviation*standard_deviation;
1269 return(y == (ssize_t) image->rows ? MagickTrue : MagickFalse);
1273 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1277 % G e t I m a g e C h a n n e l R a n g e %
1281 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1283 % GetImageChannelRange() returns the range of one or more image channels.
1285 % The format of the GetImageChannelRange method is:
1287 % MagickBooleanType GetImageChannelRange(const Image *image,
1288 % const ChannelType channel,double *minima,double *maxima,
1289 % ExceptionInfo *exception)
1291 % A description of each parameter follows:
1293 % o image: the image.
1295 % o channel: the channel.
1297 % o minima: the minimum value in the channel.
1299 % o maxima: the maximum value in the channel.
1301 % o exception: return any errors or warnings in this structure.
1305 MagickExport MagickBooleanType GetImageRange(const Image *image,
1306 double *minima,double *maxima,ExceptionInfo *exception)
1308 return(GetImageChannelRange(image,AllChannels,minima,maxima,exception));
1311 MagickExport MagickBooleanType GetImageChannelRange(const Image *image,
1312 const ChannelType channel,double *minima,double *maxima,
1313 ExceptionInfo *exception)
1321 assert(image != (Image *) NULL);
1322 assert(image->signature == MagickSignature);
1323 if (image->debug != MagickFalse)
1324 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1327 GetMagickPixelPacket(image,&pixel);
1328 for (y=0; y < (ssize_t) image->rows; y++)
1330 register const IndexPacket
1333 register const PixelPacket
1339 p=GetVirtualPixels(image,0,y,image->columns,1,exception);
1340 if (p == (const PixelPacket *) NULL)
1342 indexes=GetVirtualIndexQueue(image);
1343 for (x=0; x < (ssize_t) image->columns; x++)
1345 SetMagickPixelPacket(image,p,indexes+x,&pixel);
1346 if ((channel & RedChannel) != 0)
1348 if (pixel.red < *minima)
1349 *minima=(double) pixel.red;
1350 if (pixel.red > *maxima)
1351 *maxima=(double) pixel.red;
1353 if ((channel & GreenChannel) != 0)
1355 if (pixel.green < *minima)
1356 *minima=(double) pixel.green;
1357 if (pixel.green > *maxima)
1358 *maxima=(double) pixel.green;
1360 if ((channel & BlueChannel) != 0)
1362 if (pixel.blue < *minima)
1363 *minima=(double) pixel.blue;
1364 if (pixel.blue > *maxima)
1365 *maxima=(double) pixel.blue;
1367 if ((channel & OpacityChannel) != 0)
1369 if (pixel.opacity < *minima)
1370 *minima=(double) pixel.opacity;
1371 if (pixel.opacity > *maxima)
1372 *maxima=(double) pixel.opacity;
1374 if (((channel & IndexChannel) != 0) &&
1375 (image->colorspace == CMYKColorspace))
1377 if ((double) indexes[x] < *minima)
1378 *minima=(double) indexes[x];
1379 if ((double) indexes[x] > *maxima)
1380 *maxima=(double) indexes[x];
1385 return(y == (ssize_t) image->rows ? MagickTrue : MagickFalse);
1389 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1393 % 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 %
1397 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1399 % GetImageChannelStatistics() returns statistics for each channel in the
1400 % image. The statistics include the channel depth, its minima, maxima, mean,
1401 % standard deviation, kurtosis and skewness. You can access the red channel
1402 % mean, for example, like this:
1404 % channel_statistics=GetImageChannelStatistics(image,excepton);
1405 % red_mean=channel_statistics[RedChannel].mean;
1407 % Use MagickRelinquishMemory() to free the statistics buffer.
1409 % The format of the GetImageChannelStatistics method is:
1411 % ChannelStatistics *GetImageChannelStatistics(const Image *image,
1412 % ExceptionInfo *exception)
1414 % A description of each parameter follows:
1416 % o image: the image.
1418 % o exception: return any errors or warnings in this structure.
1421 MagickExport ChannelStatistics *GetImageChannelStatistics(const Image *image,
1422 ExceptionInfo *exception)
1425 *channel_statistics;
1449 assert(image != (Image *) NULL);
1450 assert(image->signature == MagickSignature);
1451 if (image->debug != MagickFalse)
1452 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1453 length=AllChannels+1UL;
1454 channel_statistics=(ChannelStatistics *) AcquireQuantumMemory(length,
1455 sizeof(*channel_statistics));
1456 if (channel_statistics == (ChannelStatistics *) NULL)
1457 ThrowFatalException(ResourceLimitFatalError,"MemoryAllocationFailed");
1458 (void) ResetMagickMemory(channel_statistics,0,length*
1459 sizeof(*channel_statistics));
1460 for (i=0; i <= AllChannels; i++)
1462 channel_statistics[i].depth=1;
1463 channel_statistics[i].maxima=(-1.0E-37);
1464 channel_statistics[i].minima=1.0E+37;
1466 for (y=0; y < (ssize_t) image->rows; y++)
1468 register const IndexPacket
1471 register const PixelPacket
1477 p=GetVirtualPixels(image,0,y,image->columns,1,exception);
1478 if (p == (const PixelPacket *) NULL)
1480 indexes=GetVirtualIndexQueue(image);
1481 for (x=0; x < (ssize_t) image->columns; )
1483 if (channel_statistics[RedChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1485 depth=channel_statistics[RedChannel].depth;
1486 range=GetQuantumRange(depth);
1487 status=p->red != ScaleAnyToQuantum(ScaleQuantumToAny(p->red,range),
1488 range) ? MagickTrue : MagickFalse;
1489 if (status != MagickFalse)
1491 channel_statistics[RedChannel].depth++;
1495 if (channel_statistics[GreenChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1497 depth=channel_statistics[GreenChannel].depth;
1498 range=GetQuantumRange(depth);
1499 status=p->green != ScaleAnyToQuantum(ScaleQuantumToAny(p->green,
1500 range),range) ? MagickTrue : MagickFalse;
1501 if (status != MagickFalse)
1503 channel_statistics[GreenChannel].depth++;
1507 if (channel_statistics[BlueChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1509 depth=channel_statistics[BlueChannel].depth;
1510 range=GetQuantumRange(depth);
1511 status=p->blue != ScaleAnyToQuantum(ScaleQuantumToAny(p->blue,
1512 range),range) ? MagickTrue : MagickFalse;
1513 if (status != MagickFalse)
1515 channel_statistics[BlueChannel].depth++;
1519 if (image->matte != MagickFalse)
1521 if (channel_statistics[OpacityChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1523 depth=channel_statistics[OpacityChannel].depth;
1524 range=GetQuantumRange(depth);
1525 status=p->opacity != ScaleAnyToQuantum(ScaleQuantumToAny(
1526 p->opacity,range),range) ? MagickTrue : MagickFalse;
1527 if (status != MagickFalse)
1529 channel_statistics[OpacityChannel].depth++;
1534 if (image->colorspace == CMYKColorspace)
1536 if (channel_statistics[BlackChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1538 depth=channel_statistics[BlackChannel].depth;
1539 range=GetQuantumRange(depth);
1540 status=indexes[x] != ScaleAnyToQuantum(ScaleQuantumToAny(
1541 indexes[x],range),range) ? MagickTrue : MagickFalse;
1542 if (status != MagickFalse)
1544 channel_statistics[BlackChannel].depth++;
1549 if ((double) p->red < channel_statistics[RedChannel].minima)
1550 channel_statistics[RedChannel].minima=(double) GetRedPixelComponent(p);
1551 if ((double) p->red > channel_statistics[RedChannel].maxima)
1552 channel_statistics[RedChannel].maxima=(double) GetRedPixelComponent(p);
1553 channel_statistics[RedChannel].sum+=GetRedPixelComponent(p);
1554 channel_statistics[RedChannel].sum_squared+=(double) p->red*
1555 GetRedPixelComponent(p);
1556 channel_statistics[RedChannel].kurtosis+=(double) p->red*p->red*
1557 p->red*GetRedPixelComponent(p);
1558 channel_statistics[RedChannel].skewness+=(double) p->red*p->red*
1559 GetRedPixelComponent(p);
1560 if ((double) p->green < channel_statistics[GreenChannel].minima)
1561 channel_statistics[GreenChannel].minima=(double)
1562 GetGreenPixelComponent(p);
1563 if ((double) p->green > channel_statistics[GreenChannel].maxima)
1564 channel_statistics[GreenChannel].maxima=(double)
1565 GetGreenPixelComponent(p);
1566 channel_statistics[GreenChannel].sum+=GetGreenPixelComponent(p);
1567 channel_statistics[GreenChannel].sum_squared+=(double) p->green*
1568 GetGreenPixelComponent(p);
1569 channel_statistics[GreenChannel].kurtosis+=(double) p->green*p->green*
1570 p->green*GetGreenPixelComponent(p);
1571 channel_statistics[GreenChannel].skewness+=(double) p->green*p->green*
1572 GetGreenPixelComponent(p);
1573 if ((double) p->blue < channel_statistics[BlueChannel].minima)
1574 channel_statistics[BlueChannel].minima=(double)
1575 GetBluePixelComponent(p);
1576 if ((double) p->blue > channel_statistics[BlueChannel].maxima)
1577 channel_statistics[BlueChannel].maxima=(double)
1578 GetBluePixelComponent(p);
1579 channel_statistics[BlueChannel].sum+=GetBluePixelComponent(p);
1580 channel_statistics[BlueChannel].sum_squared+=(double) p->blue*
1581 GetBluePixelComponent(p);
1582 channel_statistics[BlueChannel].kurtosis+=(double) p->blue*p->blue*
1583 p->blue*GetBluePixelComponent(p);
1584 channel_statistics[BlueChannel].skewness+=(double) p->blue*p->blue*
1585 GetBluePixelComponent(p);
1586 if (image->matte != MagickFalse)
1588 if ((double) p->opacity < channel_statistics[OpacityChannel].minima)
1589 channel_statistics[OpacityChannel].minima=(double)
1590 GetOpacityPixelComponent(p);
1591 if ((double) p->opacity > channel_statistics[OpacityChannel].maxima)
1592 channel_statistics[OpacityChannel].maxima=(double)
1593 GetOpacityPixelComponent(p);
1594 channel_statistics[OpacityChannel].sum+=GetOpacityPixelComponent(p);
1595 channel_statistics[OpacityChannel].sum_squared+=(double)
1596 p->opacity*GetOpacityPixelComponent(p);
1597 channel_statistics[OpacityChannel].kurtosis+=(double) p->opacity*
1598 p->opacity*p->opacity*GetOpacityPixelComponent(p);
1599 channel_statistics[OpacityChannel].skewness+=(double) p->opacity*
1600 p->opacity*GetOpacityPixelComponent(p);
1602 if (image->colorspace == CMYKColorspace)
1604 if ((double) indexes[x] < channel_statistics[BlackChannel].minima)
1605 channel_statistics[BlackChannel].minima=(double) indexes[x];
1606 if ((double) indexes[x] > channel_statistics[BlackChannel].maxima)
1607 channel_statistics[BlackChannel].maxima=(double) indexes[x];
1608 channel_statistics[BlackChannel].sum+=indexes[x];
1609 channel_statistics[BlackChannel].sum_squared+=(double)
1610 indexes[x]*indexes[x];
1611 channel_statistics[BlackChannel].kurtosis+=(double) indexes[x]*
1612 indexes[x]*indexes[x]*indexes[x];
1613 channel_statistics[BlackChannel].skewness+=(double) indexes[x]*
1614 indexes[x]*indexes[x];
1620 area=(double) image->columns*image->rows;
1621 for (i=0; i < AllChannels; i++)
1623 channel_statistics[i].sum/=area;
1624 channel_statistics[i].sum_squared/=area;
1625 channel_statistics[i].mean=channel_statistics[i].sum;
1626 channel_statistics[i].variance=channel_statistics[i].sum_squared;
1627 channel_statistics[i].standard_deviation=sqrt(
1628 channel_statistics[i].variance-(channel_statistics[i].mean*
1629 channel_statistics[i].mean));
1630 channel_statistics[i].kurtosis/=area;
1631 channel_statistics[i].skewness/=area;
1633 for (i=0; i < AllChannels; i++)
1635 channel_statistics[AllChannels].depth=(size_t) MagickMax((double)
1636 channel_statistics[AllChannels].depth,(double)
1637 channel_statistics[i].depth);
1638 channel_statistics[AllChannels].minima=MagickMin(
1639 channel_statistics[AllChannels].minima,channel_statistics[i].minima);
1640 channel_statistics[AllChannels].maxima=MagickMax(
1641 channel_statistics[AllChannels].maxima,channel_statistics[i].maxima);
1642 channel_statistics[AllChannels].sum+=channel_statistics[i].sum;
1643 channel_statistics[AllChannels].sum_squared+=
1644 channel_statistics[i].sum_squared;
1645 channel_statistics[AllChannels].mean+=channel_statistics[i].mean;
1646 channel_statistics[AllChannels].variance+=channel_statistics[i].variance-
1647 channel_statistics[i].mean*channel_statistics[i].mean;
1648 channel_statistics[AllChannels].standard_deviation+=
1649 channel_statistics[i].variance-channel_statistics[i].mean*
1650 channel_statistics[i].mean;
1651 channel_statistics[AllChannels].kurtosis+=channel_statistics[i].kurtosis;
1652 channel_statistics[AllChannels].skewness+=channel_statistics[i].skewness;
1655 if (image->matte != MagickFalse)
1657 if (image->colorspace == CMYKColorspace)
1659 channel_statistics[AllChannels].sum/=channels;
1660 channel_statistics[AllChannels].sum_squared/=channels;
1661 channel_statistics[AllChannels].mean/=channels;
1662 channel_statistics[AllChannels].variance/=channels;
1663 channel_statistics[AllChannels].standard_deviation=
1664 sqrt(channel_statistics[AllChannels].standard_deviation/channels);
1665 channel_statistics[AllChannels].kurtosis/=channels;
1666 channel_statistics[AllChannels].skewness/=channels;
1667 for (i=0; i <= AllChannels; i++)
1669 if (channel_statistics[i].standard_deviation == 0.0)
1671 channel_statistics[i].kurtosis=0.0;
1672 channel_statistics[i].skewness=0.0;
1676 channel_statistics[i].skewness=(channel_statistics[i].skewness-
1677 3.0*channel_statistics[i].mean*channel_statistics[i].sum_squared+
1678 2.0*channel_statistics[i].mean*channel_statistics[i].mean*
1679 channel_statistics[i].mean)/(channel_statistics[i].standard_deviation*
1680 channel_statistics[i].standard_deviation*
1681 channel_statistics[i].standard_deviation);
1682 channel_statistics[i].kurtosis=(channel_statistics[i].kurtosis-
1683 4.0*channel_statistics[i].mean*channel_statistics[i].skewness+
1684 6.0*channel_statistics[i].mean*channel_statistics[i].mean*
1685 channel_statistics[i].sum_squared-3.0*channel_statistics[i].mean*
1686 channel_statistics[i].mean*1.0*channel_statistics[i].mean*
1687 channel_statistics[i].mean)/(channel_statistics[i].standard_deviation*
1688 channel_statistics[i].standard_deviation*
1689 channel_statistics[i].standard_deviation*
1690 channel_statistics[i].standard_deviation)-3.0;
1693 return(channel_statistics);