% SSSSS T A A T IIIII SSSSS T IIIII CCCC %
% %
% %
-% MagickCore Image Methods %
+% MagickCore Image Statistical Methods %
% %
% Software Design %
% John Cristy %
#include "magick/profile.h"
#include "magick/quantize.h"
#include "magick/random_.h"
+#include "magick/random-private.h"
#include "magick/segment.h"
#include "magick/semaphore.h"
#include "magick/signature-private.h"
% %
% %
% %
-% A v e r a g e I m a g e s %
+% E v a l u a t e I m a g e %
% %
% %
% %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
-% AverageImages() takes a set of images and averages them together. Each
-% image in the set must have the same width and height. AverageImages()
-% returns a single image with each corresponding pixel component of each
-% image averaged. On failure, a NULL image is returned and exception
-% describes the reason for the failure.
+% EvaluateImage() applies a value to the image with an arithmetic, relational,
+% or logical operator to an image. Use these operations to lighten or darken
+% an image, to increase or decrease contrast in an image, or to produce the
+% "negative" of an image.
%
-% The format of the AverageImages method is:
+% The format of the EvaluateImageChannel method is:
%
-% Image *AverageImages(Image *image,ExceptionInfo *exception)
+% MagickBooleanType EvaluateImage(Image *image,
+% const MagickEvaluateOperator op,const double value,
+% ExceptionInfo *exception)
+% MagickBooleanType EvaluateImages(Image *images,
+% const MagickEvaluateOperator op,const double value,
+% ExceptionInfo *exception)
+% MagickBooleanType EvaluateImageChannel(Image *image,
+% const ChannelType channel,const MagickEvaluateOperator op,
+% const double value,ExceptionInfo *exception)
%
% A description of each parameter follows:
%
-% o image: the image sequence.
+% o image: the image.
+%
+% o channel: the channel.
+%
+% o op: A channel op.
+%
+% o value: A value value.
%
% o exception: return any errors or warnings in this structure.
%
static MagickPixelPacket **DestroyPixelThreadSet(MagickPixelPacket **pixels)
{
- register long
+ register ssize_t
i;
assert(pixels != (MagickPixelPacket **) NULL);
- for (i=0; i < (long) GetOpenMPMaximumThreads(); i++)
+ for (i=0; i < (ssize_t) GetOpenMPMaximumThreads(); i++)
if (pixels[i] != (MagickPixelPacket *) NULL)
pixels[i]=(MagickPixelPacket *) RelinquishMagickMemory(pixels[i]);
pixels=(MagickPixelPacket **) RelinquishAlignedMemory(pixels);
static MagickPixelPacket **AcquirePixelThreadSet(const Image *image)
{
- register long
+ register ssize_t
i,
j;
MagickPixelPacket
**pixels;
- unsigned long
+ size_t
number_threads;
number_threads=GetOpenMPMaximumThreads();
if (pixels == (MagickPixelPacket **) NULL)
return((MagickPixelPacket **) NULL);
(void) ResetMagickMemory(pixels,0,number_threads*sizeof(*pixels));
- for (i=0; i < (long) number_threads; i++)
+ for (i=0; i < (ssize_t) number_threads; i++)
{
pixels[i]=(MagickPixelPacket *) AcquireQuantumMemory(image->columns,
sizeof(**pixels));
if (pixels[i] == (MagickPixelPacket *) NULL)
return(DestroyPixelThreadSet(pixels));
- for (j=0; j < (long) image->columns; j++)
+ for (j=0; j < (ssize_t) image->columns; j++)
GetMagickPixelPacket(image,&pixels[i][j]);
}
return(pixels);
}
-MagickExport Image *AverageImages(const Image *image,ExceptionInfo *exception)
+static inline double MagickMax(const double x,const double y)
+{
+ if (x > y)
+ return(x);
+ return(y);
+}
+
+static inline double MagickMin(const double x,const double y)
+{
+ if (x < y)
+ return(x);
+ return(y);
+}
+
+static MagickRealType ApplyEvaluateOperator(RandomInfo *random_info,
+ Quantum pixel,const MagickEvaluateOperator op,const MagickRealType value)
+{
+ MagickRealType
+ result;
+
+ result=0.0;
+ switch (op)
+ {
+ case UndefinedEvaluateOperator:
+ break;
+ case AbsEvaluateOperator:
+ {
+ result=(MagickRealType) fabs((double) (pixel+value));
+ break;
+ }
+ case AddEvaluateOperator:
+ {
+ result=(MagickRealType) (pixel+value);
+ break;
+ }
+ case AddModulusEvaluateOperator:
+ {
+ /*
+ This returns a 'floored modulus' of the addition which is a
+ positive result. It differs from % or fmod() which returns a
+ 'truncated modulus' result, where floor() is replaced by trunc()
+ and could return a negative result (which is clipped).
+ */
+ result=pixel+value;
+ result-=(QuantumRange+1.0)*floor((double) result/(QuantumRange+1.0));
+ break;
+ }
+ case AndEvaluateOperator:
+ {
+ result=(MagickRealType) ((size_t) pixel & (size_t)
+ (value+0.5));
+ break;
+ }
+ case CosineEvaluateOperator:
+ {
+ result=(MagickRealType) (QuantumRange*(0.5*cos((double) (2.0*MagickPI*
+ QuantumScale*pixel*value))+0.5));
+ break;
+ }
+ case DivideEvaluateOperator:
+ {
+ result=pixel/(value == 0.0 ? 1.0 : value);
+ break;
+ }
+ case GaussianNoiseEvaluateOperator:
+ {
+ result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
+ GaussianNoise,value);
+ break;
+ }
+ case ImpulseNoiseEvaluateOperator:
+ {
+ result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
+ ImpulseNoise,value);
+ break;
+ }
+ case LaplacianNoiseEvaluateOperator:
+ {
+ result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
+ LaplacianNoise,value);
+ break;
+ }
+ case LeftShiftEvaluateOperator:
+ {
+ result=(MagickRealType) ((size_t) pixel << (size_t)
+ (value+0.5));
+ break;
+ }
+ case LogEvaluateOperator:
+ {
+ result=(MagickRealType) (QuantumRange*log((double) (QuantumScale*value*
+ pixel+1.0))/log((double) (value+1.0)));
+ break;
+ }
+ case MaxEvaluateOperator:
+ {
+ result=(MagickRealType) MagickMax((double) pixel,value);
+ break;
+ }
+ case MeanEvaluateOperator:
+ {
+ result=(MagickRealType) (pixel+value);
+ break;
+ }
+ case MinEvaluateOperator:
+ {
+ result=(MagickRealType) MagickMin((double) pixel,value);
+ break;
+ }
+ case MultiplicativeNoiseEvaluateOperator:
+ {
+ result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
+ MultiplicativeGaussianNoise,value);
+ break;
+ }
+ case MultiplyEvaluateOperator:
+ {
+ result=(MagickRealType) (value*pixel);
+ break;
+ }
+ case OrEvaluateOperator:
+ {
+ result=(MagickRealType) ((size_t) pixel | (size_t)
+ (value+0.5));
+ break;
+ }
+ case PoissonNoiseEvaluateOperator:
+ {
+ result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
+ PoissonNoise,value);
+ break;
+ }
+ case PowEvaluateOperator:
+ {
+ result=(MagickRealType) (QuantumRange*pow((double) (QuantumScale*pixel),
+ (double) value));
+ break;
+ }
+ case RightShiftEvaluateOperator:
+ {
+ result=(MagickRealType) ((size_t) pixel >> (size_t)
+ (value+0.5));
+ break;
+ }
+ case SetEvaluateOperator:
+ {
+ result=value;
+ break;
+ }
+ case SineEvaluateOperator:
+ {
+ result=(MagickRealType) (QuantumRange*(0.5*sin((double) (2.0*MagickPI*
+ QuantumScale*pixel*value))+0.5));
+ break;
+ }
+ case SubtractEvaluateOperator:
+ {
+ result=(MagickRealType) (pixel-value);
+ break;
+ }
+ case ThresholdEvaluateOperator:
+ {
+ result=(MagickRealType) (((MagickRealType) pixel <= value) ? 0 :
+ QuantumRange);
+ break;
+ }
+ case ThresholdBlackEvaluateOperator:
+ {
+ result=(MagickRealType) (((MagickRealType) pixel <= value) ? 0 : pixel);
+ break;
+ }
+ case ThresholdWhiteEvaluateOperator:
+ {
+ result=(MagickRealType) (((MagickRealType) pixel > value) ? QuantumRange :
+ pixel);
+ break;
+ }
+ case UniformNoiseEvaluateOperator:
+ {
+ result=(MagickRealType) GenerateDifferentialNoise(random_info,pixel,
+ UniformNoise,value);
+ break;
+ }
+ case XorEvaluateOperator:
+ {
+ result=(MagickRealType) ((size_t) pixel ^ (size_t)
+ (value+0.5));
+ break;
+ }
+ }
+ return(result);
+}
+
+MagickExport MagickBooleanType EvaluateImage(Image *image,
+ const MagickEvaluateOperator op,const double value,ExceptionInfo *exception)
+{
+ MagickBooleanType
+ status;
+
+ status=EvaluateImageChannel(image,AllChannels,op,value,exception);
+ return(status);
+}
+
+MagickExport Image *EvaluateImages(const Image *images,
+ const MagickEvaluateOperator op,ExceptionInfo *exception)
{
-#define AverageImageTag "Average/Image"
+#define EvaluateImageTag "Evaluate/Image"
CacheView
- *average_view;
+ *evaluate_view;
const Image
*next;
Image
- *average_image;
-
- long
- progress,
- y;
+ *evaluate_image;
MagickBooleanType
status;
+ MagickOffsetType
+ progress;
+
MagickPixelPacket
- **restrict average_pixels,
+ **restrict evaluate_pixels,
zero;
- unsigned long
+ RandomInfo
+ **restrict random_info;
+
+ size_t
number_images;
+ ssize_t
+ y;
+
/*
Ensure the image are the same size.
*/
- assert(image != (Image *) NULL);
- assert(image->signature == MagickSignature);
- if (image->debug != MagickFalse)
- (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
+ assert(images != (Image *) NULL);
+ assert(images->signature == MagickSignature);
+ if (images->debug != MagickFalse)
+ (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",images->filename);
assert(exception != (ExceptionInfo *) NULL);
assert(exception->signature == MagickSignature);
- for (next=image; next != (Image *) NULL; next=GetNextImageInList(next))
- if ((next->columns != image->columns) || (next->rows != image->rows))
- ThrowImageException(OptionError,"ImageWidthsOrHeightsDiffer");
+ for (next=images; next != (Image *) NULL; next=GetNextImageInList(next))
+ if ((next->columns != images->columns) || (next->rows != images->rows))
+ {
+ (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
+ "ImageWidthsOrHeightsDiffer","`%s'",images->filename);
+ return((Image *) NULL);
+ }
/*
- Initialize average next attributes.
+ Initialize evaluate next attributes.
*/
- average_image=CloneImage(image,image->columns,image->rows,MagickTrue,
+ evaluate_image=CloneImage(images,images->columns,images->rows,MagickTrue,
exception);
- if (average_image == (Image *) NULL)
+ if (evaluate_image == (Image *) NULL)
return((Image *) NULL);
- if (SetImageStorageClass(average_image,DirectClass) == MagickFalse)
+ if (SetImageStorageClass(evaluate_image,DirectClass) == MagickFalse)
{
- InheritException(exception,&average_image->exception);
- average_image=DestroyImage(average_image);
+ InheritException(exception,&evaluate_image->exception);
+ evaluate_image=DestroyImage(evaluate_image);
return((Image *) NULL);
}
- average_pixels=AcquirePixelThreadSet(image);
- if (average_pixels == (MagickPixelPacket **) NULL)
+ evaluate_pixels=AcquirePixelThreadSet(images);
+ if (evaluate_pixels == (MagickPixelPacket **) NULL)
{
- average_image=DestroyImage(average_image);
- ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
+ evaluate_image=DestroyImage(evaluate_image);
+ (void) ThrowMagickException(exception,GetMagickModule(),
+ ResourceLimitError,"MemoryAllocationFailed","`%s'",images->filename);
+ return((Image *) NULL);
}
/*
- Average image pixels.
+ Evaluate image pixels.
*/
status=MagickTrue;
progress=0;
- GetMagickPixelPacket(image,&zero);
- number_images=GetImageListLength(image);
- average_view=AcquireCacheView(average_image);
+ GetMagickPixelPacket(images,&zero);
+ random_info=AcquireRandomInfoThreadSet();
+ number_images=GetImageListLength(images);
+ evaluate_view=AcquireCacheView(evaluate_image);
#if defined(MAGICKCORE_OPENMP_SUPPORT)
#pragma omp parallel for schedule(dynamic) shared(progress,status)
#endif
- for (y=0; y < (long) average_image->rows; y++)
+ for (y=0; y < (ssize_t) evaluate_image->rows; y++)
{
CacheView
*image_view;
const Image
*next;
+ int
+ id;
+
MagickPixelPacket
pixel;
register IndexPacket
- *restrict average_indexes;
+ *restrict evaluate_indexes;
- register long
+ register ssize_t
i,
- id,
x;
register MagickPixelPacket
- *average_pixel;
+ *evaluate_pixel;
register PixelPacket
*restrict q;
if (status == MagickFalse)
continue;
- q=QueueCacheViewAuthenticPixels(average_view,0,y,average_image->columns,1,
+ q=QueueCacheViewAuthenticPixels(evaluate_view,0,y,evaluate_image->columns,1,
exception);
if (q == (PixelPacket *) NULL)
{
status=MagickFalse;
continue;
}
- average_indexes=GetCacheViewAuthenticIndexQueue(average_view);
+ evaluate_indexes=GetCacheViewAuthenticIndexQueue(evaluate_view);
pixel=zero;
id=GetOpenMPThreadId();
- average_pixel=average_pixels[id];
- for (x=0; x < (long) average_image->columns; x++)
- average_pixel[x]=zero;
- next=image;
- for (i=0; i < (long) number_images; i++)
+ evaluate_pixel=evaluate_pixels[id];
+ for (x=0; x < (ssize_t) evaluate_image->columns; x++)
+ evaluate_pixel[x]=zero;
+ next=images;
+ for (i=0; i < (ssize_t) number_images; i++)
{
register const IndexPacket
*indexes;
break;
}
indexes=GetCacheViewVirtualIndexQueue(image_view);
- for (x=0; x < (long) next->columns; x++)
+ for (x=0; x < (ssize_t) next->columns; x++)
{
- SetMagickPixelPacket(next,p,indexes+x,&pixel);
- average_pixel[x].red+=QuantumScale*pixel.red;
- average_pixel[x].green+=QuantumScale*pixel.green;
- average_pixel[x].blue+=QuantumScale*pixel.blue;
- average_pixel[x].opacity+=QuantumScale*pixel.opacity;
- if (average_image->colorspace == CMYKColorspace)
- average_pixel[x].index+=QuantumScale*pixel.index;
+ evaluate_pixel[x].red=ApplyEvaluateOperator(random_info[id],p->red,
+ i == 0 ? AddEvaluateOperator : op,evaluate_pixel[x].red);
+ evaluate_pixel[x].green=ApplyEvaluateOperator(random_info[id],p->green,
+ i == 0 ? AddEvaluateOperator : op,evaluate_pixel[x].green);
+ evaluate_pixel[x].blue=ApplyEvaluateOperator(random_info[id],p->blue,
+ i == 0 ? AddEvaluateOperator : op,evaluate_pixel[x].blue);
+ evaluate_pixel[x].opacity=ApplyEvaluateOperator(random_info[id],
+ p->opacity,i == 0 ? AddEvaluateOperator : op,
+ evaluate_pixel[x].opacity);
+ if (evaluate_image->colorspace == CMYKColorspace)
+ evaluate_pixel[x].index=ApplyEvaluateOperator(random_info[id],
+ indexes[x],i == 0 ? AddEvaluateOperator : op,
+ evaluate_pixel[x].index);
p++;
}
image_view=DestroyCacheView(image_view);
next=GetNextImageInList(next);
}
- for (x=0; x < (long) average_image->columns; x++)
+ if (op == MeanEvaluateOperator)
+ for (x=0; x < (ssize_t) evaluate_image->columns; x++)
+ {
+ evaluate_pixel[x].red/=number_images;
+ evaluate_pixel[x].green/=number_images;
+ evaluate_pixel[x].blue/=number_images;
+ evaluate_pixel[x].opacity/=number_images;
+ evaluate_pixel[x].index/=number_images;
+ }
+ for (x=0; x < (ssize_t) evaluate_image->columns; x++)
{
- average_pixel[x].red=(MagickRealType) (QuantumRange*
- average_pixel[x].red/number_images);
- average_pixel[x].green=(MagickRealType) (QuantumRange*
- average_pixel[x].green/number_images);
- average_pixel[x].blue=(MagickRealType) (QuantumRange*
- average_pixel[x].blue/number_images);
- average_pixel[x].opacity=(MagickRealType) (QuantumRange*
- average_pixel[x].opacity/number_images);
- if (average_image->colorspace == CMYKColorspace)
- average_pixel[x].index=(MagickRealType) (QuantumRange*
- average_pixel[x].index/number_images);
- SetPixelPacket(average_image,&average_pixel[x],q,average_indexes+x);
+ q->red=ClampToQuantum(evaluate_pixel[x].red);
+ q->green=ClampToQuantum(evaluate_pixel[x].green);
+ q->blue=ClampToQuantum(evaluate_pixel[x].blue);
+ if (evaluate_image->matte == MagickFalse)
+ q->opacity=ClampToQuantum(evaluate_pixel[x].opacity);
+ else
+ q->opacity=ClampToQuantum(QuantumRange-evaluate_pixel[x].opacity);
+ if (evaluate_image->colorspace == CMYKColorspace)
+ evaluate_indexes[x]=ClampToQuantum(evaluate_pixel[x].index);
q++;
}
- if (SyncCacheViewAuthenticPixels(average_view,exception) == MagickFalse)
+ if (SyncCacheViewAuthenticPixels(evaluate_view,exception) == MagickFalse)
status=MagickFalse;
- if (image->progress_monitor != (MagickProgressMonitor) NULL)
+ if (images->progress_monitor != (MagickProgressMonitor) NULL)
{
MagickBooleanType
proceed;
#if defined(MAGICKCORE_OPENMP_SUPPORT)
- #pragma omp critical (MagickCore_AverageImages)
+ #pragma omp critical (MagickCore_EvaluateImages)
#endif
- proceed=SetImageProgress(image,AverageImageTag,progress++,
- average_image->rows);
+ proceed=SetImageProgress(images,EvaluateImageTag,progress++,
+ evaluate_image->rows);
if (proceed == MagickFalse)
status=MagickFalse;
}
}
- average_view=DestroyCacheView(average_view);
- average_pixels=DestroyPixelThreadSet(average_pixels);
+ evaluate_view=DestroyCacheView(evaluate_view);
+ evaluate_pixels=DestroyPixelThreadSet(evaluate_pixels);
+ random_info=DestroyRandomInfoThreadSet(random_info);
if (status == MagickFalse)
- average_image=DestroyImage(average_image);
- return(average_image);
+ evaluate_image=DestroyImage(evaluate_image);
+ return(evaluate_image);
+}
+
+MagickExport MagickBooleanType EvaluateImageChannel(Image *image,
+ const ChannelType channel,const MagickEvaluateOperator op,const double value,
+ ExceptionInfo *exception)
+{
+ CacheView
+ *image_view;
+
+ MagickBooleanType
+ status;
+
+ MagickOffsetType
+ progress;
+
+ RandomInfo
+ **restrict random_info;
+
+ ssize_t
+ y;
+
+ assert(image != (Image *) NULL);
+ assert(image->signature == MagickSignature);
+ if (image->debug != MagickFalse)
+ (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
+ assert(exception != (ExceptionInfo *) NULL);
+ assert(exception->signature == MagickSignature);
+ if (SetImageStorageClass(image,DirectClass) == MagickFalse)
+ {
+ InheritException(exception,&image->exception);
+ return(MagickFalse);
+ }
+ status=MagickTrue;
+ progress=0;
+ random_info=AcquireRandomInfoThreadSet();
+ image_view=AcquireCacheView(image);
+#if defined(MAGICKCORE_OPENMP_SUPPORT)
+ #pragma omp parallel for schedule(dynamic,4) shared(progress,status)
+#endif
+ for (y=0; y < (ssize_t) image->rows; y++)
+ {
+ int
+ id;
+
+ register IndexPacket
+ *restrict indexes;
+
+ register ssize_t
+ x;
+
+ register PixelPacket
+ *restrict q;
+
+ if (status == MagickFalse)
+ continue;
+ q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
+ if (q == (PixelPacket *) NULL)
+ {
+ status=MagickFalse;
+ continue;
+ }
+ indexes=GetCacheViewAuthenticIndexQueue(image_view);
+ id=GetOpenMPThreadId();
+ for (x=0; x < (ssize_t) image->columns; x++)
+ {
+ if ((channel & RedChannel) != 0)
+ q->red=ClampToQuantum(ApplyEvaluateOperator(random_info[id],q->red,op,
+ value));
+ if ((channel & GreenChannel) != 0)
+ q->green=ClampToQuantum(ApplyEvaluateOperator(random_info[id],q->green,
+ op,value));
+ if ((channel & BlueChannel) != 0)
+ q->blue=ClampToQuantum(ApplyEvaluateOperator(random_info[id],q->blue,op,
+ value));
+ if ((channel & OpacityChannel) != 0)
+ {
+ if (image->matte == MagickFalse)
+ q->opacity=ClampToQuantum(ApplyEvaluateOperator(random_info[id],
+ q->opacity,op,value));
+ else
+ q->opacity=ClampToQuantum(QuantumRange-ApplyEvaluateOperator(
+ random_info[id],(Quantum) GetAlphaPixelComponent(q),op,value));
+ }
+ if (((channel & IndexChannel) != 0) && (indexes != (IndexPacket *) NULL))
+ indexes[x]=(IndexPacket) ClampToQuantum(ApplyEvaluateOperator(
+ random_info[id],indexes[x],op,value));
+ q++;
+ }
+ if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
+ status=MagickFalse;
+ if (image->progress_monitor != (MagickProgressMonitor) NULL)
+ {
+ MagickBooleanType
+ proceed;
+
+#if defined(MAGICKCORE_OPENMP_SUPPORT)
+ #pragma omp critical (MagickCore_EvaluateImageChannel)
+#endif
+ proceed=SetImageProgress(image,EvaluateImageTag,progress++,image->rows);
+ if (proceed == MagickFalse)
+ status=MagickFalse;
+ }
+ }
+ image_view=DestroyCacheView(image_view);
+ random_info=DestroyRandomInfoThreadSet(random_info);
+ return(status);
+}
+\f
+/*
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% %
+% %
+% %
+% F u n c t i o n I m a g e %
+% %
+% %
+% %
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%
+% FunctionImage() applies a value to the image with an arithmetic, relational,
+% or logical operator to an image. Use these operations to lighten or darken
+% an image, to increase or decrease contrast in an image, or to produce the
+% "negative" of an image.
+%
+% The format of the FunctionImageChannel method is:
+%
+% MagickBooleanType FunctionImage(Image *image,
+% const MagickFunction function,const ssize_t number_parameters,
+% const double *parameters,ExceptionInfo *exception)
+% MagickBooleanType FunctionImageChannel(Image *image,
+% const ChannelType channel,const MagickFunction function,
+% const ssize_t number_parameters,const double *argument,
+% ExceptionInfo *exception)
+%
+% A description of each parameter follows:
+%
+% o image: the image.
+%
+% o channel: the channel.
+%
+% o function: A channel function.
+%
+% o parameters: one or more parameters.
+%
+% o exception: return any errors or warnings in this structure.
+%
+*/
+
+static Quantum ApplyFunction(Quantum pixel,const MagickFunction function,
+ const size_t number_parameters,const double *parameters,
+ ExceptionInfo *exception)
+{
+ MagickRealType
+ result;
+
+ register ssize_t
+ i;
+
+ (void) exception;
+ result=0.0;
+ switch (function)
+ {
+ case PolynomialFunction:
+ {
+ /*
+ * Polynomial
+ * Parameters: polynomial constants, highest to lowest order
+ * For example: c0*x^3 + c1*x^2 + c2*x + c3
+ */
+ result=0.0;
+ for (i=0; i < (ssize_t) number_parameters; i++)
+ result = result*QuantumScale*pixel + parameters[i];
+ result *= QuantumRange;
+ break;
+ }
+ case SinusoidFunction:
+ {
+ /* Sinusoid Function
+ * Parameters: Freq, Phase, Ampl, bias
+ */
+ double freq,phase,ampl,bias;
+ freq = ( number_parameters >= 1 ) ? parameters[0] : 1.0;
+ phase = ( number_parameters >= 2 ) ? parameters[1] : 0.0;
+ ampl = ( number_parameters >= 3 ) ? parameters[2] : 0.5;
+ bias = ( number_parameters >= 4 ) ? parameters[3] : 0.5;
+ result=(MagickRealType) (QuantumRange*(ampl*sin((double) (2.0*MagickPI*
+ (freq*QuantumScale*pixel + phase/360.0) )) + bias ) );
+ break;
+ }
+ case ArcsinFunction:
+ {
+ /* Arcsin Function (peged at range limits for invalid results)
+ * Parameters: Width, Center, Range, Bias
+ */
+ double width,range,center,bias;
+ width = ( number_parameters >= 1 ) ? parameters[0] : 1.0;
+ center = ( number_parameters >= 2 ) ? parameters[1] : 0.5;
+ range = ( number_parameters >= 3 ) ? parameters[2] : 1.0;
+ bias = ( number_parameters >= 4 ) ? parameters[3] : 0.5;
+ result = 2.0/width*(QuantumScale*pixel - center);
+ if ( result <= -1.0 )
+ result = bias - range/2.0;
+ else if ( result >= 1.0 )
+ result = bias + range/2.0;
+ else
+ result=range/MagickPI*asin((double)result) + bias;
+ result *= QuantumRange;
+ break;
+ }
+ case ArctanFunction:
+ {
+ /* Arctan Function
+ * Parameters: Slope, Center, Range, Bias
+ */
+ double slope,range,center,bias;
+ slope = ( number_parameters >= 1 ) ? parameters[0] : 1.0;
+ center = ( number_parameters >= 2 ) ? parameters[1] : 0.5;
+ range = ( number_parameters >= 3 ) ? parameters[2] : 1.0;
+ bias = ( number_parameters >= 4 ) ? parameters[3] : 0.5;
+ result = MagickPI*slope*(QuantumScale*pixel - center);
+ result=(MagickRealType) (QuantumRange*(range/MagickPI*atan((double)
+ result) + bias ) );
+ break;
+ }
+ case UndefinedFunction:
+ break;
+ }
+ return(ClampToQuantum(result));
+}
+
+MagickExport MagickBooleanType FunctionImage(Image *image,
+ const MagickFunction function,const size_t number_parameters,
+ const double *parameters,ExceptionInfo *exception)
+{
+ MagickBooleanType
+ status;
+
+ status=FunctionImageChannel(image,AllChannels,function,number_parameters,
+ parameters,exception);
+ return(status);
+}
+
+MagickExport MagickBooleanType FunctionImageChannel(Image *image,
+ const ChannelType channel,const MagickFunction function,
+ const size_t number_parameters,const double *parameters,
+ ExceptionInfo *exception)
+{
+#define FunctionImageTag "Function/Image "
+
+ CacheView
+ *image_view;
+
+ MagickBooleanType
+ status;
+
+ MagickOffsetType
+ progress;
+
+ ssize_t
+ y;
+
+ assert(image != (Image *) NULL);
+ assert(image->signature == MagickSignature);
+ if (image->debug != MagickFalse)
+ (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
+ assert(exception != (ExceptionInfo *) NULL);
+ assert(exception->signature == MagickSignature);
+ if (SetImageStorageClass(image,DirectClass) == MagickFalse)
+ {
+ InheritException(exception,&image->exception);
+ return(MagickFalse);
+ }
+ status=MagickTrue;
+ progress=0;
+ image_view=AcquireCacheView(image);
+#if defined(MAGICKCORE_OPENMP_SUPPORT)
+ #pragma omp parallel for schedule(dynamic,4) shared(progress,status)
+#endif
+ for (y=0; y < (ssize_t) image->rows; y++)
+ {
+ register IndexPacket
+ *restrict indexes;
+
+ register ssize_t
+ x;
+
+ register PixelPacket
+ *restrict q;
+
+ if (status == MagickFalse)
+ continue;
+ q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
+ if (q == (PixelPacket *) NULL)
+ {
+ status=MagickFalse;
+ continue;
+ }
+ indexes=GetCacheViewAuthenticIndexQueue(image_view);
+ for (x=0; x < (ssize_t) image->columns; x++)
+ {
+ if ((channel & RedChannel) != 0)
+ q->red=ApplyFunction(q->red,function,number_parameters,parameters,
+ exception);
+ if ((channel & GreenChannel) != 0)
+ q->green=ApplyFunction(q->green,function,number_parameters,parameters,
+ exception);
+ if ((channel & BlueChannel) != 0)
+ q->blue=ApplyFunction(q->blue,function,number_parameters,parameters,
+ exception);
+ if ((channel & OpacityChannel) != 0)
+ {
+ if (image->matte == MagickFalse)
+ q->opacity=ApplyFunction(q->opacity,function,number_parameters,
+ parameters,exception);
+ else
+ q->opacity=(Quantum) QuantumRange-ApplyFunction((Quantum)
+ GetAlphaPixelComponent(q),function,number_parameters,parameters,
+ exception);
+ }
+ if (((channel & IndexChannel) != 0) && (indexes != (IndexPacket *) NULL))
+ indexes[x]=(IndexPacket) ApplyFunction(indexes[x],function,
+ number_parameters,parameters,exception);
+ q++;
+ }
+ if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
+ status=MagickFalse;
+ if (image->progress_monitor != (MagickProgressMonitor) NULL)
+ {
+ MagickBooleanType
+ proceed;
+
+#if defined(MAGICKCORE_OPENMP_SUPPORT)
+ #pragma omp critical (MagickCore_FunctionImageChannel)
+#endif
+ proceed=SetImageProgress(image,FunctionImageTag,progress++,image->rows);
+ if (proceed == MagickFalse)
+ status=MagickFalse;
+ }
+ }
+ image_view=DestroyCacheView(image_view);
+ return(status);
}
\f
/*
% The format of the GetImageChannelExtrema method is:
%
% MagickBooleanType GetImageChannelExtrema(const Image *image,
-% const ChannelType channel,unsigned long *minima,unsigned long *maxima,
+% const ChannelType channel,size_t *minima,size_t *maxima,
% ExceptionInfo *exception)
%
% A description of each parameter follows:
*/
MagickExport MagickBooleanType GetImageExtrema(const Image *image,
- unsigned long *minima,unsigned long *maxima,ExceptionInfo *exception)
+ size_t *minima,size_t *maxima,ExceptionInfo *exception)
{
return(GetImageChannelExtrema(image,AllChannels,minima,maxima,exception));
}
MagickExport MagickBooleanType GetImageChannelExtrema(const Image *image,
- const ChannelType channel,unsigned long *minima,unsigned long *maxima,
+ const ChannelType channel,size_t *minima,size_t *maxima,
ExceptionInfo *exception)
{
double
if (image->debug != MagickFalse)
(void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
status=GetImageChannelRange(image,channel,&min,&max,exception);
- *minima=(unsigned long) (min+0.5);
- *maxima=(unsigned long) (max+0.5);
+ *minima=(size_t) ceil(min-0.5);
+ *maxima=(size_t) floor(max+0.5);
return(status);
}
\f
const ChannelType channel,double *mean,double *standard_deviation,
ExceptionInfo *exception)
{
- double
- area;
+ ChannelStatistics
+ *channel_statistics;
- long
- y;
+ size_t
+ channels;
assert(image != (Image *) NULL);
assert(image->signature == MagickSignature);
if (image->debug != MagickFalse)
(void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
- *mean=0.0;
- *standard_deviation=0.0;
- area=0.0;
- for (y=0; y < (long) image->rows; y++)
- {
- register const IndexPacket
- *restrict indexes;
-
- register const PixelPacket
- *restrict p;
-
- register long
- x;
-
- p=GetVirtualPixels(image,0,y,image->columns,1,exception);
- if (p == (const PixelPacket *) NULL)
- break;
- indexes=GetVirtualIndexQueue(image);
- for (x=0; x < (long) image->columns; x++)
+ channel_statistics=GetImageChannelStatistics(image,exception);
+ if (channel_statistics == (ChannelStatistics *) NULL)
+ return(MagickFalse);
+ channels=0;
+ channel_statistics[AllChannels].mean=0.0;
+ channel_statistics[AllChannels].standard_deviation=0.0;
+ if ((channel & RedChannel) != 0)
{
- if ((channel & RedChannel) != 0)
- {
- *mean+=GetRedSample(p);
- *standard_deviation+=(double) p->red*GetRedSample(p);
- area++;
- }
- if ((channel & GreenChannel) != 0)
- {
- *mean+=GetGreenSample(p);
- *standard_deviation+=(double) p->green*GetGreenSample(p);
- area++;
- }
- if ((channel & BlueChannel) != 0)
- {
- *mean+=GetBlueSample(p);
- *standard_deviation+=(double) p->blue*GetBlueSample(p);
- area++;
- }
- if ((channel & OpacityChannel) != 0)
- {
- *mean+=GetOpacitySample(p);
- *standard_deviation+=(double) p->opacity*GetOpacitySample(p);
- area++;
- }
- if (((channel & IndexChannel) != 0) &&
- (image->colorspace == CMYKColorspace))
- {
- *mean+=indexes[x];
- *standard_deviation+=(double) indexes[x]*indexes[x];
- area++;
- }
- p++;
+ channel_statistics[AllChannels].mean+=
+ channel_statistics[RedChannel].mean;
+ channel_statistics[AllChannels].standard_deviation+=
+ channel_statistics[RedChannel].variance-
+ channel_statistics[RedChannel].mean*
+ channel_statistics[RedChannel].mean;
+ channels++;
}
- }
- if (y < (long) image->rows)
- return(MagickFalse);
- if (area != 0)
+ if ((channel & GreenChannel) != 0)
+ {
+ channel_statistics[AllChannels].mean+=
+ channel_statistics[GreenChannel].mean;
+ channel_statistics[AllChannels].standard_deviation+=
+ channel_statistics[GreenChannel].variance-
+ channel_statistics[GreenChannel].mean*
+ channel_statistics[GreenChannel].mean;
+ channels++;
+ }
+ if ((channel & BlueChannel) != 0)
+ {
+ channel_statistics[AllChannels].mean+=
+ channel_statistics[BlueChannel].mean;
+ channel_statistics[AllChannels].standard_deviation+=
+ channel_statistics[BlueChannel].variance-
+ channel_statistics[BlueChannel].mean*
+ channel_statistics[BlueChannel].mean;
+ channels++;
+ }
+ if (((channel & OpacityChannel) != 0) &&
+ (image->matte != MagickFalse))
+ {
+ channel_statistics[AllChannels].mean+=
+ channel_statistics[OpacityChannel].mean;
+ channel_statistics[AllChannels].standard_deviation+=
+ channel_statistics[OpacityChannel].variance-
+ channel_statistics[OpacityChannel].mean*
+ channel_statistics[OpacityChannel].mean;
+ channels++;
+ }
+ if (((channel & IndexChannel) != 0) &&
+ (image->colorspace == CMYKColorspace))
{
- *mean/=area;
- *standard_deviation/=area;
+ channel_statistics[AllChannels].mean+=
+ channel_statistics[BlackChannel].mean;
+ channel_statistics[AllChannels].standard_deviation+=
+ channel_statistics[BlackChannel].variance-
+ channel_statistics[BlackChannel].mean*
+ channel_statistics[BlackChannel].mean;
+ channels++;
}
- *standard_deviation=sqrt(*standard_deviation-(*mean*(*mean)));
- return(y == (long) image->rows ? MagickTrue : MagickFalse);
+ channel_statistics[AllChannels].mean/=channels;
+ channel_statistics[AllChannels].standard_deviation=
+ sqrt(channel_statistics[AllChannels].standard_deviation/channels);
+ *mean=channel_statistics[AllChannels].mean;
+ *standard_deviation=channel_statistics[AllChannels].standard_deviation;
+ channel_statistics=(ChannelStatistics *) RelinquishMagickMemory(
+ channel_statistics);
+ return(MagickTrue);
}
\f
/*
sum_cubes,
sum_fourth_power;
- long
+ ssize_t
y;
assert(image != (Image *) NULL);
sum_squares=0.0;
sum_cubes=0.0;
sum_fourth_power=0.0;
- for (y=0; y < (long) image->rows; y++)
+ for (y=0; y < (ssize_t) image->rows; y++)
{
register const IndexPacket
*restrict indexes;
register const PixelPacket
*restrict p;
- register long
+ register ssize_t
x;
p=GetVirtualPixels(image,0,y,image->columns,1,exception);
if (p == (const PixelPacket *) NULL)
break;
indexes=GetVirtualIndexQueue(image);
- for (x=0; x < (long) image->columns; x++)
+ for (x=0; x < (ssize_t) image->columns; x++)
{
if ((channel & RedChannel) != 0)
{
- mean+=GetRedSample(p);
- sum_squares+=(double) p->red*GetRedSample(p);
- sum_cubes+=(double) p->red*p->red*GetRedSample(p);
- sum_fourth_power+=(double) p->red*p->red*p->red*GetRedSample(p);
+ mean+=GetRedPixelComponent(p);
+ sum_squares+=(double) p->red*GetRedPixelComponent(p);
+ sum_cubes+=(double) p->red*p->red*GetRedPixelComponent(p);
+ sum_fourth_power+=(double) p->red*p->red*p->red*
+ GetRedPixelComponent(p);
area++;
}
if ((channel & GreenChannel) != 0)
{
- mean+=GetGreenSample(p);
- sum_squares+=(double) p->green*GetGreenSample(p);
- sum_cubes+=(double) p->green*p->green*GetGreenSample(p);
- sum_fourth_power+=(double) p->green*p->green*p->green*GetGreenSample(p);
+ mean+=GetGreenPixelComponent(p);
+ sum_squares+=(double) p->green*GetGreenPixelComponent(p);
+ sum_cubes+=(double) p->green*p->green*GetGreenPixelComponent(p);
+ sum_fourth_power+=(double) p->green*p->green*p->green*
+ GetGreenPixelComponent(p);
area++;
}
if ((channel & BlueChannel) != 0)
{
- mean+=GetBlueSample(p);
- sum_squares+=(double) p->blue*GetBlueSample(p);
- sum_cubes+=(double) p->blue*p->blue*GetBlueSample(p);
- sum_fourth_power+=(double) p->blue*p->blue*p->blue*GetBlueSample(p);
+ mean+=GetBluePixelComponent(p);
+ sum_squares+=(double) p->blue*GetBluePixelComponent(p);
+ sum_cubes+=(double) p->blue*p->blue*GetBluePixelComponent(p);
+ sum_fourth_power+=(double) p->blue*p->blue*p->blue*
+ GetBluePixelComponent(p);
area++;
}
if ((channel & OpacityChannel) != 0)
{
- mean+=GetOpacitySample(p);
- sum_squares+=(double) p->opacity*GetOpacitySample(p);
- sum_cubes+=(double) p->opacity*p->opacity*GetOpacitySample(p);
+ mean+=GetOpacityPixelComponent(p);
+ sum_squares+=(double) p->opacity*GetOpacityPixelComponent(p);
+ sum_cubes+=(double) p->opacity*p->opacity*GetOpacityPixelComponent(p);
sum_fourth_power+=(double) p->opacity*p->opacity*p->opacity*
- GetOpacitySample(p);
+ GetOpacityPixelComponent(p);
area++;
}
if (((channel & IndexChannel) != 0) &&
p++;
}
}
- if (y < (long) image->rows)
+ if (y < (ssize_t) image->rows)
return(MagickFalse);
if (area != 0.0)
{
*skewness=sum_cubes-3.0*mean*sum_squares+2.0*mean*mean*mean;
*skewness/=standard_deviation*standard_deviation*standard_deviation;
}
- return(y == (long) image->rows ? MagickTrue : MagickFalse);
+ return(y == (ssize_t) image->rows ? MagickTrue : MagickFalse);
}
-
+\f
/*
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% %
const ChannelType channel,double *minima,double *maxima,
ExceptionInfo *exception)
{
- long
+ ssize_t
y;
MagickPixelPacket
*maxima=(-1.0E-37);
*minima=1.0E+37;
GetMagickPixelPacket(image,&pixel);
- for (y=0; y < (long) image->rows; y++)
+ for (y=0; y < (ssize_t) image->rows; y++)
{
register const IndexPacket
*restrict indexes;
register const PixelPacket
*restrict p;
- register long
+ register ssize_t
x;
p=GetVirtualPixels(image,0,y,image->columns,1,exception);
if (p == (const PixelPacket *) NULL)
break;
indexes=GetVirtualIndexQueue(image);
- for (x=0; x < (long) image->columns; x++)
+ for (x=0; x < (ssize_t) image->columns; x++)
{
SetMagickPixelPacket(image,p,indexes+x,&pixel);
if ((channel & RedChannel) != 0)
p++;
}
}
- return(y == (long) image->rows ? MagickTrue : MagickFalse);
+ return(y == (ssize_t) image->rows ? MagickTrue : MagickFalse);
}
\f
/*
% o exception: return any errors or warnings in this structure.
%
*/
-
-static inline double MagickMax(const double x,const double y)
-{
- if (x > y)
- return(x);
- return(y);
-}
-
-static inline double MagickMin(const double x,const double y)
-{
- if (x < y)
- return(x);
- return(y);
-}
-
MagickExport ChannelStatistics *GetImageChannelStatistics(const Image *image,
ExceptionInfo *exception)
{
*channel_statistics;
double
- area,
- sum_squares,
- sum_cubes;
+ area;
- long
+ ssize_t
y;
MagickStatusType
QuantumAny
range;
- register long
+ register ssize_t
i;
size_t
length;
- unsigned long
+ size_t
channels,
depth;
channel_statistics[i].depth=1;
channel_statistics[i].maxima=(-1.0E-37);
channel_statistics[i].minima=1.0E+37;
- channel_statistics[i].mean=0.0;
- channel_statistics[i].standard_deviation=0.0;
- channel_statistics[i].kurtosis=0.0;
- channel_statistics[i].skewness=0.0;
}
- for (y=0; y < (long) image->rows; y++)
+ for (y=0; y < (ssize_t) image->rows; y++)
{
register const IndexPacket
*restrict indexes;
register const PixelPacket
*restrict p;
- register long
+ register ssize_t
x;
p=GetVirtualPixels(image,0,y,image->columns,1,exception);
if (p == (const PixelPacket *) NULL)
break;
indexes=GetVirtualIndexQueue(image);
- for (x=0; x < (long) image->columns; )
+ for (x=0; x < (ssize_t) image->columns; )
{
if (channel_statistics[RedChannel].depth != MAGICKCORE_QUANTUM_DEPTH)
{
}
}
if ((double) p->red < channel_statistics[RedChannel].minima)
- channel_statistics[RedChannel].minima=(double) GetRedSample(p);
+ channel_statistics[RedChannel].minima=(double) GetRedPixelComponent(p);
if ((double) p->red > channel_statistics[RedChannel].maxima)
- channel_statistics[RedChannel].maxima=(double) GetRedSample(p);
- channel_statistics[RedChannel].mean+=GetRedSample(p);
- channel_statistics[RedChannel].standard_deviation+=(double) p->red*GetRedSample(p);
- channel_statistics[RedChannel].kurtosis+=(double) p->red*p->red*
- p->red*GetRedSample(p);
- channel_statistics[RedChannel].skewness+=(double) p->red*p->red*GetRedSample(p);
+ channel_statistics[RedChannel].maxima=(double) GetRedPixelComponent(p);
+ channel_statistics[RedChannel].sum+=GetRedPixelComponent(p);
+ channel_statistics[RedChannel].sum_squared+=(double) p->red*
+ GetRedPixelComponent(p);
+ channel_statistics[RedChannel].sum_cubed+=(double) p->red*p->red*
+ GetRedPixelComponent(p);
+ channel_statistics[RedChannel].sum_fourth_power+=(double) p->red*p->red*
+ p->red*GetRedPixelComponent(p);
if ((double) p->green < channel_statistics[GreenChannel].minima)
- channel_statistics[GreenChannel].minima=(double) GetGreenSample(p);
+ channel_statistics[GreenChannel].minima=(double)
+ GetGreenPixelComponent(p);
if ((double) p->green > channel_statistics[GreenChannel].maxima)
- channel_statistics[GreenChannel].maxima=(double) GetGreenSample(p);
- channel_statistics[GreenChannel].mean+=GetGreenSample(p);
- channel_statistics[GreenChannel].standard_deviation+=(double) p->green*
- GetGreenSample(p);
- channel_statistics[GreenChannel].kurtosis+=(double) p->green*p->green*
- p->green*GetGreenSample(p);
- channel_statistics[GreenChannel].skewness+=(double) p->green*p->green*
- GetGreenSample(p);
+ channel_statistics[GreenChannel].maxima=(double)
+ GetGreenPixelComponent(p);
+ channel_statistics[GreenChannel].sum+=GetGreenPixelComponent(p);
+ channel_statistics[GreenChannel].sum_squared+=(double) p->green*
+ GetGreenPixelComponent(p);
+ channel_statistics[GreenChannel].sum_cubed+=(double) p->green*p->green*
+ GetGreenPixelComponent(p);
+ channel_statistics[GreenChannel].sum_fourth_power+=(double) p->green*
+ p->green*p->green*GetGreenPixelComponent(p);
if ((double) p->blue < channel_statistics[BlueChannel].minima)
- channel_statistics[BlueChannel].minima=(double) GetBlueSample(p);
+ channel_statistics[BlueChannel].minima=(double)
+ GetBluePixelComponent(p);
if ((double) p->blue > channel_statistics[BlueChannel].maxima)
- channel_statistics[BlueChannel].maxima=(double) GetBlueSample(p);
- channel_statistics[BlueChannel].mean+=GetBlueSample(p);
- channel_statistics[BlueChannel].standard_deviation+=(double) p->blue*
- GetBlueSample(p);
- channel_statistics[BlueChannel].kurtosis+=(double) p->blue*p->blue*
- p->blue*GetBlueSample(p);
- channel_statistics[BlueChannel].skewness+=(double) p->blue*p->blue*
- GetBlueSample(p);
+ channel_statistics[BlueChannel].maxima=(double)
+ GetBluePixelComponent(p);
+ channel_statistics[BlueChannel].sum+=GetBluePixelComponent(p);
+ channel_statistics[BlueChannel].sum_squared+=(double) p->blue*
+ GetBluePixelComponent(p);
+ channel_statistics[BlueChannel].sum_cubed+=(double) p->blue*p->blue*
+ GetBluePixelComponent(p);
+ channel_statistics[BlueChannel].sum_fourth_power+=(double) p->blue*
+ p->blue*p->blue*GetBluePixelComponent(p);
if (image->matte != MagickFalse)
{
if ((double) p->opacity < channel_statistics[OpacityChannel].minima)
- channel_statistics[OpacityChannel].minima=(double) GetOpacitySample(p);
+ channel_statistics[OpacityChannel].minima=(double)
+ GetOpacityPixelComponent(p);
if ((double) p->opacity > channel_statistics[OpacityChannel].maxima)
- channel_statistics[OpacityChannel].maxima=(double) GetOpacitySample(p);
- channel_statistics[OpacityChannel].mean+=GetOpacitySample(p);
- channel_statistics[OpacityChannel].standard_deviation+=(double)
- p->opacity*GetOpacitySample(p);
- channel_statistics[OpacityChannel].kurtosis+=(double) p->opacity*
- p->opacity*p->opacity*GetOpacitySample(p);
- channel_statistics[OpacityChannel].skewness+=(double) p->opacity*
- p->opacity*GetOpacitySample(p);
+ channel_statistics[OpacityChannel].maxima=(double)
+ GetOpacityPixelComponent(p);
+ channel_statistics[OpacityChannel].sum+=GetOpacityPixelComponent(p);
+ channel_statistics[OpacityChannel].sum_squared+=(double)
+ p->opacity*GetOpacityPixelComponent(p);
+ channel_statistics[OpacityChannel].sum_cubed+=(double) p->opacity*
+ p->opacity*GetOpacityPixelComponent(p);
+ channel_statistics[OpacityChannel].sum_fourth_power+=(double)
+ p->opacity*p->opacity*p->opacity*GetOpacityPixelComponent(p);
}
if (image->colorspace == CMYKColorspace)
{
channel_statistics[BlackChannel].minima=(double) indexes[x];
if ((double) indexes[x] > channel_statistics[BlackChannel].maxima)
channel_statistics[BlackChannel].maxima=(double) indexes[x];
- channel_statistics[BlackChannel].mean+=indexes[x];
- channel_statistics[BlackChannel].standard_deviation+=(double)
+ channel_statistics[BlackChannel].sum+=indexes[x];
+ channel_statistics[BlackChannel].sum_squared+=(double)
indexes[x]*indexes[x];
- channel_statistics[BlackChannel].kurtosis+=(double) indexes[x]*
- indexes[x]*indexes[x]*indexes[x];
- channel_statistics[BlackChannel].skewness+=(double) indexes[x]*
+ channel_statistics[BlackChannel].sum_cubed+=(double) indexes[x]*
indexes[x]*indexes[x];
+ channel_statistics[BlackChannel].sum_fourth_power+=(double)
+ indexes[x]*indexes[x]*indexes[x]*indexes[x];
}
x++;
p++;
area=(double) image->columns*image->rows;
for (i=0; i < AllChannels; i++)
{
- channel_statistics[i].mean/=area;
- channel_statistics[i].standard_deviation/=area;
- channel_statistics[i].kurtosis/=area;
- channel_statistics[i].skewness/=area;
+ channel_statistics[i].sum/=area;
+ channel_statistics[i].sum_squared/=area;
+ channel_statistics[i].sum_cubed/=area;
+ channel_statistics[i].sum_fourth_power/=area;
+ channel_statistics[i].mean=channel_statistics[i].sum;
+ channel_statistics[i].variance=channel_statistics[i].sum_squared;
+ channel_statistics[i].standard_deviation=sqrt(
+ channel_statistics[i].variance-(channel_statistics[i].mean*
+ channel_statistics[i].mean));
}
for (i=0; i < AllChannels; i++)
{
- channel_statistics[AllChannels].depth=(unsigned long) MagickMax((double)
+ channel_statistics[AllChannels].depth=(size_t) MagickMax((double)
channel_statistics[AllChannels].depth,(double)
channel_statistics[i].depth);
channel_statistics[AllChannels].minima=MagickMin(
channel_statistics[AllChannels].minima,channel_statistics[i].minima);
channel_statistics[AllChannels].maxima=MagickMax(
channel_statistics[AllChannels].maxima,channel_statistics[i].maxima);
+ channel_statistics[AllChannels].sum+=channel_statistics[i].sum;
+ channel_statistics[AllChannels].sum_squared+=
+ channel_statistics[i].sum_squared;
+ channel_statistics[AllChannels].sum_cubed+=channel_statistics[i].sum_cubed;
+ channel_statistics[AllChannels].sum_fourth_power+=
+ channel_statistics[i].sum_fourth_power;
channel_statistics[AllChannels].mean+=channel_statistics[i].mean;
+ channel_statistics[AllChannels].variance+=channel_statistics[i].variance-
+ channel_statistics[i].mean*channel_statistics[i].mean;
channel_statistics[AllChannels].standard_deviation+=
- channel_statistics[i].standard_deviation;
- channel_statistics[AllChannels].kurtosis+=channel_statistics[i].kurtosis;
- channel_statistics[AllChannels].skewness+=channel_statistics[i].skewness;
+ channel_statistics[i].variance-channel_statistics[i].mean*
+ channel_statistics[i].mean;
}
- channels=4;
+ channels=3;
+ if (image->matte != MagickFalse)
+ channels++;
if (image->colorspace == CMYKColorspace)
channels++;
+ channel_statistics[AllChannels].sum/=channels;
+ channel_statistics[AllChannels].sum_squared/=channels;
+ channel_statistics[AllChannels].sum_cubed/=channels;
+ channel_statistics[AllChannels].sum_fourth_power/=channels;
channel_statistics[AllChannels].mean/=channels;
- channel_statistics[AllChannels].standard_deviation/=channels;
+ channel_statistics[AllChannels].variance/=channels;
+ channel_statistics[AllChannels].standard_deviation=
+ sqrt(channel_statistics[AllChannels].standard_deviation/channels);
channel_statistics[AllChannels].kurtosis/=channels;
channel_statistics[AllChannels].skewness/=channels;
for (i=0; i <= AllChannels; i++)
{
- sum_squares=0.0;
- sum_squares=channel_statistics[i].standard_deviation;
- sum_cubes=0.0;
- sum_cubes=channel_statistics[i].skewness;
- channel_statistics[i].standard_deviation=sqrt(
- channel_statistics[i].standard_deviation-
- (channel_statistics[i].mean*channel_statistics[i].mean));
if (channel_statistics[i].standard_deviation == 0.0)
- {
- channel_statistics[i].kurtosis=0.0;
- channel_statistics[i].skewness=0.0;
- }
- else
- {
- channel_statistics[i].skewness=(channel_statistics[i].skewness-
- 3.0*channel_statistics[i].mean*sum_squares+
- 2.0*channel_statistics[i].mean*channel_statistics[i].mean*
- channel_statistics[i].mean)/
- (channel_statistics[i].standard_deviation*
- channel_statistics[i].standard_deviation*
- channel_statistics[i].standard_deviation);
- channel_statistics[i].kurtosis=(channel_statistics[i].kurtosis-
- 4.0*channel_statistics[i].mean*sum_cubes+
- 6.0*channel_statistics[i].mean*channel_statistics[i].mean*sum_squares-
- 3.0*channel_statistics[i].mean*channel_statistics[i].mean*
- 1.0*channel_statistics[i].mean*channel_statistics[i].mean)/
- (channel_statistics[i].standard_deviation*
- channel_statistics[i].standard_deviation*
- channel_statistics[i].standard_deviation*
- channel_statistics[i].standard_deviation)-3.0;
- }
+ continue;
+ channel_statistics[i].skewness=(channel_statistics[i].sum_cubed-
+ 3.0*channel_statistics[i].mean*channel_statistics[i].sum_squared+
+ 2.0*channel_statistics[i].mean*channel_statistics[i].mean*
+ channel_statistics[i].mean)/(channel_statistics[i].standard_deviation*
+ channel_statistics[i].standard_deviation*
+ channel_statistics[i].standard_deviation);
+ channel_statistics[i].kurtosis=(channel_statistics[i].sum_fourth_power-
+ 4.0*channel_statistics[i].mean*channel_statistics[i].sum_cubed+
+ 6.0*channel_statistics[i].mean*channel_statistics[i].mean*
+ channel_statistics[i].sum_squared-3.0*channel_statistics[i].mean*
+ channel_statistics[i].mean*1.0*channel_statistics[i].mean*
+ channel_statistics[i].mean)/(channel_statistics[i].standard_deviation*
+ channel_statistics[i].standard_deviation*
+ channel_statistics[i].standard_deviation*
+ channel_statistics[i].standard_deviation)-3.0;
}
return(channel_statistics);
}