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