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 < (long) 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 < (long) 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 < (long) 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)*floor(result/(QuantumRange+1));
220 case AndEvaluateOperator:
222 result=(MagickRealType) ((unsigned long) pixel & (unsigned long)
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) ((unsigned long) pixel << (unsigned long)
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) ((unsigned long) pixel | (unsigned long)
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) ((unsigned long) pixel >> (unsigned long)
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) ((unsigned long) pixel ^ (unsigned long)
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"
398 **restrict evaluate_pixels,
402 **restrict random_info;
408 Ensure the image are the same size.
410 assert(images != (Image *) NULL);
411 assert(images->signature == MagickSignature);
412 if (images->debug != MagickFalse)
413 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",images->filename);
414 assert(exception != (ExceptionInfo *) NULL);
415 assert(exception->signature == MagickSignature);
416 for (next=images; next != (Image *) NULL; next=GetNextImageInList(next))
417 if ((next->columns != images->columns) || (next->rows != images->rows))
419 (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
420 "ImageWidthsOrHeightsDiffer","`%s'",images->filename);
421 return((Image *) NULL);
424 Initialize evaluate next attributes.
426 evaluate_image=CloneImage(images,images->columns,images->rows,MagickTrue,
428 if (evaluate_image == (Image *) NULL)
429 return((Image *) NULL);
430 if (SetImageStorageClass(evaluate_image,DirectClass) == MagickFalse)
432 InheritException(exception,&evaluate_image->exception);
433 evaluate_image=DestroyImage(evaluate_image);
434 return((Image *) NULL);
436 evaluate_pixels=AcquirePixelThreadSet(images);
437 if (evaluate_pixels == (MagickPixelPacket **) NULL)
439 evaluate_image=DestroyImage(evaluate_image);
440 (void) ThrowMagickException(exception,GetMagickModule(),
441 ResourceLimitError,"MemoryAllocationFailed","`%s'",images->filename);
442 return((Image *) NULL);
445 Evaluate image pixels.
449 GetMagickPixelPacket(images,&zero);
450 random_info=AcquireRandomInfoThreadSet();
451 number_images=GetImageListLength(images);
452 evaluate_view=AcquireCacheView(evaluate_image);
453 #if defined(MAGICKCORE_OPENMP_SUPPORT)
454 #pragma omp parallel for schedule(dynamic) shared(progress,status)
456 for (y=0; y < (long) evaluate_image->rows; y++)
468 *restrict evaluate_indexes;
475 register MagickPixelPacket
481 if (status == MagickFalse)
483 q=QueueCacheViewAuthenticPixels(evaluate_view,0,y,evaluate_image->columns,1,
485 if (q == (PixelPacket *) NULL)
490 evaluate_indexes=GetCacheViewAuthenticIndexQueue(evaluate_view);
492 id=GetOpenMPThreadId();
493 evaluate_pixel=evaluate_pixels[id];
494 for (x=0; x < (long) evaluate_image->columns; x++)
495 evaluate_pixel[x]=zero;
497 for (i=0; i < (long) number_images; i++)
499 register const IndexPacket
502 register const PixelPacket
505 image_view=AcquireCacheView(next);
506 p=GetCacheViewVirtualPixels(image_view,0,y,next->columns,1,exception);
507 if (p == (const PixelPacket *) NULL)
509 image_view=DestroyCacheView(image_view);
512 indexes=GetCacheViewVirtualIndexQueue(image_view);
513 for (x=0; x < (long) next->columns; x++)
515 evaluate_pixel[x].red=ApplyEvaluateOperator(random_info[id],p->red,
516 i == 0 ? AddEvaluateOperator : op,evaluate_pixel[x].red);
517 evaluate_pixel[x].green=ApplyEvaluateOperator(random_info[id],p->green,
518 i == 0 ? AddEvaluateOperator : op,evaluate_pixel[x].green);
519 evaluate_pixel[x].blue=ApplyEvaluateOperator(random_info[id],p->blue,
520 i == 0 ? AddEvaluateOperator : op,evaluate_pixel[x].blue);
521 evaluate_pixel[x].opacity=ApplyEvaluateOperator(random_info[id],
522 p->opacity,i == 0 ? AddEvaluateOperator : op,
523 evaluate_pixel[x].opacity);
524 if (evaluate_image->colorspace == CMYKColorspace)
525 evaluate_pixel[x].index=ApplyEvaluateOperator(random_info[id],
526 indexes[x],i == 0 ? AddEvaluateOperator : op,
527 evaluate_pixel[x].index);
530 image_view=DestroyCacheView(image_view);
531 next=GetNextImageInList(next);
533 if (op == MeanEvaluateOperator)
534 for (x=0; x < (long) evaluate_image->columns; x++)
536 evaluate_pixel[x].red/=number_images;
537 evaluate_pixel[x].green/=number_images;
538 evaluate_pixel[x].blue/=number_images;
539 evaluate_pixel[x].opacity/=number_images;
540 evaluate_pixel[x].index/=number_images;
542 for (x=0; x < (long) evaluate_image->columns; x++)
544 q->red=ClampToQuantum(evaluate_pixel[x].red);
545 q->green=ClampToQuantum(evaluate_pixel[x].green);
546 q->blue=ClampToQuantum(evaluate_pixel[x].blue);
547 if (evaluate_image->matte == MagickFalse)
548 q->opacity=ClampToQuantum(evaluate_pixel[x].opacity);
550 q->opacity=ClampToQuantum(QuantumRange-evaluate_pixel[x].opacity);
551 if (evaluate_image->colorspace == CMYKColorspace)
552 evaluate_indexes[x]=ClampToQuantum(evaluate_pixel[x].index);
555 if (SyncCacheViewAuthenticPixels(evaluate_view,exception) == MagickFalse)
557 if (images->progress_monitor != (MagickProgressMonitor) NULL)
562 #if defined(MAGICKCORE_OPENMP_SUPPORT)
563 #pragma omp critical (MagickCore_EvaluateImages)
565 proceed=SetImageProgress(images,EvaluateImageTag,progress++,
566 evaluate_image->rows);
567 if (proceed == MagickFalse)
571 evaluate_view=DestroyCacheView(evaluate_view);
572 evaluate_pixels=DestroyPixelThreadSet(evaluate_pixels);
573 random_info=DestroyRandomInfoThreadSet(random_info);
574 if (status == MagickFalse)
575 evaluate_image=DestroyImage(evaluate_image);
576 return(evaluate_image);
579 MagickExport MagickBooleanType EvaluateImageChannel(Image *image,
580 const ChannelType channel,const MagickEvaluateOperator op,const double value,
581 ExceptionInfo *exception)
594 **restrict random_info;
596 assert(image != (Image *) NULL);
597 assert(image->signature == MagickSignature);
598 if (image->debug != MagickFalse)
599 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
600 assert(exception != (ExceptionInfo *) NULL);
601 assert(exception->signature == MagickSignature);
602 if (SetImageStorageClass(image,DirectClass) == MagickFalse)
604 InheritException(exception,&image->exception);
609 random_info=AcquireRandomInfoThreadSet();
610 image_view=AcquireCacheView(image);
611 #if defined(MAGICKCORE_OPENMP_SUPPORT)
612 #pragma omp parallel for schedule(dynamic,4) shared(progress,status)
614 for (y=0; y < (long) image->rows; y++)
626 if (status == MagickFalse)
628 q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
629 if (q == (PixelPacket *) NULL)
634 indexes=GetCacheViewAuthenticIndexQueue(image_view);
635 id=GetOpenMPThreadId();
636 for (x=0; x < (long) image->columns; x++)
638 if ((channel & RedChannel) != 0)
639 q->red=ClampToQuantum(ApplyEvaluateOperator(random_info[id],q->red,op,
641 if ((channel & GreenChannel) != 0)
642 q->green=ClampToQuantum(ApplyEvaluateOperator(random_info[id],q->green,
644 if ((channel & BlueChannel) != 0)
645 q->blue=ClampToQuantum(ApplyEvaluateOperator(random_info[id],q->blue,op,
647 if ((channel & OpacityChannel) != 0)
649 if (image->matte == MagickFalse)
650 q->opacity=ClampToQuantum(ApplyEvaluateOperator(random_info[id],
651 q->opacity,op,value));
653 q->opacity=ClampToQuantum(QuantumRange-ApplyEvaluateOperator(
654 random_info[id],(Quantum) GetAlphaPixelComponent(q),op,value));
656 if (((channel & IndexChannel) != 0) && (indexes != (IndexPacket *) NULL))
657 indexes[x]=(IndexPacket) ClampToQuantum(ApplyEvaluateOperator(
658 random_info[id],indexes[x],op,value));
661 if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
663 if (image->progress_monitor != (MagickProgressMonitor) NULL)
668 #if defined(MAGICKCORE_OPENMP_SUPPORT)
669 #pragma omp critical (MagickCore_EvaluateImageChannel)
671 proceed=SetImageProgress(image,EvaluateImageTag,progress++,image->rows);
672 if (proceed == MagickFalse)
676 image_view=DestroyCacheView(image_view);
677 random_info=DestroyRandomInfoThreadSet(random_info);
682 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
686 % F u n c t i o n I m a g e %
690 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
692 % FunctionImage() applies a value to the image with an arithmetic, relational,
693 % or logical operator to an image. Use these operations to lighten or darken
694 % an image, to increase or decrease contrast in an image, or to produce the
695 % "negative" of an image.
697 % The format of the FunctionImageChannel method is:
699 % MagickBooleanType FunctionImage(Image *image,
700 % const MagickFunction function,const long number_parameters,
701 % const double *parameters,ExceptionInfo *exception)
702 % MagickBooleanType FunctionImageChannel(Image *image,
703 % const ChannelType channel,const MagickFunction function,
704 % const long number_parameters,const double *argument,
705 % ExceptionInfo *exception)
707 % A description of each parameter follows:
709 % o image: the image.
711 % o channel: the channel.
713 % o function: A channel function.
715 % o parameters: one or more parameters.
717 % o exception: return any errors or warnings in this structure.
721 static Quantum ApplyFunction(Quantum pixel,const MagickFunction function,
722 const unsigned long number_parameters,const double *parameters,
723 ExceptionInfo *exception)
735 case PolynomialFunction:
739 * Parameters: polynomial constants, highest to lowest order
740 * For example: c0*x^3 + c1*x^2 + c2*x + c3
743 for (i=0; i < (long) number_parameters; i++)
744 result = result*QuantumScale*pixel + parameters[i];
745 result *= QuantumRange;
748 case SinusoidFunction:
751 * Parameters: Freq, Phase, Ampl, bias
753 double freq,phase,ampl,bias;
754 freq = ( number_parameters >= 1 ) ? parameters[0] : 1.0;
755 phase = ( number_parameters >= 2 ) ? parameters[1] : 0.0;
756 ampl = ( number_parameters >= 3 ) ? parameters[2] : 0.5;
757 bias = ( number_parameters >= 4 ) ? parameters[3] : 0.5;
758 result=(MagickRealType) (QuantumRange*(ampl*sin((double) (2.0*MagickPI*
759 (freq*QuantumScale*pixel + phase/360.0) )) + bias ) );
764 /* Arcsin Function (peged at range limits for invalid results)
765 * Parameters: Width, Center, Range, Bias
767 double width,range,center,bias;
768 width = ( number_parameters >= 1 ) ? parameters[0] : 1.0;
769 center = ( number_parameters >= 2 ) ? parameters[1] : 0.5;
770 range = ( number_parameters >= 3 ) ? parameters[2] : 1.0;
771 bias = ( number_parameters >= 4 ) ? parameters[3] : 0.5;
772 result = 2.0/width*(QuantumScale*pixel - center);
773 if ( result <= -1.0 )
774 result = bias - range/2.0;
775 else if ( result >= 1.0 )
776 result = bias + range/2.0;
778 result=range/MagickPI*asin((double)result) + bias;
779 result *= QuantumRange;
785 * Parameters: Slope, Center, Range, Bias
787 double slope,range,center,bias;
788 slope = ( number_parameters >= 1 ) ? parameters[0] : 1.0;
789 center = ( number_parameters >= 2 ) ? parameters[1] : 0.5;
790 range = ( number_parameters >= 3 ) ? parameters[2] : 1.0;
791 bias = ( number_parameters >= 4 ) ? parameters[3] : 0.5;
792 result = MagickPI*slope*(QuantumScale*pixel - center);
793 result=(MagickRealType) (QuantumRange*(range/MagickPI*atan((double)
797 case UndefinedFunction:
800 return(ClampToQuantum(result));
803 MagickExport MagickBooleanType FunctionImage(Image *image,
804 const MagickFunction function,const unsigned long number_parameters,
805 const double *parameters,ExceptionInfo *exception)
810 status=FunctionImageChannel(image,AllChannels,function,number_parameters,
811 parameters,exception);
815 MagickExport MagickBooleanType FunctionImageChannel(Image *image,
816 const ChannelType channel,const MagickFunction function,
817 const unsigned long number_parameters,const double *parameters,
818 ExceptionInfo *exception)
820 #define FunctionImageTag "Function/Image "
832 assert(image != (Image *) NULL);
833 assert(image->signature == MagickSignature);
834 if (image->debug != MagickFalse)
835 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
836 assert(exception != (ExceptionInfo *) NULL);
837 assert(exception->signature == MagickSignature);
838 if (SetImageStorageClass(image,DirectClass) == MagickFalse)
840 InheritException(exception,&image->exception);
845 image_view=AcquireCacheView(image);
846 #if defined(MAGICKCORE_OPENMP_SUPPORT)
847 #pragma omp parallel for schedule(dynamic,4) shared(progress,status)
849 for (y=0; y < (long) image->rows; y++)
860 if (status == MagickFalse)
862 q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
863 if (q == (PixelPacket *) NULL)
868 indexes=GetCacheViewAuthenticIndexQueue(image_view);
869 for (x=0; x < (long) image->columns; x++)
871 if ((channel & RedChannel) != 0)
872 q->red=ApplyFunction(q->red,function,number_parameters,parameters,
874 if ((channel & GreenChannel) != 0)
875 q->green=ApplyFunction(q->green,function,number_parameters,parameters,
877 if ((channel & BlueChannel) != 0)
878 q->blue=ApplyFunction(q->blue,function,number_parameters,parameters,
880 if ((channel & OpacityChannel) != 0)
882 if (image->matte == MagickFalse)
883 q->opacity=ApplyFunction(q->opacity,function,number_parameters,
884 parameters,exception);
886 q->opacity=(Quantum) QuantumRange-ApplyFunction((Quantum)
887 GetAlphaPixelComponent(q),function,number_parameters,parameters,
890 if (((channel & IndexChannel) != 0) && (indexes != (IndexPacket *) NULL))
891 indexes[x]=(IndexPacket) ApplyFunction(indexes[x],function,
892 number_parameters,parameters,exception);
895 if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
897 if (image->progress_monitor != (MagickProgressMonitor) NULL)
902 #if defined(MAGICKCORE_OPENMP_SUPPORT)
903 #pragma omp critical (MagickCore_FunctionImageChannel)
905 proceed=SetImageProgress(image,FunctionImageTag,progress++,image->rows);
906 if (proceed == MagickFalse)
910 image_view=DestroyCacheView(image_view);
915 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
919 + G e t I m a g e C h a n n e l E x t r e m a %
923 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
925 % GetImageChannelExtrema() returns the extrema of one or more image channels.
927 % The format of the GetImageChannelExtrema method is:
929 % MagickBooleanType GetImageChannelExtrema(const Image *image,
930 % const ChannelType channel,unsigned long *minima,unsigned long *maxima,
931 % ExceptionInfo *exception)
933 % A description of each parameter follows:
935 % o image: the image.
937 % o channel: the channel.
939 % o minima: the minimum value in the channel.
941 % o maxima: the maximum value in the channel.
943 % o exception: return any errors or warnings in this structure.
947 MagickExport MagickBooleanType GetImageExtrema(const Image *image,
948 unsigned long *minima,unsigned long *maxima,ExceptionInfo *exception)
950 return(GetImageChannelExtrema(image,AllChannels,minima,maxima,exception));
953 MagickExport MagickBooleanType GetImageChannelExtrema(const Image *image,
954 const ChannelType channel,unsigned long *minima,unsigned long *maxima,
955 ExceptionInfo *exception)
964 assert(image != (Image *) NULL);
965 assert(image->signature == MagickSignature);
966 if (image->debug != MagickFalse)
967 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
968 status=GetImageChannelRange(image,channel,&min,&max,exception);
969 *minima=(unsigned long) ceil(min-0.5);
970 *maxima=(unsigned long) floor(max+0.5);
975 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
979 % G e t I m a g e C h a n n e l M e a n %
983 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
985 % GetImageChannelMean() returns the mean and standard deviation of one or more
988 % The format of the GetImageChannelMean method is:
990 % MagickBooleanType GetImageChannelMean(const Image *image,
991 % const ChannelType channel,double *mean,double *standard_deviation,
992 % ExceptionInfo *exception)
994 % A description of each parameter follows:
996 % o image: the image.
998 % o channel: the channel.
1000 % o mean: the average value in the channel.
1002 % o standard_deviation: the standard deviation of the channel.
1004 % o exception: return any errors or warnings in this structure.
1008 MagickExport MagickBooleanType GetImageMean(const Image *image,double *mean,
1009 double *standard_deviation,ExceptionInfo *exception)
1014 status=GetImageChannelMean(image,AllChannels,mean,standard_deviation,
1019 MagickExport MagickBooleanType GetImageChannelMean(const Image *image,
1020 const ChannelType channel,double *mean,double *standard_deviation,
1021 ExceptionInfo *exception)
1029 assert(image != (Image *) NULL);
1030 assert(image->signature == MagickSignature);
1031 if (image->debug != MagickFalse)
1032 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1034 *standard_deviation=0.0;
1036 for (y=0; y < (long) image->rows; y++)
1038 register const IndexPacket
1041 register const PixelPacket
1047 p=GetVirtualPixels(image,0,y,image->columns,1,exception);
1048 if (p == (const PixelPacket *) NULL)
1050 indexes=GetVirtualIndexQueue(image);
1051 for (x=0; x < (long) image->columns; x++)
1053 if ((channel & RedChannel) != 0)
1055 *mean+=GetRedPixelComponent(p);
1056 *standard_deviation+=(double) p->red*GetRedPixelComponent(p);
1059 if ((channel & GreenChannel) != 0)
1061 *mean+=GetGreenPixelComponent(p);
1062 *standard_deviation+=(double) p->green*GetGreenPixelComponent(p);
1065 if ((channel & BlueChannel) != 0)
1067 *mean+=GetBluePixelComponent(p);
1068 *standard_deviation+=(double) p->blue*GetBluePixelComponent(p);
1071 if ((channel & OpacityChannel) != 0)
1073 *mean+=GetOpacityPixelComponent(p);
1074 *standard_deviation+=(double) p->opacity*GetOpacityPixelComponent(p);
1077 if (((channel & IndexChannel) != 0) &&
1078 (image->colorspace == CMYKColorspace))
1081 *standard_deviation+=(double) indexes[x]*indexes[x];
1087 if (y < (long) image->rows)
1088 return(MagickFalse);
1092 *standard_deviation/=area;
1094 *standard_deviation=sqrt(*standard_deviation-(*mean*(*mean)));
1095 return(y == (long) image->rows ? MagickTrue : MagickFalse);
1099 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1103 % G e t I m a g e C h a n n e l K u r t o s i s %
1107 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1109 % GetImageChannelKurtosis() returns the kurtosis and skewness of one or more
1112 % The format of the GetImageChannelKurtosis method is:
1114 % MagickBooleanType GetImageChannelKurtosis(const Image *image,
1115 % const ChannelType channel,double *kurtosis,double *skewness,
1116 % ExceptionInfo *exception)
1118 % A description of each parameter follows:
1120 % o image: the image.
1122 % o channel: the channel.
1124 % o kurtosis: the kurtosis of the channel.
1126 % o skewness: the skewness of the channel.
1128 % o exception: return any errors or warnings in this structure.
1132 MagickExport MagickBooleanType GetImageKurtosis(const Image *image,
1133 double *kurtosis,double *skewness,ExceptionInfo *exception)
1138 status=GetImageChannelKurtosis(image,AllChannels,kurtosis,skewness,
1143 MagickExport MagickBooleanType GetImageChannelKurtosis(const Image *image,
1144 const ChannelType channel,double *kurtosis,double *skewness,
1145 ExceptionInfo *exception)
1158 assert(image != (Image *) NULL);
1159 assert(image->signature == MagickSignature);
1160 if (image->debug != MagickFalse)
1161 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1166 standard_deviation=0.0;
1169 sum_fourth_power=0.0;
1170 for (y=0; y < (long) image->rows; y++)
1172 register const IndexPacket
1175 register const PixelPacket
1181 p=GetVirtualPixels(image,0,y,image->columns,1,exception);
1182 if (p == (const PixelPacket *) NULL)
1184 indexes=GetVirtualIndexQueue(image);
1185 for (x=0; x < (long) image->columns; x++)
1187 if ((channel & RedChannel) != 0)
1189 mean+=GetRedPixelComponent(p);
1190 sum_squares+=(double) p->red*GetRedPixelComponent(p);
1191 sum_cubes+=(double) p->red*p->red*GetRedPixelComponent(p);
1192 sum_fourth_power+=(double) p->red*p->red*p->red*
1193 GetRedPixelComponent(p);
1196 if ((channel & GreenChannel) != 0)
1198 mean+=GetGreenPixelComponent(p);
1199 sum_squares+=(double) p->green*GetGreenPixelComponent(p);
1200 sum_cubes+=(double) p->green*p->green*GetGreenPixelComponent(p);
1201 sum_fourth_power+=(double) p->green*p->green*p->green*
1202 GetGreenPixelComponent(p);
1205 if ((channel & BlueChannel) != 0)
1207 mean+=GetBluePixelComponent(p);
1208 sum_squares+=(double) p->blue*GetBluePixelComponent(p);
1209 sum_cubes+=(double) p->blue*p->blue*GetBluePixelComponent(p);
1210 sum_fourth_power+=(double) p->blue*p->blue*p->blue*
1211 GetBluePixelComponent(p);
1214 if ((channel & OpacityChannel) != 0)
1216 mean+=GetOpacityPixelComponent(p);
1217 sum_squares+=(double) p->opacity*GetOpacityPixelComponent(p);
1218 sum_cubes+=(double) p->opacity*p->opacity*GetOpacityPixelComponent(p);
1219 sum_fourth_power+=(double) p->opacity*p->opacity*p->opacity*
1220 GetOpacityPixelComponent(p);
1223 if (((channel & IndexChannel) != 0) &&
1224 (image->colorspace == CMYKColorspace))
1227 sum_squares+=(double) indexes[x]*indexes[x];
1228 sum_cubes+=(double) indexes[x]*indexes[x]*indexes[x];
1229 sum_fourth_power+=(double) indexes[x]*indexes[x]*indexes[x]*
1236 if (y < (long) image->rows)
1237 return(MagickFalse);
1243 sum_fourth_power/=area;
1245 standard_deviation=sqrt(sum_squares-(mean*mean));
1246 if (standard_deviation != 0.0)
1248 *kurtosis=sum_fourth_power-4.0*mean*sum_cubes+6.0*mean*mean*sum_squares-
1249 3.0*mean*mean*mean*mean;
1250 *kurtosis/=standard_deviation*standard_deviation*standard_deviation*
1253 *skewness=sum_cubes-3.0*mean*sum_squares+2.0*mean*mean*mean;
1254 *skewness/=standard_deviation*standard_deviation*standard_deviation;
1256 return(y == (long) image->rows ? MagickTrue : MagickFalse);
1260 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1264 % G e t I m a g e C h a n n e l R a n g e %
1268 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1270 % GetImageChannelRange() returns the range of one or more image channels.
1272 % The format of the GetImageChannelRange method is:
1274 % MagickBooleanType GetImageChannelRange(const Image *image,
1275 % const ChannelType channel,double *minima,double *maxima,
1276 % ExceptionInfo *exception)
1278 % A description of each parameter follows:
1280 % o image: the image.
1282 % o channel: the channel.
1284 % o minima: the minimum value in the channel.
1286 % o maxima: the maximum value in the channel.
1288 % o exception: return any errors or warnings in this structure.
1292 MagickExport MagickBooleanType GetImageRange(const Image *image,
1293 double *minima,double *maxima,ExceptionInfo *exception)
1295 return(GetImageChannelRange(image,AllChannels,minima,maxima,exception));
1298 MagickExport MagickBooleanType GetImageChannelRange(const Image *image,
1299 const ChannelType channel,double *minima,double *maxima,
1300 ExceptionInfo *exception)
1308 assert(image != (Image *) NULL);
1309 assert(image->signature == MagickSignature);
1310 if (image->debug != MagickFalse)
1311 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1314 GetMagickPixelPacket(image,&pixel);
1315 for (y=0; y < (long) image->rows; y++)
1317 register const IndexPacket
1320 register const PixelPacket
1326 p=GetVirtualPixels(image,0,y,image->columns,1,exception);
1327 if (p == (const PixelPacket *) NULL)
1329 indexes=GetVirtualIndexQueue(image);
1330 for (x=0; x < (long) image->columns; x++)
1332 SetMagickPixelPacket(image,p,indexes+x,&pixel);
1333 if ((channel & RedChannel) != 0)
1335 if (pixel.red < *minima)
1336 *minima=(double) pixel.red;
1337 if (pixel.red > *maxima)
1338 *maxima=(double) pixel.red;
1340 if ((channel & GreenChannel) != 0)
1342 if (pixel.green < *minima)
1343 *minima=(double) pixel.green;
1344 if (pixel.green > *maxima)
1345 *maxima=(double) pixel.green;
1347 if ((channel & BlueChannel) != 0)
1349 if (pixel.blue < *minima)
1350 *minima=(double) pixel.blue;
1351 if (pixel.blue > *maxima)
1352 *maxima=(double) pixel.blue;
1354 if ((channel & OpacityChannel) != 0)
1356 if (pixel.opacity < *minima)
1357 *minima=(double) pixel.opacity;
1358 if (pixel.opacity > *maxima)
1359 *maxima=(double) pixel.opacity;
1361 if (((channel & IndexChannel) != 0) &&
1362 (image->colorspace == CMYKColorspace))
1364 if ((double) indexes[x] < *minima)
1365 *minima=(double) indexes[x];
1366 if ((double) indexes[x] > *maxima)
1367 *maxima=(double) indexes[x];
1372 return(y == (long) image->rows ? MagickTrue : MagickFalse);
1376 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1380 % 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 %
1384 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1386 % GetImageChannelStatistics() returns statistics for each channel in the
1387 % image. The statistics include the channel depth, its minima, maxima, mean,
1388 % standard deviation, kurtosis and skewness. You can access the red channel
1389 % mean, for example, like this:
1391 % channel_statistics=GetImageChannelStatistics(image,excepton);
1392 % red_mean=channel_statistics[RedChannel].mean;
1394 % Use MagickRelinquishMemory() to free the statistics buffer.
1396 % The format of the GetImageChannelStatistics method is:
1398 % ChannelStatistics *GetImageChannelStatistics(const Image *image,
1399 % ExceptionInfo *exception)
1401 % A description of each parameter follows:
1403 % o image: the image.
1405 % o exception: return any errors or warnings in this structure.
1408 MagickExport ChannelStatistics *GetImageChannelStatistics(const Image *image,
1409 ExceptionInfo *exception)
1412 *channel_statistics;
1438 assert(image != (Image *) NULL);
1439 assert(image->signature == MagickSignature);
1440 if (image->debug != MagickFalse)
1441 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1442 length=AllChannels+1UL;
1443 channel_statistics=(ChannelStatistics *) AcquireQuantumMemory(length,
1444 sizeof(*channel_statistics));
1445 if (channel_statistics == (ChannelStatistics *) NULL)
1446 ThrowFatalException(ResourceLimitFatalError,"MemoryAllocationFailed");
1447 (void) ResetMagickMemory(channel_statistics,0,length*
1448 sizeof(*channel_statistics));
1449 for (i=0; i <= AllChannels; i++)
1451 channel_statistics[i].depth=1;
1452 channel_statistics[i].maxima=(-1.0E-37);
1453 channel_statistics[i].minima=1.0E+37;
1454 channel_statistics[i].mean=0.0;
1455 channel_statistics[i].standard_deviation=0.0;
1456 channel_statistics[i].kurtosis=0.0;
1457 channel_statistics[i].skewness=0.0;
1459 for (y=0; y < (long) image->rows; y++)
1461 register const IndexPacket
1464 register const PixelPacket
1470 p=GetVirtualPixels(image,0,y,image->columns,1,exception);
1471 if (p == (const PixelPacket *) NULL)
1473 indexes=GetVirtualIndexQueue(image);
1474 for (x=0; x < (long) image->columns; )
1476 if (channel_statistics[RedChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1478 depth=channel_statistics[RedChannel].depth;
1479 range=GetQuantumRange(depth);
1480 status=p->red != ScaleAnyToQuantum(ScaleQuantumToAny(p->red,range),
1481 range) ? MagickTrue : MagickFalse;
1482 if (status != MagickFalse)
1484 channel_statistics[RedChannel].depth++;
1488 if (channel_statistics[GreenChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1490 depth=channel_statistics[GreenChannel].depth;
1491 range=GetQuantumRange(depth);
1492 status=p->green != ScaleAnyToQuantum(ScaleQuantumToAny(p->green,
1493 range),range) ? MagickTrue : MagickFalse;
1494 if (status != MagickFalse)
1496 channel_statistics[GreenChannel].depth++;
1500 if (channel_statistics[BlueChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1502 depth=channel_statistics[BlueChannel].depth;
1503 range=GetQuantumRange(depth);
1504 status=p->blue != ScaleAnyToQuantum(ScaleQuantumToAny(p->blue,
1505 range),range) ? MagickTrue : MagickFalse;
1506 if (status != MagickFalse)
1508 channel_statistics[BlueChannel].depth++;
1512 if (image->matte != MagickFalse)
1514 if (channel_statistics[OpacityChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1516 depth=channel_statistics[OpacityChannel].depth;
1517 range=GetQuantumRange(depth);
1518 status=p->opacity != ScaleAnyToQuantum(ScaleQuantumToAny(
1519 p->opacity,range),range) ? MagickTrue : MagickFalse;
1520 if (status != MagickFalse)
1522 channel_statistics[OpacityChannel].depth++;
1527 if (image->colorspace == CMYKColorspace)
1529 if (channel_statistics[BlackChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1531 depth=channel_statistics[BlackChannel].depth;
1532 range=GetQuantumRange(depth);
1533 status=indexes[x] != ScaleAnyToQuantum(ScaleQuantumToAny(
1534 indexes[x],range),range) ? MagickTrue : MagickFalse;
1535 if (status != MagickFalse)
1537 channel_statistics[BlackChannel].depth++;
1542 if ((double) p->red < channel_statistics[RedChannel].minima)
1543 channel_statistics[RedChannel].minima=(double) GetRedPixelComponent(p);
1544 if ((double) p->red > channel_statistics[RedChannel].maxima)
1545 channel_statistics[RedChannel].maxima=(double) GetRedPixelComponent(p);
1546 channel_statistics[RedChannel].mean+=GetRedPixelComponent(p);
1547 channel_statistics[RedChannel].standard_deviation+=(double) p->red*
1548 GetRedPixelComponent(p);
1549 channel_statistics[RedChannel].kurtosis+=(double) p->red*p->red*
1550 p->red*GetRedPixelComponent(p);
1551 channel_statistics[RedChannel].skewness+=(double) p->red*p->red*
1552 GetRedPixelComponent(p);
1553 if ((double) p->green < channel_statistics[GreenChannel].minima)
1554 channel_statistics[GreenChannel].minima=(double)
1555 GetGreenPixelComponent(p);
1556 if ((double) p->green > channel_statistics[GreenChannel].maxima)
1557 channel_statistics[GreenChannel].maxima=(double)
1558 GetGreenPixelComponent(p);
1559 channel_statistics[GreenChannel].mean+=GetGreenPixelComponent(p);
1560 channel_statistics[GreenChannel].standard_deviation+=(double) p->green*
1561 GetGreenPixelComponent(p);
1562 channel_statistics[GreenChannel].kurtosis+=(double) p->green*p->green*
1563 p->green*GetGreenPixelComponent(p);
1564 channel_statistics[GreenChannel].skewness+=(double) p->green*p->green*
1565 GetGreenPixelComponent(p);
1566 if ((double) p->blue < channel_statistics[BlueChannel].minima)
1567 channel_statistics[BlueChannel].minima=(double)
1568 GetBluePixelComponent(p);
1569 if ((double) p->blue > channel_statistics[BlueChannel].maxima)
1570 channel_statistics[BlueChannel].maxima=(double)
1571 GetBluePixelComponent(p);
1572 channel_statistics[BlueChannel].mean+=GetBluePixelComponent(p);
1573 channel_statistics[BlueChannel].standard_deviation+=(double) p->blue*
1574 GetBluePixelComponent(p);
1575 channel_statistics[BlueChannel].kurtosis+=(double) p->blue*p->blue*
1576 p->blue*GetBluePixelComponent(p);
1577 channel_statistics[BlueChannel].skewness+=(double) p->blue*p->blue*
1578 GetBluePixelComponent(p);
1579 if (image->matte != MagickFalse)
1581 if ((double) p->opacity < channel_statistics[OpacityChannel].minima)
1582 channel_statistics[OpacityChannel].minima=(double)
1583 GetOpacityPixelComponent(p);
1584 if ((double) p->opacity > channel_statistics[OpacityChannel].maxima)
1585 channel_statistics[OpacityChannel].maxima=(double)
1586 GetOpacityPixelComponent(p);
1587 channel_statistics[OpacityChannel].mean+=GetOpacityPixelComponent(p);
1588 channel_statistics[OpacityChannel].standard_deviation+=(double)
1589 p->opacity*GetOpacityPixelComponent(p);
1590 channel_statistics[OpacityChannel].kurtosis+=(double) p->opacity*
1591 p->opacity*p->opacity*GetOpacityPixelComponent(p);
1592 channel_statistics[OpacityChannel].skewness+=(double) p->opacity*
1593 p->opacity*GetOpacityPixelComponent(p);
1595 if (image->colorspace == CMYKColorspace)
1597 if ((double) indexes[x] < channel_statistics[BlackChannel].minima)
1598 channel_statistics[BlackChannel].minima=(double) indexes[x];
1599 if ((double) indexes[x] > channel_statistics[BlackChannel].maxima)
1600 channel_statistics[BlackChannel].maxima=(double) indexes[x];
1601 channel_statistics[BlackChannel].mean+=indexes[x];
1602 channel_statistics[BlackChannel].standard_deviation+=(double)
1603 indexes[x]*indexes[x];
1604 channel_statistics[BlackChannel].kurtosis+=(double) indexes[x]*
1605 indexes[x]*indexes[x]*indexes[x];
1606 channel_statistics[BlackChannel].skewness+=(double) indexes[x]*
1607 indexes[x]*indexes[x];
1613 area=(double) image->columns*image->rows;
1614 for (i=0; i < AllChannels; i++)
1616 channel_statistics[i].mean/=area;
1617 channel_statistics[i].standard_deviation/=area;
1618 channel_statistics[i].kurtosis/=area;
1619 channel_statistics[i].skewness/=area;
1621 for (i=0; i < AllChannels; i++)
1623 channel_statistics[AllChannels].depth=(unsigned long) MagickMax((double)
1624 channel_statistics[AllChannels].depth,(double)
1625 channel_statistics[i].depth);
1626 channel_statistics[AllChannels].minima=MagickMin(
1627 channel_statistics[AllChannels].minima,channel_statistics[i].minima);
1628 channel_statistics[AllChannels].maxima=MagickMax(
1629 channel_statistics[AllChannels].maxima,channel_statistics[i].maxima);
1630 channel_statistics[AllChannels].mean+=channel_statistics[i].mean;
1631 channel_statistics[AllChannels].standard_deviation+=
1632 channel_statistics[i].standard_deviation;
1633 channel_statistics[AllChannels].kurtosis+=channel_statistics[i].kurtosis;
1634 channel_statistics[AllChannels].skewness+=channel_statistics[i].skewness;
1637 if (image->colorspace == CMYKColorspace)
1639 channel_statistics[AllChannels].mean/=channels;
1640 channel_statistics[AllChannels].standard_deviation/=channels;
1641 channel_statistics[AllChannels].kurtosis/=channels;
1642 channel_statistics[AllChannels].skewness/=channels;
1643 for (i=0; i <= AllChannels; i++)
1646 sum_squares=channel_statistics[i].standard_deviation;
1648 sum_cubes=channel_statistics[i].skewness;
1649 channel_statistics[i].standard_deviation=sqrt(
1650 channel_statistics[i].standard_deviation-
1651 (channel_statistics[i].mean*channel_statistics[i].mean));
1652 if (channel_statistics[i].standard_deviation == 0.0)
1654 channel_statistics[i].kurtosis=0.0;
1655 channel_statistics[i].skewness=0.0;
1659 channel_statistics[i].skewness=(channel_statistics[i].skewness-
1660 3.0*channel_statistics[i].mean*sum_squares+
1661 2.0*channel_statistics[i].mean*channel_statistics[i].mean*
1662 channel_statistics[i].mean)/
1663 (channel_statistics[i].standard_deviation*
1664 channel_statistics[i].standard_deviation*
1665 channel_statistics[i].standard_deviation);
1666 channel_statistics[i].kurtosis=(channel_statistics[i].kurtosis-
1667 4.0*channel_statistics[i].mean*sum_cubes+
1668 6.0*channel_statistics[i].mean*channel_statistics[i].mean*sum_squares-
1669 3.0*channel_statistics[i].mean*channel_statistics[i].mean*
1670 1.0*channel_statistics[i].mean*channel_statistics[i].mean)/
1671 (channel_statistics[i].standard_deviation*
1672 channel_statistics[i].standard_deviation*
1673 channel_statistics[i].standard_deviation*
1674 channel_statistics[i].standard_deviation)-3.0;
1677 return(channel_statistics);