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