/*

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%                                                                             %

%                                                                             %

%                                                                             %

%    M   M    OOO    RRRR   PPPP   H   H   OOO   L       OOO    GGGG  Y   Y   %

%    MM MM   O   O   R   R  P   P  H   H  O   O  L      O   O  G       Y Y    %

%    M M M   O   O   RRRR   PPPP   HHHHH  O   O  L      O   O  G GGG    Y     %

%    M   M   O   O   R R    P      H   H  O   O  L      O   O  G   G    Y     %

%    M   M    OOO    R  R   P      H   H   OOO   LLLLL   OOO    GGG     Y     %

%                                                                             %

%                                                                             %

%                        MagickCore Morphology Methods                        %

%                                                                             %

%                              Software Design                                %

%                              Anthony Thyssen                                %

%                               January 2010                                  %

%                                                                             %

%                                                                             %

%  Copyright @ 1999 ImageMagick Studio LLC, a non-profit organization         %

%  dedicated to making software imaging solutions freely available.           %

%                                                                             %

%  You may not use this file except in compliance with the License.  You may  %

%  obtain a copy of the License at                                            %

%                                                                             %

%    https://imagemagick.org/script/license.php                               %

%                                                                             %

%  Unless required by applicable law or agreed to in writing, software        %

%  distributed under the License is distributed on an "AS IS" BASIS,          %

%  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.   %

%  See the License for the specific language governing permissions and        %

%  limitations under the License.                                             %

%                                                                             %

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%

% Morphology is the application of various kernels, of any size or shape, to an

% image in various ways (typically binary, but not always).

%

% Convolution (weighted sum or average) is just one specific type of

% morphology. Just one that is very common for image blurring and sharpening

% effects.  Not only 2D Gaussian blurring, but also 2-pass 1D Blurring.

%

% This module provides not only a general morphology function, and the ability

% to apply more advanced or iterative morphologies, but also functions for the

% generation of many different types of kernel arrays from user supplied

% arguments. Prehaps even the generation of a kernel from a small image.

*/

/*

  Include declarations.

*/

#include "MagickCore/studio.h"

#include "MagickCore/artifact.h"

#include "MagickCore/cache-view.h"

#include "MagickCore/channel.h"

#include "MagickCore/color-private.h"

#include "MagickCore/enhance.h"

#include "MagickCore/exception.h"

#include "MagickCore/exception-private.h"

#include "MagickCore/gem.h"

#include "MagickCore/gem-private.h"

#include "MagickCore/image.h"

#include "MagickCore/image-private.h"

#include "MagickCore/linked-list.h"

#include "MagickCore/list.h"

#include "MagickCore/magick.h"

#include "MagickCore/memory_.h"

#include "MagickCore/memory-private.h"

#include "MagickCore/monitor-private.h"

#include "MagickCore/morphology.h"

#include "MagickCore/morphology-private.h"

#include "MagickCore/option.h"

#include "MagickCore/pixel-accessor.h"

#include "MagickCore/prepress.h"

#include "MagickCore/quantize.h"

#include "MagickCore/resource_.h"

#include "MagickCore/registry.h"

#include "MagickCore/semaphore.h"

#include "MagickCore/splay-tree.h"

#include "MagickCore/statistic.h"

#include "MagickCore/string_.h"

#include "MagickCore/string-private.h"

#include "MagickCore/thread-private.h"

#include "MagickCore/token.h"

#include "MagickCore/utility.h"

#include "MagickCore/utility-private.h"

/*

  Other global definitions used by module.

*/

#define Minimize(assign,value) assign=MagickMin(assign,value)

#define Maximize(assign,value) assign=MagickMax(assign,value)

/* Integer Factorial Function - for a Binomial kernel */

#if 1

static inline size_t fact(size_t n)

{

  size_t f,l;

  for(f=1, l=2; l <= n; f=f*l, l++);

  return(f);

}

#elif 1 /* glibc floating point alternatives */

#define fact(n) ((size_t)tgamma((double)n+1))

#else

#define fact(n) ((size_t)lgamma((double)n+1))

#endif

/* Currently these are only internal to this module */

static void

  CalcKernelMetaData(KernelInfo *),

  ExpandMirrorKernelInfo(KernelInfo *),

  ExpandRotateKernelInfo(KernelInfo *, const double),

  RotateKernelInfo(KernelInfo *, double);

/* Quick function to find last kernel in a kernel list */

static inline KernelInfo *LastKernelInfo(KernelInfo *kernel)

{

  while (kernel->next != (KernelInfo *) NULL)

    kernel=kernel->next;

  return(kernel);

}

/*

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%                                                                             %

%                                                                             %

%                                                                             %

%     A c q u i r e K e r n e l I n f o                                       %

%                                                                             %

%                                                                             %

%                                                                             %

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%

%  AcquireKernelInfo() takes the given string (generally supplied by the

%  user) and converts it into a Morphology/Convolution Kernel.  This allows

%  users to specify a kernel from a number of pre-defined kernels, or to fully

%  specify their own kernel for a specific Convolution or Morphology

%  Operation.

%

%  The kernel so generated can be any rectangular array of floating point

%  values (doubles) with the 'control point' or 'pixel being affected'

%  anywhere within that array of values.

%

%  Previously IM was restricted to a square of odd size using the exact

%  center as origin, this is no longer the case, and any rectangular kernel

%  with any value being declared the origin. This in turn allows the use of

%  highly asymmetrical kernels.

%

%  The floating point values in the kernel can also include a special value

%  known as 'nan' or 'not a number' to indicate that this value is not part

%  of the kernel array. This allows you to shaped the kernel within its

%  rectangular area. That is 'nan' values provide a 'mask' for the kernel

%  shape.  However at least one non-nan value must be provided for correct

%  working of a kernel.

%

%  The returned kernel should be freed using the DestroyKernelInfo() when you

%  are finished with it.  Do not free this memory yourself.

%

%  Input kernel definition strings can consist of any of three types.

%

%    "name:args[[@><]"

%         Select from one of the built in kernels, using the name and

%         geometry arguments supplied.  See AcquireKernelBuiltIn()

%

%    "WxH[+X+Y][@><]:num, num, num ..."

%         a kernel of size W by H, with W*H floating point numbers following.

%         the 'center' can be optionally be defined at +X+Y (such that +0+0

%         is top left corner). If not defined the pixel in the center, for

%         odd sizes, or to the immediate top or left of center for even sizes

%         is automatically selected.

%

%    "num, num, num, num, ..."

%         list of floating point numbers defining an 'old style' odd sized

%         square kernel.  At least 9 values should be provided for a 3x3

%         square kernel, 25 for a 5x5 square kernel, 49 for 7x7, etc.

%         Values can be space or comma separated.  This is not recommended.

%

%  You can define a 'list of kernels' which can be used by some morphology

%  operators A list is defined as a semi-colon separated list kernels.

%

%     " kernel ; kernel ; kernel ; "

%

%  Any extra ';' characters, at start, end or between kernel definitions are

%  simply ignored.

%

%  The special flags will expand a single kernel, into a list of rotated

%  kernels. A '@' flag will expand a 3x3 kernel into a list of 45-degree

%  cyclic rotations, while a '>' will generate a list of 90-degree rotations.

%  The '<' also expands using 90-degree rotates, but giving a 180-degree

%  reflected kernel before the +/- 90-degree rotations, which can be important

%  for Thinning operations.

%

%  Note that 'name' kernels will start with an alphabetic character while the

%  new kernel specification has a ':' character in its specification string.

%  If neither is the case, it is assumed an old style of a simple list of

%  numbers generating a odd-sized square kernel has been given.

%

%  The format of the AcquireKernel method is:

%

%      KernelInfo *AcquireKernelInfo(const char *kernel_string)

%

%  A description of each parameter follows:

%

%    o kernel_string: the Morphology/Convolution kernel wanted.

%

*/

/* This was separated so that it could be used as a separate

** array input handling function, such as for -color-matrix

*/

static KernelInfo *ParseKernelArray(const char *kernel_string)

{

  KernelInfo

    *kernel;

  char

    token[MagickPathExtent];

  const char

    *p,

    *end;

  ssize_t

    i;

  double

    nan = sqrt(-1.0);  /* Special Value : Not A Number */

  MagickStatusType

    flags;

  GeometryInfo

    args;

  kernel=(KernelInfo *) AcquireMagickMemory(sizeof(*kernel));

  if (kernel == (KernelInfo *) NULL)

    return(kernel);

  (void) memset(kernel,0,sizeof(*kernel));

  kernel->minimum = kernel->maximum = kernel->angle = 0.0;

  kernel->negative_range = kernel->positive_range = 0.0;

  kernel->type = UserDefinedKernel;

  kernel->next = (KernelInfo *) NULL;

  kernel->signature=MagickCoreSignature;

  if (kernel_string == (const char *) NULL)

    return(kernel);

  /* find end of this specific kernel definition string */

  end = strchr(kernel_string, ';');

  if ( end == (char *) NULL )

    end = strchr(kernel_string, '\0');

  /* clear flags - for Expanding kernel lists through rotations */

   flags = NoValue;

  /* Has a ':' in argument - New user kernel specification

     FUTURE: this split on ':' could be done by StringToken()

   */

  p = strchr(kernel_string, ':');

  if ( p != (char *) NULL && p < end)

    {

      /* ParseGeometry() needs the geometry separated! -- Arrgghh */

      (void) memcpy(token, kernel_string, (size_t) (p-kernel_string));

      token[p-kernel_string] = '\0';

      SetGeometryInfo(&args);

      flags = ParseGeometry(token, &args);

      /* Size handling and checks of geometry settings */

      if ( (flags & WidthValue) == 0 ) /* if no width then */

        args.rho = args.sigma;         /* then  width = height */

      if ( args.rho < 1.0 )            /* if width too small */

         args.rho = 1.0;               /* then  width = 1 */

      if ( args.sigma < 1.0 )          /* if height too small */

        args.sigma = args.rho;         /* then  height = width */

      kernel->width = (size_t)args.rho;

      kernel->height = (size_t)args.sigma;

      /* Offset Handling and Checks */

      if ( args.xi  < 0.0 || args.psi < 0.0 )

        return(DestroyKernelInfo(kernel));

      kernel->x = ((flags & XValue)!=0) ? (ssize_t)args.xi

                                        : (ssize_t) (kernel->width-1)/2;

      kernel->y = ((flags & YValue)!=0) ? (ssize_t)args.psi

                                        : (ssize_t) (kernel->height-1)/2;

      if ( kernel->x >= (ssize_t) kernel->width ||

           kernel->y >= (ssize_t) kernel->height )

        return(DestroyKernelInfo(kernel));

      p++; /* advance beyond the ':' */

    }

  else

    { /* ELSE - Old old specification, forming odd-square kernel */

      /* count up number of values given */

      p=(const char *) kernel_string;

      while ((isspace((int) ((unsigned char) *p)) != 0) || (*p == '\''))

        p++;  /* ignore "'" chars for convolve filter usage - Cristy */

      for (i=0; p < end; i++)

      {

        (void) GetNextToken(p,&p,MagickPathExtent,token);

        if (*token == ',')

          (void) GetNextToken(p,&p,MagickPathExtent,token);

      }

      /* set the size of the kernel - old sized square */

      kernel->width = kernel->height= (size_t) sqrt((double) i+1.0);

      kernel->x = kernel->y = (ssize_t) (kernel->width-1)/2;

      p=(const char *) kernel_string;

      while ((isspace((int) ((unsigned char) *p)) != 0) || (*p == '\''))

        p++;  /* ignore "'" chars for convolve filter usage - Cristy */

    }

  /* Read in the kernel values from rest of input string argument */

  kernel->values=(MagickRealType *) MagickAssumeAligned(AcquireAlignedMemory(

    kernel->width,kernel->height*sizeof(*kernel->values)));

  if (kernel->values == (MagickRealType *) NULL)

    return(DestroyKernelInfo(kernel));

  kernel->minimum=MagickMaximumValue;

  kernel->maximum=(-MagickMaximumValue);

  kernel->negative_range = kernel->positive_range = 0.0;

  for (i=0; (i < (ssize_t) (kernel->width*kernel->height)) && (p < end); i++)

  {

    (void) GetNextToken(p,&p,MagickPathExtent,token);

    if (*token == ',')

      (void) GetNextToken(p,&p,MagickPathExtent,token);

    if (    LocaleCompare("nan",token) == 0

        || LocaleCompare("-",token) == 0 ) {

      kernel->values[i] = nan; /* this value is not part of neighbourhood */

    }

    else {

      kernel->values[i] = StringToDouble(token,(char **) NULL);

      ( kernel->values[i] < 0)

          ?  ( kernel->negative_range += kernel->values[i] )

          :  ( kernel->positive_range += kernel->values[i] );

      Minimize(kernel->minimum, kernel->values[i]);

      Maximize(kernel->maximum, kernel->values[i]);

    }

  }

  /* sanity check -- no more values in kernel definition */

  (void) GetNextToken(p,&p,MagickPathExtent,token);

  if ( *token != '\0' && *token != ';' && *token != '\'' )

    return(DestroyKernelInfo(kernel));

#if 0

  /* this was the old method of handling a incomplete kernel */

  if ( i < (ssize_t) (kernel->width*kernel->height) ) {

    Minimize(kernel->minimum, kernel->values[i]);

    Maximize(kernel->maximum, kernel->values[i]);

    for ( ; i < (ssize_t) (kernel->width*kernel->height); i++)

      kernel->values[i]=0.0;

  }

#else

  /* Number of values for kernel was not enough - Report Error */

  if ( i < (ssize_t) (kernel->width*kernel->height) )

    return(DestroyKernelInfo(kernel));

#endif

  /* check that we received at least one real (non-nan) value! */

  if (kernel->minimum == MagickMaximumValue)

    return(DestroyKernelInfo(kernel));

  if ( (flags & AreaValue) != 0 )         /* '@' symbol in kernel size */

    ExpandRotateKernelInfo(kernel, 45.0); /* cyclic rotate 3x3 kernels */

  else if ( (flags & GreaterValue) != 0 ) /* '>' symbol in kernel args */

    ExpandRotateKernelInfo(kernel, 90.0); /* 90 degree rotate of kernel */

  else if ( (flags & LessValue) != 0 )    /* '<' symbol in kernel args */

    ExpandMirrorKernelInfo(kernel);       /* 90 degree mirror rotate */

  return(kernel);

}

static KernelInfo *ParseKernelName(const char *kernel_string,

  ExceptionInfo *exception)

{

  char

    token[MagickPathExtent] = "";

  const char

    *p,

    *end;

  GeometryInfo

    args;

  KernelInfo

    *kernel;

  MagickStatusType

    flags;

  ssize_t

    type;

  /* Parse special 'named' kernel */

  (void) GetNextToken(kernel_string,&p,MagickPathExtent,token);

  type=ParseCommandOption(MagickKernelOptions,MagickFalse,token);

  if ( type < 0 || type == UserDefinedKernel )

    return((KernelInfo *) NULL);  /* not a valid named kernel */

  while (((isspace((int) ((unsigned char) *p)) != 0) ||

          (*p == ',') || (*p == ':' )) && (*p != '\0') && (*p != ';'))

    p++;

  end = strchr(p, ';'); /* end of this kernel definition */

  if ( end == (char *) NULL )

    end = strchr(p, '\0');

  /* ParseGeometry() needs the geometry separated! -- Arrgghh */

  (void) memcpy(token, p, (size_t) (end-p));

  token[end-p] = '\0';

  SetGeometryInfo(&args);

  flags = ParseGeometry(token, &args);

#if 0

  /* For Debugging Geometry Input */

  (void) FormatLocaleFile(stderr, "Geometry = 0x%04X : %lg x %lg %+lg %+lg\n",

    flags, args.rho, args.sigma, args.xi, args.psi );

#endif

  /* special handling of missing values in input string */

  switch( type ) {

    /* Shape Kernel Defaults */

    case UnityKernel:

      if ( (flags & WidthValue) == 0 )

        args.rho = 1.0;    /* Default scale = 1.0, zero is valid */

      break;

    case SquareKernel:

    case DiamondKernel:

    case OctagonKernel:

    case DiskKernel:

    case PlusKernel:

    case CrossKernel:

      if ( (flags & HeightValue) == 0 )

        args.sigma = 1.0;    /* Default scale = 1.0, zero is valid */

      break;

    case RingKernel:

      if ( (flags & XValue) == 0 )

        args.xi = 1.0;       /* Default scale = 1.0, zero is valid */

      break;

    case RectangleKernel:    /* Rectangle - set size defaults */

      if ( (flags & WidthValue) == 0 ) /* if no width then */

        args.rho = args.sigma;         /* then  width = height */

      if ( args.rho < 1.0 )            /* if width too small */

          args.rho = 3;                /* then  width = 3 */

      if ( args.sigma < 1.0 )          /* if height too small */

        args.sigma = args.rho;         /* then  height = width */

      if ( (flags & XValue) == 0 )     /* center offset if not defined */

        args.xi = (double)(((ssize_t)args.rho-1)/2);

      if ( (flags & YValue) == 0 )

        args.psi = (double)(((ssize_t)args.sigma-1)/2);

      break;

    /* Distance Kernel Defaults */

    case ChebyshevKernel:

    case ManhattanKernel:

    case OctagonalKernel:

    case EuclideanKernel:

      if ( (flags & HeightValue) == 0 )           /* no distance scale */

        args.sigma = 100.0;                       /* default distance scaling */

      else if ( (flags & AspectValue ) != 0 )     /* '!' flag */

        args.sigma = (double) QuantumRange/(args.sigma+1); /* maximum pixel distance */

      else if ( (flags & PercentValue ) != 0 )    /* '%' flag */

        args.sigma *= (double) QuantumRange/100.0;         /* percentage of color range */

      break;

    default:

      break;

  }

  kernel = AcquireKernelBuiltIn((KernelInfoType)type, &args, exception);

  if ( kernel == (KernelInfo *) NULL )

    return(kernel);

  /* global expand to rotated kernel list - only for single kernels */

  if ( kernel->next == (KernelInfo *) NULL ) {

    if ( (flags & AreaValue) != 0 )         /* '@' symbol in kernel args */

      ExpandRotateKernelInfo(kernel, 45.0);

    else if ( (flags & GreaterValue) != 0 ) /* '>' symbol in kernel args */

      ExpandRotateKernelInfo(kernel, 90.0);

    else if ( (flags & LessValue) != 0 )    /* '<' symbol in kernel args */

      ExpandMirrorKernelInfo(kernel);

  }

  return(kernel);

}

MagickExport KernelInfo *AcquireKernelInfo(const char *kernel_string,

  ExceptionInfo *exception)

{

  KernelInfo

    *kernel,

    *new_kernel;

  char

    *kernel_cache,

    token[MagickPathExtent];

  const char

    *p;

  if (kernel_string == (const char *) NULL)

    return(ParseKernelArray(kernel_string));

  p=kernel_string;

  kernel_cache=(char *) NULL;

  if (*kernel_string == '@')

    {

      kernel_cache=FileToString(kernel_string,~0UL,exception);

      if (kernel_cache == (char *) NULL)

        return((KernelInfo *) NULL);

      p=(const char *) kernel_cache;

    }

  kernel=NULL;

  while (GetNextToken(p,(const char **) NULL,MagickPathExtent,token), *token != '\0')

  {

    /* ignore extra or multiple ';' kernel separators */

    if (*token != ';')

      {

        /* tokens starting with alpha is a Named kernel */

        if (isalpha((int) ((unsigned char) *token)) != 0)

          new_kernel=ParseKernelName(p,exception);

        else /* otherwise a user defined kernel array */

          new_kernel=ParseKernelArray(p);

        /* Error handling -- this is not proper error handling! */

        if (new_kernel == (KernelInfo *) NULL)

          {

            if (kernel != (KernelInfo *) NULL)

              kernel=DestroyKernelInfo(kernel);

            return((KernelInfo *) NULL);

          }

        /* initialise or append the kernel list */

        if (kernel == (KernelInfo *) NULL)

          kernel=new_kernel;

        else

          LastKernelInfo(kernel)->next=new_kernel;

      }

    /* look for the next kernel in list */

    p=strchr(p,';');

    if (p == (char *) NULL)

      break;

    p++;

  }

  if (kernel_cache != (char *) NULL)

    kernel_cache=DestroyString(kernel_cache);

  return(kernel);

}

/*

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%                                                                             %

%                                                                             %

%                                                                             %

%     A c q u i r e K e r n e l B u i l t I n                                 %

%                                                                             %

%                                                                             %

%                                                                             %

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%

%  AcquireKernelBuiltIn() returned one of the 'named' built-in types of

%  kernels used for special purposes such as gaussian blurring, skeleton

%  pruning, and edge distance determination.

%

%  They take a KernelType, and a set of geometry style arguments, which were

%  typically decoded from a user supplied string, or from a more complex

%  Morphology Method that was requested.

%

%  The format of the AcquireKernelBuiltIn method is:

%

%      KernelInfo *AcquireKernelBuiltIn(const KernelInfoType type,

%           const GeometryInfo args)

%

%  A description of each parameter follows:

%

%    o type: the pre-defined type of kernel wanted

%

%    o args: arguments defining or modifying the kernel

%

%  Convolution Kernels

%

%    Unity

%       The a No-Op or Scaling single element kernel.

%

%    Gaussian:{radius},{sigma}

%       Generate a two-dimensional gaussian kernel, as used by -gaussian.

%       The sigma for the curve is required.  The resulting kernel is

%       normalized,

%

%       If 'sigma' is zero, you get a single pixel on a field of zeros.

%

%       NOTE: that the 'radius' is optional, but if provided can limit (clip)

%       the final size of the resulting kernel to a square 2*radius+1 in size.

%       The radius should be at least 2 times that of the sigma value, or

%       sever clipping and aliasing may result.  If not given or set to 0 the

%       radius will be determined so as to produce the best minimal error

%       result, which is usually much larger than is normally needed.

%

%    LoG:{radius},{sigma}

%        "Laplacian of a Gaussian" or "Mexican Hat" Kernel.

%        The supposed ideal edge detection, zero-summing kernel.

%

%        An alternative to this kernel is to use a "DoG" with a sigma ratio of

%        approx 1.6 (according to wikipedia).

%

%    DoG:{radius},{sigma1},{sigma2}

%        "Difference of Gaussians" Kernel.

%        As "Gaussian" but with a gaussian produced by 'sigma2' subtracted

%        from the gaussian produced by 'sigma1'. Typically sigma2 > sigma1.

%        The result is a zero-summing kernel.

%

%    Blur:{radius},{sigma}[,{angle}]

%       Generates a 1 dimensional or linear gaussian blur, at the angle given

%       (current restricted to orthogonal angles).  If a 'radius' is given the

%       kernel is clipped to a width of 2*radius+1.  Kernel can be rotated

%       by a 90 degree angle.

%

%       If 'sigma' is zero, you get a single pixel on a field of zeros.

%

%       Note that two convolutions with two "Blur" kernels perpendicular to

%       each other, is equivalent to a far larger "Gaussian" kernel with the

%       same sigma value, However it is much faster to apply. This is how the

%       "-blur" operator actually works.

%

%    Comet:{width},{sigma},{angle}

%       Blur in one direction only, much like how a bright object leaves

%       a comet like trail.  The Kernel is actually half a gaussian curve,

%       Adding two such blurs in opposite directions produces a Blur Kernel.

%       Angle can be rotated in multiples of 90 degrees.

%

%       Note that the first argument is the width of the kernel and not the

%       radius of the kernel.

%

%    Binomial:[{radius}]

%       Generate a discrete kernel using a 2 dimensional Pascal's Triangle

%       of values. Used for special forma of image filters.

%

%    # Still to be implemented...

%    #

%    # Filter2D

%    # Filter1D

%    #    Set kernel values using a resize filter, and given scale (sigma)

%    #    Cylindrical or Linear.   Is this possible with an image?

%    #

%

%  Named Constant Convolution Kernels

%

%  All these are unscaled, zero-summing kernels by default. As such for

%  non-HDRI version of ImageMagick some form of normalization, user scaling,

%  and biasing the results is recommended, to prevent the resulting image

%  being 'clipped'.

%

%  The 3x3 kernels (most of these) can be circularly rotated in multiples of

%  45 degrees to generate the 8 angled variants of each of the kernels.

%

%    Laplacian:{type}

%      Discrete Laplacian Kernels, (without normalization)

%        Type 0 :  3x3 with center:8 surrounded by -1  (8 neighbourhood)

%        Type 1 :  3x3 with center:4 edge:-1 corner:0 (4 neighbourhood)

%        Type 2 :  3x3 with center:4 edge:1 corner:-2

%        Type 3 :  3x3 with center:4 edge:-2 corner:1

%        Type 5 :  5x5 laplacian

%        Type 7 :  7x7 laplacian

%        Type 15 : 5x5 LoG (sigma approx 1.4)

%        Type 19 : 9x9 LoG (sigma approx 1.4)

%

%    Sobel:{angle}

%      Sobel 'Edge' convolution kernel (3x3)

%          | -1, 0, 1 |

%          | -2, 0,-2 |

%          | -1, 0, 1 |

%

%    Roberts:{angle}

%      Roberts convolution kernel (3x3)

%          |  0, 0, 0 |

%          | -1, 1, 0 |

%          |  0, 0, 0 |

%

%    Prewitt:{angle}

%      Prewitt Edge convolution kernel (3x3)

%          | -1, 0, 1 |

%          | -1, 0, 1 |

%          | -1, 0, 1 |

%

%    Compass:{angle}

%      Prewitt's "Compass" convolution kernel (3x3)

%          | -1, 1, 1 |

%          | -1,-2, 1 |

%          | -1, 1, 1 |

%

%    Kirsch:{angle}

%      Kirsch's "Compass" convolution kernel (3x3)

%          | -3,-3, 5 |

%          | -3, 0, 5 |

%          | -3,-3, 5 |

%

%    FreiChen:{angle}

%      Frei-Chen Edge Detector is based on a kernel that is similar to

%      the Sobel Kernel, but is designed to be isotropic. That is it takes

%      into account the distance of the diagonal in the kernel.

%

%          |   1,     0,   -1     |

%          | sqrt(2), 0, -sqrt(2) |

%          |   1,     0,   -1     |

%

%    FreiChen:{type},{angle}

%

%      Frei-Chen Pre-weighted kernels...

%

%        Type 0:  default un-normalized version shown above.

%

%        Type 1: Orthogonal Kernel (same as type 11 below)

%          |   1,     0,   -1     |

%          | sqrt(2), 0, -sqrt(2) | / 2*sqrt(2)

%          |   1,     0,   -1     |

%

%        Type 2: Diagonal form of Kernel...

%          |   1,     sqrt(2),    0     |

%          | sqrt(2),   0,     -sqrt(2) | / 2*sqrt(2)

%          |   0,    -sqrt(2)    -1     |

%

%      However this kernel is als at the heart of the FreiChen Edge Detection

%      Process which uses a set of 9 specially weighted kernel.  These 9

%      kernels not be normalized, but directly applied to the image. The

%      results is then added together, to produce the intensity of an edge in

%      a specific direction.  The square root of the pixel value can then be

%      taken as the cosine of the edge, and at least 2 such runs at 90 degrees

%      from each other, both the direction and the strength of the edge can be

%      determined.

%

%        Type 10: All 9 of the following pre-weighted kernels...

%

%        Type 11: |   1,     0,   -1     |

%                 | sqrt(2), 0, -sqrt(2) | / 2*sqrt(2)

%                 |   1,     0,   -1     |

%

%        Type 12: | 1, sqrt(2), 1 |

%                 | 0,   0,     0 | / 2*sqrt(2)

%                 | 1, sqrt(2), 1 |

%

%        Type 13: | sqrt(2), -1,    0     |

%                 |  -1,      0,    1     | / 2*sqrt(2)

%                 |   0,      1, -sqrt(2) |

%

%        Type 14: |    0,     1, -sqrt(2) |

%                 |   -1,     0,     1    | / 2*sqrt(2)

%                 | sqrt(2), -1,     0    |

%

%        Type 15: | 0, -1, 0 |

%                 | 1,  0, 1 | / 2

%                 | 0, -1, 0 |

%

%        Type 16: |  1, 0, -1 |

%                 |  0, 0,  0 | / 2

%                 | -1, 0,  1 |

%

%        Type 17: |  1, -2,  1 |

%                 | -2,  4, -2 | / 6

%                 | -1, -2,  1 |

%

%        Type 18: | -2, 1, -2 |

%                 |  1, 4,  1 | / 6

%                 | -2, 1, -2 |

%

%        Type 19: | 1, 1, 1 |

%                 | 1, 1, 1 | / 3

%                 | 1, 1, 1 |

%

%      The first 4 are for edge detection, the next 4 are for line detection

%      and the last is to add a average component to the results.

%

%      Using a special type of '-1' will return all 9 pre-weighted kernels

%      as a multi-kernel list, so that you can use them directly (without

%      normalization) with the special "-set option:morphology:compose Plus"

%      setting to apply the full FreiChen Edge Detection Technique.

%

%      If 'type' is large it will be taken to be an actual rotation angle for

%      the default FreiChen (type 0) kernel.  As such  FreiChen:45  will look

%      like a  Sobel:45  but with 'sqrt(2)' instead of '2' values.

%

%      WARNING: The above was layed out as per

%          http://www.math.tau.ac.il/~turkel/notes/edge_detectors.pdf

%      But rotated 90 degrees so direction is from left rather than the top.

%      I have yet to find any secondary confirmation of the above. The only

%      other source found was actual source code at

%          http://ltswww.epfl.ch/~courstiv/exos_labos/sol3.pdf

%      Neither paper defines the kernels in a way that looks logical or

%      correct when taken as a whole.

%

%  Boolean Kernels

%

%    Diamond:[{radius}[,{scale}]]

%       Generate a diamond shaped kernel with given radius to the points.

%       Kernel size will again be radius*2+1 square and defaults to radius 1,

%       generating a 3x3 kernel that is slightly larger than a square.

%

%    Square:[{radius}[,{scale}]]

%       Generate a square shaped kernel of size radius*2+1, and defaulting

%       to a 3x3 (radius 1).

%

%    Octagon:[{radius}[,{scale}]]

%       Generate octagonal shaped kernel of given radius and constant scale.

%       Default radius is 3 producing a 7x7 kernel. A radius of 1 will result

%       in "Diamond" kernel.

%

%    Disk:[{radius}[,{scale}]]

%       Generate a binary disk, thresholded at the radius given, the radius

%       may be a float-point value. Final Kernel size is floor(radius)*2+1

%       square. A radius of 5.3 is the default.

%

%       NOTE: That a low radii Disk kernels produce the same results as

%       many of the previously defined kernels, but differ greatly at larger

%       radii.  Here is a table of equivalences...

%          "Disk:1"    => "Diamond", "Octagon:1", or "Cross:1"

%          "Disk:1.5"  => "Square"

%          "Disk:2"    => "Diamond:2"

%          "Disk:2.5"  => "Octagon"

%          "Disk:2.9"  => "Square:2"

%          "Disk:3.5"  => "Octagon:3"

%          "Disk:4.5"  => "Octagon:4"

%          "Disk:5.4"  => "Octagon:5"

%          "Disk:6.4"  => "Octagon:6"

%       All other Disk shapes are unique to this kernel, but because a "Disk"

%       is more circular when using a larger radius, using a larger radius is

%       preferred over iterating the morphological operation.

%

%    Rectangle:{geometry}

%       Simply generate a rectangle of 1's with the size given. You can also

%       specify the location of the 'control point', otherwise the closest

%       pixel to the center of the rectangle is selected.

%

%       Properly centered and odd sized rectangles work the best.

%

%  Symbol Dilation Kernels

%

%    These kernel is not a good general morphological kernel, but is used

%    more for highlighting and marking any single pixels in an image using,

%    a "Dilate" method as appropriate.

%

%    For the same reasons iterating these kernels does not produce the

%    same result as using a larger radius for the symbol.

%

%    Plus:[{radius}[,{scale}]]

%    Cross:[{radius}[,{scale}]]

%       Generate a kernel in the shape of a 'plus' or a 'cross' with

%       a each arm the length of the given radius (default 2).

%

%       NOTE: "plus:1" is equivalent to a "Diamond" kernel.

%

%    Ring:{radius1},{radius2}[,{scale}]

%       A ring of the values given that falls between the two radii.

%       Defaults to a ring of approximately 3 radius in a 7x7 kernel.

%       This is the 'edge' pixels of the default "Disk" kernel,

%       More specifically, "Ring" -> "Ring:2.5,3.5,1.0"

%

%  Hit and Miss Kernels

%

%    Peak:radius1,radius2

%       Find any peak larger than the pixels the fall between the two radii.

%       The default ring of pixels is as per "Ring".

%    Edges

%       Find flat orthogonal edges of a binary shape

%    Corners

%       Find 90 degree corners of a binary shape

%    Diagonals:type

%       A special kernel to thin the 'outside' of diagonals

%    LineEnds:type

%       Find end points of lines (for pruning a skeleton)

%       Two types of lines ends (default to both) can be searched for

%         Type 0: All line ends

%         Type 1: single kernel for 4-connected line ends

%         Type 2: single kernel for simple line ends

%    LineJunctions

%       Find three line junctions (within a skeleton)

%         Type 0: all line junctions

%         Type 1: Y Junction kernel

%         Type 2: Diagonal T Junction kernel

%         Type 3: Orthogonal T Junction kernel

%         Type 4: Diagonal X Junction kernel

%         Type 5: Orthogonal + Junction kernel

%    Ridges:type

%       Find single pixel ridges or thin lines

%         Type 1: Fine single pixel thick lines and ridges

%         Type 2: Find two pixel thick lines and ridges

%    ConvexHull

%       Octagonal Thickening Kernel, to generate convex hulls of 45 degrees

%    Skeleton:type

%       Traditional skeleton generating kernels.

%         Type 1: Traditional Skeleton kernel (4 connected skeleton)

%         Type 2: HIPR2 Skeleton kernel (8 connected skeleton)

%         Type 3: Thinning skeleton based on a research paper by

%                 Dan S. Bloomberg (Default Type)

%    ThinSE:type

%       A huge variety of Thinning Kernels designed to preserve connectivity.

%       many other kernel sets use these kernels as source definitions.

%       Type numbers are 41-49, 81-89, 481, and 482 which are based on

%       the super and sub notations used in the source research paper.

%

%  Distance Measuring Kernels

%

%    Different types of distance measuring methods, which are used with the

%    a 'Distance' morphology method for generating a gradient based on

%    distance from an edge of a binary shape, though there is a technique

%    for handling a anti-aliased shape.

%

%    See the 'Distance' Morphological Method, for information of how it is

%    applied.

%

%    Chebyshev:[{radius}][x{scale}[%!]]

%       Chebyshev Distance (also known as Tchebychev or Chessboard distance)

%       is a value of one to any neighbour, orthogonal or diagonal. One why

%       of thinking of it is the number of squares a 'King' or 'Queen' in

%       chess needs to traverse reach any other position on a chess board.

%       It results in a 'square' like distance function, but one where

%       diagonals are given a value that is closer than expected.

%

%    Manhattan:[{radius}][x{scale}[%!]]

%       Manhattan Distance (also known as Rectilinear, City Block, or the Taxi

%       Cab distance metric), it is the distance needed when you can only

%       travel in horizontal or vertical directions only.  It is the

%       distance a 'Rook' in chess would have to travel, and results in a

%       diamond like distances, where diagonals are further than expected.

%

%    Octagonal:[{radius}][x{scale}[%!]]

%       An interleaving of Manhattan and Chebyshev metrics producing an

%       increasing octagonally shaped distance.  Distances matches those of

%       the "Octagon" shaped kernel of the same radius.  The minimum radius

%       and default is 2, producing a 5x5 kernel.

%

%    Euclidean:[{radius}][x{scale}[%!]]

%       Euclidean distance is the 'direct' or 'as the crow flys' distance.

%       However by default the kernel size only has a radius of 1, which

%       limits the distance to 'Knight' like moves, with only orthogonal and

%       diagonal measurements being correct.  As such for the default kernel

%       you will get octagonal like distance function.

%

%       However using a larger radius such as "Euclidean:4" you will get a

%       much smoother distance gradient from the edge of the shape. Especially

%       if the image is pre-processed to include any anti-aliasing pixels.

%       Of course a larger kernel is slower to use, and not always needed.

%

%    The first three Distance Measuring Kernels will only generate distances

%    of exact multiples of {scale} in binary images. As such you can use a

%    scale of 1 without loosing any information.  However you also need some

%    scaling when handling non-binary anti-aliased shapes.

%

%    The "Euclidean" Distance Kernel however does generate a non-integer

%    fractional results, and as such scaling is vital even for binary shapes.

%

*/

MagickExport KernelInfo *AcquireKernelBuiltIn(const KernelInfoType type,

  const GeometryInfo *args,ExceptionInfo *exception)

{

  KernelInfo

    *kernel;

  ssize_t

    i;

  ssize_t

    u,

    v;

  double

    nan = sqrt(-1.0);  /* Special Value : Not A Number */

  /* Generate a new empty kernel if needed */

  kernel=(KernelInfo *) NULL;

  switch(type) {

    case UndefinedKernel:    /* These should not call this function */

    case UserDefinedKernel:

      (void) ThrowMagickException(exception,GetMagickModule(),OptionWarning,

        "InvalidOption","`%s'","Should not call this function");

      return((KernelInfo *) NULL);

    case LaplacianKernel:   /* Named Discrete Convolution Kernels */

    case SobelKernel:       /* these are defined using other kernels */

    case RobertsKernel:

    case PrewittKernel:

    case CompassKernel:

    case KirschKernel:

    case FreiChenKernel:

    case EdgesKernel:       /* Hit and Miss kernels */

    case CornersKernel:

    case DiagonalsKernel:

    case LineEndsKernel:

    case LineJunctionsKernel:

    case RidgesKernel:

    case ConvexHullKernel:

    case SkeletonKernel:

    case ThinSEKernel:

      break;               /* A pre-generated kernel is not needed */

#if 0

    /* set to 1 to do a compile-time check that we haven't missed anything */

    case UnityKernel:

    case GaussianKernel:

    case DoGKernel:

    case LoGKernel:

    case BlurKernel:

    case CometKernel:

    case BinomialKernel:

    case DiamondKernel:

    case SquareKernel: