MagickCore 6.9.13-26
Convert, Edit, Or Compose Bitmap Images
Loading...
Searching...
No Matches
compare.c
1/*
2%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3% %
4% %
5% %
6% CCCC OOO M M PPPP AAA RRRR EEEEE %
7% C O O MM MM P P A A R R E %
8% C O O M M M PPPP AAAAA RRRR EEE %
9% C O O M M P A A R R E %
10% CCCC OOO M M P A A R R EEEEE %
11% %
12% %
13% MagickCore Image Comparison Methods %
14% %
15% Software Design %
16% Cristy %
17% December 2003 %
18% %
19% %
20% Copyright 1999 ImageMagick Studio LLC, a non-profit organization dedicated %
21% 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 "magick/studio.h"
44#include "magick/artifact.h"
45#include "magick/attribute.h"
46#include "magick/cache-view.h"
47#include "magick/channel.h"
48#include "magick/client.h"
49#include "magick/color.h"
50#include "magick/color-private.h"
51#include "magick/colorspace.h"
52#include "magick/colorspace-private.h"
53#include "magick/compare.h"
54#include "magick/composite-private.h"
55#include "magick/constitute.h"
56#include "magick/exception-private.h"
57#include "magick/geometry.h"
58#include "magick/image-private.h"
59#include "magick/list.h"
60#include "magick/log.h"
61#include "magick/memory_.h"
62#include "magick/monitor.h"
63#include "magick/monitor-private.h"
64#include "magick/option.h"
65#include "magick/pixel-private.h"
66#include "magick/property.h"
67#include "magick/resource_.h"
68#include "magick/statistic-private.h"
69#include "magick/string_.h"
70#include "magick/string-private.h"
71#include "magick/statistic.h"
72#include "magick/thread-private.h"
73#include "magick/transform.h"
74#include "magick/utility.h"
75#include "magick/version.h"
76
77/*
78%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
79% %
80% %
81% %
82% C o m p a r e I m a g e C h a n n e l s %
83% %
84% %
85% %
86%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
87%
88% CompareImageChannels() compares one or more image channels of an image
89% to a reconstructed image and returns the difference image.
90%
91% The format of the CompareImageChannels method is:
92%
93% Image *CompareImageChannels(const Image *image,
94% const Image *reconstruct_image,const ChannelType channel,
95% const MetricType metric,double *distortion,ExceptionInfo *exception)
96%
97% A description of each parameter follows:
98%
99% o image: the image.
100%
101% o reconstruct_image: the reconstruct image.
102%
103% o channel: the channel.
104%
105% o metric: the metric.
106%
107% o distortion: the computed distortion between the images.
108%
109% o exception: return any errors or warnings in this structure.
110%
111*/
112
113MagickExport Image *CompareImages(Image *image,const Image *reconstruct_image,
114 const MetricType metric,double *distortion,ExceptionInfo *exception)
115{
116 Image
117 *highlight_image;
118
119 highlight_image=CompareImageChannels(image,reconstruct_image,
120 CompositeChannels,metric,distortion,exception);
121 return(highlight_image);
122}
123
124static size_t GetNumberChannels(const Image *image,const ChannelType channel)
125{
126 size_t
127 channels;
128
129 channels=0;
130 if ((channel & RedChannel) != 0)
131 channels++;
132 if ((channel & GreenChannel) != 0)
133 channels++;
134 if ((channel & BlueChannel) != 0)
135 channels++;
136 if (((channel & OpacityChannel) != 0) && (image->matte != MagickFalse))
137 channels++;
138 if (((channel & IndexChannel) != 0) && (image->colorspace == CMYKColorspace))
139 channels++;
140 return(channels == 0 ? 1UL : channels);
141}
142
143static void SetImageDistortionBounds(const Image *image,
144 const Image *reconstruct_image,size_t *columns,size_t *rows)
145{
146 const char
147 *artifact;
148
149 *columns=MagickMax(image->columns,reconstruct_image->columns);
150 *rows=MagickMax(image->rows,reconstruct_image->rows);
151 artifact=GetImageArtifact(image,"compare:virtual-pixels");
152 if ((artifact != (const char *) NULL) &&
153 (IsStringTrue(artifact) == MagickFalse))
154 {
155 *columns=MagickMin(image->columns,reconstruct_image->columns);
156 *rows=MagickMin(image->rows,reconstruct_image->rows);
157 }
158}
159
160static inline MagickBooleanType ValidateImageMorphology(
161 const Image *magick_restrict image,
162 const Image *magick_restrict reconstruct_image)
163{
164 /*
165 Does the image match the reconstructed image morphology?
166 */
167 if (GetNumberChannels(image,DefaultChannels) !=
168 GetNumberChannels(reconstruct_image,DefaultChannels))
169 return(MagickFalse);
170 return(MagickTrue);
171}
172
173MagickExport Image *CompareImageChannels(Image *image,
174 const Image *reconstruct_image,const ChannelType channel,
175 const MetricType metric,double *distortion,ExceptionInfo *exception)
176{
177 CacheView
178 *highlight_view,
179 *image_view,
180 *reconstruct_view;
181
182 const char
183 *artifact;
184
185 double
186 fuzz;
187
188 Image
189 *clone_image,
190 *difference_image,
191 *highlight_image;
192
193 MagickBooleanType
194 status;
195
196 MagickPixelPacket
197 highlight,
198 lowlight,
199 zero;
200
201 size_t
202 columns,
203 rows;
204
205 ssize_t
206 y;
207
208 assert(image != (Image *) NULL);
209 assert(image->signature == MagickCoreSignature);
210 assert(reconstruct_image != (const Image *) NULL);
211 assert(reconstruct_image->signature == MagickCoreSignature);
212 assert(distortion != (double *) NULL);
213 if (IsEventLogging() != MagickFalse)
214 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
215 *distortion=0.0;
216 if (metric != PerceptualHashErrorMetric)
217 if (ValidateImageMorphology(image,reconstruct_image) == MagickFalse)
218 ThrowImageException(ImageError,"ImageMorphologyDiffers");
219 status=GetImageChannelDistortion(image,reconstruct_image,channel,metric,
220 distortion,exception);
221 if (status == MagickFalse)
222 return((Image *) NULL);
223 clone_image=CloneImage(image,0,0,MagickTrue,exception);
224 if (clone_image == (Image *) NULL)
225 return((Image *) NULL);
226 (void) SetImageMask(clone_image,(Image *) NULL);
227 difference_image=CloneImage(clone_image,0,0,MagickTrue,exception);
228 clone_image=DestroyImage(clone_image);
229 if (difference_image == (Image *) NULL)
230 return((Image *) NULL);
231 (void) SetImageAlphaChannel(difference_image,OpaqueAlphaChannel);
232 SetImageDistortionBounds(image,reconstruct_image,&columns,&rows);
233 highlight_image=CloneImage(image,columns,rows,MagickTrue,exception);
234 if (highlight_image == (Image *) NULL)
235 {
236 difference_image=DestroyImage(difference_image);
237 return((Image *) NULL);
238 }
239 if (SetImageStorageClass(highlight_image,DirectClass) == MagickFalse)
240 {
241 InheritException(exception,&highlight_image->exception);
242 difference_image=DestroyImage(difference_image);
243 highlight_image=DestroyImage(highlight_image);
244 return((Image *) NULL);
245 }
246 (void) SetImageMask(highlight_image,(Image *) NULL);
247 (void) SetImageAlphaChannel(highlight_image,OpaqueAlphaChannel);
248 (void) QueryMagickColor("#f1001ecc",&highlight,exception);
249 artifact=GetImageArtifact(image,"compare:highlight-color");
250 if (artifact != (const char *) NULL)
251 (void) QueryMagickColor(artifact,&highlight,exception);
252 (void) QueryMagickColor("#ffffffcc",&lowlight,exception);
253 artifact=GetImageArtifact(image,"compare:lowlight-color");
254 if (artifact != (const char *) NULL)
255 (void) QueryMagickColor(artifact,&lowlight,exception);
256 if (highlight_image->colorspace == CMYKColorspace)
257 {
258 ConvertRGBToCMYK(&highlight);
259 ConvertRGBToCMYK(&lowlight);
260 }
261 /*
262 Generate difference image.
263 */
264 status=MagickTrue;
265 fuzz=GetFuzzyColorDistance(image,reconstruct_image);
266 GetMagickPixelPacket(image,&zero);
267 image_view=AcquireVirtualCacheView(image,exception);
268 reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
269 highlight_view=AcquireAuthenticCacheView(highlight_image,exception);
270#if defined(MAGICKCORE_OPENMP_SUPPORT)
271 #pragma omp parallel for schedule(static) shared(status) \
272 magick_number_threads(image,highlight_image,rows,1)
273#endif
274 for (y=0; y < (ssize_t) rows; y++)
275 {
276 MagickBooleanType
277 sync;
278
279 MagickPixelPacket
280 pixel,
281 reconstruct_pixel;
282
283 const IndexPacket
284 *magick_restrict indexes,
285 *magick_restrict reconstruct_indexes;
286
287 const PixelPacket
288 *magick_restrict p,
289 *magick_restrict q;
290
291 IndexPacket
292 *magick_restrict highlight_indexes;
293
294 PixelPacket
295 *magick_restrict r;
296
297 ssize_t
298 x;
299
300 if (status == MagickFalse)
301 continue;
302 p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
303 q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
304 r=QueueCacheViewAuthenticPixels(highlight_view,0,y,columns,1,exception);
305 if ((p == (const PixelPacket *) NULL) ||
306 (q == (const PixelPacket *) NULL) || (r == (PixelPacket *) NULL))
307 {
308 status=MagickFalse;
309 continue;
310 }
311 indexes=GetCacheViewVirtualIndexQueue(image_view);
312 reconstruct_indexes=GetCacheViewVirtualIndexQueue(reconstruct_view);
313 highlight_indexes=GetCacheViewAuthenticIndexQueue(highlight_view);
314 pixel=zero;
315 reconstruct_pixel=zero;
316 for (x=0; x < (ssize_t) columns; x++)
317 {
318 MagickStatusType
319 difference;
320
321 SetMagickPixelPacket(image,p,indexes == (IndexPacket *) NULL ? NULL :
322 indexes+x,&pixel);
323 SetMagickPixelPacket(reconstruct_image,q,reconstruct_indexes ==
324 (IndexPacket *) NULL ? NULL : reconstruct_indexes+x,&reconstruct_pixel);
325 difference=MagickFalse;
326 if (channel == CompositeChannels)
327 {
328 if (IsMagickColorSimilar(&pixel,&reconstruct_pixel) == MagickFalse)
329 difference=MagickTrue;
330 }
331 else
332 {
333 double
334 Da,
335 distance,
336 pixel,
337 Sa;
338
339 Sa=QuantumScale*(image->matte != MagickFalse ? (double)
340 GetPixelAlpha(p) : ((double) QuantumRange-(double) OpaqueOpacity));
341 Da=QuantumScale*(image->matte != MagickFalse ? (double)
342 GetPixelAlpha(q) : ((double) QuantumRange-(double) OpaqueOpacity));
343 if ((channel & RedChannel) != 0)
344 {
345 pixel=Sa*(double) GetPixelRed(p)-Da*(double) GetPixelRed(q);
346 distance=pixel*pixel;
347 if (distance >= fuzz)
348 difference=MagickTrue;
349 }
350 if ((channel & GreenChannel) != 0)
351 {
352 pixel=Sa*(double) GetPixelGreen(p)-Da*(double) GetPixelGreen(q);
353 distance=pixel*pixel;
354 if (distance >= fuzz)
355 difference=MagickTrue;
356 }
357 if ((channel & BlueChannel) != 0)
358 {
359 pixel=Sa*(double) GetPixelBlue(p)-Da*(double) GetPixelBlue(q);
360 distance=pixel*pixel;
361 if (distance >= fuzz)
362 difference=MagickTrue;
363 }
364 if (((channel & OpacityChannel) != 0) &&
365 (image->matte != MagickFalse))
366 {
367 pixel=(double) GetPixelOpacity(p)-(double) GetPixelOpacity(q);
368 distance=pixel*pixel;
369 if (distance >= fuzz)
370 difference=MagickTrue;
371 }
372 if (((channel & IndexChannel) != 0) &&
373 (image->colorspace == CMYKColorspace))
374 {
375 pixel=Sa*(double) indexes[x]-Da*(double) reconstruct_indexes[x];
376 distance=pixel*pixel;
377 if (distance >= fuzz)
378 difference=MagickTrue;
379 }
380 }
381 if (difference != MagickFalse)
382 SetPixelPacket(highlight_image,&highlight,r,highlight_indexes ==
383 (IndexPacket *) NULL ? NULL : highlight_indexes+x);
384 else
385 SetPixelPacket(highlight_image,&lowlight,r,highlight_indexes ==
386 (IndexPacket *) NULL ? NULL : highlight_indexes+x);
387 p++;
388 q++;
389 r++;
390 }
391 sync=SyncCacheViewAuthenticPixels(highlight_view,exception);
392 if (sync == MagickFalse)
393 status=MagickFalse;
394 }
395 highlight_view=DestroyCacheView(highlight_view);
396 reconstruct_view=DestroyCacheView(reconstruct_view);
397 image_view=DestroyCacheView(image_view);
398 (void) CompositeImage(difference_image,image->compose,highlight_image,0,0);
399 highlight_image=DestroyImage(highlight_image);
400 if (status == MagickFalse)
401 difference_image=DestroyImage(difference_image);
402 return(difference_image);
403}
404
405/*
406%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
407% %
408% %
409% %
410% G e t I m a g e C h a n n e l D i s t o r t i o n %
411% %
412% %
413% %
414%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
415%
416% GetImageChannelDistortion() compares one or more image channels of an image
417% to a reconstructed image and returns the specified distortion metric.
418%
419% The format of the GetImageChannelDistortion method is:
420%
421% MagickBooleanType GetImageChannelDistortion(const Image *image,
422% const Image *reconstruct_image,const ChannelType channel,
423% const MetricType metric,double *distortion,ExceptionInfo *exception)
424%
425% A description of each parameter follows:
426%
427% o image: the image.
428%
429% o reconstruct_image: the reconstruct image.
430%
431% o channel: the channel.
432%
433% o metric: the metric.
434%
435% o distortion: the computed distortion between the images.
436%
437% o exception: return any errors or warnings in this structure.
438%
439*/
440
441MagickExport MagickBooleanType GetImageDistortion(Image *image,
442 const Image *reconstruct_image,const MetricType metric,double *distortion,
443 ExceptionInfo *exception)
444{
445 MagickBooleanType
446 status;
447
448 status=GetImageChannelDistortion(image,reconstruct_image,CompositeChannels,
449 metric,distortion,exception);
450 return(status);
451}
452
453static MagickBooleanType GetAbsoluteDistortion(const Image *image,
454 const Image *reconstruct_image,const ChannelType channel,double *distortion,
455 ExceptionInfo *exception)
456{
457 CacheView
458 *image_view,
459 *reconstruct_view;
460
461 double
462 fuzz;
463
464 MagickBooleanType
465 status;
466
467 size_t
468 columns,
469 rows;
470
471 ssize_t
472 y;
473
474 /*
475 Compute the absolute difference in pixels between two images.
476 */
477 status=MagickTrue;
478 fuzz=GetFuzzyColorDistance(image,reconstruct_image);
479 SetImageDistortionBounds(image,reconstruct_image,&columns,&rows);
480 image_view=AcquireVirtualCacheView(image,exception);
481 reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
482#if defined(MAGICKCORE_OPENMP_SUPPORT)
483 #pragma omp parallel for schedule(static) shared(status) \
484 magick_number_threads(image,image,rows,1)
485#endif
486 for (y=0; y < (ssize_t) rows; y++)
487 {
488 double
489 channel_distortion[CompositeChannels+1];
490
491 const IndexPacket
492 *magick_restrict indexes,
493 *magick_restrict reconstruct_indexes;
494
495 const PixelPacket
496 *magick_restrict p,
497 *magick_restrict q;
498
499 ssize_t
500 i,
501 x;
502
503 if (status == MagickFalse)
504 continue;
505 p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
506 q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
507 if ((p == (const PixelPacket *) NULL) || (q == (const PixelPacket *) NULL))
508 {
509 status=MagickFalse;
510 continue;
511 }
512 indexes=GetCacheViewVirtualIndexQueue(image_view);
513 reconstruct_indexes=GetCacheViewVirtualIndexQueue(reconstruct_view);
514 (void) memset(channel_distortion,0,sizeof(channel_distortion));
515 for (x=0; x < (ssize_t) columns; x++)
516 {
517 double
518 Da,
519 distance,
520 pixel,
521 Sa;
522
523 MagickBooleanType
524 difference;
525
526 difference=MagickFalse;
527 Sa=QuantumScale*(image->matte != MagickFalse ? (double) GetPixelAlpha(p) :
528 ((double) QuantumRange-(double) OpaqueOpacity));
529 Da=QuantumScale*(image->matte != MagickFalse ? (double) GetPixelAlpha(q) :
530 ((double) QuantumRange-(double) OpaqueOpacity));
531 if ((channel & RedChannel) != 0)
532 {
533 pixel=Sa*(double) GetPixelRed(p)-Da*(double) GetPixelRed(q);
534 distance=pixel*pixel;
535 if (distance >= fuzz)
536 {
537 channel_distortion[RedChannel]++;
538 difference=MagickTrue;
539 }
540 }
541 if ((channel & GreenChannel) != 0)
542 {
543 pixel=Sa*(double) GetPixelGreen(p)-Da*(double) GetPixelGreen(q);
544 distance=pixel*pixel;
545 if (distance >= fuzz)
546 {
547 channel_distortion[GreenChannel]++;
548 difference=MagickTrue;
549 }
550 }
551 if ((channel & BlueChannel) != 0)
552 {
553 pixel=Sa*(double) GetPixelBlue(p)-Da*(double) GetPixelBlue(q);
554 distance=pixel*pixel;
555 if (distance >= fuzz)
556 {
557 channel_distortion[BlueChannel]++;
558 difference=MagickTrue;
559 }
560 }
561 if (((channel & OpacityChannel) != 0) &&
562 (image->matte != MagickFalse))
563 {
564 pixel=(double) GetPixelOpacity(p)-(double) GetPixelOpacity(q);
565 distance=pixel*pixel;
566 if (distance >= fuzz)
567 {
568 channel_distortion[OpacityChannel]++;
569 difference=MagickTrue;
570 }
571 }
572 if (((channel & IndexChannel) != 0) &&
573 (image->colorspace == CMYKColorspace))
574 {
575 pixel=Sa*(double) indexes[x]-Da*(double) reconstruct_indexes[x];
576 distance=pixel*pixel;
577 if (distance >= fuzz)
578 {
579 channel_distortion[BlackChannel]++;
580 difference=MagickTrue;
581 }
582 }
583 if (difference != MagickFalse)
584 channel_distortion[CompositeChannels]++;
585 p++;
586 q++;
587 }
588#if defined(MAGICKCORE_OPENMP_SUPPORT)
589 #pragma omp critical (MagickCore_GetAbsoluteDistortion)
590#endif
591 for (i=0; i <= (ssize_t) CompositeChannels; i++)
592 distortion[i]+=channel_distortion[i];
593 }
594 reconstruct_view=DestroyCacheView(reconstruct_view);
595 image_view=DestroyCacheView(image_view);
596 return(status);
597}
598
599static MagickBooleanType GetFuzzDistortion(const Image *image,
600 const Image *reconstruct_image,const ChannelType channel,
601 double *distortion,ExceptionInfo *exception)
602{
603 CacheView
604 *image_view,
605 *reconstruct_view;
606
607 MagickBooleanType
608 status;
609
610 ssize_t
611 i;
612
613 size_t
614 columns,
615 rows;
616
617 ssize_t
618 y;
619
620 status=MagickTrue;
621 SetImageDistortionBounds(image,reconstruct_image,&columns,&rows);
622 image_view=AcquireVirtualCacheView(image,exception);
623 reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
624#if defined(MAGICKCORE_OPENMP_SUPPORT)
625 #pragma omp parallel for schedule(static) shared(status) \
626 magick_number_threads(image,image,rows,1)
627#endif
628 for (y=0; y < (ssize_t) rows; y++)
629 {
630 double
631 channel_distortion[CompositeChannels+1];
632
633 const IndexPacket
634 *magick_restrict indexes,
635 *magick_restrict reconstruct_indexes;
636
637 const PixelPacket
638 *magick_restrict p,
639 *magick_restrict q;
640
641 ssize_t
642 i,
643 x;
644
645 if (status == MagickFalse)
646 continue;
647 p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
648 q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
649 if ((p == (const PixelPacket *) NULL) || (q == (const PixelPacket *) NULL))
650 {
651 status=MagickFalse;
652 continue;
653 }
654 indexes=GetCacheViewVirtualIndexQueue(image_view);
655 reconstruct_indexes=GetCacheViewVirtualIndexQueue(reconstruct_view);
656 (void) memset(channel_distortion,0,sizeof(channel_distortion));
657 for (x=0; x < (ssize_t) columns; x++)
658 {
659 MagickRealType
660 distance,
661 Da,
662 Sa;
663
664 Sa=QuantumScale*(image->matte != MagickFalse ? (double) GetPixelAlpha(p) :
665 ((double) QuantumRange-(double) OpaqueOpacity));
666 Da=QuantumScale*(reconstruct_image->matte != MagickFalse ?
667 (double) GetPixelAlpha(q) : ((double) QuantumRange-(double)
668 OpaqueOpacity));
669 if ((channel & RedChannel) != 0)
670 {
671 distance=QuantumScale*(Sa*(double) GetPixelRed(p)-Da*(double)
672 GetPixelRed(q));
673 channel_distortion[RedChannel]+=distance*distance;
674 channel_distortion[CompositeChannels]+=distance*distance;
675 }
676 if ((channel & GreenChannel) != 0)
677 {
678 distance=QuantumScale*(Sa*(double) GetPixelGreen(p)-Da*(double)
679 GetPixelGreen(q));
680 channel_distortion[GreenChannel]+=distance*distance;
681 channel_distortion[CompositeChannels]+=distance*distance;
682 }
683 if ((channel & BlueChannel) != 0)
684 {
685 distance=QuantumScale*(Sa*(double) GetPixelBlue(p)-Da*(double)
686 GetPixelBlue(q));
687 channel_distortion[BlueChannel]+=distance*distance;
688 channel_distortion[CompositeChannels]+=distance*distance;
689 }
690 if (((channel & OpacityChannel) != 0) && ((image->matte != MagickFalse) ||
691 (reconstruct_image->matte != MagickFalse)))
692 {
693 distance=QuantumScale*((image->matte != MagickFalse ? (double)
694 GetPixelOpacity(p) : (double) OpaqueOpacity)-
695 (reconstruct_image->matte != MagickFalse ?
696 (double) GetPixelOpacity(q): (double) OpaqueOpacity));
697 channel_distortion[OpacityChannel]+=distance*distance;
698 channel_distortion[CompositeChannels]+=distance*distance;
699 }
700 if (((channel & IndexChannel) != 0) &&
701 (image->colorspace == CMYKColorspace) &&
702 (reconstruct_image->colorspace == CMYKColorspace))
703 {
704 distance=QuantumScale*(Sa*(double) GetPixelIndex(indexes+x)-
705 Da*(double) GetPixelIndex(reconstruct_indexes+x));
706 channel_distortion[BlackChannel]+=distance*distance;
707 channel_distortion[CompositeChannels]+=distance*distance;
708 }
709 p++;
710 q++;
711 }
712#if defined(MAGICKCORE_OPENMP_SUPPORT)
713 #pragma omp critical (MagickCore_GetFuzzDistortion)
714#endif
715 for (i=0; i <= (ssize_t) CompositeChannels; i++)
716 distortion[i]+=channel_distortion[i];
717 }
718 reconstruct_view=DestroyCacheView(reconstruct_view);
719 image_view=DestroyCacheView(image_view);
720 for (i=0; i <= (ssize_t) CompositeChannels; i++)
721 distortion[i]/=((double) columns*rows);
722 distortion[CompositeChannels]/=(double) GetNumberChannels(image,channel);
723 distortion[CompositeChannels]=sqrt(distortion[CompositeChannels]);
724 return(status);
725}
726
727static MagickBooleanType GetMeanAbsoluteDistortion(const Image *image,
728 const Image *reconstruct_image,const ChannelType channel,
729 double *distortion,ExceptionInfo *exception)
730{
731 CacheView
732 *image_view,
733 *reconstruct_view;
734
735 MagickBooleanType
736 status;
737
738 size_t
739 columns,
740 rows;
741
742 ssize_t
743 i,
744 y;
745
746 status=MagickTrue;
747 SetImageDistortionBounds(image,reconstruct_image,&columns,&rows);
748 image_view=AcquireVirtualCacheView(image,exception);
749 reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
750#if defined(MAGICKCORE_OPENMP_SUPPORT)
751 #pragma omp parallel for schedule(static) shared(status) \
752 magick_number_threads(image,image,rows,1)
753#endif
754 for (y=0; y < (ssize_t) rows; y++)
755 {
756 double
757 channel_distortion[CompositeChannels+1];
758
759 const IndexPacket
760 *magick_restrict indexes,
761 *magick_restrict reconstruct_indexes;
762
763 const PixelPacket
764 *magick_restrict p,
765 *magick_restrict q;
766
767 ssize_t
768 i,
769 x;
770
771 if (status == MagickFalse)
772 continue;
773 p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
774 q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
775 if ((p == (const PixelPacket *) NULL) || (q == (const PixelPacket *) NULL))
776 {
777 status=MagickFalse;
778 continue;
779 }
780 indexes=GetCacheViewVirtualIndexQueue(image_view);
781 reconstruct_indexes=GetCacheViewVirtualIndexQueue(reconstruct_view);
782 (void) memset(channel_distortion,0,sizeof(channel_distortion));
783 for (x=0; x < (ssize_t) columns; x++)
784 {
785 MagickRealType
786 distance,
787 Da,
788 Sa;
789
790 Sa=QuantumScale*(image->matte != MagickFalse ? (double) GetPixelAlpha(p) :
791 ((double) QuantumRange-(double) OpaqueOpacity));
792 Da=QuantumScale*(reconstruct_image->matte != MagickFalse ?
793 (double) GetPixelAlpha(q) : ((double) QuantumRange-(double)
794 OpaqueOpacity));
795 if ((channel & RedChannel) != 0)
796 {
797 distance=QuantumScale*fabs(Sa*(double) GetPixelRed(p)-Da*
798 (double) GetPixelRed(q));
799 channel_distortion[RedChannel]+=distance;
800 channel_distortion[CompositeChannels]+=distance;
801 }
802 if ((channel & GreenChannel) != 0)
803 {
804 distance=QuantumScale*fabs(Sa*(double) GetPixelGreen(p)-Da*
805 (double) GetPixelGreen(q));
806 channel_distortion[GreenChannel]+=distance;
807 channel_distortion[CompositeChannels]+=distance;
808 }
809 if ((channel & BlueChannel) != 0)
810 {
811 distance=QuantumScale*fabs(Sa*(double) GetPixelBlue(p)-Da*
812 (double) GetPixelBlue(q));
813 channel_distortion[BlueChannel]+=distance;
814 channel_distortion[CompositeChannels]+=distance;
815 }
816 if (((channel & OpacityChannel) != 0) &&
817 (image->matte != MagickFalse))
818 {
819 distance=QuantumScale*fabs((double) GetPixelOpacity(p)-(double)
820 GetPixelOpacity(q));
821 channel_distortion[OpacityChannel]+=distance;
822 channel_distortion[CompositeChannels]+=distance;
823 }
824 if (((channel & IndexChannel) != 0) &&
825 (image->colorspace == CMYKColorspace))
826 {
827 distance=QuantumScale*fabs(Sa*(double) GetPixelIndex(indexes+x)-Da*
828 (double) GetPixelIndex(reconstruct_indexes+x));
829 channel_distortion[BlackChannel]+=distance;
830 channel_distortion[CompositeChannels]+=distance;
831 }
832 p++;
833 q++;
834 }
835#if defined(MAGICKCORE_OPENMP_SUPPORT)
836 #pragma omp critical (MagickCore_GetMeanAbsoluteError)
837#endif
838 for (i=0; i <= (ssize_t) CompositeChannels; i++)
839 distortion[i]+=channel_distortion[i];
840 }
841 reconstruct_view=DestroyCacheView(reconstruct_view);
842 image_view=DestroyCacheView(image_view);
843 for (i=0; i <= (ssize_t) CompositeChannels; i++)
844 distortion[i]/=((double) columns*rows);
845 distortion[CompositeChannels]/=(double) GetNumberChannels(image,channel);
846 return(status);
847}
848
849static MagickBooleanType GetMeanErrorPerPixel(Image *image,
850 const Image *reconstruct_image,const ChannelType channel,double *distortion,
851 ExceptionInfo *exception)
852{
853 CacheView
854 *image_view,
855 *reconstruct_view;
856
857 double
858 maximum_error = MagickMinimumValue;
859
860 MagickBooleanType
861 status;
862
863 size_t
864 columns,
865 rows;
866
867 ssize_t
868 i,
869 y;
870
871 status=MagickTrue;
872 SetImageDistortionBounds(image,reconstruct_image,&columns,&rows);
873 image_view=AcquireVirtualCacheView(image,exception);
874 reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
875#if defined(MAGICKCORE_OPENMP_SUPPORT)
876 #pragma omp parallel for schedule(static) shared(status) \
877 magick_number_threads(image,image,rows,1)
878#endif
879 for (y=0; y < (ssize_t) rows; y++)
880 {
881 double
882 channel_distortion[CompositeChannels+1],
883 local_maximum = MinimumValue;
884
885 const IndexPacket
886 *magick_restrict indexes,
887 *magick_restrict reconstruct_indexes;
888
889 const PixelPacket
890 *magick_restrict p,
891 *magick_restrict q;
892
893 ssize_t
894 i,
895 x;
896
897 if (status == MagickFalse)
898 continue;
899 p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
900 q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
901 if ((p == (const PixelPacket *) NULL) || (q == (const PixelPacket *) NULL))
902 {
903 status=MagickFalse;
904 continue;
905 }
906 indexes=GetCacheViewVirtualIndexQueue(image_view);
907 reconstruct_indexes=GetCacheViewVirtualIndexQueue(reconstruct_view);
908 (void) memset(channel_distortion,0,sizeof(channel_distortion));
909 for (x=0; x < (ssize_t) columns; x++)
910 {
911 MagickRealType
912 distance,
913 Da,
914 Sa;
915
916 Sa=QuantumScale*(image->matte != MagickFalse ? (double) GetPixelAlpha(p) :
917 ((double) QuantumRange-(double) OpaqueOpacity));
918 Da=QuantumScale*(reconstruct_image->matte != MagickFalse ?
919 (double) GetPixelAlpha(q) : ((double) QuantumRange-(double)
920 OpaqueOpacity));
921 if ((channel & RedChannel) != 0)
922 {
923 distance=QuantumScale*fabs(Sa*(double) GetPixelRed(p)-Da*
924 (double) GetPixelRed(q));
925 channel_distortion[RedChannel]+=distance;
926 channel_distortion[CompositeChannels]+=distance;
927 if (distance > local_maximum)
928 local_maximum=distance;
929 }
930 if ((channel & GreenChannel) != 0)
931 {
932 distance=QuantumScale*fabs(Sa*(double) GetPixelGreen(p)-Da*
933 (double) GetPixelGreen(q));
934 channel_distortion[GreenChannel]+=distance;
935 channel_distortion[CompositeChannels]+=distance;
936 if (distance > local_maximum)
937 local_maximum=distance;
938 }
939 if ((channel & BlueChannel) != 0)
940 {
941 distance=QuantumScale*fabs(Sa*(double) GetPixelBlue(p)-Da*
942 (double) GetPixelBlue(q));
943 channel_distortion[BlueChannel]+=distance;
944 channel_distortion[CompositeChannels]+=distance;
945 if (distance > local_maximum)
946 local_maximum=distance;
947 }
948 if (((channel & OpacityChannel) != 0) &&
949 (image->matte != MagickFalse))
950 {
951 distance=QuantumScale*fabs((double) GetPixelOpacity(p)-(double)
952 GetPixelOpacity(q));
953 channel_distortion[OpacityChannel]+=distance;
954 channel_distortion[CompositeChannels]+=distance;
955 if (distance > local_maximum)
956 local_maximum=distance;
957 }
958 if (((channel & IndexChannel) != 0) &&
959 (image->colorspace == CMYKColorspace))
960 {
961 distance=QuantumScale*fabs(Sa*(double) GetPixelIndex(indexes+x)-Da*
962 (double) GetPixelIndex(reconstruct_indexes+x));
963 channel_distortion[BlackChannel]+=distance;
964 channel_distortion[CompositeChannels]+=distance;
965 if (distance > local_maximum)
966 local_maximum=distance;
967 }
968 p++;
969 q++;
970 }
971#if defined(MAGICKCORE_OPENMP_SUPPORT)
972 #pragma omp critical (MagickCore_GetMeanAbsoluteError)
973#endif
974 {
975 for (i=0; i <= (ssize_t) CompositeChannels; i++)
976 distortion[i]+=channel_distortion[i];
977 if (local_maximum > maximum_error)
978 maximum_error=local_maximum;
979 }
980 }
981 reconstruct_view=DestroyCacheView(reconstruct_view);
982 image_view=DestroyCacheView(image_view);
983 distortion[CompositeChannels]/=(double) GetNumberChannels(image,channel);
984 for (i=0; i <= (ssize_t) CompositeChannels; i++)
985 distortion[i]/=((double) columns*rows);
986 image->error.mean_error_per_pixel=distortion[CompositeChannels];
987 image->error.normalized_mean_error=distortion[CompositeChannels];
988 image->error.normalized_maximum_error=QuantumScale*maximum_error;
989 return(status);
990}
991
992static MagickBooleanType GetMeanSquaredDistortion(const Image *image,
993 const Image *reconstruct_image,const ChannelType channel,
994 double *distortion,ExceptionInfo *exception)
995{
996 CacheView
997 *image_view,
998 *reconstruct_view;
999
1000 MagickBooleanType
1001 status;
1002
1003 ssize_t
1004 i;
1005
1006 size_t
1007 columns,
1008 rows;
1009
1010 ssize_t
1011 y;
1012
1013 status=MagickTrue;
1014 SetImageDistortionBounds(image,reconstruct_image,&columns,&rows);
1015 image_view=AcquireVirtualCacheView(image,exception);
1016 reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
1017#if defined(MAGICKCORE_OPENMP_SUPPORT)
1018 #pragma omp parallel for schedule(static) shared(status) \
1019 magick_number_threads(image,image,rows,1)
1020#endif
1021 for (y=0; y < (ssize_t) rows; y++)
1022 {
1023 double
1024 channel_distortion[CompositeChannels+1];
1025
1026 const IndexPacket
1027 *magick_restrict indexes,
1028 *magick_restrict reconstruct_indexes;
1029
1030 const PixelPacket
1031 *magick_restrict p,
1032 *magick_restrict q;
1033
1034 ssize_t
1035 i,
1036 x;
1037
1038 if (status == MagickFalse)
1039 continue;
1040 p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
1041 q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
1042 if ((p == (const PixelPacket *) NULL) || (q == (const PixelPacket *) NULL))
1043 {
1044 status=MagickFalse;
1045 continue;
1046 }
1047 indexes=GetCacheViewVirtualIndexQueue(image_view);
1048 reconstruct_indexes=GetCacheViewVirtualIndexQueue(reconstruct_view);
1049 (void) memset(channel_distortion,0,sizeof(channel_distortion));
1050 for (x=0; x < (ssize_t) columns; x++)
1051 {
1052 MagickRealType
1053 distance,
1054 Da,
1055 Sa;
1056
1057 Sa=QuantumScale*(image->matte != MagickFalse ? (double) GetPixelAlpha(p) :
1058 ((double) QuantumRange-(double) OpaqueOpacity));
1059 Da=QuantumScale*(reconstruct_image->matte != MagickFalse ?
1060 (double) GetPixelAlpha(q) : ((double) QuantumRange-(double)
1061 OpaqueOpacity));
1062 if ((channel & RedChannel) != 0)
1063 {
1064 distance=QuantumScale*(Sa*(double) GetPixelRed(p)-Da*(double)
1065 GetPixelRed(q));
1066 channel_distortion[RedChannel]+=distance*distance;
1067 channel_distortion[CompositeChannels]+=distance*distance;
1068 }
1069 if ((channel & GreenChannel) != 0)
1070 {
1071 distance=QuantumScale*(Sa*(double) GetPixelGreen(p)-Da*(double)
1072 GetPixelGreen(q));
1073 channel_distortion[GreenChannel]+=distance*distance;
1074 channel_distortion[CompositeChannels]+=distance*distance;
1075 }
1076 if ((channel & BlueChannel) != 0)
1077 {
1078 distance=QuantumScale*(Sa*(double) GetPixelBlue(p)-Da*(double)
1079 GetPixelBlue(q));
1080 channel_distortion[BlueChannel]+=distance*distance;
1081 channel_distortion[CompositeChannels]+=distance*distance;
1082 }
1083 if (((channel & OpacityChannel) != 0) &&
1084 (image->matte != MagickFalse))
1085 {
1086 distance=QuantumScale*((double) GetPixelOpacity(p)-(double)
1087 GetPixelOpacity(q));
1088 channel_distortion[OpacityChannel]+=distance*distance;
1089 channel_distortion[CompositeChannels]+=distance*distance;
1090 }
1091 if (((channel & IndexChannel) != 0) &&
1092 (image->colorspace == CMYKColorspace) &&
1093 (reconstruct_image->colorspace == CMYKColorspace))
1094 {
1095 distance=QuantumScale*(Sa*(double) GetPixelIndex(indexes+x)-Da*
1096 (double) GetPixelIndex(reconstruct_indexes+x));
1097 channel_distortion[BlackChannel]+=distance*distance;
1098 channel_distortion[CompositeChannels]+=distance*distance;
1099 }
1100 p++;
1101 q++;
1102 }
1103#if defined(MAGICKCORE_OPENMP_SUPPORT)
1104 #pragma omp critical (MagickCore_GetMeanSquaredError)
1105#endif
1106 for (i=0; i <= (ssize_t) CompositeChannels; i++)
1107 distortion[i]+=channel_distortion[i];
1108 }
1109 reconstruct_view=DestroyCacheView(reconstruct_view);
1110 image_view=DestroyCacheView(image_view);
1111 for (i=0; i <= (ssize_t) CompositeChannels; i++)
1112 distortion[i]/=((double) columns*rows);
1113 distortion[CompositeChannels]/=(double) GetNumberChannels(image,channel);
1114 return(status);
1115}
1116
1117static MagickBooleanType GetNormalizedCrossCorrelationDistortion(
1118 const Image *image,const Image *reconstruct_image,const ChannelType channel,
1119 double *distortion,ExceptionInfo *exception)
1120{
1121#define SimilarityImageTag "Similarity/Image"
1122
1123 CacheView
1124 *image_view,
1125 *reconstruct_view;
1126
1127 ChannelStatistics
1128 *image_statistics,
1129 *reconstruct_statistics;
1130
1131 double
1132 alpha_variance[CompositeChannels+1],
1133 beta_variance[CompositeChannels+1];
1134
1135 MagickBooleanType
1136 status;
1137
1138 MagickOffsetType
1139 progress;
1140
1141 ssize_t
1142 i;
1143
1144 size_t
1145 columns,
1146 rows;
1147
1148 ssize_t
1149 y;
1150
1151 /*
1152 Normalize to account for variation due to lighting and exposure condition.
1153 */
1154 image_statistics=GetImageChannelStatistics(image,exception);
1155 reconstruct_statistics=GetImageChannelStatistics(reconstruct_image,exception);
1156 if ((image_statistics == (ChannelStatistics *) NULL) ||
1157 (reconstruct_statistics == (ChannelStatistics *) NULL))
1158 {
1159 if (image_statistics != (ChannelStatistics *) NULL)
1160 image_statistics=(ChannelStatistics *) RelinquishMagickMemory(
1161 image_statistics);
1162 if (reconstruct_statistics != (ChannelStatistics *) NULL)
1163 reconstruct_statistics=(ChannelStatistics *) RelinquishMagickMemory(
1164 reconstruct_statistics);
1165 return(MagickFalse);
1166 }
1167 (void) memset(distortion,0,(CompositeChannels+1)*sizeof(*distortion));
1168 (void) memset(alpha_variance,0,(CompositeChannels+1)*sizeof(*alpha_variance));
1169 (void) memset(beta_variance,0,(CompositeChannels+1)*sizeof(*beta_variance));
1170 status=MagickTrue;
1171 progress=0;
1172 SetImageDistortionBounds(image,reconstruct_image,&columns,&rows);
1173 image_view=AcquireVirtualCacheView(image,exception);
1174 reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
1175 for (y=0; y < (ssize_t) rows; y++)
1176 {
1177 const IndexPacket
1178 *magick_restrict indexes,
1179 *magick_restrict reconstruct_indexes;
1180
1181 const PixelPacket
1182 *magick_restrict p,
1183 *magick_restrict q;
1184
1185 ssize_t
1186 x;
1187
1188 if (status == MagickFalse)
1189 continue;
1190 p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
1191 q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
1192 if ((p == (const PixelPacket *) NULL) || (q == (const PixelPacket *) NULL))
1193 {
1194 status=MagickFalse;
1195 continue;
1196 }
1197 indexes=GetCacheViewVirtualIndexQueue(image_view);
1198 reconstruct_indexes=GetCacheViewVirtualIndexQueue(reconstruct_view);
1199 for (x=0; x < (ssize_t) columns; x++)
1200 {
1201 MagickRealType
1202 alpha,
1203 beta,
1204 Da,
1205 Sa;
1206
1207 Sa=QuantumScale*(image->matte != MagickFalse ? (double) GetPixelAlpha(p) :
1208 (double) QuantumRange);
1209 Da=QuantumScale*(reconstruct_image->matte != MagickFalse ?
1210 (double) GetPixelAlpha(q) : (double) QuantumRange);
1211 if ((channel & RedChannel) != 0)
1212 {
1213 alpha=QuantumScale*fabs(Sa*(double) GetPixelRed(p)-
1214 image_statistics[RedChannel].mean);
1215 beta=QuantumScale*fabs(Da*(double) GetPixelRed(q)-
1216 reconstruct_statistics[RedChannel].mean);
1217 distortion[RedChannel]+=alpha*beta;
1218 alpha_variance[RedChannel]+=alpha*alpha;
1219 beta_variance[RedChannel]+=beta*beta;
1220 }
1221 if ((channel & GreenChannel) != 0)
1222 {
1223 alpha=QuantumScale*fabs(Sa*(double) GetPixelGreen(p)-
1224 image_statistics[GreenChannel].mean);
1225 beta=QuantumScale*fabs(Da*(double) GetPixelGreen(q)-
1226 reconstruct_statistics[GreenChannel].mean);
1227 distortion[GreenChannel]+=alpha*beta;
1228 alpha_variance[GreenChannel]+=alpha*alpha;
1229 beta_variance[GreenChannel]+=beta*beta;
1230 }
1231 if ((channel & BlueChannel) != 0)
1232 {
1233 alpha=QuantumScale*fabs(Sa*(double) GetPixelBlue(p)-
1234 image_statistics[BlueChannel].mean);
1235 beta=QuantumScale*fabs(Da*(double) GetPixelBlue(q)-
1236 reconstruct_statistics[BlueChannel].mean);
1237 distortion[BlueChannel]+=alpha*beta;
1238 alpha_variance[BlueChannel]+=alpha*alpha;
1239 beta_variance[BlueChannel]+=beta*beta;
1240 }
1241 if (((channel & OpacityChannel) != 0) && (image->matte != MagickFalse))
1242 {
1243 alpha=QuantumScale*fabs((double) GetPixelAlpha(p)-
1244 image_statistics[AlphaChannel].mean);
1245 beta=QuantumScale*fabs((double) GetPixelAlpha(q)-
1246 reconstruct_statistics[AlphaChannel].mean);
1247 distortion[OpacityChannel]+=alpha*beta;
1248 alpha_variance[OpacityChannel]+=alpha*alpha;
1249 beta_variance[OpacityChannel]+=beta*beta;
1250 }
1251 if (((channel & IndexChannel) != 0) &&
1252 (image->colorspace == CMYKColorspace) &&
1253 (reconstruct_image->colorspace == CMYKColorspace))
1254 {
1255 alpha=QuantumScale*fabs(Sa*(double) GetPixelIndex(indexes+x)-
1256 image_statistics[BlackChannel].mean);
1257 beta=QuantumScale*fabs(Da*(double) GetPixelIndex(reconstruct_indexes+
1258 x)-reconstruct_statistics[BlackChannel].mean);
1259 distortion[BlackChannel]+=alpha*beta;
1260 alpha_variance[BlackChannel]+=alpha*alpha;
1261 beta_variance[BlackChannel]+=beta*beta;
1262 }
1263 p++;
1264 q++;
1265 }
1266 if (image->progress_monitor != (MagickProgressMonitor) NULL)
1267 {
1268 MagickBooleanType
1269 proceed;
1270
1271#if defined(MAGICKCORE_OPENMP_SUPPORT)
1272 #pragma omp atomic
1273#endif
1274 progress++;
1275 proceed=SetImageProgress(image,SimilarityImageTag,progress,rows);
1276 if (proceed == MagickFalse)
1277 status=MagickFalse;
1278 }
1279 }
1280 reconstruct_view=DestroyCacheView(reconstruct_view);
1281 image_view=DestroyCacheView(image_view);
1282 /*
1283 Divide by the standard deviation.
1284 */
1285 for (i=0; i < (ssize_t) CompositeChannels; i++)
1286 {
1287 distortion[i]/=sqrt(alpha_variance[i]*beta_variance[i]);
1288 if (fabs(distortion[i]) >= MagickEpsilon)
1289 distortion[CompositeChannels]+=distortion[i];
1290 }
1291 distortion[CompositeChannels]=distortion[CompositeChannels]/
1292 GetNumberChannels(image,channel);
1293 /*
1294 Free resources.
1295 */
1296 reconstruct_statistics=(ChannelStatistics *) RelinquishMagickMemory(
1297 reconstruct_statistics);
1298 image_statistics=(ChannelStatistics *) RelinquishMagickMemory(
1299 image_statistics);
1300 return(status);
1301}
1302
1303static MagickBooleanType GetPeakAbsoluteDistortion(const Image *image,
1304 const Image *reconstruct_image,const ChannelType channel,
1305 double *distortion,ExceptionInfo *exception)
1306{
1307 CacheView
1308 *image_view,
1309 *reconstruct_view;
1310
1311 MagickBooleanType
1312 status;
1313
1314 size_t
1315 columns,
1316 rows;
1317
1318 ssize_t
1319 y;
1320
1321 status=MagickTrue;
1322 (void) memset(distortion,0,(CompositeChannels+1)*sizeof(*distortion));
1323 SetImageDistortionBounds(image,reconstruct_image,&columns,&rows);
1324 image_view=AcquireVirtualCacheView(image,exception);
1325 reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
1326#if defined(MAGICKCORE_OPENMP_SUPPORT)
1327 #pragma omp parallel for schedule(static) shared(status) \
1328 magick_number_threads(image,image,rows,1)
1329#endif
1330 for (y=0; y < (ssize_t) rows; y++)
1331 {
1332 double
1333 channel_distortion[CompositeChannels+1];
1334
1335 const IndexPacket
1336 *magick_restrict indexes,
1337 *magick_restrict reconstruct_indexes;
1338
1339 const PixelPacket
1340 *magick_restrict p,
1341 *magick_restrict q;
1342
1343 ssize_t
1344 i,
1345 x;
1346
1347 if (status == MagickFalse)
1348 continue;
1349 p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
1350 q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
1351 if ((p == (const PixelPacket *) NULL) || (q == (const PixelPacket *) NULL))
1352 {
1353 status=MagickFalse;
1354 continue;
1355 }
1356 indexes=GetCacheViewVirtualIndexQueue(image_view);
1357 reconstruct_indexes=GetCacheViewVirtualIndexQueue(reconstruct_view);
1358 (void) memset(channel_distortion,0,(CompositeChannels+1)*
1359 sizeof(*channel_distortion));
1360 for (x=0; x < (ssize_t) columns; x++)
1361 {
1362 MagickRealType
1363 distance,
1364 Da,
1365 Sa;
1366
1367 Sa=QuantumScale*(image->matte != MagickFalse ? (double) GetPixelAlpha(p) :
1368 ((double) QuantumRange-(double) OpaqueOpacity));
1369 Da=QuantumScale*(reconstruct_image->matte != MagickFalse ?
1370 (double) GetPixelAlpha(q) : ((double) QuantumRange-(double)
1371 OpaqueOpacity));
1372 if ((channel & RedChannel) != 0)
1373 {
1374 distance=QuantumScale*fabs(Sa*(double) GetPixelRed(p)-Da*
1375 (double) GetPixelRed(q));
1376 if (distance > channel_distortion[RedChannel])
1377 channel_distortion[RedChannel]=distance;
1378 if (distance > channel_distortion[CompositeChannels])
1379 channel_distortion[CompositeChannels]=distance;
1380 }
1381 if ((channel & GreenChannel) != 0)
1382 {
1383 distance=QuantumScale*fabs(Sa*(double) GetPixelGreen(p)-Da*
1384 (double) GetPixelGreen(q));
1385 if (distance > channel_distortion[GreenChannel])
1386 channel_distortion[GreenChannel]=distance;
1387 if (distance > channel_distortion[CompositeChannels])
1388 channel_distortion[CompositeChannels]=distance;
1389 }
1390 if ((channel & BlueChannel) != 0)
1391 {
1392 distance=QuantumScale*fabs(Sa*(double) GetPixelBlue(p)-Da*
1393 (double) GetPixelBlue(q));
1394 if (distance > channel_distortion[BlueChannel])
1395 channel_distortion[BlueChannel]=distance;
1396 if (distance > channel_distortion[CompositeChannels])
1397 channel_distortion[CompositeChannels]=distance;
1398 }
1399 if (((channel & OpacityChannel) != 0) &&
1400 (image->matte != MagickFalse))
1401 {
1402 distance=QuantumScale*fabs((double) GetPixelOpacity(p)-(double)
1403 GetPixelOpacity(q));
1404 if (distance > channel_distortion[OpacityChannel])
1405 channel_distortion[OpacityChannel]=distance;
1406 if (distance > channel_distortion[CompositeChannels])
1407 channel_distortion[CompositeChannels]=distance;
1408 }
1409 if (((channel & IndexChannel) != 0) &&
1410 (image->colorspace == CMYKColorspace) &&
1411 (reconstruct_image->colorspace == CMYKColorspace))
1412 {
1413 distance=QuantumScale*fabs(Sa*(double) GetPixelIndex(indexes+x)-Da*
1414 (double) GetPixelIndex(reconstruct_indexes+x));
1415 if (distance > channel_distortion[BlackChannel])
1416 channel_distortion[BlackChannel]=distance;
1417 if (distance > channel_distortion[CompositeChannels])
1418 channel_distortion[CompositeChannels]=distance;
1419 }
1420 p++;
1421 q++;
1422 }
1423#if defined(MAGICKCORE_OPENMP_SUPPORT)
1424 #pragma omp critical (MagickCore_GetPeakAbsoluteError)
1425#endif
1426 for (i=0; i <= (ssize_t) CompositeChannels; i++)
1427 if (channel_distortion[i] > distortion[i])
1428 distortion[i]=channel_distortion[i];
1429 }
1430 reconstruct_view=DestroyCacheView(reconstruct_view);
1431 image_view=DestroyCacheView(image_view);
1432 return(status);
1433}
1434
1435static MagickBooleanType GetPeakSignalToNoiseRatio(const Image *image,
1436 const Image *reconstruct_image,const ChannelType channel,
1437 double *distortion,ExceptionInfo *exception)
1438{
1439 MagickBooleanType
1440 status;
1441
1442 status=GetMeanSquaredDistortion(image,reconstruct_image,channel,distortion,
1443 exception);
1444 if ((channel & RedChannel) != 0)
1445 {
1446 if (fabs(distortion[RedChannel]) >= MagickEpsilon)
1447 distortion[RedChannel]=fabs(-10.0*MagickLog10(PerceptibleReciprocal(
1448 distortion[RedChannel])))/48.1647;
1449 }
1450 if ((channel & GreenChannel) != 0)
1451 {
1452 if (fabs(distortion[GreenChannel]) >= MagickEpsilon)
1453 distortion[GreenChannel]=fabs(-10.0*MagickLog10(PerceptibleReciprocal(
1454 distortion[GreenChannel])))/48.1647;
1455 }
1456 if ((channel & BlueChannel) != 0)
1457 {
1458 if (fabs(distortion[BlueChannel]) >= MagickEpsilon)
1459 distortion[BlueChannel]=fabs(-10.0*MagickLog10(PerceptibleReciprocal(
1460 distortion[BlueChannel])))/48.1647;
1461 }
1462 if (((channel & OpacityChannel) != 0) && (image->matte != MagickFalse))
1463 {
1464 if (fabs(distortion[OpacityChannel]) >= MagickEpsilon)
1465 distortion[OpacityChannel]=fabs(-10.0*MagickLog10(PerceptibleReciprocal(
1466 distortion[OpacityChannel])))/48.1647;
1467 }
1468 if (((channel & IndexChannel) != 0) && (image->colorspace == CMYKColorspace))
1469 {
1470 if (fabs(distortion[BlackChannel]) >= MagickEpsilon)
1471 distortion[BlackChannel]=fabs(-10.0*MagickLog10(PerceptibleReciprocal(
1472 distortion[BlackChannel])))/48.1647;
1473 }
1474 if (fabs(distortion[CompositeChannels]) >= MagickEpsilon)
1475 distortion[CompositeChannels]=fabs(-10.0*MagickLog10(PerceptibleReciprocal(
1476 distortion[CompositeChannels])))/48.1647;
1477 return(status);
1478}
1479
1480static MagickBooleanType GetPerceptualHashDistortion(const Image *image,
1481 const Image *reconstruct_image,const ChannelType channel,double *distortion,
1482 ExceptionInfo *exception)
1483{
1484 ChannelPerceptualHash
1485 *image_phash,
1486 *reconstruct_phash;
1487
1488 double
1489 difference;
1490
1491 ssize_t
1492 i;
1493
1494 /*
1495 Compute perceptual hash in the sRGB colorspace.
1496 */
1497 image_phash=GetImageChannelPerceptualHash(image,exception);
1498 if (image_phash == (ChannelPerceptualHash *) NULL)
1499 return(MagickFalse);
1500 reconstruct_phash=GetImageChannelPerceptualHash(reconstruct_image,exception);
1501 if (reconstruct_phash == (ChannelPerceptualHash *) NULL)
1502 {
1503 image_phash=(ChannelPerceptualHash *) RelinquishMagickMemory(image_phash);
1504 return(MagickFalse);
1505 }
1506 for (i=0; i < MaximumNumberOfImageMoments; i++)
1507 {
1508 /*
1509 Compute sum of moment differences squared.
1510 */
1511 if ((channel & RedChannel) != 0)
1512 {
1513 difference=reconstruct_phash[RedChannel].P[i]-
1514 image_phash[RedChannel].P[i];
1515 distortion[RedChannel]+=difference*difference;
1516 distortion[CompositeChannels]+=difference*difference;
1517 }
1518 if ((channel & GreenChannel) != 0)
1519 {
1520 difference=reconstruct_phash[GreenChannel].P[i]-
1521 image_phash[GreenChannel].P[i];
1522 distortion[GreenChannel]+=difference*difference;
1523 distortion[CompositeChannels]+=difference*difference;
1524 }
1525 if ((channel & BlueChannel) != 0)
1526 {
1527 difference=reconstruct_phash[BlueChannel].P[i]-
1528 image_phash[BlueChannel].P[i];
1529 distortion[BlueChannel]+=difference*difference;
1530 distortion[CompositeChannels]+=difference*difference;
1531 }
1532 if (((channel & OpacityChannel) != 0) && (image->matte != MagickFalse) &&
1533 (reconstruct_image->matte != MagickFalse))
1534 {
1535 difference=reconstruct_phash[OpacityChannel].P[i]-
1536 image_phash[OpacityChannel].P[i];
1537 distortion[OpacityChannel]+=difference*difference;
1538 distortion[CompositeChannels]+=difference*difference;
1539 }
1540 if (((channel & IndexChannel) != 0) &&
1541 (image->colorspace == CMYKColorspace) &&
1542 (reconstruct_image->colorspace == CMYKColorspace))
1543 {
1544 difference=reconstruct_phash[IndexChannel].P[i]-
1545 image_phash[IndexChannel].P[i];
1546 distortion[IndexChannel]+=difference*difference;
1547 distortion[CompositeChannels]+=difference*difference;
1548 }
1549 }
1550 /*
1551 Compute perceptual hash in the HCLP colorspace.
1552 */
1553 for (i=0; i < MaximumNumberOfImageMoments; i++)
1554 {
1555 /*
1556 Compute sum of moment differences squared.
1557 */
1558 if ((channel & RedChannel) != 0)
1559 {
1560 difference=reconstruct_phash[RedChannel].Q[i]-
1561 image_phash[RedChannel].Q[i];
1562 distortion[RedChannel]+=difference*difference;
1563 distortion[CompositeChannels]+=difference*difference;
1564 }
1565 if ((channel & GreenChannel) != 0)
1566 {
1567 difference=reconstruct_phash[GreenChannel].Q[i]-
1568 image_phash[GreenChannel].Q[i];
1569 distortion[GreenChannel]+=difference*difference;
1570 distortion[CompositeChannels]+=difference*difference;
1571 }
1572 if ((channel & BlueChannel) != 0)
1573 {
1574 difference=reconstruct_phash[BlueChannel].Q[i]-
1575 image_phash[BlueChannel].Q[i];
1576 distortion[BlueChannel]+=difference*difference;
1577 distortion[CompositeChannels]+=difference*difference;
1578 }
1579 if (((channel & OpacityChannel) != 0) && (image->matte != MagickFalse) &&
1580 (reconstruct_image->matte != MagickFalse))
1581 {
1582 difference=reconstruct_phash[OpacityChannel].Q[i]-
1583 image_phash[OpacityChannel].Q[i];
1584 distortion[OpacityChannel]+=difference*difference;
1585 distortion[CompositeChannels]+=difference*difference;
1586 }
1587 if (((channel & IndexChannel) != 0) &&
1588 (image->colorspace == CMYKColorspace) &&
1589 (reconstruct_image->colorspace == CMYKColorspace))
1590 {
1591 difference=reconstruct_phash[IndexChannel].Q[i]-
1592 image_phash[IndexChannel].Q[i];
1593 distortion[IndexChannel]+=difference*difference;
1594 distortion[CompositeChannels]+=difference*difference;
1595 }
1596 }
1597 /*
1598 Free resources.
1599 */
1600 reconstruct_phash=(ChannelPerceptualHash *) RelinquishMagickMemory(
1601 reconstruct_phash);
1602 image_phash=(ChannelPerceptualHash *) RelinquishMagickMemory(image_phash);
1603 return(MagickTrue);
1604}
1605
1606static MagickBooleanType GetRootMeanSquaredDistortion(const Image *image,
1607 const Image *reconstruct_image,const ChannelType channel,double *distortion,
1608 ExceptionInfo *exception)
1609{
1610 MagickBooleanType
1611 status;
1612
1613 status=GetMeanSquaredDistortion(image,reconstruct_image,channel,distortion,
1614 exception);
1615 if ((channel & RedChannel) != 0)
1616 distortion[RedChannel]=sqrt(distortion[RedChannel]);
1617 if ((channel & GreenChannel) != 0)
1618 distortion[GreenChannel]=sqrt(distortion[GreenChannel]);
1619 if ((channel & BlueChannel) != 0)
1620 distortion[BlueChannel]=sqrt(distortion[BlueChannel]);
1621 if (((channel & OpacityChannel) != 0) &&
1622 (image->matte != MagickFalse))
1623 distortion[OpacityChannel]=sqrt(distortion[OpacityChannel]);
1624 if (((channel & IndexChannel) != 0) &&
1625 (image->colorspace == CMYKColorspace))
1626 distortion[BlackChannel]=sqrt(distortion[BlackChannel]);
1627 distortion[CompositeChannels]=sqrt(distortion[CompositeChannels]);
1628 return(status);
1629}
1630
1631MagickExport MagickBooleanType GetImageChannelDistortion(Image *image,
1632 const Image *reconstruct_image,const ChannelType channel,
1633 const MetricType metric,double *distortion,ExceptionInfo *exception)
1634{
1635 double
1636 *channel_distortion;
1637
1638 MagickBooleanType
1639 status;
1640
1641 size_t
1642 length;
1643
1644 assert(image != (Image *) NULL);
1645 assert(image->signature == MagickCoreSignature);
1646 assert(reconstruct_image != (const Image *) NULL);
1647 assert(reconstruct_image->signature == MagickCoreSignature);
1648 assert(distortion != (double *) NULL);
1649 if (IsEventLogging() != MagickFalse)
1650 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1651 *distortion=0.0;
1652 if (metric != PerceptualHashErrorMetric)
1653 if (ValidateImageMorphology(image,reconstruct_image) == MagickFalse)
1654 ThrowBinaryException(ImageError,"ImageMorphologyDiffers",image->filename);
1655 /*
1656 Get image distortion.
1657 */
1658 length=CompositeChannels+1UL;
1659 channel_distortion=(double *) AcquireQuantumMemory(length,
1660 sizeof(*channel_distortion));
1661 if (channel_distortion == (double *) NULL)
1662 ThrowFatalException(ResourceLimitFatalError,"MemoryAllocationFailed");
1663 (void) memset(channel_distortion,0,length*sizeof(*channel_distortion));
1664 switch (metric)
1665 {
1666 case AbsoluteErrorMetric:
1667 {
1668 status=GetAbsoluteDistortion(image,reconstruct_image,channel,
1669 channel_distortion,exception);
1670 break;
1671 }
1672 case FuzzErrorMetric:
1673 {
1674 status=GetFuzzDistortion(image,reconstruct_image,channel,
1675 channel_distortion,exception);
1676 break;
1677 }
1678 case MeanAbsoluteErrorMetric:
1679 {
1680 status=GetMeanAbsoluteDistortion(image,reconstruct_image,channel,
1681 channel_distortion,exception);
1682 break;
1683 }
1684 case MeanErrorPerPixelMetric:
1685 {
1686 status=GetMeanErrorPerPixel(image,reconstruct_image,channel,
1687 channel_distortion,exception);
1688 break;
1689 }
1690 case MeanSquaredErrorMetric:
1691 {
1692 status=GetMeanSquaredDistortion(image,reconstruct_image,channel,
1693 channel_distortion,exception);
1694 break;
1695 }
1696 case NormalizedCrossCorrelationErrorMetric:
1697 default:
1698 {
1699 status=GetNormalizedCrossCorrelationDistortion(image,reconstruct_image,
1700 channel,channel_distortion,exception);
1701 break;
1702 }
1703 case PeakAbsoluteErrorMetric:
1704 {
1705 status=GetPeakAbsoluteDistortion(image,reconstruct_image,channel,
1706 channel_distortion,exception);
1707 break;
1708 }
1709 case PeakSignalToNoiseRatioMetric:
1710 {
1711 status=GetPeakSignalToNoiseRatio(image,reconstruct_image,channel,
1712 channel_distortion,exception);
1713 break;
1714 }
1715 case PerceptualHashErrorMetric:
1716 {
1717 status=GetPerceptualHashDistortion(image,reconstruct_image,channel,
1718 channel_distortion,exception);
1719 break;
1720 }
1721 case RootMeanSquaredErrorMetric:
1722 {
1723 status=GetRootMeanSquaredDistortion(image,reconstruct_image,channel,
1724 channel_distortion,exception);
1725 break;
1726 }
1727 }
1728 *distortion=channel_distortion[CompositeChannels];
1729 channel_distortion=(double *) RelinquishMagickMemory(channel_distortion);
1730 (void) FormatImageProperty(image,"distortion","%.*g",GetMagickPrecision(),
1731 *distortion);
1732 return(status);
1733}
1734
1735/*
1736%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1737% %
1738% %
1739% %
1740% G e t I m a g e C h a n n e l D i s t o r t i o n s %
1741% %
1742% %
1743% %
1744%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1745%
1746% GetImageChannelDistortions() compares the image channels of an image to a
1747% reconstructed image and returns the specified distortion metric for each
1748% channel.
1749%
1750% The format of the GetImageChannelDistortions method is:
1751%
1752% double *GetImageChannelDistortions(const Image *image,
1753% const Image *reconstruct_image,const MetricType metric,
1754% ExceptionInfo *exception)
1755%
1756% A description of each parameter follows:
1757%
1758% o image: the image.
1759%
1760% o reconstruct_image: the reconstruct image.
1761%
1762% o metric: the metric.
1763%
1764% o exception: return any errors or warnings in this structure.
1765%
1766*/
1767MagickExport double *GetImageChannelDistortions(Image *image,
1768 const Image *reconstruct_image,const MetricType metric,
1769 ExceptionInfo *exception)
1770{
1771 double
1772 *channel_distortion;
1773
1774 MagickBooleanType
1775 status;
1776
1777 size_t
1778 length;
1779
1780 assert(image != (Image *) NULL);
1781 assert(image->signature == MagickCoreSignature);
1782 assert(reconstruct_image != (const Image *) NULL);
1783 assert(reconstruct_image->signature == MagickCoreSignature);
1784 if (IsEventLogging() != MagickFalse)
1785 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1786 if (metric != PerceptualHashErrorMetric)
1787 if (ValidateImageMorphology(image,reconstruct_image) == MagickFalse)
1788 {
1789 (void) ThrowMagickException(&image->exception,GetMagickModule(),
1790 ImageError,"ImageMorphologyDiffers","`%s'",image->filename);
1791 return((double *) NULL);
1792 }
1793 /*
1794 Get image distortion.
1795 */
1796 length=CompositeChannels+1UL;
1797 channel_distortion=(double *) AcquireQuantumMemory(length,
1798 sizeof(*channel_distortion));
1799 if (channel_distortion == (double *) NULL)
1800 ThrowFatalException(ResourceLimitFatalError,"MemoryAllocationFailed");
1801 (void) memset(channel_distortion,0,length*
1802 sizeof(*channel_distortion));
1803 status=MagickTrue;
1804 switch (metric)
1805 {
1806 case AbsoluteErrorMetric:
1807 {
1808 status=GetAbsoluteDistortion(image,reconstruct_image,CompositeChannels,
1809 channel_distortion,exception);
1810 break;
1811 }
1812 case FuzzErrorMetric:
1813 {
1814 status=GetFuzzDistortion(image,reconstruct_image,CompositeChannels,
1815 channel_distortion,exception);
1816 break;
1817 }
1818 case MeanAbsoluteErrorMetric:
1819 {
1820 status=GetMeanAbsoluteDistortion(image,reconstruct_image,
1821 CompositeChannels,channel_distortion,exception);
1822 break;
1823 }
1824 case MeanErrorPerPixelMetric:
1825 {
1826 status=GetMeanErrorPerPixel(image,reconstruct_image,CompositeChannels,
1827 channel_distortion,exception);
1828 break;
1829 }
1830 case MeanSquaredErrorMetric:
1831 {
1832 status=GetMeanSquaredDistortion(image,reconstruct_image,CompositeChannels,
1833 channel_distortion,exception);
1834 break;
1835 }
1836 case NormalizedCrossCorrelationErrorMetric:
1837 default:
1838 {
1839 status=GetNormalizedCrossCorrelationDistortion(image,reconstruct_image,
1840 CompositeChannels,channel_distortion,exception);
1841 break;
1842 }
1843 case PeakAbsoluteErrorMetric:
1844 {
1845 status=GetPeakAbsoluteDistortion(image,reconstruct_image,
1846 CompositeChannels,channel_distortion,exception);
1847 break;
1848 }
1849 case PeakSignalToNoiseRatioMetric:
1850 {
1851 status=GetPeakSignalToNoiseRatio(image,reconstruct_image,
1852 CompositeChannels,channel_distortion,exception);
1853 break;
1854 }
1855 case PerceptualHashErrorMetric:
1856 {
1857 status=GetPerceptualHashDistortion(image,reconstruct_image,
1858 CompositeChannels,channel_distortion,exception);
1859 break;
1860 }
1861 case RootMeanSquaredErrorMetric:
1862 {
1863 status=GetRootMeanSquaredDistortion(image,reconstruct_image,
1864 CompositeChannels,channel_distortion,exception);
1865 break;
1866 }
1867 }
1868 if (status == MagickFalse)
1869 {
1870 channel_distortion=(double *) RelinquishMagickMemory(channel_distortion);
1871 return((double *) NULL);
1872 }
1873 return(channel_distortion);
1874}
1875
1876/*
1877%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1878% %
1879% %
1880% %
1881% I s I m a g e s E q u a l %
1882% %
1883% %
1884% %
1885%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1886%
1887% IsImagesEqual() measures the difference between colors at each pixel
1888% location of two images. A value other than 0 means the colors match
1889% exactly. Otherwise an error measure is computed by summing over all
1890% pixels in an image the distance squared in RGB space between each image
1891% pixel and its corresponding pixel in the reconstruct image. The error
1892% measure is assigned to these image members:
1893%
1894% o mean_error_per_pixel: The mean error for any single pixel in
1895% the image.
1896%
1897% o normalized_mean_error: The normalized mean quantization error for
1898% any single pixel in the image. This distance measure is normalized to
1899% a range between 0 and 1. It is independent of the range of red, green,
1900% and blue values in the image.
1901%
1902% o normalized_maximum_error: The normalized maximum quantization
1903% error for any single pixel in the image. This distance measure is
1904% normalized to a range between 0 and 1. It is independent of the range
1905% of red, green, and blue values in your image.
1906%
1907% A small normalized mean square error, accessed as
1908% image->normalized_mean_error, suggests the images are very similar in
1909% spatial layout and color.
1910%
1911% The format of the IsImagesEqual method is:
1912%
1913% MagickBooleanType IsImagesEqual(Image *image,
1914% const Image *reconstruct_image)
1915%
1916% A description of each parameter follows.
1917%
1918% o image: the image.
1919%
1920% o reconstruct_image: the reconstruct image.
1921%
1922*/
1923MagickExport MagickBooleanType IsImagesEqual(Image *image,
1924 const Image *reconstruct_image)
1925{
1926 CacheView
1927 *image_view,
1928 *reconstruct_view;
1929
1930 ExceptionInfo
1931 *exception;
1932
1933 MagickBooleanType
1934 status;
1935
1936 MagickRealType
1937 area,
1938 gamma,
1939 maximum_error,
1940 mean_error,
1941 mean_error_per_pixel;
1942
1943 size_t
1944 columns,
1945 rows;
1946
1947 ssize_t
1948 y;
1949
1950 assert(image != (Image *) NULL);
1951 assert(image->signature == MagickCoreSignature);
1952 assert(reconstruct_image != (const Image *) NULL);
1953 assert(reconstruct_image->signature == MagickCoreSignature);
1954 exception=(&image->exception);
1955 if (ValidateImageMorphology(image,reconstruct_image) == MagickFalse)
1956 ThrowBinaryException(ImageError,"ImageMorphologyDiffers",image->filename);
1957 area=0.0;
1958 maximum_error=0.0;
1959 mean_error_per_pixel=0.0;
1960 mean_error=0.0;
1961 SetImageDistortionBounds(image,reconstruct_image,&columns,&rows);
1962 image_view=AcquireVirtualCacheView(image,exception);
1963 reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
1964 for (y=0; y < (ssize_t) rows; y++)
1965 {
1966 const IndexPacket
1967 *magick_restrict indexes,
1968 *magick_restrict reconstruct_indexes;
1969
1970 const PixelPacket
1971 *magick_restrict p,
1972 *magick_restrict q;
1973
1974 ssize_t
1975 x;
1976
1977 p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
1978 q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
1979 if ((p == (const PixelPacket *) NULL) || (q == (const PixelPacket *) NULL))
1980 break;
1981 indexes=GetCacheViewVirtualIndexQueue(image_view);
1982 reconstruct_indexes=GetCacheViewVirtualIndexQueue(reconstruct_view);
1983 for (x=0; x < (ssize_t) columns; x++)
1984 {
1985 MagickRealType
1986 distance;
1987
1988 distance=fabs((double) GetPixelRed(p)-(double) GetPixelRed(q));
1989 mean_error_per_pixel+=distance;
1990 mean_error+=distance*distance;
1991 if (distance > maximum_error)
1992 maximum_error=distance;
1993 area++;
1994 distance=fabs((double) GetPixelGreen(p)-(double) GetPixelGreen(q));
1995 mean_error_per_pixel+=distance;
1996 mean_error+=distance*distance;
1997 if (distance > maximum_error)
1998 maximum_error=distance;
1999 area++;
2000 distance=fabs((double) GetPixelBlue(p)-(double) GetPixelBlue(q));
2001 mean_error_per_pixel+=distance;
2002 mean_error+=distance*distance;
2003 if (distance > maximum_error)
2004 maximum_error=distance;
2005 area++;
2006 if (image->matte != MagickFalse)
2007 {
2008 distance=fabs((double) GetPixelOpacity(p)-(double)
2009 GetPixelOpacity(q));
2010 mean_error_per_pixel+=distance;
2011 mean_error+=distance*distance;
2012 if (distance > maximum_error)
2013 maximum_error=distance;
2014 area++;
2015 }
2016 if ((image->colorspace == CMYKColorspace) &&
2017 (reconstruct_image->colorspace == CMYKColorspace))
2018 {
2019 distance=fabs((double) GetPixelIndex(indexes+x)-(double)
2020 GetPixelIndex(reconstruct_indexes+x));
2021 mean_error_per_pixel+=distance;
2022 mean_error+=distance*distance;
2023 if (distance > maximum_error)
2024 maximum_error=distance;
2025 area++;
2026 }
2027 p++;
2028 q++;
2029 }
2030 }
2031 reconstruct_view=DestroyCacheView(reconstruct_view);
2032 image_view=DestroyCacheView(image_view);
2033 gamma=PerceptibleReciprocal(area);
2034 image->error.mean_error_per_pixel=gamma*mean_error_per_pixel;
2035 image->error.normalized_mean_error=gamma*QuantumScale*QuantumScale*mean_error;
2036 image->error.normalized_maximum_error=QuantumScale*maximum_error;
2037 status=image->error.mean_error_per_pixel == 0.0 ? MagickTrue : MagickFalse;
2038 return(status);
2039}
2040
2041/*
2042%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2043% %
2044% %
2045% %
2046% S i m i l a r i t y I m a g e %
2047% %
2048% %
2049% %
2050%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2051%
2052% SimilarityImage() compares the reference image of the image and returns the
2053% best match offset. In addition, it returns a similarity image such that an
2054% exact match location is completely white and if none of the pixels match,
2055% black, otherwise some gray level in-between.
2056%
2057% The format of the SimilarityImageImage method is:
2058%
2059% Image *SimilarityImage(const Image *image,const Image *reference,
2060% RectangleInfo *offset,double *similarity,ExceptionInfo *exception)
2061%
2062% A description of each parameter follows:
2063%
2064% o image: the image.
2065%
2066% o reference: find an area of the image that closely resembles this image.
2067%
2068% o the best match offset of the reference image within the image.
2069%
2070% o similarity: the computed similarity between the images.
2071%
2072% o exception: return any errors or warnings in this structure.
2073%
2074*/
2075
2076static double GetSimilarityMetric(const Image *image,const Image *reference,
2077 const MetricType metric,const ssize_t x_offset,const ssize_t y_offset,
2078 ExceptionInfo *exception)
2079{
2080 double
2081 distortion;
2082
2083 Image
2084 *similarity_image;
2085
2086 MagickBooleanType
2087 status;
2088
2089 RectangleInfo
2090 geometry;
2091
2092 SetGeometry(reference,&geometry);
2093 geometry.x=x_offset;
2094 geometry.y=y_offset;
2095 similarity_image=CropImage(image,&geometry,exception);
2096 if (similarity_image == (Image *) NULL)
2097 return(0.0);
2098 distortion=0.0;
2099 status=GetImageDistortion(similarity_image,reference,metric,&distortion,
2100 exception);
2101 (void) status;
2102 similarity_image=DestroyImage(similarity_image);
2103 return(distortion);
2104}
2105
2106MagickExport Image *SimilarityImage(Image *image,const Image *reference,
2107 RectangleInfo *offset,double *similarity_metric,ExceptionInfo *exception)
2108{
2109 Image
2110 *similarity_image;
2111
2112 similarity_image=SimilarityMetricImage(image,reference,
2113 RootMeanSquaredErrorMetric,offset,similarity_metric,exception);
2114 return(similarity_image);
2115}
2116
2117MagickExport Image *SimilarityMetricImage(Image *image,const Image *reference,
2118 const MetricType metric,RectangleInfo *offset,double *similarity_metric,
2119 ExceptionInfo *exception)
2120{
2121#define SimilarityImageTag "Similarity/Image"
2122
2123 CacheView
2124 *similarity_view;
2125
2126 const char
2127 *artifact;
2128
2129 double
2130 similarity_threshold;
2131
2132 Image
2133 *similarity_image = (Image *) NULL;
2134
2135 MagickBooleanType
2136 status;
2137
2138 MagickOffsetType
2139 progress;
2140
2141 size_t
2142 rows;
2143
2144 ssize_t
2145 y;
2146
2147 assert(image != (const Image *) NULL);
2148 assert(image->signature == MagickCoreSignature);
2149 assert(exception != (ExceptionInfo *) NULL);
2150 assert(exception->signature == MagickCoreSignature);
2151 assert(offset != (RectangleInfo *) NULL);
2152 if (IsEventLogging() != MagickFalse)
2153 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2154 SetGeometry(reference,offset);
2155 *similarity_metric=MagickMaximumValue;
2156 if (ValidateImageMorphology(image,reference) == MagickFalse)
2157 ThrowImageException(ImageError,"ImageMorphologyDiffers");
2158 if ((image->columns < reference->columns) || (image->rows < reference->rows))
2159 {
2160 (void) ThrowMagickException(&image->exception,GetMagickModule(),
2161 OptionWarning,"GeometryDoesNotContainImage","`%s'",image->filename);
2162 return((Image *) NULL);
2163 }
2164 if (metric == PeakAbsoluteErrorMetric)
2165 {
2166 (void) ThrowMagickException(&image->exception,GetMagickModule(),
2167 OptionError,"InvalidUseOfOption","`%s'",image->filename);
2168 return((Image *) NULL);
2169 }
2170 similarity_image=CloneImage(image,image->columns,image->rows,MagickTrue,
2171 exception);
2172 if (similarity_image == (Image *) NULL)
2173 return((Image *) NULL);
2174 similarity_image->depth=MAGICKCORE_QUANTUM_DEPTH;
2175 similarity_image->matte=MagickFalse;
2176 status=SetImageStorageClass(similarity_image,DirectClass);
2177 if (status == MagickFalse)
2178 {
2179 InheritException(exception,&similarity_image->exception);
2180 return(DestroyImage(similarity_image));
2181 }
2182 /*
2183 Measure similarity of reference image against image.
2184 */
2185 similarity_threshold=(-1.0);
2186 artifact=GetImageArtifact(image,"compare:similarity-threshold");
2187 if (artifact != (const char *) NULL)
2188 similarity_threshold=StringToDouble(artifact,(char **) NULL);
2189 status=MagickTrue;
2190 progress=0;
2191 similarity_view=AcquireVirtualCacheView(similarity_image,exception);
2192 rows=similarity_image->rows;
2193 for (y=0; y < (ssize_t) rows; y++)
2194 {
2195 double
2196 similarity;
2197
2198 ssize_t
2199 x;
2200
2201 PixelPacket
2202 *magick_restrict q;
2203
2204 if (status == MagickFalse)
2205 continue;
2206 if (*similarity_metric <= similarity_threshold)
2207 continue;
2208 q=GetCacheViewAuthenticPixels(similarity_view,0,y,similarity_image->columns,
2209 1,exception);
2210 if (q == (const PixelPacket *) NULL)
2211 {
2212 status=MagickFalse;
2213 continue;
2214 }
2215 for (x=0; x < (ssize_t) similarity_image->columns; x++)
2216 {
2217 if (*similarity_metric <= similarity_threshold)
2218 break;
2219 similarity=GetSimilarityMetric(image,reference,metric,x,y,exception);
2220 switch (metric)
2221 {
2222 case NormalizedCrossCorrelationErrorMetric:
2223 case UndefinedErrorMetric:
2224 {
2225 similarity=1.0-similarity;
2226 break;
2227 }
2228 case PerceptualHashErrorMetric:
2229 {
2230 similarity=1.0-MagickMin(0.01*similarity,1.0);
2231 break;
2232 }
2233 default:
2234 break;
2235 }
2236 if (similarity < *similarity_metric)
2237 {
2238 *similarity_metric=similarity;
2239 offset->x=x;
2240 offset->y=y;
2241 }
2242 if ((metric == PeakSignalToNoiseRatioMetric) &&
2243 (fabs(similarity) < MagickEpsilon))
2244 similarity=1.0-similarity;
2245 switch (metric)
2246 {
2247 case AbsoluteErrorMetric:
2248 case FuzzErrorMetric:
2249 case MeanAbsoluteErrorMetric:
2250 case MeanErrorPerPixelMetric:
2251 case MeanSquaredErrorMetric:
2252 case NormalizedCrossCorrelationErrorMetric:
2253 case PeakAbsoluteErrorMetric:
2254 case PerceptualHashErrorMetric:
2255 case RootMeanSquaredErrorMetric:
2256 {
2257 SetPixelRed(q,ClampToQuantum((double) QuantumRange-QuantumRange*
2258 similarity));
2259 break;
2260 }
2261 default:
2262 {
2263 SetPixelRed(q,ClampToQuantum((double) QuantumRange*similarity));
2264 break;
2265 }
2266 }
2267 SetPixelGreen(q,GetPixelRed(q));
2268 SetPixelBlue(q,GetPixelRed(q));
2269 q++;
2270 }
2271 if (SyncCacheViewAuthenticPixels(similarity_view,exception) == MagickFalse)
2272 status=MagickFalse;
2273 if (image->progress_monitor != (MagickProgressMonitor) NULL)
2274 {
2275 MagickBooleanType
2276 proceed;
2277
2278 progress++;
2279 proceed=SetImageProgress(image,SimilarityImageTag,progress,image->rows);
2280 if (proceed == MagickFalse)
2281 status=MagickFalse;
2282 }
2283 }
2284 similarity_view=DestroyCacheView(similarity_view);
2285 (void) SetImageType(similarity_image,GrayscaleType);
2286 if (status == MagickFalse)
2287 similarity_image=DestroyImage(similarity_image);
2288 return(similarity_image);
2289}