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