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-2011 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 "MagickCore/studio.h"
44 #include "MagickCore/property.h"
45 #include "MagickCore/animate.h"
46 #include "MagickCore/blob.h"
47 #include "MagickCore/blob-private.h"
48 #include "MagickCore/cache.h"
49 #include "MagickCore/cache-private.h"
50 #include "MagickCore/cache-view.h"
51 #include "MagickCore/client.h"
52 #include "MagickCore/color.h"
53 #include "MagickCore/color-private.h"
54 #include "MagickCore/colorspace.h"
55 #include "MagickCore/colorspace-private.h"
56 #include "MagickCore/composite.h"
57 #include "MagickCore/composite-private.h"
58 #include "MagickCore/compress.h"
59 #include "MagickCore/constitute.h"
60 #include "MagickCore/display.h"
61 #include "MagickCore/draw.h"
62 #include "MagickCore/enhance.h"
63 #include "MagickCore/exception.h"
64 #include "MagickCore/exception-private.h"
65 #include "MagickCore/gem.h"
66 #include "MagickCore/gem-private.h"
67 #include "MagickCore/geometry.h"
68 #include "MagickCore/list.h"
69 #include "MagickCore/image-private.h"
70 #include "MagickCore/magic.h"
71 #include "MagickCore/magick.h"
72 #include "MagickCore/memory_.h"
73 #include "MagickCore/module.h"
74 #include "MagickCore/monitor.h"
75 #include "MagickCore/monitor-private.h"
76 #include "MagickCore/option.h"
77 #include "MagickCore/paint.h"
78 #include "MagickCore/pixel-accessor.h"
79 #include "MagickCore/profile.h"
80 #include "MagickCore/quantize.h"
81 #include "MagickCore/quantum-private.h"
82 #include "MagickCore/random_.h"
83 #include "MagickCore/random-private.h"
84 #include "MagickCore/segment.h"
85 #include "MagickCore/semaphore.h"
86 #include "MagickCore/signature-private.h"
87 #include "MagickCore/statistic.h"
88 #include "MagickCore/string_.h"
89 #include "MagickCore/thread-private.h"
90 #include "MagickCore/timer.h"
91 #include "MagickCore/utility.h"
92 #include "MagickCore/version.h"
95 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
99 % E v a l u a t e I m a g e %
103 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
105 % EvaluateImage() applies a value to the image with an arithmetic, relational,
106 % or logical operator to an image. Use these operations to lighten or darken
107 % an image, to increase or decrease contrast in an image, or to produce the
108 % "negative" of an image.
110 % The format of the EvaluateImage method is:
112 % MagickBooleanType EvaluateImage(Image *image,
113 % const MagickEvaluateOperator op,const double value,
114 % ExceptionInfo *exception)
115 % MagickBooleanType EvaluateImages(Image *images,
116 % const MagickEvaluateOperator op,const double value,
117 % ExceptionInfo *exception)
119 % A description of each parameter follows:
121 % o image: the image.
123 % o op: A channel op.
125 % o value: A value value.
127 % o exception: return any errors or warnings in this structure.
131 static PixelInfo **DestroyPixelThreadSet(PixelInfo **pixels)
136 assert(pixels != (PixelInfo **) NULL);
137 for (i=0; i < (ssize_t) GetOpenMPMaximumThreads(); i++)
138 if (pixels[i] != (PixelInfo *) NULL)
139 pixels[i]=(PixelInfo *) RelinquishMagickMemory(pixels[i]);
140 pixels=(PixelInfo **) RelinquishMagickMemory(pixels);
144 static PixelInfo **AcquirePixelThreadSet(const Image *image,
145 const size_t number_images)
158 number_threads=GetOpenMPMaximumThreads();
159 pixels=(PixelInfo **) AcquireQuantumMemory(number_threads,
161 if (pixels == (PixelInfo **) NULL)
162 return((PixelInfo **) NULL);
163 (void) ResetMagickMemory(pixels,0,number_threads*sizeof(*pixels));
164 for (i=0; i < (ssize_t) number_threads; i++)
166 length=image->columns;
167 if (length < number_images)
168 length=number_images;
169 pixels[i]=(PixelInfo *) AcquireQuantumMemory(length,
171 if (pixels[i] == (PixelInfo *) NULL)
172 return(DestroyPixelThreadSet(pixels));
173 for (j=0; j < (ssize_t) length; j++)
174 GetPixelInfo(image,&pixels[i][j]);
179 static inline double MagickMax(const double x,const double y)
186 #if defined(__cplusplus) || defined(c_plusplus)
190 static int IntensityCompare(const void *x,const void *y)
199 color_1=(const PixelInfo *) x;
200 color_2=(const PixelInfo *) y;
201 intensity=(int) GetPixelInfoIntensity(color_2)-(int)
202 GetPixelInfoIntensity(color_1);
206 #if defined(__cplusplus) || defined(c_plusplus)
210 static inline double MagickMin(const double x,const double y)
217 static MagickRealType ApplyEvaluateOperator(RandomInfo *random_info,
218 Quantum pixel,const MagickEvaluateOperator op,const MagickRealType value)
226 case UndefinedEvaluateOperator:
228 case AbsEvaluateOperator:
230 result=(MagickRealType) fabs((double) (pixel+value));
233 case AddEvaluateOperator:
235 result=(MagickRealType) (pixel+value);
238 case AddModulusEvaluateOperator:
241 This returns a 'floored modulus' of the addition which is a
242 positive result. It differs from % or fmod() which returns a
243 'truncated modulus' result, where floor() is replaced by trunc()
244 and could return a negative result (which is clipped).
247 result-=(QuantumRange+1.0)*floor((double) result/(QuantumRange+1.0));
250 case AndEvaluateOperator:
252 result=(MagickRealType) ((size_t) pixel & (size_t) (value+0.5));
255 case CosineEvaluateOperator:
257 result=(MagickRealType) (QuantumRange*(0.5*cos((double) (2.0*MagickPI*
258 QuantumScale*pixel*value))+0.5));
261 case DivideEvaluateOperator:
263 result=pixel/(value == 0.0 ? 1.0 : value);
266 case ExponentialEvaluateOperator:
268 result=(MagickRealType) (QuantumRange*exp((double) (value*QuantumScale*
272 case GaussianNoiseEvaluateOperator:
274 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
275 GaussianNoise,value);
278 case ImpulseNoiseEvaluateOperator:
280 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
284 case LaplacianNoiseEvaluateOperator:
286 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
287 LaplacianNoise,value);
290 case LeftShiftEvaluateOperator:
292 result=(MagickRealType) ((size_t) pixel << (size_t) (value+0.5));
295 case LogEvaluateOperator:
297 result=(MagickRealType) (QuantumRange*log((double) (QuantumScale*value*
298 pixel+1.0))/log((double) (value+1.0)));
301 case MaxEvaluateOperator:
303 result=(MagickRealType) MagickMax((double) pixel,value);
306 case MeanEvaluateOperator:
308 result=(MagickRealType) (pixel+value);
311 case MedianEvaluateOperator:
313 result=(MagickRealType) (pixel+value);
316 case MinEvaluateOperator:
318 result=(MagickRealType) MagickMin((double) pixel,value);
321 case MultiplicativeNoiseEvaluateOperator:
323 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
324 MultiplicativeGaussianNoise,value);
327 case MultiplyEvaluateOperator:
329 result=(MagickRealType) (value*pixel);
332 case OrEvaluateOperator:
334 result=(MagickRealType) ((size_t) pixel | (size_t) (value+0.5));
337 case PoissonNoiseEvaluateOperator:
339 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
343 case PowEvaluateOperator:
345 result=(MagickRealType) (QuantumRange*pow((double) (QuantumScale*pixel),
349 case RightShiftEvaluateOperator:
351 result=(MagickRealType) ((size_t) pixel >> (size_t) (value+0.5));
354 case SetEvaluateOperator:
359 case SineEvaluateOperator:
361 result=(MagickRealType) (QuantumRange*(0.5*sin((double) (2.0*MagickPI*
362 QuantumScale*pixel*value))+0.5));
365 case SubtractEvaluateOperator:
367 result=(MagickRealType) (pixel-value);
370 case ThresholdEvaluateOperator:
372 result=(MagickRealType) (((MagickRealType) pixel <= value) ? 0 :
376 case ThresholdBlackEvaluateOperator:
378 result=(MagickRealType) (((MagickRealType) pixel <= value) ? 0 : pixel);
381 case ThresholdWhiteEvaluateOperator:
383 result=(MagickRealType) (((MagickRealType) pixel > value) ? QuantumRange :
387 case UniformNoiseEvaluateOperator:
389 result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
393 case XorEvaluateOperator:
395 result=(MagickRealType) ((size_t) pixel ^ (size_t) (value+0.5));
402 MagickExport Image *EvaluateImages(const Image *images,
403 const MagickEvaluateOperator op,ExceptionInfo *exception)
405 #define EvaluateImageTag "Evaluate/Image"
423 **restrict evaluate_pixels,
427 **restrict random_info;
436 Ensure the image are the same size.
438 assert(images != (Image *) NULL);
439 assert(images->signature == MagickSignature);
440 if (images->debug != MagickFalse)
441 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",images->filename);
442 assert(exception != (ExceptionInfo *) NULL);
443 assert(exception->signature == MagickSignature);
444 for (next=images; next != (Image *) NULL; next=GetNextImageInList(next))
445 if ((next->columns != images->columns) || (next->rows != images->rows))
447 (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
448 "ImageWidthsOrHeightsDiffer","`%s'",images->filename);
449 return((Image *) NULL);
452 Initialize evaluate next attributes.
454 evaluate_image=CloneImage(images,images->columns,images->rows,MagickTrue,
456 if (evaluate_image == (Image *) NULL)
457 return((Image *) NULL);
458 if (SetImageStorageClass(evaluate_image,DirectClass,exception) == MagickFalse)
460 evaluate_image=DestroyImage(evaluate_image);
461 return((Image *) NULL);
463 number_images=GetImageListLength(images);
464 evaluate_pixels=AcquirePixelThreadSet(images,number_images);
465 if (evaluate_pixels == (PixelInfo **) NULL)
467 evaluate_image=DestroyImage(evaluate_image);
468 (void) ThrowMagickException(exception,GetMagickModule(),
469 ResourceLimitError,"MemoryAllocationFailed","`%s'",images->filename);
470 return((Image *) NULL);
473 Evaluate image pixels.
477 GetPixelInfo(images,&zero);
478 random_info=AcquireRandomInfoThreadSet();
479 evaluate_view=AcquireCacheView(evaluate_image);
480 if (op == MedianEvaluateOperator)
481 #if defined(MAGICKCORE_OPENMP_SUPPORT)
482 #pragma omp parallel for schedule(dynamic) shared(progress,status)
484 for (y=0; y < (ssize_t) evaluate_image->rows; y++)
493 id = GetOpenMPThreadId();
504 if (status == MagickFalse)
506 q=QueueCacheViewAuthenticPixels(evaluate_view,0,y,evaluate_image->columns,
508 if (q == (Quantum *) NULL)
513 evaluate_pixel=evaluate_pixels[id];
514 for (x=0; x < (ssize_t) evaluate_image->columns; x++)
519 for (i=0; i < (ssize_t) number_images; i++)
520 evaluate_pixel[i]=zero;
522 for (i=0; i < (ssize_t) number_images; i++)
524 register const Quantum
527 image_view=AcquireCacheView(next);
528 p=GetCacheViewVirtualPixels(image_view,x,y,1,1,exception);
529 if (p == (const Quantum *) NULL)
531 image_view=DestroyCacheView(image_view);
534 evaluate_pixel[i].red=ApplyEvaluateOperator(random_info[id],
535 GetPixelRed(next,p),op,evaluate_pixel[i].red);
536 evaluate_pixel[i].green=ApplyEvaluateOperator(random_info[id],
537 GetPixelGreen(next,p),op,evaluate_pixel[i].green);
538 evaluate_pixel[i].blue=ApplyEvaluateOperator(random_info[id],
539 GetPixelBlue(next,p),op,evaluate_pixel[i].blue);
540 if (evaluate_image->colorspace == CMYKColorspace)
541 evaluate_pixel[i].black=ApplyEvaluateOperator(random_info[id],
542 GetPixelBlack(next,p),op,evaluate_pixel[i].black);
543 evaluate_pixel[i].alpha=ApplyEvaluateOperator(random_info[id],
544 GetPixelAlpha(next,p),op,evaluate_pixel[i].alpha);
545 image_view=DestroyCacheView(image_view);
546 next=GetNextImageInList(next);
548 qsort((void *) evaluate_pixel,number_images,sizeof(*evaluate_pixel),
550 SetPixelRed(evaluate_image,
551 ClampToQuantum(evaluate_pixel[i/2].red),q);
552 SetPixelGreen(evaluate_image,
553 ClampToQuantum(evaluate_pixel[i/2].green),q);
554 SetPixelBlue(evaluate_image,
555 ClampToQuantum(evaluate_pixel[i/2].blue),q);
556 if (evaluate_image->colorspace == CMYKColorspace)
557 SetPixelBlack(evaluate_image,
558 ClampToQuantum(evaluate_pixel[i/2].black),q);
559 if (evaluate_image->matte == MagickFalse)
560 SetPixelAlpha(evaluate_image,
561 ClampToQuantum(evaluate_pixel[i/2].alpha),q);
563 SetPixelAlpha(evaluate_image,
564 ClampToQuantum(evaluate_pixel[i/2].alpha),q);
565 q+=GetPixelChannels(evaluate_image);
567 if (SyncCacheViewAuthenticPixels(evaluate_view,exception) == MagickFalse)
569 if (images->progress_monitor != (MagickProgressMonitor) NULL)
574 #if defined(MAGICKCORE_OPENMP_SUPPORT)
575 #pragma omp critical (MagickCore_EvaluateImages)
577 proceed=SetImageProgress(images,EvaluateImageTag,progress++,
578 evaluate_image->rows);
579 if (proceed == MagickFalse)
584 #if defined(MAGICKCORE_OPENMP_SUPPORT)
585 #pragma omp parallel for schedule(dynamic) shared(progress,status)
587 for (y=0; y < (ssize_t) evaluate_image->rows; y++)
596 id = GetOpenMPThreadId();
608 if (status == MagickFalse)
610 q=QueueCacheViewAuthenticPixels(evaluate_view,0,y,evaluate_image->columns,
612 if (q == (Quantum *) NULL)
617 evaluate_pixel=evaluate_pixels[id];
618 for (x=0; x < (ssize_t) evaluate_image->columns; x++)
619 evaluate_pixel[x]=zero;
621 for (i=0; i < (ssize_t) number_images; i++)
623 register const Quantum
626 image_view=AcquireCacheView(next);
627 p=GetCacheViewVirtualPixels(image_view,0,y,next->columns,1,exception);
628 if (p == (const Quantum *) NULL)
630 image_view=DestroyCacheView(image_view);
633 for (x=0; x < (ssize_t) next->columns; x++)
635 evaluate_pixel[x].red=ApplyEvaluateOperator(random_info[id],
636 GetPixelRed(next,p),i == 0 ? AddEvaluateOperator : op,
637 evaluate_pixel[x].red);
638 evaluate_pixel[x].green=ApplyEvaluateOperator(random_info[id],
639 GetPixelGreen(next,p),i == 0 ? AddEvaluateOperator : op,
640 evaluate_pixel[x].green);
641 evaluate_pixel[x].blue=ApplyEvaluateOperator(random_info[id],
642 GetPixelBlue(next,p),i == 0 ? AddEvaluateOperator : op,
643 evaluate_pixel[x].blue);
644 if (evaluate_image->colorspace == CMYKColorspace)
645 evaluate_pixel[x].black=ApplyEvaluateOperator(random_info[id],
646 GetPixelBlack(next,p),i == 0 ? AddEvaluateOperator : op,
647 evaluate_pixel[x].black);
648 evaluate_pixel[x].alpha=ApplyEvaluateOperator(random_info[id],
649 GetPixelAlpha(next,p),i == 0 ? AddEvaluateOperator : op,
650 evaluate_pixel[x].alpha);
651 p+=GetPixelChannels(next);
653 image_view=DestroyCacheView(image_view);
654 next=GetNextImageInList(next);
656 if (op == MeanEvaluateOperator)
657 for (x=0; x < (ssize_t) evaluate_image->columns; x++)
659 evaluate_pixel[x].red/=number_images;
660 evaluate_pixel[x].green/=number_images;
661 evaluate_pixel[x].blue/=number_images;
662 evaluate_pixel[x].black/=number_images;
663 evaluate_pixel[x].alpha/=number_images;
665 for (x=0; x < (ssize_t) evaluate_image->columns; x++)
667 SetPixelRed(evaluate_image,ClampToQuantum(evaluate_pixel[x].red),q);
668 SetPixelGreen(evaluate_image,ClampToQuantum(evaluate_pixel[x].green),q);
669 SetPixelBlue(evaluate_image,ClampToQuantum(evaluate_pixel[x].blue),q);
670 if (evaluate_image->colorspace == CMYKColorspace)
671 SetPixelBlack(evaluate_image,ClampToQuantum(evaluate_pixel[x].black),
673 if (evaluate_image->matte == MagickFalse)
674 SetPixelAlpha(evaluate_image,ClampToQuantum(evaluate_pixel[x].alpha),
677 SetPixelAlpha(evaluate_image,ClampToQuantum(evaluate_pixel[x].alpha),
679 q+=GetPixelChannels(evaluate_image);
681 if (SyncCacheViewAuthenticPixels(evaluate_view,exception) == MagickFalse)
683 if (images->progress_monitor != (MagickProgressMonitor) NULL)
688 #if defined(MAGICKCORE_OPENMP_SUPPORT)
689 #pragma omp critical (MagickCore_EvaluateImages)
691 proceed=SetImageProgress(images,EvaluateImageTag,progress++,
692 evaluate_image->rows);
693 if (proceed == MagickFalse)
697 evaluate_view=DestroyCacheView(evaluate_view);
698 evaluate_pixels=DestroyPixelThreadSet(evaluate_pixels);
699 random_info=DestroyRandomInfoThreadSet(random_info);
700 if (status == MagickFalse)
701 evaluate_image=DestroyImage(evaluate_image);
702 return(evaluate_image);
705 MagickExport MagickBooleanType EvaluateImage(Image *image,
706 const MagickEvaluateOperator op,const double value,ExceptionInfo *exception)
718 **restrict random_info;
723 assert(image != (Image *) NULL);
724 assert(image->signature == MagickSignature);
725 if (image->debug != MagickFalse)
726 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
727 assert(exception != (ExceptionInfo *) NULL);
728 assert(exception->signature == MagickSignature);
729 if (SetImageStorageClass(image,DirectClass,exception) == MagickFalse)
733 random_info=AcquireRandomInfoThreadSet();
734 image_view=AcquireCacheView(image);
735 #if defined(MAGICKCORE_OPENMP_SUPPORT)
736 #pragma omp parallel for schedule(dynamic,4) shared(progress,status)
738 for (y=0; y < (ssize_t) image->rows; y++)
741 id = GetOpenMPThreadId();
749 if (status == MagickFalse)
751 q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
752 if (q == (Quantum *) NULL)
757 for (x=0; x < (ssize_t) image->columns; x++)
759 if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
760 SetPixelRed(image,ClampToQuantum(ApplyEvaluateOperator(
761 random_info[id],GetPixelRed(image,q),op,value)),q);
762 if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
763 SetPixelGreen(image,ClampToQuantum(ApplyEvaluateOperator(
764 random_info[id],GetPixelGreen(image,q),op,value)),q);
765 if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
766 SetPixelBlue(image,ClampToQuantum(ApplyEvaluateOperator(
767 random_info[id],GetPixelBlue(image,q),op,value)),q);
768 if (((GetPixelBlackTraits(image) & UpdatePixelTrait) != 0) &&
769 (image->colorspace == CMYKColorspace))
770 SetPixelBlack(image,ClampToQuantum(ApplyEvaluateOperator(
771 random_info[id],GetPixelBlack(image,q),op,value)),q);
772 if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
774 if (image->matte == MagickFalse)
775 SetPixelAlpha(image,ClampToQuantum(ApplyEvaluateOperator(
776 random_info[id],GetPixelAlpha(image,q),op,value)),q);
778 SetPixelAlpha(image,ClampToQuantum(ApplyEvaluateOperator(
779 random_info[id],GetPixelAlpha(image,q),op,value)),q);
781 q+=GetPixelChannels(image);
783 if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
785 if (image->progress_monitor != (MagickProgressMonitor) NULL)
790 #if defined(MAGICKCORE_OPENMP_SUPPORT)
791 #pragma omp critical (MagickCore_EvaluateImage)
793 proceed=SetImageProgress(image,EvaluateImageTag,progress++,image->rows);
794 if (proceed == MagickFalse)
798 image_view=DestroyCacheView(image_view);
799 random_info=DestroyRandomInfoThreadSet(random_info);
804 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
808 % F u n c t i o n I m a g e %
812 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
814 % FunctionImage() applies a value to the image with an arithmetic, relational,
815 % or logical operator to an image. Use these operations to lighten or darken
816 % an image, to increase or decrease contrast in an image, or to produce the
817 % "negative" of an image.
819 % The format of the FunctionImage method is:
821 % MagickBooleanType FunctionImage(Image *image,
822 % const MagickFunction function,const ssize_t number_parameters,
823 % const double *parameters,ExceptionInfo *exception)
825 % A description of each parameter follows:
827 % o image: the image.
829 % o function: A channel function.
831 % o parameters: one or more parameters.
833 % o exception: return any errors or warnings in this structure.
837 static Quantum ApplyFunction(Quantum pixel,const MagickFunction function,
838 const size_t number_parameters,const double *parameters,
839 ExceptionInfo *exception)
851 case PolynomialFunction:
855 * Parameters: polynomial constants, highest to lowest order
856 * For example: c0*x^3 + c1*x^2 + c2*x + c3
859 for (i=0; i < (ssize_t) number_parameters; i++)
860 result = result*QuantumScale*pixel + parameters[i];
861 result *= QuantumRange;
864 case SinusoidFunction:
867 * Parameters: Freq, Phase, Ampl, bias
869 double freq,phase,ampl,bias;
870 freq = ( number_parameters >= 1 ) ? parameters[0] : 1.0;
871 phase = ( number_parameters >= 2 ) ? parameters[1] : 0.0;
872 ampl = ( number_parameters >= 3 ) ? parameters[2] : 0.5;
873 bias = ( number_parameters >= 4 ) ? parameters[3] : 0.5;
874 result=(MagickRealType) (QuantumRange*(ampl*sin((double) (2.0*MagickPI*
875 (freq*QuantumScale*pixel + phase/360.0) )) + bias ) );
880 /* Arcsin Function (peged at range limits for invalid results)
881 * Parameters: Width, Center, Range, Bias
883 double width,range,center,bias;
884 width = ( number_parameters >= 1 ) ? parameters[0] : 1.0;
885 center = ( number_parameters >= 2 ) ? parameters[1] : 0.5;
886 range = ( number_parameters >= 3 ) ? parameters[2] : 1.0;
887 bias = ( number_parameters >= 4 ) ? parameters[3] : 0.5;
888 result = 2.0/width*(QuantumScale*pixel - center);
889 if ( result <= -1.0 )
890 result = bias - range/2.0;
891 else if ( result >= 1.0 )
892 result = bias + range/2.0;
894 result=(MagickRealType) (range/MagickPI*asin((double) result)+bias);
895 result *= QuantumRange;
901 * Parameters: Slope, Center, Range, Bias
903 double slope,range,center,bias;
904 slope = ( number_parameters >= 1 ) ? parameters[0] : 1.0;
905 center = ( number_parameters >= 2 ) ? parameters[1] : 0.5;
906 range = ( number_parameters >= 3 ) ? parameters[2] : 1.0;
907 bias = ( number_parameters >= 4 ) ? parameters[3] : 0.5;
908 result=(MagickRealType) (MagickPI*slope*(QuantumScale*pixel-center));
909 result=(MagickRealType) (QuantumRange*(range/MagickPI*atan((double)
913 case UndefinedFunction:
916 return(ClampToQuantum(result));
919 MagickExport MagickBooleanType FunctionImage(Image *image,
920 const MagickFunction function,const size_t number_parameters,
921 const double *parameters,ExceptionInfo *exception)
923 #define FunctionImageTag "Function/Image "
937 assert(image != (Image *) NULL);
938 assert(image->signature == MagickSignature);
939 if (image->debug != MagickFalse)
940 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
941 assert(exception != (ExceptionInfo *) NULL);
942 assert(exception->signature == MagickSignature);
943 if (SetImageStorageClass(image,DirectClass,exception) == MagickFalse)
947 image_view=AcquireCacheView(image);
948 #if defined(MAGICKCORE_OPENMP_SUPPORT)
949 #pragma omp parallel for schedule(dynamic,4) shared(progress,status)
951 for (y=0; y < (ssize_t) image->rows; y++)
959 if (status == MagickFalse)
961 q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
962 if (q == (Quantum *) NULL)
967 for (x=0; x < (ssize_t) image->columns; x++)
969 if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
970 SetPixelRed(image,ApplyFunction(GetPixelRed(image,q),function,
971 number_parameters,parameters,exception),q);
972 if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
973 SetPixelGreen(image,ApplyFunction(GetPixelGreen(image,q),function,
974 number_parameters,parameters,exception),q);
975 if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
976 SetPixelBlue(image,ApplyFunction(GetPixelBlue(image,q),function,
977 number_parameters,parameters,exception),q);
978 if (((GetPixelBlackTraits(image) & UpdatePixelTrait) != 0) &&
979 (image->colorspace == CMYKColorspace))
980 SetPixelBlack(image,ApplyFunction(GetPixelBlack(image,q),function,
981 number_parameters,parameters,exception),q);
982 if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
984 if (image->matte == MagickFalse)
985 SetPixelAlpha(image,ApplyFunction(GetPixelAlpha(image,q),function,
986 number_parameters,parameters,exception),q);
988 SetPixelAlpha(image,ApplyFunction(GetPixelAlpha(image,q),function,
989 number_parameters,parameters,exception),q);
991 q+=GetPixelChannels(image);
993 if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
995 if (image->progress_monitor != (MagickProgressMonitor) NULL)
1000 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1001 #pragma omp critical (MagickCore_FunctionImage)
1003 proceed=SetImageProgress(image,FunctionImageTag,progress++,image->rows);
1004 if (proceed == MagickFalse)
1008 image_view=DestroyCacheView(image_view);
1013 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1017 % G e t I m a g e E x t r e m a %
1021 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1023 % GetImageExtrema() returns the extrema of one or more image channels.
1025 % The format of the GetImageExtrema method is:
1027 % MagickBooleanType GetImageExtrema(const Image *image,size_t *minima,
1028 % size_t *maxima,ExceptionInfo *exception)
1030 % A description of each parameter follows:
1032 % o image: the image.
1034 % o minima: the minimum value in the channel.
1036 % o maxima: the maximum value in the channel.
1038 % o exception: return any errors or warnings in this structure.
1041 MagickExport MagickBooleanType GetImageExtrema(const Image *image,
1042 size_t *minima,size_t *maxima,ExceptionInfo *exception)
1051 assert(image != (Image *) NULL);
1052 assert(image->signature == MagickSignature);
1053 if (image->debug != MagickFalse)
1054 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1055 status=GetImageRange(image,&min,&max,exception);
1056 *minima=(size_t) ceil(min-0.5);
1057 *maxima=(size_t) floor(max+0.5);
1062 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1066 % G e t I m a g e M e a n %
1070 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1072 % GetImageMean() returns the mean and standard deviation of one or more
1075 % The format of the GetImageMean method is:
1077 % MagickBooleanType GetImageMean(const Image *image,double *mean,
1078 % double *standard_deviation,ExceptionInfo *exception)
1080 % A description of each parameter follows:
1082 % o image: the image.
1084 % o mean: the average value in the channel.
1086 % o standard_deviation: the standard deviation of the channel.
1088 % o exception: return any errors or warnings in this structure.
1091 MagickExport MagickBooleanType GetImageMean(const Image *image,double *mean,
1092 double *standard_deviation,ExceptionInfo *exception)
1095 *channel_statistics;
1100 assert(image != (Image *) NULL);
1101 assert(image->signature == MagickSignature);
1102 if (image->debug != MagickFalse)
1103 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1104 channel_statistics=GetImageStatistics(image,exception);
1105 if (channel_statistics == (ChannelStatistics *) NULL)
1106 return(MagickFalse);
1108 channel_statistics[CompositeChannels].mean=0.0;
1109 channel_statistics[CompositeChannels].standard_deviation=0.0;
1110 if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
1112 channel_statistics[CompositeChannels].mean+=
1113 channel_statistics[RedChannel].mean;
1114 channel_statistics[CompositeChannels].standard_deviation+=
1115 channel_statistics[RedChannel].variance-
1116 channel_statistics[RedChannel].mean*
1117 channel_statistics[RedChannel].mean;
1120 if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
1122 channel_statistics[CompositeChannels].mean+=
1123 channel_statistics[GreenChannel].mean;
1124 channel_statistics[CompositeChannels].standard_deviation+=
1125 channel_statistics[GreenChannel].variance-
1126 channel_statistics[GreenChannel].mean*
1127 channel_statistics[GreenChannel].mean;
1130 if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
1132 channel_statistics[CompositeChannels].mean+=
1133 channel_statistics[BlueChannel].mean;
1134 channel_statistics[CompositeChannels].standard_deviation+=
1135 channel_statistics[BlueChannel].variance-
1136 channel_statistics[BlueChannel].mean*
1137 channel_statistics[BlueChannel].mean;
1140 if (((GetPixelBlackTraits(image) & UpdatePixelTrait) != 0) &&
1141 (image->colorspace == CMYKColorspace))
1143 channel_statistics[CompositeChannels].mean+=
1144 channel_statistics[BlackChannel].mean;
1145 channel_statistics[CompositeChannels].standard_deviation+=
1146 channel_statistics[BlackChannel].variance-
1147 channel_statistics[BlackChannel].mean*
1148 channel_statistics[BlackChannel].mean;
1151 if (((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0) &&
1152 (image->matte != MagickFalse))
1154 channel_statistics[CompositeChannels].mean+=
1155 channel_statistics[AlphaChannel].mean;
1156 channel_statistics[CompositeChannels].standard_deviation+=
1157 channel_statistics[AlphaChannel].variance-
1158 channel_statistics[AlphaChannel].mean*
1159 channel_statistics[AlphaChannel].mean;
1162 channel_statistics[CompositeChannels].mean/=channels;
1163 channel_statistics[CompositeChannels].standard_deviation=
1164 sqrt(channel_statistics[CompositeChannels].standard_deviation/channels);
1165 *mean=channel_statistics[CompositeChannels].mean;
1166 *standard_deviation=channel_statistics[CompositeChannels].standard_deviation;
1167 channel_statistics=(ChannelStatistics *) RelinquishMagickMemory(
1168 channel_statistics);
1173 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1177 % G e t I m a g e K u r t o s i s %
1181 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1183 % GetImageKurtosis() returns the kurtosis and skewness of one or more
1186 % The format of the GetImageKurtosis method is:
1188 % MagickBooleanType GetImageKurtosis(const Image *image,double *kurtosis,
1189 % double *skewness,ExceptionInfo *exception)
1191 % A description of each parameter follows:
1193 % o image: the image.
1195 % o kurtosis: the kurtosis of the channel.
1197 % o skewness: the skewness of the channel.
1199 % o exception: return any errors or warnings in this structure.
1202 MagickExport MagickBooleanType GetImageKurtosis(const Image *image,
1203 double *kurtosis,double *skewness,ExceptionInfo *exception)
1216 assert(image != (Image *) NULL);
1217 assert(image->signature == MagickSignature);
1218 if (image->debug != MagickFalse)
1219 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1224 standard_deviation=0.0;
1227 sum_fourth_power=0.0;
1228 for (y=0; y < (ssize_t) image->rows; y++)
1230 register const Quantum
1236 p=GetVirtualPixels(image,0,y,image->columns,1,exception);
1237 if (p == (const Quantum *) NULL)
1239 for (x=0; x < (ssize_t) image->columns; x++)
1241 if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
1243 mean+=GetPixelRed(image,p);
1244 sum_squares+=(double) GetPixelRed(image,p)*GetPixelRed(image,p);
1245 sum_cubes+=(double) GetPixelRed(image,p)*GetPixelRed(image,p)*
1246 GetPixelRed(image,p);
1247 sum_fourth_power+=(double) GetPixelRed(image,p)*
1248 GetPixelRed(image,p)*GetPixelRed(image,p)*GetPixelRed(image,p);
1251 if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
1253 mean+=GetPixelGreen(image,p);
1254 sum_squares+=(double) GetPixelGreen(image,p)*GetPixelGreen(image,p);
1255 sum_cubes+=(double) GetPixelGreen(image,p)*GetPixelGreen(image,p)*
1256 GetPixelGreen(image,p);
1257 sum_fourth_power+=(double) GetPixelGreen(image,p)*
1258 GetPixelGreen(image,p)*GetPixelGreen(image,p)*
1259 GetPixelGreen(image,p);
1262 if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
1264 mean+=GetPixelBlue(image,p);
1265 sum_squares+=(double) GetPixelBlue(image,p)*GetPixelBlue(image,p);
1266 sum_cubes+=(double) GetPixelBlue(image,p)*GetPixelBlue(image,p)*
1267 GetPixelBlue(image,p);
1268 sum_fourth_power+=(double) GetPixelBlue(image,p)*
1269 GetPixelBlue(image,p)*GetPixelBlue(image,p)*
1270 GetPixelBlue(image,p);
1273 if (((GetPixelBlackTraits(image) & UpdatePixelTrait) != 0) &&
1274 (image->colorspace == CMYKColorspace))
1276 mean+=GetPixelBlack(image,p);
1277 sum_squares+=(double) GetPixelBlack(image,p)*GetPixelBlack(image,p);
1278 sum_cubes+=(double) GetPixelBlack(image,p)*GetPixelBlack(image,p)*
1279 GetPixelBlack(image,p);
1280 sum_fourth_power+=(double) GetPixelBlack(image,p)*
1281 GetPixelBlack(image,p)*GetPixelBlack(image,p)*
1282 GetPixelBlack(image,p);
1285 if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
1287 mean+=GetPixelAlpha(image,p);
1288 sum_squares+=(double) GetPixelAlpha(image,p)*GetPixelAlpha(image,p);
1289 sum_cubes+=(double) GetPixelAlpha(image,p)*GetPixelAlpha(image,p)*
1290 GetPixelAlpha(image,p);
1291 sum_fourth_power+=(double) GetPixelAlpha(image,p)*
1292 GetPixelAlpha(image,p)*GetPixelAlpha(image,p)*
1293 GetPixelAlpha(image,p);
1296 p+=GetPixelChannels(image);
1299 if (y < (ssize_t) image->rows)
1300 return(MagickFalse);
1306 sum_fourth_power/=area;
1308 standard_deviation=sqrt(sum_squares-(mean*mean));
1309 if (standard_deviation != 0.0)
1311 *kurtosis=sum_fourth_power-4.0*mean*sum_cubes+6.0*mean*mean*sum_squares-
1312 3.0*mean*mean*mean*mean;
1313 *kurtosis/=standard_deviation*standard_deviation*standard_deviation*
1316 *skewness=sum_cubes-3.0*mean*sum_squares+2.0*mean*mean*mean;
1317 *skewness/=standard_deviation*standard_deviation*standard_deviation;
1319 return(y == (ssize_t) image->rows ? MagickTrue : MagickFalse);
1323 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1327 % G e t I m a g e R a n g e %
1331 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1333 % GetImageRange() returns the range of one or more image channels.
1335 % The format of the GetImageRange method is:
1337 % MagickBooleanType GetImageRange(const Image *image,double *minima,
1338 % double *maxima,ExceptionInfo *exception)
1340 % A description of each parameter follows:
1342 % o image: the image.
1344 % o minima: the minimum value in the channel.
1346 % o maxima: the maximum value in the channel.
1348 % o exception: return any errors or warnings in this structure.
1351 MagickExport MagickBooleanType GetImageRange(const Image *image,double *minima,
1352 double *maxima,ExceptionInfo *exception)
1360 assert(image != (Image *) NULL);
1361 assert(image->signature == MagickSignature);
1362 if (image->debug != MagickFalse)
1363 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1366 GetPixelInfo(image,&pixel);
1367 for (y=0; y < (ssize_t) image->rows; y++)
1369 register const Quantum
1375 p=GetVirtualPixels(image,0,y,image->columns,1,exception);
1376 if (p == (const Quantum *) NULL)
1378 for (x=0; x < (ssize_t) image->columns; x++)
1380 SetPixelInfo(image,p,&pixel);
1381 if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
1383 if (pixel.red < *minima)
1384 *minima=(double) pixel.red;
1385 if (pixel.red > *maxima)
1386 *maxima=(double) pixel.red;
1388 if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
1390 if (pixel.green < *minima)
1391 *minima=(double) pixel.green;
1392 if (pixel.green > *maxima)
1393 *maxima=(double) pixel.green;
1395 if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
1397 if (pixel.blue < *minima)
1398 *minima=(double) pixel.blue;
1399 if (pixel.blue > *maxima)
1400 *maxima=(double) pixel.blue;
1402 if (((GetPixelBlackTraits(image) & UpdatePixelTrait) != 0) &&
1403 (image->colorspace == CMYKColorspace))
1405 if (pixel.black < *minima)
1406 *minima=(double) pixel.black;
1407 if (pixel.black > *maxima)
1408 *maxima=(double) pixel.black;
1410 if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
1412 if (pixel.alpha < *minima)
1413 *minima=(double) pixel.alpha;
1414 if (pixel.alpha > *maxima)
1415 *maxima=(double) pixel.alpha;
1417 p+=GetPixelChannels(image);
1420 return(y == (ssize_t) image->rows ? MagickTrue : MagickFalse);
1424 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1428 % G e t I m a g e S t a t i s t i c s %
1432 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1434 % GetImageStatistics() returns statistics for each channel in the
1435 % image. The statistics include the channel depth, its minima, maxima, mean,
1436 % standard deviation, kurtosis and skewness. You can access the red channel
1437 % mean, for example, like this:
1439 % channel_statistics=GetImageStatistics(image,exception);
1440 % red_mean=channel_statistics[RedChannel].mean;
1442 % Use MagickRelinquishMemory() to free the statistics buffer.
1444 % The format of the GetImageStatistics method is:
1446 % ChannelStatistics *GetImageStatistics(const Image *image,
1447 % ExceptionInfo *exception)
1449 % A description of each parameter follows:
1451 % o image: the image.
1453 % o exception: return any errors or warnings in this structure.
1456 MagickExport ChannelStatistics *GetImageStatistics(const Image *image,
1457 ExceptionInfo *exception)
1460 *channel_statistics;
1482 assert(image != (Image *) NULL);
1483 assert(image->signature == MagickSignature);
1484 if (image->debug != MagickFalse)
1485 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1486 length=CompositeChannels+1UL;
1487 channel_statistics=(ChannelStatistics *) AcquireQuantumMemory(length,
1488 sizeof(*channel_statistics));
1489 if (channel_statistics == (ChannelStatistics *) NULL)
1490 ThrowFatalException(ResourceLimitFatalError,"MemoryAllocationFailed");
1491 (void) ResetMagickMemory(channel_statistics,0,length*
1492 sizeof(*channel_statistics));
1493 for (i=0; i <= (ssize_t) CompositeChannels; i++)
1495 channel_statistics[i].depth=1;
1496 channel_statistics[i].maxima=(-1.0E-37);
1497 channel_statistics[i].minima=1.0E+37;
1499 for (y=0; y < (ssize_t) image->rows; y++)
1501 register const Quantum
1507 p=GetVirtualPixels(image,0,y,image->columns,1,exception);
1508 if (p == (const Quantum *) NULL)
1510 for (x=0; x < (ssize_t) image->columns; )
1512 if (channel_statistics[RedChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1514 depth=channel_statistics[RedChannel].depth;
1515 range=GetQuantumRange(depth);
1516 status=GetPixelRed(image,p) != ScaleAnyToQuantum(ScaleQuantumToAny(
1517 GetPixelRed(image,p),range),range) ? MagickTrue : MagickFalse;
1518 if (status != MagickFalse)
1520 channel_statistics[RedChannel].depth++;
1524 if (channel_statistics[GreenChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1526 depth=channel_statistics[GreenChannel].depth;
1527 range=GetQuantumRange(depth);
1528 status=GetPixelGreen(image,p) != ScaleAnyToQuantum(ScaleQuantumToAny(
1529 GetPixelGreen(image,p),range),range) ? MagickTrue : MagickFalse;
1530 if (status != MagickFalse)
1532 channel_statistics[GreenChannel].depth++;
1536 if (channel_statistics[BlueChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1538 depth=channel_statistics[BlueChannel].depth;
1539 range=GetQuantumRange(depth);
1540 status=GetPixelBlue(image,p) != ScaleAnyToQuantum(ScaleQuantumToAny(
1541 GetPixelBlue(image,p),range),range) ? MagickTrue : MagickFalse;
1542 if (status != MagickFalse)
1544 channel_statistics[BlueChannel].depth++;
1548 if (image->matte != MagickFalse)
1550 if (channel_statistics[AlphaChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1552 depth=channel_statistics[AlphaChannel].depth;
1553 range=GetQuantumRange(depth);
1554 status=GetPixelAlpha(image,p) != ScaleAnyToQuantum(
1555 ScaleQuantumToAny(GetPixelAlpha(image,p),range),range) ?
1556 MagickTrue : MagickFalse;
1557 if (status != MagickFalse)
1559 channel_statistics[AlphaChannel].depth++;
1564 if (image->colorspace == CMYKColorspace)
1566 if (channel_statistics[BlackChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
1568 depth=channel_statistics[BlackChannel].depth;
1569 range=GetQuantumRange(depth);
1570 status=GetPixelBlack(image,p) != ScaleAnyToQuantum(
1571 ScaleQuantumToAny(GetPixelBlack(image,p),range),range) ?
1572 MagickTrue : MagickFalse;
1573 if (status != MagickFalse)
1575 channel_statistics[BlackChannel].depth++;
1580 if ((double) GetPixelRed(image,p) < channel_statistics[RedChannel].minima)
1581 channel_statistics[RedChannel].minima=(double) GetPixelRed(image,p);
1582 if ((double) GetPixelRed(image,p) > channel_statistics[RedChannel].maxima)
1583 channel_statistics[RedChannel].maxima=(double) GetPixelRed(image,p);
1584 channel_statistics[RedChannel].sum+=GetPixelRed(image,p);
1585 channel_statistics[RedChannel].sum_squared+=(double)
1586 GetPixelRed(image,p)*GetPixelRed(image,p);
1587 channel_statistics[RedChannel].sum_cubed+=(double)
1588 GetPixelRed(image,p)*GetPixelRed(image,p)*GetPixelRed(image,p);
1589 channel_statistics[RedChannel].sum_fourth_power+=(double)
1590 GetPixelRed(image,p)*GetPixelRed(image,p)*GetPixelRed(image,p)*
1591 GetPixelRed(image,p);
1592 if ((double) GetPixelGreen(image,p) < channel_statistics[GreenChannel].minima)
1593 channel_statistics[GreenChannel].minima=(double)
1594 GetPixelGreen(image,p);
1595 if ((double) GetPixelGreen(image,p) > channel_statistics[GreenChannel].maxima)
1596 channel_statistics[GreenChannel].maxima=(double) GetPixelGreen(image,p);
1597 channel_statistics[GreenChannel].sum+=GetPixelGreen(image,p);
1598 channel_statistics[GreenChannel].sum_squared+=(double)
1599 GetPixelGreen(image,p)*GetPixelGreen(image,p);
1600 channel_statistics[GreenChannel].sum_cubed+=(double)
1601 GetPixelGreen(image,p)*GetPixelGreen(image,p)*GetPixelGreen(image,p);
1602 channel_statistics[GreenChannel].sum_fourth_power+=(double)
1603 GetPixelGreen(image,p)*GetPixelGreen(image,p)*GetPixelGreen(image,p)*
1604 GetPixelGreen(image,p);
1605 if ((double) GetPixelBlue(image,p) < channel_statistics[BlueChannel].minima)
1606 channel_statistics[BlueChannel].minima=(double) GetPixelBlue(image,p);
1607 if ((double) GetPixelBlue(image,p) > channel_statistics[BlueChannel].maxima)
1608 channel_statistics[BlueChannel].maxima=(double) GetPixelBlue(image,p);
1609 channel_statistics[BlueChannel].sum+=GetPixelBlue(image,p);
1610 channel_statistics[BlueChannel].sum_squared+=(double)
1611 GetPixelBlue(image,p)*GetPixelBlue(image,p);
1612 channel_statistics[BlueChannel].sum_cubed+=(double)
1613 GetPixelBlue(image,p)*GetPixelBlue(image,p)*GetPixelBlue(image,p);
1614 channel_statistics[BlueChannel].sum_fourth_power+=(double)
1615 GetPixelBlue(image,p)*GetPixelBlue(image,p)*GetPixelBlue(image,p)*
1616 GetPixelBlue(image,p);
1617 if (image->colorspace == CMYKColorspace)
1619 if ((double) GetPixelBlack(image,p) < channel_statistics[BlackChannel].minima)
1620 channel_statistics[BlackChannel].minima=(double)
1621 GetPixelBlack(image,p);
1622 if ((double) GetPixelBlack(image,p) > channel_statistics[BlackChannel].maxima)
1623 channel_statistics[BlackChannel].maxima=(double)
1624 GetPixelBlack(image,p);
1625 channel_statistics[BlackChannel].sum+=GetPixelBlack(image,p);
1626 channel_statistics[BlackChannel].sum_squared+=(double)
1627 GetPixelBlack(image,p)*GetPixelBlack(image,p);
1628 channel_statistics[BlackChannel].sum_cubed+=(double)
1629 GetPixelBlack(image,p)*GetPixelBlack(image,p)*
1630 GetPixelBlack(image,p);
1631 channel_statistics[BlackChannel].sum_fourth_power+=(double)
1632 GetPixelBlack(image,p)*GetPixelBlack(image,p)*
1633 GetPixelBlack(image,p)*GetPixelBlack(image,p);
1635 if (image->matte != MagickFalse)
1637 if ((double) GetPixelAlpha(image,p) < channel_statistics[AlphaChannel].minima)
1638 channel_statistics[AlphaChannel].minima=(double)
1639 GetPixelAlpha(image,p);
1640 if ((double) GetPixelAlpha(image,p) > channel_statistics[AlphaChannel].maxima)
1641 channel_statistics[AlphaChannel].maxima=(double)
1642 GetPixelAlpha(image,p);
1643 channel_statistics[AlphaChannel].sum+=GetPixelAlpha(image,p);
1644 channel_statistics[AlphaChannel].sum_squared+=(double)
1645 GetPixelAlpha(image,p)*GetPixelAlpha(image,p);
1646 channel_statistics[AlphaChannel].sum_cubed+=(double)
1647 GetPixelAlpha(image,p)*GetPixelAlpha(image,p)*
1648 GetPixelAlpha(image,p);
1649 channel_statistics[AlphaChannel].sum_fourth_power+=(double)
1650 GetPixelAlpha(image,p)*GetPixelAlpha(image,p)*
1651 GetPixelAlpha(image,p)*GetPixelAlpha(image,p);
1654 p+=GetPixelChannels(image);
1657 area=(double) image->columns*image->rows;
1658 for (i=0; i < (ssize_t) CompositeChannels; i++)
1660 channel_statistics[i].sum/=area;
1661 channel_statistics[i].sum_squared/=area;
1662 channel_statistics[i].sum_cubed/=area;
1663 channel_statistics[i].sum_fourth_power/=area;
1664 channel_statistics[i].mean=channel_statistics[i].sum;
1665 channel_statistics[i].variance=channel_statistics[i].sum_squared;
1666 channel_statistics[i].standard_deviation=sqrt(
1667 channel_statistics[i].variance-(channel_statistics[i].mean*
1668 channel_statistics[i].mean));
1670 for (i=0; i < (ssize_t) CompositeChannels; i++)
1672 channel_statistics[CompositeChannels].depth=(size_t) MagickMax((double)
1673 channel_statistics[CompositeChannels].depth,(double)
1674 channel_statistics[i].depth);
1675 channel_statistics[CompositeChannels].minima=MagickMin(
1676 channel_statistics[CompositeChannels].minima,
1677 channel_statistics[i].minima);
1678 channel_statistics[CompositeChannels].maxima=MagickMax(
1679 channel_statistics[CompositeChannels].maxima,
1680 channel_statistics[i].maxima);
1681 channel_statistics[CompositeChannels].sum+=channel_statistics[i].sum;
1682 channel_statistics[CompositeChannels].sum_squared+=
1683 channel_statistics[i].sum_squared;
1684 channel_statistics[CompositeChannels].sum_cubed+=
1685 channel_statistics[i].sum_cubed;
1686 channel_statistics[CompositeChannels].sum_fourth_power+=
1687 channel_statistics[i].sum_fourth_power;
1688 channel_statistics[CompositeChannels].mean+=channel_statistics[i].mean;
1689 channel_statistics[CompositeChannels].variance+=
1690 channel_statistics[i].variance-channel_statistics[i].mean*
1691 channel_statistics[i].mean;
1692 channel_statistics[CompositeChannels].standard_deviation+=
1693 channel_statistics[i].variance-channel_statistics[i].mean*
1694 channel_statistics[i].mean;
1697 if (image->matte != MagickFalse)
1699 if (image->colorspace == CMYKColorspace)
1701 channel_statistics[CompositeChannels].sum/=channels;
1702 channel_statistics[CompositeChannels].sum_squared/=channels;
1703 channel_statistics[CompositeChannels].sum_cubed/=channels;
1704 channel_statistics[CompositeChannels].sum_fourth_power/=channels;
1705 channel_statistics[CompositeChannels].mean/=channels;
1706 channel_statistics[CompositeChannels].variance/=channels;
1707 channel_statistics[CompositeChannels].standard_deviation=
1708 sqrt(channel_statistics[CompositeChannels].standard_deviation/channels);
1709 channel_statistics[CompositeChannels].kurtosis/=channels;
1710 channel_statistics[CompositeChannels].skewness/=channels;
1711 for (i=0; i <= (ssize_t) CompositeChannels; i++)
1713 if (channel_statistics[i].standard_deviation == 0.0)
1715 channel_statistics[i].skewness=(channel_statistics[i].sum_cubed-
1716 3.0*channel_statistics[i].mean*channel_statistics[i].sum_squared+
1717 2.0*channel_statistics[i].mean*channel_statistics[i].mean*
1718 channel_statistics[i].mean)/(channel_statistics[i].standard_deviation*
1719 channel_statistics[i].standard_deviation*
1720 channel_statistics[i].standard_deviation);
1721 channel_statistics[i].kurtosis=(channel_statistics[i].sum_fourth_power-
1722 4.0*channel_statistics[i].mean*channel_statistics[i].sum_cubed+
1723 6.0*channel_statistics[i].mean*channel_statistics[i].mean*
1724 channel_statistics[i].sum_squared-3.0*channel_statistics[i].mean*
1725 channel_statistics[i].mean*1.0*channel_statistics[i].mean*
1726 channel_statistics[i].mean)/(channel_statistics[i].standard_deviation*
1727 channel_statistics[i].standard_deviation*
1728 channel_statistics[i].standard_deviation*
1729 channel_statistics[i].standard_deviation)-3.0;
1731 return(channel_statistics);