• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
1 /*
2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3 %                                                                             %
4 %                                                                             %
5 %                                                                             %
6 %       TTTTT  H   H  RRRR   EEEEE  SSSSS  H   H   OOO   L      DDDD          %
7 %         T    H   H  R   R  E      SS     H   H  O   O  L      D   D         %
8 %         T    HHHHH  RRRR   EEE     SSS   HHHHH  O   O  L      D   D         %
9 %         T    H   H  R R    E         SS  H   H  O   O  L      D   D         %
10 %         T    H   H  R  R   EEEEE  SSSSS  H   H   OOO   LLLLL  DDDD          %
11 %                                                                             %
12 %                                                                             %
13 %                      MagickCore Image Threshold Methods                     %
14 %                                                                             %
15 %                               Software Design                               %
16 %                                    Cristy                                   %
17 %                                 October 1996                                %
18 %                                                                             %
19 %                                                                             %
20 %  Copyright 1999-2019 ImageMagick Studio LLC, a non-profit organization      %
21 %  dedicated to making software imaging solutions freely available.           %
22 %                                                                             %
23 %  You may not use this file except in compliance with the License.  You may  %
24 %  obtain a copy of the License at                                            %
25 %                                                                             %
26 %    https://imagemagick.org/script/license.php                               %
27 %                                                                             %
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.                                             %
33 %                                                                             %
34 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
35 %
36 %
37 %
38 */
39 
40 /*
41   Include declarations.
42 */
43 #include "MagickCore/studio.h"
44 #include "MagickCore/property.h"
45 #include "MagickCore/blob.h"
46 #include "MagickCore/cache-view.h"
47 #include "MagickCore/color.h"
48 #include "MagickCore/color-private.h"
49 #include "MagickCore/colormap.h"
50 #include "MagickCore/colorspace.h"
51 #include "MagickCore/colorspace-private.h"
52 #include "MagickCore/configure.h"
53 #include "MagickCore/constitute.h"
54 #include "MagickCore/decorate.h"
55 #include "MagickCore/draw.h"
56 #include "MagickCore/enhance.h"
57 #include "MagickCore/exception.h"
58 #include "MagickCore/exception-private.h"
59 #include "MagickCore/effect.h"
60 #include "MagickCore/fx.h"
61 #include "MagickCore/gem.h"
62 #include "MagickCore/geometry.h"
63 #include "MagickCore/image-private.h"
64 #include "MagickCore/list.h"
65 #include "MagickCore/log.h"
66 #include "MagickCore/memory_.h"
67 #include "MagickCore/memory-private.h"
68 #include "MagickCore/monitor.h"
69 #include "MagickCore/monitor-private.h"
70 #include "MagickCore/montage.h"
71 #include "MagickCore/option.h"
72 #include "MagickCore/pixel-accessor.h"
73 #include "MagickCore/pixel-private.h"
74 #include "MagickCore/quantize.h"
75 #include "MagickCore/quantum.h"
76 #include "MagickCore/quantum-private.h"
77 #include "MagickCore/random_.h"
78 #include "MagickCore/random-private.h"
79 #include "MagickCore/resize.h"
80 #include "MagickCore/resource_.h"
81 #include "MagickCore/segment.h"
82 #include "MagickCore/shear.h"
83 #include "MagickCore/signature-private.h"
84 #include "MagickCore/string_.h"
85 #include "MagickCore/string-private.h"
86 #include "MagickCore/thread-private.h"
87 #include "MagickCore/threshold.h"
88 #include "MagickCore/token.h"
89 #include "MagickCore/transform.h"
90 #include "MagickCore/xml-tree.h"
91 #include "MagickCore/xml-tree-private.h"
92 
93 /*
94   Define declarations.
95 */
96 #define ThresholdsFilename  "thresholds.xml"
97 
98 /*
99   Typedef declarations.
100 */
101 struct _ThresholdMap
102 {
103   char
104     *map_id,
105     *description;
106 
107   size_t
108     width,
109     height;
110 
111   ssize_t
112     divisor,
113     *levels;
114 };
115 
116 /*
117   Static declarations.
118 */
119 static const char
120   *MinimalThresholdMap =
121     "<?xml version=\"1.0\"?>"
122     "<thresholds>"
123     "  <threshold map=\"threshold\" alias=\"1x1\">"
124     "    <description>Threshold 1x1 (non-dither)</description>"
125     "    <levels width=\"1\" height=\"1\" divisor=\"2\">"
126     "        1"
127     "    </levels>"
128     "  </threshold>"
129     "  <threshold map=\"checks\" alias=\"2x1\">"
130     "    <description>Checkerboard 2x1 (dither)</description>"
131     "    <levels width=\"2\" height=\"2\" divisor=\"3\">"
132     "       1 2"
133     "       2 1"
134     "    </levels>"
135     "  </threshold>"
136     "</thresholds>";
137 
138 /*
139   Forward declarations.
140 */
141 static ThresholdMap
142   *GetThresholdMapFile(const char *,const char *,const char *,ExceptionInfo *);
143 
144 /*
145 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
146 %                                                                             %
147 %                                                                             %
148 %                                                                             %
149 %     A d a p t i v e T h r e s h o l d I m a g e                             %
150 %                                                                             %
151 %                                                                             %
152 %                                                                             %
153 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
154 %
155 %  AdaptiveThresholdImage() selects an individual threshold for each pixel
156 %  based on the range of intensity values in its local neighborhood.  This
157 %  allows for thresholding of an image whose global intensity histogram
158 %  doesn't contain distinctive peaks.
159 %
160 %  The format of the AdaptiveThresholdImage method is:
161 %
162 %      Image *AdaptiveThresholdImage(const Image *image,const size_t width,
163 %        const size_t height,const double bias,ExceptionInfo *exception)
164 %
165 %  A description of each parameter follows:
166 %
167 %    o image: the image.
168 %
169 %    o width: the width of the local neighborhood.
170 %
171 %    o height: the height of the local neighborhood.
172 %
173 %    o bias: the mean bias.
174 %
175 %    o exception: return any errors or warnings in this structure.
176 %
177 */
AdaptiveThresholdImage(const Image * image,const size_t width,const size_t height,const double bias,ExceptionInfo * exception)178 MagickExport Image *AdaptiveThresholdImage(const Image *image,
179   const size_t width,const size_t height,const double bias,
180   ExceptionInfo *exception)
181 {
182 #define AdaptiveThresholdImageTag  "AdaptiveThreshold/Image"
183 
184   CacheView
185     *image_view,
186     *threshold_view;
187 
188   Image
189     *threshold_image;
190 
191   MagickBooleanType
192     status;
193 
194   MagickOffsetType
195     progress;
196 
197   MagickSizeType
198     number_pixels;
199 
200   ssize_t
201     y;
202 
203   /*
204     Initialize threshold image attributes.
205   */
206   assert(image != (Image *) NULL);
207   assert(image->signature == MagickCoreSignature);
208   if (image->debug != MagickFalse)
209     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
210   assert(exception != (ExceptionInfo *) NULL);
211   assert(exception->signature == MagickCoreSignature);
212   threshold_image=CloneImage(image,0,0,MagickTrue,exception);
213   if (threshold_image == (Image *) NULL)
214     return((Image *) NULL);
215   status=SetImageStorageClass(threshold_image,DirectClass,exception);
216   if (status == MagickFalse)
217     {
218       threshold_image=DestroyImage(threshold_image);
219       return((Image *) NULL);
220     }
221   /*
222     Threshold image.
223   */
224   status=MagickTrue;
225   progress=0;
226   number_pixels=(MagickSizeType) width*height;
227   image_view=AcquireVirtualCacheView(image,exception);
228   threshold_view=AcquireAuthenticCacheView(threshold_image,exception);
229 #if defined(MAGICKCORE_OPENMP_SUPPORT)
230   #pragma omp parallel for schedule(static) shared(progress,status) \
231     magick_number_threads(image,threshold_image,image->rows,1)
232 #endif
233   for (y=0; y < (ssize_t) image->rows; y++)
234   {
235     double
236       channel_bias[MaxPixelChannels],
237       channel_sum[MaxPixelChannels];
238 
239     register const Quantum
240       *magick_restrict p,
241       *magick_restrict pixels;
242 
243     register Quantum
244       *magick_restrict q;
245 
246     register ssize_t
247       i,
248       x;
249 
250     ssize_t
251       center,
252       u,
253       v;
254 
255     if (status == MagickFalse)
256       continue;
257     p=GetCacheViewVirtualPixels(image_view,-((ssize_t) width/2L),y-(ssize_t)
258       (height/2L),image->columns+width,height,exception);
259     q=QueueCacheViewAuthenticPixels(threshold_view,0,y,threshold_image->columns,
260       1,exception);
261     if ((p == (const Quantum *) NULL) || (q == (Quantum *) NULL))
262       {
263         status=MagickFalse;
264         continue;
265       }
266     center=(ssize_t) GetPixelChannels(image)*(image->columns+width)*(height/2L)+
267       GetPixelChannels(image)*(width/2);
268     for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
269     {
270       PixelChannel channel = GetPixelChannelChannel(image,i);
271       PixelTrait traits = GetPixelChannelTraits(image,channel);
272       PixelTrait threshold_traits=GetPixelChannelTraits(threshold_image,
273         channel);
274       if ((traits == UndefinedPixelTrait) ||
275           (threshold_traits == UndefinedPixelTrait))
276         continue;
277       if ((threshold_traits & CopyPixelTrait) != 0)
278         {
279           SetPixelChannel(threshold_image,channel,p[center+i],q);
280           continue;
281         }
282       pixels=p;
283       channel_bias[channel]=0.0;
284       channel_sum[channel]=0.0;
285       for (v=0; v < (ssize_t) height; v++)
286       {
287         for (u=0; u < (ssize_t) width; u++)
288         {
289           if (u == (ssize_t) (width-1))
290             channel_bias[channel]+=pixels[i];
291           channel_sum[channel]+=pixels[i];
292           pixels+=GetPixelChannels(image);
293         }
294         pixels+=GetPixelChannels(image)*image->columns;
295       }
296     }
297     for (x=0; x < (ssize_t) image->columns; x++)
298     {
299       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
300       {
301         double
302           mean;
303 
304         PixelChannel channel = GetPixelChannelChannel(image,i);
305         PixelTrait traits = GetPixelChannelTraits(image,channel);
306         PixelTrait threshold_traits=GetPixelChannelTraits(threshold_image,
307           channel);
308         if ((traits == UndefinedPixelTrait) ||
309             (threshold_traits == UndefinedPixelTrait))
310           continue;
311         if ((threshold_traits & CopyPixelTrait) != 0)
312           {
313             SetPixelChannel(threshold_image,channel,p[center+i],q);
314             continue;
315           }
316         channel_sum[channel]-=channel_bias[channel];
317         channel_bias[channel]=0.0;
318         pixels=p;
319         for (v=0; v < (ssize_t) height; v++)
320         {
321           channel_bias[channel]+=pixels[i];
322           pixels+=(width-1)*GetPixelChannels(image);
323           channel_sum[channel]+=pixels[i];
324           pixels+=GetPixelChannels(image)*(image->columns+1);
325         }
326         mean=(double) (channel_sum[channel]/number_pixels+bias);
327         SetPixelChannel(threshold_image,channel,(Quantum) ((double)
328           p[center+i] <= mean ? 0 : QuantumRange),q);
329       }
330       p+=GetPixelChannels(image);
331       q+=GetPixelChannels(threshold_image);
332     }
333     if (SyncCacheViewAuthenticPixels(threshold_view,exception) == MagickFalse)
334       status=MagickFalse;
335     if (image->progress_monitor != (MagickProgressMonitor) NULL)
336       {
337         MagickBooleanType
338           proceed;
339 
340 #if defined(MAGICKCORE_OPENMP_SUPPORT)
341         #pragma omp atomic
342 #endif
343         progress++;
344         proceed=SetImageProgress(image,AdaptiveThresholdImageTag,progress,
345           image->rows);
346         if (proceed == MagickFalse)
347           status=MagickFalse;
348       }
349   }
350   threshold_image->type=image->type;
351   threshold_view=DestroyCacheView(threshold_view);
352   image_view=DestroyCacheView(image_view);
353   if (status == MagickFalse)
354     threshold_image=DestroyImage(threshold_image);
355   return(threshold_image);
356 }
357 
358 /*
359 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
360 %                                                                             %
361 %                                                                             %
362 %                                                                             %
363 %     A u t o T h r e s h o l d I m a g e                                     %
364 %                                                                             %
365 %                                                                             %
366 %                                                                             %
367 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
368 %
369 %  AutoThresholdImage() automatically selects a threshold and replaces each
370 %  pixel in the image with a black pixel if the image intentsity is less than
371 %  the selected threshold otherwise white.
372 %
373 %  The format of the AutoThresholdImage method is:
374 %
375 %      MagickBooleanType AutoThresholdImage(Image *image,
376 %        const AutoThresholdMethod method,ExceptionInfo *exception)
377 %
378 %  A description of each parameter follows:
379 %
380 %    o image: The image to auto-threshold.
381 %
382 %    o method: choose from Kapur, OTSU, or Triangle.
383 %
384 %    o exception: return any errors or warnings in this structure.
385 %
386 */
387 
KapurThreshold(const Image * image,const double * histogram,ExceptionInfo * exception)388 static double KapurThreshold(const Image *image,const double *histogram,
389   ExceptionInfo *exception)
390 {
391 #define MaxIntensity  255
392 
393   double
394     *black_entropy,
395     *cumulative_histogram,
396     entropy,
397     epsilon,
398     maximum_entropy,
399     *white_entropy;
400 
401   register ssize_t
402     i,
403     j;
404 
405   size_t
406     threshold;
407 
408   /*
409     Compute optimal threshold from the entopy of the histogram.
410   */
411   cumulative_histogram=(double *) AcquireQuantumMemory(MaxIntensity+1UL,
412     sizeof(*cumulative_histogram));
413   black_entropy=(double *) AcquireQuantumMemory(MaxIntensity+1UL,
414     sizeof(*black_entropy));
415   white_entropy=(double *) AcquireQuantumMemory(MaxIntensity+1UL,
416     sizeof(*white_entropy));
417   if ((cumulative_histogram == (double *) NULL) ||
418       (black_entropy == (double *) NULL) || (white_entropy == (double *) NULL))
419     {
420       if (white_entropy != (double *) NULL)
421         white_entropy=(double *) RelinquishMagickMemory(white_entropy);
422       if (black_entropy != (double *) NULL)
423         black_entropy=(double *) RelinquishMagickMemory(black_entropy);
424       if (cumulative_histogram != (double *) NULL)
425         cumulative_histogram=(double *)
426           RelinquishMagickMemory(cumulative_histogram);
427       (void) ThrowMagickException(exception,GetMagickModule(),
428         ResourceLimitError,"MemoryAllocationFailed","`%s'",image->filename);
429       return(-1.0);
430     }
431    /*
432      Entropy for black and white parts of the histogram.
433    */
434    cumulative_histogram[0]=histogram[0];
435    for (i=1; i <= MaxIntensity; i++)
436      cumulative_histogram[i]=cumulative_histogram[i-1]+histogram[i];
437    epsilon=MagickMinimumValue;
438    for (j=0; j <= MaxIntensity; j++)
439    {
440      /*
441        Black entropy.
442      */
443      black_entropy[j]=0.0;
444      if (cumulative_histogram[j] > epsilon)
445        {
446          entropy=0.0;
447          for (i=0; i <= j; i++)
448            if (histogram[i] > epsilon)
449              entropy-=histogram[i]/cumulative_histogram[j]*
450                log(histogram[i]/cumulative_histogram[j]);
451          black_entropy[j]=entropy;
452        }
453      /*
454        White entropy.
455      */
456      white_entropy[j]=0.0;
457      if ((1.0-cumulative_histogram[j]) > epsilon)
458        {
459          entropy=0.0;
460          for (i=j+1; i <= MaxIntensity; i++)
461            if (histogram[i] > epsilon)
462              entropy-=histogram[i]/(1.0-cumulative_histogram[j])*
463                log(histogram[i]/(1.0-cumulative_histogram[j]));
464          white_entropy[j]=entropy;
465        }
466    }
467   /*
468     Find histogram bin with maximum entropy.
469   */
470   maximum_entropy=black_entropy[0]+white_entropy[0];
471   threshold=0;
472   for (j=1; j <= MaxIntensity; j++)
473     if ((black_entropy[j]+white_entropy[j]) > maximum_entropy)
474       {
475         maximum_entropy=black_entropy[j]+white_entropy[j];
476         threshold=(size_t) j;
477       }
478   /*
479     Free resources.
480   */
481   white_entropy=(double *) RelinquishMagickMemory(white_entropy);
482   black_entropy=(double *) RelinquishMagickMemory(black_entropy);
483   cumulative_histogram=(double *) RelinquishMagickMemory(cumulative_histogram);
484   return(100.0*threshold/MaxIntensity);
485 }
486 
OTSUThreshold(const Image * image,const double * histogram,ExceptionInfo * exception)487 static double OTSUThreshold(const Image *image,const double *histogram,
488   ExceptionInfo *exception)
489 {
490   double
491     max_sigma,
492     *myu,
493     *omega,
494     *probability,
495     *sigma,
496     threshold;
497 
498   register ssize_t
499     i;
500 
501   /*
502     Compute optimal threshold from maximization of inter-class variance.
503   */
504   myu=(double *) AcquireQuantumMemory(MaxIntensity+1UL,sizeof(*myu));
505   omega=(double *) AcquireQuantumMemory(MaxIntensity+1UL,sizeof(*omega));
506   probability=(double *) AcquireQuantumMemory(MaxIntensity+1UL,
507     sizeof(*probability));
508   sigma=(double *) AcquireQuantumMemory(MaxIntensity+1UL,sizeof(*sigma));
509   if ((myu == (double *) NULL) || (omega == (double *) NULL) ||
510       (probability == (double *) NULL) || (sigma == (double *) NULL))
511     {
512       if (sigma != (double *) NULL)
513         sigma=(double *) RelinquishMagickMemory(sigma);
514       if (probability != (double *) NULL)
515         probability=(double *) RelinquishMagickMemory(probability);
516       if (omega != (double *) NULL)
517         omega=(double *) RelinquishMagickMemory(omega);
518       if (myu != (double *) NULL)
519         myu=(double *) RelinquishMagickMemory(myu);
520       (void) ThrowMagickException(exception,GetMagickModule(),
521         ResourceLimitError,"MemoryAllocationFailed","`%s'",image->filename);
522       return(-1.0);
523     }
524   /*
525     Calculate probability density.
526   */
527   for (i=0; i <= (ssize_t) MaxIntensity; i++)
528     probability[i]=histogram[i];
529   /*
530     Generate probability of graylevels and mean value for separation.
531   */
532   omega[0]=probability[0];
533   myu[0]=0.0;
534   for (i=1; i <= (ssize_t) MaxIntensity; i++)
535   {
536     omega[i]=omega[i-1]+probability[i];
537     myu[i]=myu[i-1]+i*probability[i];
538   }
539   /*
540     Sigma maximization: inter-class variance and compute optimal threshold.
541   */
542   threshold=0;
543   max_sigma=0.0;
544   for (i=0; i < (ssize_t) MaxIntensity; i++)
545   {
546     sigma[i]=0.0;
547     if ((omega[i] != 0.0) && (omega[i] != 1.0))
548       sigma[i]=pow(myu[MaxIntensity]*omega[i]-myu[i],2.0)/(omega[i]*(1.0-
549         omega[i]));
550     if (sigma[i] > max_sigma)
551       {
552         max_sigma=sigma[i];
553         threshold=(double) i;
554       }
555   }
556   /*
557     Free resources.
558   */
559   myu=(double *) RelinquishMagickMemory(myu);
560   omega=(double *) RelinquishMagickMemory(omega);
561   probability=(double *) RelinquishMagickMemory(probability);
562   sigma=(double *) RelinquishMagickMemory(sigma);
563   return(100.0*threshold/MaxIntensity);
564 }
565 
TriangleThreshold(const double * histogram,ExceptionInfo * exception)566 static double TriangleThreshold(const double *histogram,
567   ExceptionInfo *exception)
568 {
569   double
570     a,
571     b,
572     c,
573     count,
574     distance,
575     inverse_ratio,
576     max_distance,
577     segment,
578     x1,
579     x2,
580     y1,
581     y2;
582 
583   register ssize_t
584     i;
585 
586   ssize_t
587     end,
588     max,
589     start,
590     threshold;
591 
592   /*
593     Compute optimal threshold with triangle algorithm.
594   */
595   (void) exception;
596   start=0;  /* find start bin, first bin not zero count */
597   for (i=0; i <= (ssize_t) MaxIntensity; i++)
598     if (histogram[i] > 0.0)
599       {
600         start=i;
601         break;
602       }
603   end=0;  /* find end bin, last bin not zero count */
604   for (i=(ssize_t) MaxIntensity; i >= 0; i--)
605     if (histogram[i] > 0.0)
606       {
607         end=i;
608         break;
609       }
610   max=0;  /* find max bin, bin with largest count */
611   count=0.0;
612   for (i=0; i <= (ssize_t) MaxIntensity; i++)
613     if (histogram[i] > count)
614       {
615         max=i;
616         count=histogram[i];
617       }
618   /*
619     Compute threshold at split point.
620   */
621   x1=(double) max;
622   y1=histogram[max];
623   x2=(double) end;
624   if ((max-start) >= (end-max))
625     x2=(double) start;
626   y2=0.0;
627   a=y1-y2;
628   b=x2-x1;
629   c=(-1.0)*(a*x1+b*y1);
630   inverse_ratio=1.0/sqrt(a*a+b*b+c*c);
631   threshold=0;
632   max_distance=0.0;
633   if (x2 == (double) start)
634     for (i=start; i < max; i++)
635     {
636       segment=inverse_ratio*(a*i+b*histogram[i]+c);
637       distance=sqrt(segment*segment);
638       if ((distance > max_distance) && (segment > 0.0))
639         {
640           threshold=i;
641           max_distance=distance;
642         }
643     }
644   else
645     for (i=end; i > max; i--)
646     {
647       segment=inverse_ratio*(a*i+b*histogram[i]+c);
648       distance=sqrt(segment*segment);
649       if ((distance > max_distance) && (segment < 0.0))
650         {
651           threshold=i;
652           max_distance=distance;
653         }
654     }
655   return(100.0*threshold/MaxIntensity);
656 }
657 
AutoThresholdImage(Image * image,const AutoThresholdMethod method,ExceptionInfo * exception)658 MagickExport MagickBooleanType AutoThresholdImage(Image *image,
659   const AutoThresholdMethod method,ExceptionInfo *exception)
660 {
661   CacheView
662     *image_view;
663 
664   char
665     property[MagickPathExtent];
666 
667   double
668     gamma,
669     *histogram,
670     sum,
671     threshold;
672 
673   MagickBooleanType
674     status;
675 
676   register ssize_t
677     i;
678 
679   ssize_t
680     y;
681 
682   /*
683     Form histogram.
684   */
685   assert(image != (Image *) NULL);
686   assert(image->signature == MagickCoreSignature);
687   if (image->debug != MagickFalse)
688     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
689   histogram=(double *) AcquireQuantumMemory(MaxIntensity+1UL,
690     sizeof(*histogram));
691   if (histogram == (double *) NULL)
692     ThrowBinaryException(ResourceLimitError,"MemoryAllocationFailed",
693       image->filename);
694   status=MagickTrue;
695   (void) memset(histogram,0,(MaxIntensity+1UL)*sizeof(*histogram));
696   image_view=AcquireVirtualCacheView(image,exception);
697   for (y=0; y < (ssize_t) image->rows; y++)
698   {
699     register const Quantum
700       *magick_restrict p;
701 
702     register ssize_t
703       x;
704 
705     p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
706     if (p == (const Quantum *) NULL)
707       break;
708     for (x=0; x < (ssize_t) image->columns; x++)
709     {
710       double intensity = GetPixelIntensity(image,p);
711       histogram[ScaleQuantumToChar(ClampToQuantum(intensity))]++;
712       p+=GetPixelChannels(image);
713     }
714   }
715   image_view=DestroyCacheView(image_view);
716   /*
717     Normalize histogram.
718   */
719   sum=0.0;
720   for (i=0; i <= (ssize_t) MaxIntensity; i++)
721     sum+=histogram[i];
722   gamma=PerceptibleReciprocal(sum);
723   for (i=0; i <= (ssize_t) MaxIntensity; i++)
724     histogram[i]=gamma*histogram[i];
725   /*
726     Discover threshold from histogram.
727   */
728   switch (method)
729   {
730     case KapurThresholdMethod:
731     {
732       threshold=KapurThreshold(image,histogram,exception);
733       break;
734     }
735     case OTSUThresholdMethod:
736     default:
737     {
738       threshold=OTSUThreshold(image,histogram,exception);
739       break;
740     }
741     case TriangleThresholdMethod:
742     {
743       threshold=TriangleThreshold(histogram,exception);
744       break;
745     }
746   }
747   histogram=(double *) RelinquishMagickMemory(histogram);
748   if (threshold < 0.0)
749     status=MagickFalse;
750   if (status == MagickFalse)
751     return(MagickFalse);
752   /*
753     Threshold image.
754   */
755   (void) FormatLocaleString(property,MagickPathExtent,"%g%%",threshold);
756   (void) SetImageProperty(image,"auto-threshold:threshold",property,exception);
757   return(BilevelImage(image,QuantumRange*threshold/100.0,exception));
758 }
759 
760 /*
761 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
762 %                                                                             %
763 %                                                                             %
764 %                                                                             %
765 %     B i l e v e l I m a g e                                                 %
766 %                                                                             %
767 %                                                                             %
768 %                                                                             %
769 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
770 %
771 %  BilevelImage() changes the value of individual pixels based on the
772 %  intensity of each pixel channel.  The result is a high-contrast image.
773 %
774 %  More precisely each channel value of the image is 'thresholded' so that if
775 %  it is equal to or less than the given value it is set to zero, while any
776 %  value greater than that give is set to it maximum or QuantumRange.
777 %
778 %  This function is what is used to implement the "-threshold" operator for
779 %  the command line API.
780 %
781 %  If the default channel setting is given the image is thresholded using just
782 %  the gray 'intensity' of the image, rather than the individual channels.
783 %
784 %  The format of the BilevelImage method is:
785 %
786 %      MagickBooleanType BilevelImage(Image *image,const double threshold,
787 %        ExceptionInfo *exception)
788 %
789 %  A description of each parameter follows:
790 %
791 %    o image: the image.
792 %
793 %    o threshold: define the threshold values.
794 %
795 %    o exception: return any errors or warnings in this structure.
796 %
797 %  Aside: You can get the same results as operator using LevelImages()
798 %  with the 'threshold' value for both the black_point and the white_point.
799 %
800 */
BilevelImage(Image * image,const double threshold,ExceptionInfo * exception)801 MagickExport MagickBooleanType BilevelImage(Image *image,const double threshold,
802   ExceptionInfo *exception)
803 {
804 #define ThresholdImageTag  "Threshold/Image"
805 
806   CacheView
807     *image_view;
808 
809   MagickBooleanType
810     status;
811 
812   MagickOffsetType
813     progress;
814 
815   ssize_t
816     y;
817 
818   assert(image != (Image *) NULL);
819   assert(image->signature == MagickCoreSignature);
820   if (image->debug != MagickFalse)
821     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
822   if (SetImageStorageClass(image,DirectClass,exception) == MagickFalse)
823     return(MagickFalse);
824   if (IsGrayColorspace(image->colorspace) != MagickFalse)
825     (void) SetImageColorspace(image,sRGBColorspace,exception);
826   /*
827     Bilevel threshold image.
828   */
829   status=MagickTrue;
830   progress=0;
831   image_view=AcquireAuthenticCacheView(image,exception);
832 #if defined(MAGICKCORE_OPENMP_SUPPORT)
833   #pragma omp parallel for schedule(static) shared(progress,status) \
834     magick_number_threads(image,image,image->rows,1)
835 #endif
836   for (y=0; y < (ssize_t) image->rows; y++)
837   {
838     register ssize_t
839       x;
840 
841     register Quantum
842       *magick_restrict q;
843 
844     if (status == MagickFalse)
845       continue;
846     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
847     if (q == (Quantum *) NULL)
848       {
849         status=MagickFalse;
850         continue;
851       }
852     for (x=0; x < (ssize_t) image->columns; x++)
853     {
854       double
855         pixel;
856 
857       register ssize_t
858         i;
859 
860       pixel=GetPixelIntensity(image,q);
861       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
862       {
863         PixelChannel channel = GetPixelChannelChannel(image,i);
864         PixelTrait traits = GetPixelChannelTraits(image,channel);
865         if ((traits & UpdatePixelTrait) == 0)
866           continue;
867         if (image->channel_mask != DefaultChannels)
868           pixel=(double) q[i];
869         q[i]=(Quantum) (pixel <= threshold ? 0 : QuantumRange);
870       }
871       q+=GetPixelChannels(image);
872     }
873     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
874       status=MagickFalse;
875     if (image->progress_monitor != (MagickProgressMonitor) NULL)
876       {
877         MagickBooleanType
878           proceed;
879 
880 #if defined(MAGICKCORE_OPENMP_SUPPORT)
881         #pragma omp atomic
882 #endif
883         progress++;
884         proceed=SetImageProgress(image,ThresholdImageTag,progress++,
885           image->rows);
886         if (proceed == MagickFalse)
887           status=MagickFalse;
888       }
889   }
890   image_view=DestroyCacheView(image_view);
891   return(status);
892 }
893 
894 /*
895 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
896 %                                                                             %
897 %                                                                             %
898 %                                                                             %
899 %     B l a c k T h r e s h o l d I m a g e                                   %
900 %                                                                             %
901 %                                                                             %
902 %                                                                             %
903 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
904 %
905 %  BlackThresholdImage() is like ThresholdImage() but forces all pixels below
906 %  the threshold into black while leaving all pixels at or above the threshold
907 %  unchanged.
908 %
909 %  The format of the BlackThresholdImage method is:
910 %
911 %      MagickBooleanType BlackThresholdImage(Image *image,
912 %        const char *threshold,ExceptionInfo *exception)
913 %
914 %  A description of each parameter follows:
915 %
916 %    o image: the image.
917 %
918 %    o threshold: define the threshold value.
919 %
920 %    o exception: return any errors or warnings in this structure.
921 %
922 */
BlackThresholdImage(Image * image,const char * thresholds,ExceptionInfo * exception)923 MagickExport MagickBooleanType BlackThresholdImage(Image *image,
924   const char *thresholds,ExceptionInfo *exception)
925 {
926 #define ThresholdImageTag  "Threshold/Image"
927 
928   CacheView
929     *image_view;
930 
931   GeometryInfo
932     geometry_info;
933 
934   MagickBooleanType
935     status;
936 
937   MagickOffsetType
938     progress;
939 
940   PixelInfo
941     threshold;
942 
943   MagickStatusType
944     flags;
945 
946   ssize_t
947     y;
948 
949   assert(image != (Image *) NULL);
950   assert(image->signature == MagickCoreSignature);
951   if (image->debug != MagickFalse)
952     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
953   if (thresholds == (const char *) NULL)
954     return(MagickTrue);
955   if (SetImageStorageClass(image,DirectClass,exception) == MagickFalse)
956     return(MagickFalse);
957   if (IsGrayColorspace(image->colorspace) != MagickFalse)
958     (void) SetImageColorspace(image,sRGBColorspace,exception);
959   GetPixelInfo(image,&threshold);
960   flags=ParseGeometry(thresholds,&geometry_info);
961   threshold.red=geometry_info.rho;
962   threshold.green=geometry_info.rho;
963   threshold.blue=geometry_info.rho;
964   threshold.black=geometry_info.rho;
965   threshold.alpha=100.0;
966   if ((flags & SigmaValue) != 0)
967     threshold.green=geometry_info.sigma;
968   if ((flags & XiValue) != 0)
969     threshold.blue=geometry_info.xi;
970   if ((flags & PsiValue) != 0)
971     threshold.alpha=geometry_info.psi;
972   if (threshold.colorspace == CMYKColorspace)
973     {
974       if ((flags & PsiValue) != 0)
975         threshold.black=geometry_info.psi;
976       if ((flags & ChiValue) != 0)
977         threshold.alpha=geometry_info.chi;
978     }
979   if ((flags & PercentValue) != 0)
980     {
981       threshold.red*=(MagickRealType) (QuantumRange/100.0);
982       threshold.green*=(MagickRealType) (QuantumRange/100.0);
983       threshold.blue*=(MagickRealType) (QuantumRange/100.0);
984       threshold.black*=(MagickRealType) (QuantumRange/100.0);
985       threshold.alpha*=(MagickRealType) (QuantumRange/100.0);
986     }
987   /*
988     White threshold image.
989   */
990   status=MagickTrue;
991   progress=0;
992   image_view=AcquireAuthenticCacheView(image,exception);
993 #if defined(MAGICKCORE_OPENMP_SUPPORT)
994   #pragma omp parallel for schedule(static) shared(progress,status) \
995     magick_number_threads(image,image,image->rows,1)
996 #endif
997   for (y=0; y < (ssize_t) image->rows; y++)
998   {
999     register ssize_t
1000       x;
1001 
1002     register Quantum
1003       *magick_restrict q;
1004 
1005     if (status == MagickFalse)
1006       continue;
1007     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
1008     if (q == (Quantum *) NULL)
1009       {
1010         status=MagickFalse;
1011         continue;
1012       }
1013     for (x=0; x < (ssize_t) image->columns; x++)
1014     {
1015       double
1016         pixel;
1017 
1018       register ssize_t
1019         i;
1020 
1021       pixel=GetPixelIntensity(image,q);
1022       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
1023       {
1024         PixelChannel channel = GetPixelChannelChannel(image,i);
1025         PixelTrait traits = GetPixelChannelTraits(image,channel);
1026         if ((traits & UpdatePixelTrait) == 0)
1027           continue;
1028         if (image->channel_mask != DefaultChannels)
1029           pixel=(double) q[i];
1030         if (pixel < GetPixelInfoChannel(&threshold,channel))
1031           q[i]=(Quantum) 0;
1032       }
1033       q+=GetPixelChannels(image);
1034     }
1035     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
1036       status=MagickFalse;
1037     if (image->progress_monitor != (MagickProgressMonitor) NULL)
1038       {
1039         MagickBooleanType
1040           proceed;
1041 
1042 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1043         #pragma omp atomic
1044 #endif
1045         progress++;
1046         proceed=SetImageProgress(image,ThresholdImageTag,progress,
1047           image->rows);
1048         if (proceed == MagickFalse)
1049           status=MagickFalse;
1050       }
1051   }
1052   image_view=DestroyCacheView(image_view);
1053   return(status);
1054 }
1055 
1056 /*
1057 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1058 %                                                                             %
1059 %                                                                             %
1060 %                                                                             %
1061 %     C l a m p I m a g e                                                     %
1062 %                                                                             %
1063 %                                                                             %
1064 %                                                                             %
1065 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1066 %
1067 %  ClampImage() set each pixel whose value is below zero to zero and any the
1068 %  pixel whose value is above the quantum range to the quantum range (e.g.
1069 %  65535) otherwise the pixel value remains unchanged.
1070 %
1071 %  The format of the ClampImage method is:
1072 %
1073 %      MagickBooleanType ClampImage(Image *image,ExceptionInfo *exception)
1074 %
1075 %  A description of each parameter follows:
1076 %
1077 %    o image: the image.
1078 %
1079 %    o exception: return any errors or warnings in this structure.
1080 %
1081 */
1082 
ClampImage(Image * image,ExceptionInfo * exception)1083 MagickExport MagickBooleanType ClampImage(Image *image,ExceptionInfo *exception)
1084 {
1085 #define ClampImageTag  "Clamp/Image"
1086 
1087   CacheView
1088     *image_view;
1089 
1090   MagickBooleanType
1091     status;
1092 
1093   MagickOffsetType
1094     progress;
1095 
1096   ssize_t
1097     y;
1098 
1099   assert(image != (Image *) NULL);
1100   assert(image->signature == MagickCoreSignature);
1101   if (image->debug != MagickFalse)
1102     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1103   if (image->storage_class == PseudoClass)
1104     {
1105       register ssize_t
1106         i;
1107 
1108       register PixelInfo
1109         *magick_restrict q;
1110 
1111       q=image->colormap;
1112       for (i=0; i < (ssize_t) image->colors; i++)
1113       {
1114         q->red=(double) ClampPixel(q->red);
1115         q->green=(double) ClampPixel(q->green);
1116         q->blue=(double) ClampPixel(q->blue);
1117         q->alpha=(double) ClampPixel(q->alpha);
1118         q++;
1119       }
1120       return(SyncImage(image,exception));
1121     }
1122   /*
1123     Clamp image.
1124   */
1125   status=MagickTrue;
1126   progress=0;
1127   image_view=AcquireAuthenticCacheView(image,exception);
1128 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1129   #pragma omp parallel for schedule(static) shared(progress,status) \
1130     magick_number_threads(image,image,image->rows,1)
1131 #endif
1132   for (y=0; y < (ssize_t) image->rows; y++)
1133   {
1134     register ssize_t
1135       x;
1136 
1137     register Quantum
1138       *magick_restrict q;
1139 
1140     if (status == MagickFalse)
1141       continue;
1142     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
1143     if (q == (Quantum *) NULL)
1144       {
1145         status=MagickFalse;
1146         continue;
1147       }
1148     for (x=0; x < (ssize_t) image->columns; x++)
1149     {
1150       register ssize_t
1151         i;
1152 
1153       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
1154       {
1155         PixelChannel channel = GetPixelChannelChannel(image,i);
1156         PixelTrait traits = GetPixelChannelTraits(image,channel);
1157         if ((traits & UpdatePixelTrait) == 0)
1158           continue;
1159         q[i]=ClampPixel((MagickRealType) q[i]);
1160       }
1161       q+=GetPixelChannels(image);
1162     }
1163     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
1164       status=MagickFalse;
1165     if (image->progress_monitor != (MagickProgressMonitor) NULL)
1166       {
1167         MagickBooleanType
1168           proceed;
1169 
1170 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1171         #pragma omp atomic
1172 #endif
1173         progress++;
1174         proceed=SetImageProgress(image,ClampImageTag,progress,image->rows);
1175         if (proceed == MagickFalse)
1176           status=MagickFalse;
1177       }
1178   }
1179   image_view=DestroyCacheView(image_view);
1180   return(status);
1181 }
1182 
1183 /*
1184 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1185 %                                                                             %
1186 %                                                                             %
1187 %                                                                             %
1188 %  D e s t r o y T h r e s h o l d M a p                                      %
1189 %                                                                             %
1190 %                                                                             %
1191 %                                                                             %
1192 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1193 %
1194 %  DestroyThresholdMap() de-allocate the given ThresholdMap
1195 %
1196 %  The format of the ListThresholdMaps method is:
1197 %
1198 %      ThresholdMap *DestroyThresholdMap(Threshold *map)
1199 %
1200 %  A description of each parameter follows.
1201 %
1202 %    o map:    Pointer to the Threshold map to destroy
1203 %
1204 */
DestroyThresholdMap(ThresholdMap * map)1205 MagickExport ThresholdMap *DestroyThresholdMap(ThresholdMap *map)
1206 {
1207   assert(map != (ThresholdMap *) NULL);
1208   if (map->map_id != (char *) NULL)
1209     map->map_id=DestroyString(map->map_id);
1210   if (map->description != (char *) NULL)
1211     map->description=DestroyString(map->description);
1212   if (map->levels != (ssize_t *) NULL)
1213     map->levels=(ssize_t *) RelinquishMagickMemory(map->levels);
1214   map=(ThresholdMap *) RelinquishMagickMemory(map);
1215   return(map);
1216 }
1217 
1218 /*
1219 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1220 %                                                                             %
1221 %                                                                             %
1222 %                                                                             %
1223 %  G e t T h r e s h o l d M a p                                              %
1224 %                                                                             %
1225 %                                                                             %
1226 %                                                                             %
1227 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1228 %
1229 %  GetThresholdMap() loads and searches one or more threshold map files for the
1230 %  map matching the given name or alias.
1231 %
1232 %  The format of the GetThresholdMap method is:
1233 %
1234 %      ThresholdMap *GetThresholdMap(const char *map_id,
1235 %        ExceptionInfo *exception)
1236 %
1237 %  A description of each parameter follows.
1238 %
1239 %    o map_id:  ID of the map to look for.
1240 %
1241 %    o exception: return any errors or warnings in this structure.
1242 %
1243 */
GetThresholdMap(const char * map_id,ExceptionInfo * exception)1244 MagickExport ThresholdMap *GetThresholdMap(const char *map_id,
1245   ExceptionInfo *exception)
1246 {
1247   ThresholdMap
1248     *map;
1249 
1250   map=GetThresholdMapFile(MinimalThresholdMap,"built-in",map_id,exception);
1251   if (map != (ThresholdMap *) NULL)
1252     return(map);
1253 #if !defined(MAGICKCORE_ZERO_CONFIGURATION_SUPPORT)
1254   {
1255     const StringInfo
1256       *option;
1257 
1258     LinkedListInfo
1259       *options;
1260 
1261     options=GetConfigureOptions(ThresholdsFilename,exception);
1262     option=(const StringInfo *) GetNextValueInLinkedList(options);
1263     while (option != (const StringInfo *) NULL)
1264     {
1265       map=GetThresholdMapFile((const char *) GetStringInfoDatum(option),
1266         GetStringInfoPath(option),map_id,exception);
1267       if (map != (ThresholdMap *) NULL)
1268         break;
1269       option=(const StringInfo *) GetNextValueInLinkedList(options);
1270     }
1271     options=DestroyConfigureOptions(options);
1272   }
1273 #endif
1274   return(map);
1275 }
1276 
1277 /*
1278 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1279 %                                                                             %
1280 %                                                                             %
1281 %                                                                             %
1282 +  G e t T h r e s h o l d M a p F i l e                                      %
1283 %                                                                             %
1284 %                                                                             %
1285 %                                                                             %
1286 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1287 %
1288 %  GetThresholdMapFile() look for a given threshold map name or alias in the
1289 %  given XML file data, and return the allocated the map when found.
1290 %
1291 %  The format of the ListThresholdMaps method is:
1292 %
1293 %      ThresholdMap *GetThresholdMap(const char *xml,const char *filename,
1294 %         const char *map_id,ExceptionInfo *exception)
1295 %
1296 %  A description of each parameter follows.
1297 %
1298 %    o xml:  The threshold map list in XML format.
1299 %
1300 %    o filename:  The threshold map XML filename.
1301 %
1302 %    o map_id:  ID of the map to look for in XML list.
1303 %
1304 %    o exception: return any errors or warnings in this structure.
1305 %
1306 */
GetThresholdMapFile(const char * xml,const char * filename,const char * map_id,ExceptionInfo * exception)1307 static ThresholdMap *GetThresholdMapFile(const char *xml,const char *filename,
1308   const char *map_id,ExceptionInfo *exception)
1309 {
1310   char
1311     *p;
1312 
1313   const char
1314     *attribute,
1315     *content;
1316 
1317   double
1318     value;
1319 
1320   register ssize_t
1321     i;
1322 
1323   ThresholdMap
1324     *map;
1325 
1326   XMLTreeInfo
1327     *description,
1328     *levels,
1329     *threshold,
1330     *thresholds;
1331 
1332   (void) LogMagickEvent(ConfigureEvent,GetMagickModule(),
1333     "Loading threshold map file \"%s\" ...",filename);
1334   map=(ThresholdMap *) NULL;
1335   thresholds=NewXMLTree(xml,exception);
1336   if (thresholds == (XMLTreeInfo *) NULL)
1337     return(map);
1338   for (threshold=GetXMLTreeChild(thresholds,"threshold");
1339        threshold != (XMLTreeInfo *) NULL;
1340        threshold=GetNextXMLTreeTag(threshold))
1341   {
1342     attribute=GetXMLTreeAttribute(threshold,"map");
1343     if ((attribute != (char *) NULL) && (LocaleCompare(map_id,attribute) == 0))
1344       break;
1345     attribute=GetXMLTreeAttribute(threshold,"alias");
1346     if ((attribute != (char *) NULL) && (LocaleCompare(map_id,attribute) == 0))
1347       break;
1348   }
1349   if (threshold == (XMLTreeInfo *) NULL)
1350     {
1351       thresholds=DestroyXMLTree(thresholds);
1352       return(map);
1353     }
1354   description=GetXMLTreeChild(threshold,"description");
1355   if (description == (XMLTreeInfo *) NULL)
1356     {
1357       (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
1358         "XmlMissingElement", "<description>, map \"%s\"",map_id);
1359       thresholds=DestroyXMLTree(thresholds);
1360       return(map);
1361     }
1362   levels=GetXMLTreeChild(threshold,"levels");
1363   if (levels == (XMLTreeInfo *) NULL)
1364     {
1365       (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
1366         "XmlMissingElement", "<levels>, map \"%s\"", map_id);
1367       thresholds=DestroyXMLTree(thresholds);
1368       return(map);
1369     }
1370   map=(ThresholdMap *) AcquireCriticalMemory(sizeof(*map));
1371   map->map_id=(char *) NULL;
1372   map->description=(char *) NULL;
1373   map->levels=(ssize_t *) NULL;
1374   attribute=GetXMLTreeAttribute(threshold,"map");
1375   if (attribute != (char *) NULL)
1376     map->map_id=ConstantString(attribute);
1377   content=GetXMLTreeContent(description);
1378   if (content != (char *) NULL)
1379     map->description=ConstantString(content);
1380   attribute=GetXMLTreeAttribute(levels,"width");
1381   if (attribute == (char *) NULL)
1382     {
1383       (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
1384         "XmlMissingAttribute", "<levels width>, map \"%s\"",map_id);
1385       thresholds=DestroyXMLTree(thresholds);
1386       map=DestroyThresholdMap(map);
1387       return(map);
1388     }
1389   map->width=StringToUnsignedLong(attribute);
1390   if (map->width == 0)
1391     {
1392       (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
1393        "XmlInvalidAttribute", "<levels width>, map \"%s\"",map_id);
1394       thresholds=DestroyXMLTree(thresholds);
1395       map=DestroyThresholdMap(map);
1396       return(map);
1397     }
1398   attribute=GetXMLTreeAttribute(levels,"height");
1399   if (attribute == (char *) NULL)
1400     {
1401       (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
1402         "XmlMissingAttribute", "<levels height>, map \"%s\"",map_id);
1403       thresholds=DestroyXMLTree(thresholds);
1404       map=DestroyThresholdMap(map);
1405       return(map);
1406     }
1407   map->height=StringToUnsignedLong(attribute);
1408   if (map->height == 0)
1409     {
1410       (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
1411         "XmlInvalidAttribute", "<levels height>, map \"%s\"",map_id);
1412       thresholds=DestroyXMLTree(thresholds);
1413       map=DestroyThresholdMap(map);
1414       return(map);
1415     }
1416   attribute=GetXMLTreeAttribute(levels,"divisor");
1417   if (attribute == (char *) NULL)
1418     {
1419       (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
1420         "XmlMissingAttribute", "<levels divisor>, map \"%s\"",map_id);
1421       thresholds=DestroyXMLTree(thresholds);
1422       map=DestroyThresholdMap(map);
1423       return(map);
1424     }
1425   map->divisor=(ssize_t) StringToLong(attribute);
1426   if (map->divisor < 2)
1427     {
1428       (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
1429         "XmlInvalidAttribute", "<levels divisor>, map \"%s\"",map_id);
1430       thresholds=DestroyXMLTree(thresholds);
1431       map=DestroyThresholdMap(map);
1432       return(map);
1433     }
1434   content=GetXMLTreeContent(levels);
1435   if (content == (char *) NULL)
1436     {
1437       (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
1438         "XmlMissingContent", "<levels>, map \"%s\"",map_id);
1439       thresholds=DestroyXMLTree(thresholds);
1440       map=DestroyThresholdMap(map);
1441       return(map);
1442     }
1443   map->levels=(ssize_t *) AcquireQuantumMemory((size_t) map->width,map->height*
1444     sizeof(*map->levels));
1445   if (map->levels == (ssize_t *) NULL)
1446     ThrowFatalException(ResourceLimitFatalError,"UnableToAcquireThresholdMap");
1447   for (i=0; i < (ssize_t) (map->width*map->height); i++)
1448   {
1449     map->levels[i]=(ssize_t) strtol(content,&p,10);
1450     if (p == content)
1451       {
1452         (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
1453           "XmlInvalidContent", "<level> too few values, map \"%s\"",map_id);
1454         thresholds=DestroyXMLTree(thresholds);
1455         map=DestroyThresholdMap(map);
1456         return(map);
1457       }
1458     if ((map->levels[i] < 0) || (map->levels[i] > map->divisor))
1459       {
1460         (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
1461           "XmlInvalidContent", "<level> %.20g out of range, map \"%s\"",
1462           (double) map->levels[i],map_id);
1463         thresholds=DestroyXMLTree(thresholds);
1464         map=DestroyThresholdMap(map);
1465         return(map);
1466       }
1467     content=p;
1468   }
1469   value=(double) strtol(content,&p,10);
1470   (void) value;
1471   if (p != content)
1472     {
1473       (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
1474         "XmlInvalidContent", "<level> too many values, map \"%s\"",map_id);
1475      thresholds=DestroyXMLTree(thresholds);
1476      map=DestroyThresholdMap(map);
1477      return(map);
1478    }
1479   thresholds=DestroyXMLTree(thresholds);
1480   return(map);
1481 }
1482 
1483 /*
1484 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1485 %                                                                             %
1486 %                                                                             %
1487 %                                                                             %
1488 +  L i s t T h r e s h o l d M a p F i l e                                    %
1489 %                                                                             %
1490 %                                                                             %
1491 %                                                                             %
1492 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1493 %
1494 %  ListThresholdMapFile() lists the threshold maps and their descriptions
1495 %  in the given XML file data.
1496 %
1497 %  The format of the ListThresholdMaps method is:
1498 %
1499 %      MagickBooleanType ListThresholdMaps(FILE *file,const char*xml,
1500 %         const char *filename,ExceptionInfo *exception)
1501 %
1502 %  A description of each parameter follows.
1503 %
1504 %    o file:  An pointer to the output FILE.
1505 %
1506 %    o xml:  The threshold map list in XML format.
1507 %
1508 %    o filename:  The threshold map XML filename.
1509 %
1510 %    o exception: return any errors or warnings in this structure.
1511 %
1512 */
ListThresholdMapFile(FILE * file,const char * xml,const char * filename,ExceptionInfo * exception)1513 MagickBooleanType ListThresholdMapFile(FILE *file,const char *xml,
1514   const char *filename,ExceptionInfo *exception)
1515 {
1516   const char
1517     *alias,
1518     *content,
1519     *map;
1520 
1521   XMLTreeInfo
1522     *description,
1523     *threshold,
1524     *thresholds;
1525 
1526   assert( xml != (char *) NULL );
1527   assert( file != (FILE *) NULL );
1528   (void) LogMagickEvent(ConfigureEvent,GetMagickModule(),
1529     "Loading threshold map file \"%s\" ...",filename);
1530   thresholds=NewXMLTree(xml,exception);
1531   if ( thresholds == (XMLTreeInfo *) NULL )
1532     return(MagickFalse);
1533   (void) FormatLocaleFile(file,"%-16s %-12s %s\n","Map","Alias","Description");
1534   (void) FormatLocaleFile(file,
1535     "----------------------------------------------------\n");
1536   threshold=GetXMLTreeChild(thresholds,"threshold");
1537   for ( ; threshold != (XMLTreeInfo *) NULL;
1538           threshold=GetNextXMLTreeTag(threshold))
1539   {
1540     map=GetXMLTreeAttribute(threshold,"map");
1541     if (map == (char *) NULL)
1542       {
1543         (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
1544           "XmlMissingAttribute", "<map>");
1545         thresholds=DestroyXMLTree(thresholds);
1546         return(MagickFalse);
1547       }
1548     alias=GetXMLTreeAttribute(threshold,"alias");
1549     description=GetXMLTreeChild(threshold,"description");
1550     if (description == (XMLTreeInfo *) NULL)
1551       {
1552         (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
1553           "XmlMissingElement", "<description>, map \"%s\"",map);
1554         thresholds=DestroyXMLTree(thresholds);
1555         return(MagickFalse);
1556       }
1557     content=GetXMLTreeContent(description);
1558     if (content == (char *) NULL)
1559       {
1560         (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
1561           "XmlMissingContent", "<description>, map \"%s\"", map);
1562         thresholds=DestroyXMLTree(thresholds);
1563         return(MagickFalse);
1564       }
1565     (void) FormatLocaleFile(file,"%-16s %-12s %s\n",map,alias ? alias : "",
1566       content);
1567   }
1568   thresholds=DestroyXMLTree(thresholds);
1569   return(MagickTrue);
1570 }
1571 
1572 /*
1573 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1574 %                                                                             %
1575 %                                                                             %
1576 %                                                                             %
1577 %  L i s t T h r e s h o l d M a p s                                          %
1578 %                                                                             %
1579 %                                                                             %
1580 %                                                                             %
1581 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1582 %
1583 %  ListThresholdMaps() lists the threshold maps and their descriptions
1584 %  as defined by "threshold.xml" to a file.
1585 %
1586 %  The format of the ListThresholdMaps method is:
1587 %
1588 %      MagickBooleanType ListThresholdMaps(FILE *file,ExceptionInfo *exception)
1589 %
1590 %  A description of each parameter follows.
1591 %
1592 %    o file:  An pointer to the output FILE.
1593 %
1594 %    o exception: return any errors or warnings in this structure.
1595 %
1596 */
ListThresholdMaps(FILE * file,ExceptionInfo * exception)1597 MagickExport MagickBooleanType ListThresholdMaps(FILE *file,
1598   ExceptionInfo *exception)
1599 {
1600   const StringInfo
1601     *option;
1602 
1603   LinkedListInfo
1604     *options;
1605 
1606   MagickStatusType
1607     status;
1608 
1609   status=MagickTrue;
1610   if (file == (FILE *) NULL)
1611     file=stdout;
1612   options=GetConfigureOptions(ThresholdsFilename,exception);
1613   (void) FormatLocaleFile(file,
1614     "\n   Threshold Maps for Ordered Dither Operations\n");
1615   option=(const StringInfo *) GetNextValueInLinkedList(options);
1616   while (option != (const StringInfo *) NULL)
1617   {
1618     (void) FormatLocaleFile(file,"\nPath: %s\n\n",GetStringInfoPath(option));
1619     status&=ListThresholdMapFile(file,(const char *) GetStringInfoDatum(option),
1620       GetStringInfoPath(option),exception);
1621     option=(const StringInfo *) GetNextValueInLinkedList(options);
1622   }
1623   options=DestroyConfigureOptions(options);
1624   return(status != 0 ? MagickTrue : MagickFalse);
1625 }
1626 
1627 /*
1628 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1629 %                                                                             %
1630 %                                                                             %
1631 %                                                                             %
1632 %     O r d e r e d D i t h e r I m a g e                                     %
1633 %                                                                             %
1634 %                                                                             %
1635 %                                                                             %
1636 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1637 %
1638 %  OrderedDitherImage() will perform a ordered dither based on a number
1639 %  of pre-defined dithering threshold maps, but over multiple intensity
1640 %  levels, which can be different for different channels, according to the
1641 %  input argument.
1642 %
1643 %  The format of the OrderedDitherImage method is:
1644 %
1645 %      MagickBooleanType OrderedDitherImage(Image *image,
1646 %        const char *threshold_map,ExceptionInfo *exception)
1647 %
1648 %  A description of each parameter follows:
1649 %
1650 %    o image: the image.
1651 %
1652 %    o threshold_map: A string containing the name of the threshold dither
1653 %      map to use, followed by zero or more numbers representing the number
1654 %      of color levels tho dither between.
1655 %
1656 %      Any level number less than 2 will be equivalent to 2, and means only
1657 %      binary dithering will be applied to each color channel.
1658 %
1659 %      No numbers also means a 2 level (bitmap) dither will be applied to all
1660 %      channels, while a single number is the number of levels applied to each
1661 %      channel in sequence.  More numbers will be applied in turn to each of
1662 %      the color channels.
1663 %
1664 %      For example: "o3x3,6" will generate a 6 level posterization of the
1665 %      image with a ordered 3x3 diffused pixel dither being applied between
1666 %      each level. While checker,8,8,4 will produce a 332 colormaped image
1667 %      with only a single checkerboard hash pattern (50% grey) between each
1668 %      color level, to basically double the number of color levels with
1669 %      a bare minimim of dithering.
1670 %
1671 %    o exception: return any errors or warnings in this structure.
1672 %
1673 */
OrderedDitherImage(Image * image,const char * threshold_map,ExceptionInfo * exception)1674 MagickExport MagickBooleanType OrderedDitherImage(Image *image,
1675   const char *threshold_map,ExceptionInfo *exception)
1676 {
1677 #define DitherImageTag  "Dither/Image"
1678 
1679   CacheView
1680     *image_view;
1681 
1682   char
1683     token[MagickPathExtent];
1684 
1685   const char
1686     *p;
1687 
1688   double
1689     levels[CompositePixelChannel];
1690 
1691   MagickBooleanType
1692     status;
1693 
1694   MagickOffsetType
1695     progress;
1696 
1697   register ssize_t
1698     i;
1699 
1700   ssize_t
1701     y;
1702 
1703   ThresholdMap
1704     *map;
1705 
1706   assert(image != (Image *) NULL);
1707   assert(image->signature == MagickCoreSignature);
1708   if (image->debug != MagickFalse)
1709     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1710   assert(exception != (ExceptionInfo *) NULL);
1711   assert(exception->signature == MagickCoreSignature);
1712   if (threshold_map == (const char *) NULL)
1713     return(MagickTrue);
1714   p=(char *) threshold_map;
1715   while (((isspace((int) ((unsigned char) *p)) != 0) || (*p == ',')) &&
1716          (*p != '\0'))
1717     p++;
1718   threshold_map=p;
1719   while (((isspace((int) ((unsigned char) *p)) == 0) && (*p != ',')) &&
1720          (*p != '\0'))
1721   {
1722     if ((p-threshold_map) >= (MagickPathExtent-1))
1723       break;
1724     token[p-threshold_map]=(*p);
1725     p++;
1726   }
1727   token[p-threshold_map]='\0';
1728   map=GetThresholdMap(token,exception);
1729   if (map == (ThresholdMap *) NULL)
1730     {
1731       (void) ThrowMagickException(exception,GetMagickModule(),OptionError,
1732         "InvalidArgument","%s : '%s'","ordered-dither",threshold_map);
1733       return(MagickFalse);
1734     }
1735   for (i=0; i < MaxPixelChannels; i++)
1736     levels[i]=2.0;
1737   p=strchr((char *) threshold_map,',');
1738   if ((p != (char *) NULL) && (isdigit((int) ((unsigned char) *(++p))) != 0))
1739     {
1740       GetNextToken(p,&p,MagickPathExtent,token);
1741       for (i=0; (i < MaxPixelChannels); i++)
1742         levels[i]=StringToDouble(token,(char **) NULL);
1743       for (i=0; (*p != '\0') && (i < MaxPixelChannels); i++)
1744       {
1745         GetNextToken(p,&p,MagickPathExtent,token);
1746         if (*token == ',')
1747           GetNextToken(p,&p,MagickPathExtent,token);
1748         levels[i]=StringToDouble(token,(char **) NULL);
1749       }
1750     }
1751   for (i=0; i < MaxPixelChannels; i++)
1752     if (fabs(levels[i]) >= 1)
1753       levels[i]-=1.0;
1754   if (SetImageStorageClass(image,DirectClass,exception) == MagickFalse)
1755     return(MagickFalse);
1756   status=MagickTrue;
1757   progress=0;
1758   image_view=AcquireAuthenticCacheView(image,exception);
1759 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1760   #pragma omp parallel for schedule(static) shared(progress,status) \
1761     magick_number_threads(image,image,image->rows,1)
1762 #endif
1763   for (y=0; y < (ssize_t) image->rows; y++)
1764   {
1765     register ssize_t
1766       x;
1767 
1768     register Quantum
1769       *magick_restrict q;
1770 
1771     if (status == MagickFalse)
1772       continue;
1773     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
1774     if (q == (Quantum *) NULL)
1775       {
1776         status=MagickFalse;
1777         continue;
1778       }
1779     for (x=0; x < (ssize_t) image->columns; x++)
1780     {
1781       register ssize_t
1782         i;
1783 
1784       ssize_t
1785         n;
1786 
1787       n=0;
1788       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
1789       {
1790         ssize_t
1791           level,
1792           threshold;
1793 
1794         PixelChannel channel = GetPixelChannelChannel(image,i);
1795         PixelTrait traits = GetPixelChannelTraits(image,channel);
1796         if ((traits & UpdatePixelTrait) == 0)
1797           continue;
1798         if (fabs(levels[n]) < MagickEpsilon)
1799           {
1800             n++;
1801             continue;
1802           }
1803         threshold=(ssize_t) (QuantumScale*q[i]*(levels[n]*(map->divisor-1)+1));
1804         level=threshold/(map->divisor-1);
1805         threshold-=level*(map->divisor-1);
1806         q[i]=ClampToQuantum((double) (level+(threshold >=
1807           map->levels[(x % map->width)+map->width*(y % map->height)]))*
1808           QuantumRange/levels[n]);
1809         n++;
1810       }
1811       q+=GetPixelChannels(image);
1812     }
1813     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
1814       status=MagickFalse;
1815     if (image->progress_monitor != (MagickProgressMonitor) NULL)
1816       {
1817         MagickBooleanType
1818           proceed;
1819 
1820 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1821         #pragma omp atomic
1822 #endif
1823         progress++;
1824         proceed=SetImageProgress(image,DitherImageTag,progress,image->rows);
1825         if (proceed == MagickFalse)
1826           status=MagickFalse;
1827       }
1828   }
1829   image_view=DestroyCacheView(image_view);
1830   map=DestroyThresholdMap(map);
1831   return(MagickTrue);
1832 }
1833 
1834 /*
1835 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1836 %                                                                             %
1837 %                                                                             %
1838 %                                                                             %
1839 %     P e r c e p t i b l e I m a g e                                         %
1840 %                                                                             %
1841 %                                                                             %
1842 %                                                                             %
1843 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1844 %
1845 %  PerceptibleImage() set each pixel whose value is less than |epsilon| to
1846 %  epsilon or -epsilon (whichever is closer) otherwise the pixel value remains
1847 %  unchanged.
1848 %
1849 %  The format of the PerceptibleImage method is:
1850 %
1851 %      MagickBooleanType PerceptibleImage(Image *image,const double epsilon,
1852 %        ExceptionInfo *exception)
1853 %
1854 %  A description of each parameter follows:
1855 %
1856 %    o image: the image.
1857 %
1858 %    o epsilon: the epsilon threshold (e.g. 1.0e-9).
1859 %
1860 %    o exception: return any errors or warnings in this structure.
1861 %
1862 */
1863 
PerceptibleThreshold(const Quantum quantum,const double epsilon)1864 static inline Quantum PerceptibleThreshold(const Quantum quantum,
1865   const double epsilon)
1866 {
1867   double
1868     sign;
1869 
1870   sign=(double) quantum < 0.0 ? -1.0 : 1.0;
1871   if ((sign*quantum) >= epsilon)
1872     return(quantum);
1873   return((Quantum) (sign*epsilon));
1874 }
1875 
PerceptibleImage(Image * image,const double epsilon,ExceptionInfo * exception)1876 MagickExport MagickBooleanType PerceptibleImage(Image *image,
1877   const double epsilon,ExceptionInfo *exception)
1878 {
1879 #define PerceptibleImageTag  "Perceptible/Image"
1880 
1881   CacheView
1882     *image_view;
1883 
1884   MagickBooleanType
1885     status;
1886 
1887   MagickOffsetType
1888     progress;
1889 
1890   ssize_t
1891     y;
1892 
1893   assert(image != (Image *) NULL);
1894   assert(image->signature == MagickCoreSignature);
1895   if (image->debug != MagickFalse)
1896     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1897   if (image->storage_class == PseudoClass)
1898     {
1899       register ssize_t
1900         i;
1901 
1902       register PixelInfo
1903         *magick_restrict q;
1904 
1905       q=image->colormap;
1906       for (i=0; i < (ssize_t) image->colors; i++)
1907       {
1908         q->red=(double) PerceptibleThreshold(ClampToQuantum(q->red),
1909           epsilon);
1910         q->green=(double) PerceptibleThreshold(ClampToQuantum(q->green),
1911           epsilon);
1912         q->blue=(double) PerceptibleThreshold(ClampToQuantum(q->blue),
1913           epsilon);
1914         q->alpha=(double) PerceptibleThreshold(ClampToQuantum(q->alpha),
1915           epsilon);
1916         q++;
1917       }
1918       return(SyncImage(image,exception));
1919     }
1920   /*
1921     Perceptible image.
1922   */
1923   status=MagickTrue;
1924   progress=0;
1925   image_view=AcquireAuthenticCacheView(image,exception);
1926 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1927   #pragma omp parallel for schedule(static) shared(progress,status) \
1928     magick_number_threads(image,image,image->rows,1)
1929 #endif
1930   for (y=0; y < (ssize_t) image->rows; y++)
1931   {
1932     register ssize_t
1933       x;
1934 
1935     register Quantum
1936       *magick_restrict q;
1937 
1938     if (status == MagickFalse)
1939       continue;
1940     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
1941     if (q == (Quantum *) NULL)
1942       {
1943         status=MagickFalse;
1944         continue;
1945       }
1946     for (x=0; x < (ssize_t) image->columns; x++)
1947     {
1948       register ssize_t
1949         i;
1950 
1951       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
1952       {
1953         PixelChannel channel = GetPixelChannelChannel(image,i);
1954         PixelTrait traits = GetPixelChannelTraits(image,channel);
1955         if (traits == UndefinedPixelTrait)
1956           continue;
1957         q[i]=PerceptibleThreshold(q[i],epsilon);
1958       }
1959       q+=GetPixelChannels(image);
1960     }
1961     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
1962       status=MagickFalse;
1963     if (image->progress_monitor != (MagickProgressMonitor) NULL)
1964       {
1965         MagickBooleanType
1966           proceed;
1967 
1968 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1969         #pragma omp atomic
1970 #endif
1971         progress++;
1972         proceed=SetImageProgress(image,PerceptibleImageTag,progress,
1973           image->rows);
1974         if (proceed == MagickFalse)
1975           status=MagickFalse;
1976       }
1977   }
1978   image_view=DestroyCacheView(image_view);
1979   return(status);
1980 }
1981 
1982 /*
1983 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1984 %                                                                             %
1985 %                                                                             %
1986 %                                                                             %
1987 %     R a n d o m T h r e s h o l d I m a g e                                 %
1988 %                                                                             %
1989 %                                                                             %
1990 %                                                                             %
1991 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1992 %
1993 %  RandomThresholdImage() changes the value of individual pixels based on the
1994 %  intensity of each pixel compared to a random threshold.  The result is a
1995 %  low-contrast, two color image.
1996 %
1997 %  The format of the RandomThresholdImage method is:
1998 %
1999 %      MagickBooleanType RandomThresholdImage(Image *image,
2000 %        const char *thresholds,ExceptionInfo *exception)
2001 %
2002 %  A description of each parameter follows:
2003 %
2004 %    o image: the image.
2005 %
2006 %    o low,high: Specify the high and low thresholds. These values range from
2007 %      0 to QuantumRange.
2008 %
2009 %    o exception: return any errors or warnings in this structure.
2010 %
2011 */
RandomThresholdImage(Image * image,const double min_threshold,const double max_threshold,ExceptionInfo * exception)2012 MagickExport MagickBooleanType RandomThresholdImage(Image *image,
2013   const double min_threshold, const double max_threshold,ExceptionInfo *exception)
2014 {
2015 #define ThresholdImageTag  "Threshold/Image"
2016 
2017   CacheView
2018     *image_view;
2019 
2020   MagickBooleanType
2021     status;
2022 
2023   MagickOffsetType
2024     progress;
2025 
2026   PixelInfo
2027     threshold;
2028 
2029   RandomInfo
2030     **magick_restrict random_info;
2031 
2032   ssize_t
2033     y;
2034 
2035 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2036   unsigned long
2037     key;
2038 #endif
2039 
2040   assert(image != (Image *) NULL);
2041   assert(image->signature == MagickCoreSignature);
2042   if (image->debug != MagickFalse)
2043     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2044   assert(exception != (ExceptionInfo *) NULL);
2045   assert(exception->signature == MagickCoreSignature);
2046   if (SetImageStorageClass(image,DirectClass,exception) == MagickFalse)
2047     return(MagickFalse);
2048   GetPixelInfo(image,&threshold);
2049   /*
2050     Random threshold image.
2051   */
2052   status=MagickTrue;
2053   progress=0;
2054   random_info=AcquireRandomInfoThreadSet();
2055   image_view=AcquireAuthenticCacheView(image,exception);
2056 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2057   key=GetRandomSecretKey(random_info[0]);
2058   #pragma omp parallel for schedule(static) shared(progress,status) \
2059     magick_number_threads(image,image,image->rows,key == ~0UL)
2060 #endif
2061   for (y=0; y < (ssize_t) image->rows; y++)
2062   {
2063     const int
2064       id = GetOpenMPThreadId();
2065 
2066     register Quantum
2067       *magick_restrict q;
2068 
2069     register ssize_t
2070       x;
2071 
2072     if (status == MagickFalse)
2073       continue;
2074     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
2075     if (q == (Quantum *) NULL)
2076       {
2077         status=MagickFalse;
2078         continue;
2079       }
2080     for (x=0; x < (ssize_t) image->columns; x++)
2081     {
2082       register ssize_t
2083         i;
2084 
2085       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
2086       {
2087         double
2088           threshold;
2089 
2090         PixelChannel channel = GetPixelChannelChannel(image,i);
2091         PixelTrait traits = GetPixelChannelTraits(image,channel);
2092         if ((traits & UpdatePixelTrait) == 0)
2093           continue;
2094         if ((double) q[i] < min_threshold)
2095           threshold=min_threshold;
2096         else
2097           if ((double) q[i] > max_threshold)
2098             threshold=max_threshold;
2099           else
2100             threshold=(double) (QuantumRange*
2101               GetPseudoRandomValue(random_info[id]));
2102         q[i]=(double) q[i] <= threshold ? 0 : QuantumRange;
2103       }
2104       q+=GetPixelChannels(image);
2105     }
2106     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
2107       status=MagickFalse;
2108     if (image->progress_monitor != (MagickProgressMonitor) NULL)
2109       {
2110         MagickBooleanType
2111           proceed;
2112 
2113 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2114         #pragma omp atomic
2115 #endif
2116         progress++;
2117         proceed=SetImageProgress(image,ThresholdImageTag,progress,
2118           image->rows);
2119         if (proceed == MagickFalse)
2120           status=MagickFalse;
2121       }
2122   }
2123   image_view=DestroyCacheView(image_view);
2124   random_info=DestroyRandomInfoThreadSet(random_info);
2125   return(status);
2126 }
2127 
2128 /*
2129 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2130 %                                                                             %
2131 %                                                                             %
2132 %                                                                             %
2133 %     R a n g e T h r e s h o l d I m a g e                                   %
2134 %                                                                             %
2135 %                                                                             %
2136 %                                                                             %
2137 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2138 %
2139 %  RangeThresholdImage() applies soft and hard thresholding.
2140 %
2141 %  The format of the RangeThresholdImage method is:
2142 %
2143 %      MagickBooleanType RangeThresholdImage(Image *image,
2144 %        const double low_black,const double low_white,const double high_white,
2145 %        const double high_black,ExceptionInfo *exception)
2146 %
2147 %  A description of each parameter follows:
2148 %
2149 %    o image: the image.
2150 %
2151 %    o low_black: Define the minimum threshold value.
2152 %
2153 %    o low_white: Define the maximum threshold value.
2154 %
2155 %    o high_white: Define the minimum threshold value.
2156 %
2157 %    o low_white: Define the maximum threshold value.
2158 %
2159 %    o exception: return any errors or warnings in this structure.
2160 %
2161 */
RangeThresholdImage(Image * image,const double low_black,const double low_white,const double high_white,const double high_black,ExceptionInfo * exception)2162 MagickExport MagickBooleanType RangeThresholdImage(Image *image,
2163   const double low_black,const double low_white,const double high_white,
2164   const double high_black,ExceptionInfo *exception)
2165 {
2166 #define ThresholdImageTag  "Threshold/Image"
2167 
2168   CacheView
2169     *image_view;
2170 
2171   MagickBooleanType
2172     status;
2173 
2174   MagickOffsetType
2175     progress;
2176 
2177   ssize_t
2178     y;
2179 
2180   assert(image != (Image *) NULL);
2181   assert(image->signature == MagickCoreSignature);
2182   if (image->debug != MagickFalse)
2183     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2184   if (SetImageStorageClass(image,DirectClass,exception) == MagickFalse)
2185     return(MagickFalse);
2186   if (IsGrayColorspace(image->colorspace) != MagickFalse)
2187     (void) TransformImageColorspace(image,sRGBColorspace,exception);
2188   /*
2189     Range threshold image.
2190   */
2191   status=MagickTrue;
2192   progress=0;
2193   image_view=AcquireAuthenticCacheView(image,exception);
2194 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2195   #pragma omp parallel for schedule(static) shared(progress,status) \
2196     magick_number_threads(image,image,image->rows,1)
2197 #endif
2198   for (y=0; y < (ssize_t) image->rows; y++)
2199   {
2200     register ssize_t
2201       x;
2202 
2203     register Quantum
2204       *magick_restrict q;
2205 
2206     if (status == MagickFalse)
2207       continue;
2208     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
2209     if (q == (Quantum *) NULL)
2210       {
2211         status=MagickFalse;
2212         continue;
2213       }
2214     for (x=0; x < (ssize_t) image->columns; x++)
2215     {
2216       double
2217         pixel;
2218 
2219       register ssize_t
2220         i;
2221 
2222       pixel=GetPixelIntensity(image,q);
2223       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
2224       {
2225         PixelChannel channel = GetPixelChannelChannel(image,i);
2226         PixelTrait traits = GetPixelChannelTraits(image,channel);
2227         if ((traits & UpdatePixelTrait) == 0)
2228           continue;
2229         if (image->channel_mask != DefaultChannels)
2230           pixel=(double) q[i];
2231         if (pixel < low_black)
2232           q[i]=0;
2233         else
2234           if ((pixel >= low_black) && (pixel < low_white))
2235             q[i]=ClampToQuantum(QuantumRange*
2236               PerceptibleReciprocal(low_white-low_black)*(pixel-low_black));
2237           else
2238             if ((pixel >= low_white) && (pixel <= high_white))
2239               q[i]=QuantumRange;
2240             else
2241               if ((pixel > high_white) && (pixel <= high_black))
2242                 q[i]=ClampToQuantum(QuantumRange*PerceptibleReciprocal(
2243                   high_black-high_white)*(high_black-pixel));
2244               else
2245                 if (pixel > high_black)
2246                   q[i]=0;
2247                 else
2248                   q[i]=0;
2249       }
2250       q+=GetPixelChannels(image);
2251     }
2252     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
2253       status=MagickFalse;
2254     if (image->progress_monitor != (MagickProgressMonitor) NULL)
2255       {
2256         MagickBooleanType
2257           proceed;
2258 
2259 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2260         #pragma omp atomic
2261 #endif
2262         progress++;
2263         proceed=SetImageProgress(image,ThresholdImageTag,progress,
2264           image->rows);
2265         if (proceed == MagickFalse)
2266           status=MagickFalse;
2267       }
2268   }
2269   image_view=DestroyCacheView(image_view);
2270   return(status);
2271 }
2272 
2273 /*
2274 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2275 %                                                                             %
2276 %                                                                             %
2277 %                                                                             %
2278 %     W h i t e T h r e s h o l d I m a g e                                   %
2279 %                                                                             %
2280 %                                                                             %
2281 %                                                                             %
2282 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2283 %
2284 %  WhiteThresholdImage() is like ThresholdImage() but forces all pixels above
2285 %  the threshold into white while leaving all pixels at or below the threshold
2286 %  unchanged.
2287 %
2288 %  The format of the WhiteThresholdImage method is:
2289 %
2290 %      MagickBooleanType WhiteThresholdImage(Image *image,
2291 %        const char *threshold,ExceptionInfo *exception)
2292 %
2293 %  A description of each parameter follows:
2294 %
2295 %    o image: the image.
2296 %
2297 %    o threshold: Define the threshold value.
2298 %
2299 %    o exception: return any errors or warnings in this structure.
2300 %
2301 */
WhiteThresholdImage(Image * image,const char * thresholds,ExceptionInfo * exception)2302 MagickExport MagickBooleanType WhiteThresholdImage(Image *image,
2303   const char *thresholds,ExceptionInfo *exception)
2304 {
2305 #define ThresholdImageTag  "Threshold/Image"
2306 
2307   CacheView
2308     *image_view;
2309 
2310   GeometryInfo
2311     geometry_info;
2312 
2313   MagickBooleanType
2314     status;
2315 
2316   MagickOffsetType
2317     progress;
2318 
2319   PixelInfo
2320     threshold;
2321 
2322   MagickStatusType
2323     flags;
2324 
2325   ssize_t
2326     y;
2327 
2328   assert(image != (Image *) NULL);
2329   assert(image->signature == MagickCoreSignature);
2330   if (image->debug != MagickFalse)
2331     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2332   if (thresholds == (const char *) NULL)
2333     return(MagickTrue);
2334   if (SetImageStorageClass(image,DirectClass,exception) == MagickFalse)
2335     return(MagickFalse);
2336   if (IsGrayColorspace(image->colorspace) != MagickFalse)
2337     (void) TransformImageColorspace(image,sRGBColorspace,exception);
2338   GetPixelInfo(image,&threshold);
2339   flags=ParseGeometry(thresholds,&geometry_info);
2340   threshold.red=geometry_info.rho;
2341   threshold.green=geometry_info.rho;
2342   threshold.blue=geometry_info.rho;
2343   threshold.black=geometry_info.rho;
2344   threshold.alpha=100.0;
2345   if ((flags & SigmaValue) != 0)
2346     threshold.green=geometry_info.sigma;
2347   if ((flags & XiValue) != 0)
2348     threshold.blue=geometry_info.xi;
2349   if ((flags & PsiValue) != 0)
2350     threshold.alpha=geometry_info.psi;
2351   if (threshold.colorspace == CMYKColorspace)
2352     {
2353       if ((flags & PsiValue) != 0)
2354         threshold.black=geometry_info.psi;
2355       if ((flags & ChiValue) != 0)
2356         threshold.alpha=geometry_info.chi;
2357     }
2358   if ((flags & PercentValue) != 0)
2359     {
2360       threshold.red*=(MagickRealType) (QuantumRange/100.0);
2361       threshold.green*=(MagickRealType) (QuantumRange/100.0);
2362       threshold.blue*=(MagickRealType) (QuantumRange/100.0);
2363       threshold.black*=(MagickRealType) (QuantumRange/100.0);
2364       threshold.alpha*=(MagickRealType) (QuantumRange/100.0);
2365     }
2366   /*
2367     White threshold image.
2368   */
2369   status=MagickTrue;
2370   progress=0;
2371   image_view=AcquireAuthenticCacheView(image,exception);
2372 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2373   #pragma omp parallel for schedule(static) shared(progress,status) \
2374     magick_number_threads(image,image,image->rows,1)
2375 #endif
2376   for (y=0; y < (ssize_t) image->rows; y++)
2377   {
2378     register ssize_t
2379       x;
2380 
2381     register Quantum
2382       *magick_restrict q;
2383 
2384     if (status == MagickFalse)
2385       continue;
2386     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
2387     if (q == (Quantum *) NULL)
2388       {
2389         status=MagickFalse;
2390         continue;
2391       }
2392     for (x=0; x < (ssize_t) image->columns; x++)
2393     {
2394       double
2395         pixel;
2396 
2397       register ssize_t
2398         i;
2399 
2400       pixel=GetPixelIntensity(image,q);
2401       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
2402       {
2403         PixelChannel channel = GetPixelChannelChannel(image,i);
2404         PixelTrait traits = GetPixelChannelTraits(image,channel);
2405         if ((traits & UpdatePixelTrait) == 0)
2406           continue;
2407         if (image->channel_mask != DefaultChannels)
2408           pixel=(double) q[i];
2409         if (pixel > GetPixelInfoChannel(&threshold,channel))
2410           q[i]=QuantumRange;
2411       }
2412       q+=GetPixelChannels(image);
2413     }
2414     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
2415       status=MagickFalse;
2416     if (image->progress_monitor != (MagickProgressMonitor) NULL)
2417       {
2418         MagickBooleanType
2419           proceed;
2420 
2421 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2422         #pragma omp atomic
2423 #endif
2424         progress++;
2425         proceed=SetImageProgress(image,ThresholdImageTag,progress,image->rows);
2426         if (proceed == MagickFalse)
2427           status=MagickFalse;
2428       }
2429   }
2430   image_view=DestroyCacheView(image_view);
2431   return(status);
2432 }
2433