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