2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
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 %
13 % MagickCore Image Comparison Methods %
20 % Copyright 1999-2010 ImageMagick Studio LLC, a non-profit organization %
21 % dedicated to making software imaging solutions freely available. %
23 % You may not use this file except in compliance with the License. You may %
24 % obtain a copy of the License at %
26 % http://www.imagemagick.org/script/license.php %
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. %
34 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
43 #include "magick/studio.h"
44 #include "magick/artifact.h"
45 #include "magick/cache-view.h"
46 #include "magick/client.h"
47 #include "magick/color.h"
48 #include "magick/color-private.h"
49 #include "magick/colorspace.h"
50 #include "magick/colorspace-private.h"
51 #include "magick/compare.h"
52 #include "magick/composite-private.h"
53 #include "magick/constitute.h"
54 #include "magick/exception-private.h"
55 #include "magick/geometry.h"
56 #include "magick/image-private.h"
57 #include "magick/list.h"
58 #include "magick/log.h"
59 #include "magick/memory_.h"
60 #include "magick/monitor.h"
61 #include "magick/monitor-private.h"
62 #include "magick/option.h"
63 #include "magick/pixel-private.h"
64 #include "magick/resource_.h"
65 #include "magick/string_.h"
66 #include "magick/statistic.h"
67 #include "magick/utility.h"
68 #include "magick/version.h"
71 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
75 % C o m p a r e I m a g e C h a n n e l s %
79 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
81 % CompareImageChannels() compares one or more image channels of an image
82 % to a reconstructed image and returns the difference image.
84 % The format of the CompareImageChannels method is:
86 % Image *CompareImageChannels(const Image *image,
87 % const Image *reconstruct_image,const ChannelType channel,
88 % const MetricType metric,double *distortion,ExceptionInfo *exception)
90 % A description of each parameter follows:
94 % o reconstruct_image: the reconstruct image.
96 % o channel: the channel.
98 % o metric: the metric.
100 % o distortion: the computed distortion between the images.
102 % o exception: return any errors or warnings in this structure.
106 MagickExport Image *CompareImages(Image *image,const Image *reconstruct_image,
107 const MetricType metric,double *distortion,ExceptionInfo *exception)
112 highlight_image=CompareImageChannels(image,reconstruct_image,AllChannels,
113 metric,distortion,exception);
114 return(highlight_image);
117 MagickExport Image *CompareImageChannels(Image *image,
118 const Image *reconstruct_image,const ChannelType channel,
119 const MetricType metric,double *distortion,ExceptionInfo *exception)
144 assert(image != (Image *) NULL);
145 assert(image->signature == MagickSignature);
146 if (image->debug != MagickFalse)
147 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
148 assert(reconstruct_image != (const Image *) NULL);
149 assert(reconstruct_image->signature == MagickSignature);
150 assert(distortion != (double *) NULL);
152 if (image->debug != MagickFalse)
153 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
154 if ((reconstruct_image->columns != image->columns) ||
155 (reconstruct_image->rows != image->rows))
156 ThrowImageException(ImageError,"ImageSizeDiffers");
157 status=GetImageChannelDistortion(image,reconstruct_image,channel,metric,
158 distortion,exception);
159 if (status == MagickFalse)
160 return((Image *) NULL);
161 difference_image=CloneImage(image,0,0,MagickTrue,exception);
162 if (difference_image == (Image *) NULL)
163 return((Image *) NULL);
164 (void) SetImageAlphaChannel(difference_image,OpaqueAlphaChannel);
165 highlight_image=CloneImage(image,image->columns,image->rows,MagickTrue,
167 if (highlight_image == (Image *) NULL)
169 difference_image=DestroyImage(difference_image);
170 return((Image *) NULL);
172 if (SetImageStorageClass(highlight_image,DirectClass) == MagickFalse)
174 InheritException(exception,&highlight_image->exception);
175 difference_image=DestroyImage(difference_image);
176 highlight_image=DestroyImage(highlight_image);
177 return((Image *) NULL);
179 (void) SetImageAlphaChannel(highlight_image,OpaqueAlphaChannel);
180 (void) QueryMagickColor("#f1001ecc",&highlight,exception);
181 artifact=GetImageArtifact(image,"highlight-color");
182 if (artifact != (const char *) NULL)
183 (void) QueryMagickColor(artifact,&highlight,exception);
184 (void) QueryMagickColor("#ffffffcc",&lowlight,exception);
185 artifact=GetImageArtifact(image,"lowlight-color");
186 if (artifact != (const char *) NULL)
187 (void) QueryMagickColor(artifact,&lowlight,exception);
188 if (highlight_image->colorspace == CMYKColorspace)
190 ConvertRGBToCMYK(&highlight);
191 ConvertRGBToCMYK(&lowlight);
194 Generate difference image.
197 GetMagickPixelPacket(image,&zero);
198 image_view=AcquireCacheView(image);
199 reconstruct_view=AcquireCacheView(reconstruct_image);
200 highlight_view=AcquireCacheView(highlight_image);
201 #if defined(MAGICKCORE_OPENMP_SUPPORT)
202 #pragma omp parallel for schedule(dynamic,4) shared(status)
204 for (y=0; y < (ssize_t) image->rows; y++)
213 register const IndexPacket
215 *restrict reconstruct_indexes;
217 register const PixelPacket
222 *restrict highlight_indexes;
230 if (status == MagickFalse)
232 p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
233 q=GetCacheViewVirtualPixels(reconstruct_view,0,y,reconstruct_image->columns,
235 r=QueueCacheViewAuthenticPixels(highlight_view,0,y,highlight_image->columns,
237 if ((p == (const PixelPacket *) NULL) ||
238 (q == (const PixelPacket *) NULL) || (r == (PixelPacket *) NULL))
243 indexes=GetCacheViewVirtualIndexQueue(image_view);
244 reconstruct_indexes=GetCacheViewVirtualIndexQueue(reconstruct_view);
245 highlight_indexes=GetCacheViewAuthenticIndexQueue(highlight_view);
247 reconstruct_pixel=zero;
248 for (x=0; x < (ssize_t) image->columns; x++)
253 SetMagickPixelPacket(image,p,indexes+x,&pixel);
254 SetMagickPixelPacket(reconstruct_image,q,reconstruct_indexes+x,
256 difference=MagickFalse;
257 if (channel == AllChannels)
259 if (IsMagickColorSimilar(&pixel,&reconstruct_pixel) == MagickFalse)
260 difference=MagickTrue;
264 if (((channel & RedChannel) != 0) && (p->red != q->red))
265 difference=MagickTrue;
266 if (((channel & GreenChannel) != 0) && (p->green != q->green))
267 difference=MagickTrue;
268 if (((channel & BlueChannel) != 0) && (p->blue != q->blue))
269 difference=MagickTrue;
270 if (((channel & OpacityChannel) != 0) &&
271 (image->matte != MagickFalse) && (p->opacity != q->opacity))
272 difference=MagickTrue;
273 if ((((channel & IndexChannel) != 0) &&
274 (image->colorspace == CMYKColorspace) &&
275 (reconstruct_image->colorspace == CMYKColorspace)) &&
276 (indexes[x] != reconstruct_indexes[x]))
277 difference=MagickTrue;
279 if (difference != MagickFalse)
280 SetPixelPacket(highlight_image,&highlight,r,highlight_indexes+x);
282 SetPixelPacket(highlight_image,&lowlight,r,highlight_indexes+x);
287 sync=SyncCacheViewAuthenticPixels(highlight_view,exception);
288 if (sync == MagickFalse)
291 highlight_view=DestroyCacheView(highlight_view);
292 reconstruct_view=DestroyCacheView(reconstruct_view);
293 image_view=DestroyCacheView(image_view);
294 (void) CompositeImage(difference_image,image->compose,highlight_image,0,0);
295 highlight_image=DestroyImage(highlight_image);
296 if (status == MagickFalse)
297 difference_image=DestroyImage(difference_image);
298 return(difference_image);
302 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
306 % 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 %
310 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
312 % GetImageChannelDistortion() compares one or more image channels of an image
313 % to a reconstructed image and returns the specified distortion metric.
315 % The format of the CompareImageChannels method is:
317 % MagickBooleanType GetImageChannelDistortion(const Image *image,
318 % const Image *reconstruct_image,const ChannelType channel,
319 % const MetricType metric,double *distortion,ExceptionInfo *exception)
321 % A description of each parameter follows:
323 % o image: the image.
325 % o reconstruct_image: the reconstruct image.
327 % o channel: the channel.
329 % o metric: the metric.
331 % o distortion: the computed distortion between the images.
333 % o exception: return any errors or warnings in this structure.
337 MagickExport MagickBooleanType GetImageDistortion(Image *image,
338 const Image *reconstruct_image,const MetricType metric,double *distortion,
339 ExceptionInfo *exception)
344 status=GetImageChannelDistortion(image,reconstruct_image,AllChannels,
345 metric,distortion,exception);
349 static MagickBooleanType GetAbsoluteError(const Image *image,
350 const Image *reconstruct_image,const ChannelType channel,double *distortion,
351 ExceptionInfo *exception)
367 Compute the absolute difference in pixels between two images.
370 GetMagickPixelPacket(image,&zero);
371 image_view=AcquireCacheView(image);
372 reconstruct_view=AcquireCacheView(reconstruct_image);
373 #if defined(MAGICKCORE_OPENMP_SUPPORT)
374 #pragma omp parallel for schedule(dynamic,4) shared(status)
376 for (y=0; y < (ssize_t) image->rows; y++)
379 channel_distortion[AllChannels+1];
385 register const IndexPacket
387 *restrict reconstruct_indexes;
389 register const PixelPacket
397 if (status == MagickFalse)
399 p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
400 q=GetCacheViewVirtualPixels(reconstruct_view,0,y,reconstruct_image->columns,
402 if ((p == (const PixelPacket *) NULL) || (q == (const PixelPacket *) NULL))
407 indexes=GetCacheViewVirtualIndexQueue(image_view);
408 reconstruct_indexes=GetCacheViewVirtualIndexQueue(reconstruct_view);
410 reconstruct_pixel=pixel;
411 (void) ResetMagickMemory(channel_distortion,0,sizeof(channel_distortion));
412 for (x=0; x < (ssize_t) image->columns; x++)
414 SetMagickPixelPacket(image,p,indexes+x,&pixel);
415 SetMagickPixelPacket(reconstruct_image,q,reconstruct_indexes+x,
417 if (IsMagickColorSimilar(&pixel,&reconstruct_pixel) == MagickFalse)
419 if ((channel & RedChannel) != 0)
420 channel_distortion[RedChannel]++;
421 if ((channel & GreenChannel) != 0)
422 channel_distortion[GreenChannel]++;
423 if ((channel & BlueChannel) != 0)
424 channel_distortion[BlueChannel]++;
425 if (((channel & OpacityChannel) != 0) &&
426 (image->matte != MagickFalse))
427 channel_distortion[OpacityChannel]++;
428 if (((channel & IndexChannel) != 0) &&
429 (image->colorspace == CMYKColorspace))
430 channel_distortion[BlackChannel]++;
431 channel_distortion[AllChannels]++;
436 #if defined(MAGICKCORE_OPENMP_SUPPORT)
437 #pragma omp critical (MagickCore_GetAbsoluteError)
439 for (i=0; i <= (ssize_t) AllChannels; i++)
440 distortion[i]+=channel_distortion[i];
442 reconstruct_view=DestroyCacheView(reconstruct_view);
443 image_view=DestroyCacheView(image_view);
447 static size_t GetNumberChannels(const Image *image,
448 const ChannelType channel)
454 if ((channel & RedChannel) != 0)
456 if ((channel & GreenChannel) != 0)
458 if ((channel & BlueChannel) != 0)
460 if (((channel & OpacityChannel) != 0) &&
461 (image->matte != MagickFalse))
463 if (((channel & IndexChannel) != 0) &&
464 (image->colorspace == CMYKColorspace))
469 static MagickBooleanType GetMeanAbsoluteError(const Image *image,
470 const Image *reconstruct_image,const ChannelType channel,
471 double *distortion,ExceptionInfo *exception)
487 image_view=AcquireCacheView(image);
488 reconstruct_view=AcquireCacheView(reconstruct_image);
489 #if defined(MAGICKCORE_OPENMP_SUPPORT)
490 #pragma omp parallel for schedule(dynamic,4) shared(status)
492 for (y=0; y < (ssize_t) image->rows; y++)
495 channel_distortion[AllChannels+1];
497 register const IndexPacket
499 *restrict reconstruct_indexes;
501 register const PixelPacket
509 if (status == MagickFalse)
511 p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
512 q=GetCacheViewVirtualPixels(reconstruct_view,0,y,
513 reconstruct_image->columns,1,exception);
514 if ((p == (const PixelPacket *) NULL) || (q == (const PixelPacket *) NULL))
519 indexes=GetCacheViewVirtualIndexQueue(image_view);
520 reconstruct_indexes=GetCacheViewVirtualIndexQueue(reconstruct_view);
521 (void) ResetMagickMemory(channel_distortion,0,sizeof(channel_distortion));
522 for (x=0; x < (ssize_t) image->columns; x++)
527 if ((channel & RedChannel) != 0)
529 distance=QuantumScale*fabs(p->red-(double) q->red);
530 channel_distortion[RedChannel]+=distance;
531 channel_distortion[AllChannels]+=distance;
533 if ((channel & GreenChannel) != 0)
535 distance=QuantumScale*fabs(p->green-(double) q->green);
536 channel_distortion[GreenChannel]+=distance;
537 channel_distortion[AllChannels]+=distance;
539 if ((channel & BlueChannel) != 0)
541 distance=QuantumScale*fabs(p->blue-(double) q->blue);
542 channel_distortion[BlueChannel]+=distance;
543 channel_distortion[AllChannels]+=distance;
545 if (((channel & OpacityChannel) != 0) &&
546 (image->matte != MagickFalse))
548 distance=QuantumScale*fabs(p->opacity-(double) q->opacity);
549 channel_distortion[OpacityChannel]+=distance;
550 channel_distortion[AllChannels]+=distance;
552 if (((channel & IndexChannel) != 0) &&
553 (image->colorspace == CMYKColorspace))
555 distance=QuantumScale*fabs(indexes[x]-(double)
556 reconstruct_indexes[x]);
557 channel_distortion[BlackChannel]+=distance;
558 channel_distortion[AllChannels]+=distance;
563 #if defined(MAGICKCORE_OPENMP_SUPPORT)
564 #pragma omp critical (MagickCore_GetMeanAbsoluteError)
566 for (i=0; i <= (ssize_t) AllChannels; i++)
567 distortion[i]+=channel_distortion[i];
569 reconstruct_view=DestroyCacheView(reconstruct_view);
570 image_view=DestroyCacheView(image_view);
571 for (i=0; i <= (ssize_t) AllChannels; i++)
572 distortion[i]/=((double) image->columns*image->rows);
573 distortion[AllChannels]/=(double) GetNumberChannels(image,channel);
577 static MagickBooleanType GetMeanErrorPerPixel(Image *image,
578 const Image *reconstruct_image,const ChannelType channel,double *distortion,
579 ExceptionInfo *exception)
604 image_view=AcquireCacheView(image);
605 reconstruct_view=AcquireCacheView(reconstruct_image);
606 for (y=0; y < (ssize_t) image->rows; y++)
608 register const IndexPacket
610 *restrict reconstruct_indexes;
612 register const PixelPacket
619 p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
620 q=GetCacheViewVirtualPixels(reconstruct_view,0,y,reconstruct_image->columns,
622 if ((p == (const PixelPacket *) NULL) || (q == (const PixelPacket *) NULL))
627 indexes=GetCacheViewVirtualIndexQueue(image_view);
628 reconstruct_indexes=GetCacheViewVirtualIndexQueue(reconstruct_view);
629 for (x=0; x < (ssize_t) image->columns; x++)
634 if ((channel & OpacityChannel) != 0)
636 if (image->matte != MagickFalse)
637 alpha=(MagickRealType) (QuantumScale*(GetAlphaPixelComponent(p)));
638 if (reconstruct_image->matte != MagickFalse)
639 beta=(MagickRealType) (QuantumScale*GetAlphaPixelComponent(q));
641 if ((channel & RedChannel) != 0)
643 distance=fabs(alpha*p->red-beta*q->red);
644 distortion[RedChannel]+=distance;
645 distortion[AllChannels]+=distance;
646 mean_error+=distance*distance;
647 if (distance > maximum_error)
648 maximum_error=distance;
651 if ((channel & GreenChannel) != 0)
653 distance=fabs(alpha*p->green-beta*q->green);
654 distortion[GreenChannel]+=distance;
655 distortion[AllChannels]+=distance;
656 mean_error+=distance*distance;
657 if (distance > maximum_error)
658 maximum_error=distance;
661 if ((channel & BlueChannel) != 0)
663 distance=fabs(alpha*p->blue-beta*q->blue);
664 distortion[BlueChannel]+=distance;
665 distortion[AllChannels]+=distance;
666 mean_error+=distance*distance;
667 if (distance > maximum_error)
668 maximum_error=distance;
671 if (((channel & OpacityChannel) != 0) &&
672 (image->matte != MagickFalse))
674 distance=fabs((double) p->opacity-q->opacity);
675 distortion[OpacityChannel]+=distance;
676 distortion[AllChannels]+=distance;
677 mean_error+=distance*distance;
678 if (distance > maximum_error)
679 maximum_error=distance;
682 if (((channel & IndexChannel) != 0) &&
683 (image->colorspace == CMYKColorspace) &&
684 (reconstruct_image->colorspace == CMYKColorspace))
686 distance=fabs(alpha*indexes[x]-beta*reconstruct_indexes[x]);
687 distortion[BlackChannel]+=distance;
688 distortion[AllChannels]+=distance;
689 mean_error+=distance*distance;
690 if (distance > maximum_error)
691 maximum_error=distance;
698 reconstruct_view=DestroyCacheView(reconstruct_view);
699 image_view=DestroyCacheView(image_view);
700 image->error.mean_error_per_pixel=distortion[AllChannels]/area;
701 image->error.normalized_mean_error=QuantumScale*QuantumScale*mean_error/area;
702 image->error.normalized_maximum_error=QuantumScale*maximum_error;
706 static MagickBooleanType GetMeanSquaredError(const Image *image,
707 const Image *reconstruct_image,const ChannelType channel,
708 double *distortion,ExceptionInfo *exception)
724 image_view=AcquireCacheView(image);
725 reconstruct_view=AcquireCacheView(reconstruct_image);
726 #if defined(MAGICKCORE_OPENMP_SUPPORT)
727 #pragma omp parallel for schedule(dynamic,4) shared(status)
729 for (y=0; y < (ssize_t) image->rows; y++)
732 channel_distortion[AllChannels+1];
734 register const IndexPacket
736 *restrict reconstruct_indexes;
738 register const PixelPacket
746 if (status == MagickFalse)
748 p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
749 q=GetCacheViewVirtualPixels(reconstruct_view,0,y,
750 reconstruct_image->columns,1,exception);
751 if ((p == (const PixelPacket *) NULL) || (q == (const PixelPacket *) NULL))
756 indexes=GetCacheViewVirtualIndexQueue(image_view);
757 reconstruct_indexes=GetCacheViewVirtualIndexQueue(reconstruct_view);
758 (void) ResetMagickMemory(channel_distortion,0,sizeof(channel_distortion));
759 for (x=0; x < (ssize_t) image->columns; x++)
764 if ((channel & RedChannel) != 0)
766 distance=QuantumScale*(p->red-(MagickRealType) q->red);
767 channel_distortion[RedChannel]+=distance*distance;
768 channel_distortion[AllChannels]+=distance*distance;
770 if ((channel & GreenChannel) != 0)
772 distance=QuantumScale*(p->green-(MagickRealType) q->green);
773 channel_distortion[GreenChannel]+=distance*distance;
774 channel_distortion[AllChannels]+=distance*distance;
776 if ((channel & BlueChannel) != 0)
778 distance=QuantumScale*(p->blue-(MagickRealType) q->blue);
779 channel_distortion[BlueChannel]+=distance*distance;
780 channel_distortion[AllChannels]+=distance*distance;
782 if (((channel & OpacityChannel) != 0) &&
783 (image->matte != MagickFalse))
785 distance=QuantumScale*(p->opacity-(MagickRealType) q->opacity);
786 channel_distortion[OpacityChannel]+=distance*distance;
787 channel_distortion[AllChannels]+=distance*distance;
789 if (((channel & IndexChannel) != 0) &&
790 (image->colorspace == CMYKColorspace) &&
791 (reconstruct_image->colorspace == CMYKColorspace))
793 distance=QuantumScale*(indexes[x]-(MagickRealType)
794 reconstruct_indexes[x]);
795 channel_distortion[BlackChannel]+=distance*distance;
796 channel_distortion[AllChannels]+=distance*distance;
801 #if defined(MAGICKCORE_OPENMP_SUPPORT)
802 #pragma omp critical (MagickCore_GetMeanSquaredError)
804 for (i=0; i <= (ssize_t) AllChannels; i++)
805 distortion[i]+=channel_distortion[i];
807 reconstruct_view=DestroyCacheView(reconstruct_view);
808 image_view=DestroyCacheView(image_view);
809 for (i=0; i <= (ssize_t) AllChannels; i++)
810 distortion[i]/=((double) image->columns*image->rows);
811 distortion[AllChannels]/=(double) GetNumberChannels(image,channel);
815 static MagickBooleanType GetNormalizedCrossCorrelationError(const Image *image,
816 const Image *reconstruct_image,const ChannelType channel,
817 double *distortion,ExceptionInfo *exception)
819 #define SimilarityImageTag "Similarity/Image"
827 *correlation_statistics,
829 *reconstruct_statistics;
846 correlation_image=CloneImage(image,image->columns,image->rows,MagickTrue,
848 if (correlation_image == (Image *) NULL)
850 if (SetImageStorageClass(correlation_image,DirectClass) == MagickFalse)
852 InheritException(exception,&correlation_image->exception);
853 correlation_image=DestroyImage(correlation_image);
856 image_statistics=GetImageChannelStatistics(image,exception);
857 reconstruct_statistics=GetImageChannelStatistics(reconstruct_image,exception);
860 image_view=AcquireCacheView(image);
861 reconstruct_view=AcquireCacheView(reconstruct_image);
862 correlation_view=AcquireCacheView(correlation_image);
863 #if defined(MAGICKCORE_OPENMP_SUPPORT)
864 #pragma omp parallel for schedule(dynamic,4) shared(status,progress)
866 for (y=0; y < (ssize_t) image->rows; y++)
868 register const IndexPacket
870 *restrict reconstruct_indexes;
872 register const PixelPacket
877 *restrict correlation_indexes;
885 if (status == MagickFalse)
887 p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
888 q=GetCacheViewVirtualPixels(reconstruct_view,0,y,
889 reconstruct_image->columns,1,exception);
890 c=GetCacheViewAuthenticPixels(correlation_view,0,y,
891 correlation_image->columns,1,exception);
892 if ((p == (const PixelPacket *) NULL) ||
893 (q == (const PixelPacket *) NULL) || (c == (PixelPacket *) NULL))
898 indexes=GetCacheViewVirtualIndexQueue(image_view);
899 reconstruct_indexes=GetCacheViewVirtualIndexQueue(reconstruct_view);
900 correlation_indexes=GetCacheViewAuthenticIndexQueue(correlation_view);
901 for (x=0; x < (ssize_t) image->columns; x++)
903 if ((channel & RedChannel) != 0)
904 c->red=(p->red-image_statistics[RedChannel].mean)*
905 (q->red-reconstruct_statistics[RedChannel].mean);
906 if ((channel & GreenChannel) != 0)
907 c->green=(p->green-image_statistics[GreenChannel].mean)*
908 (q->green-reconstruct_statistics[GreenChannel].mean);
909 if ((channel & BlueChannel) != 0)
910 c->blue=(p->blue-image_statistics[BlueChannel].mean)*
911 (q->blue-reconstruct_statistics[BlueChannel].mean);
912 if (((channel & OpacityChannel) != 0) &&
913 (image->matte != MagickFalse))
914 c->opacity=(p->opacity-image_statistics[OpacityChannel].mean)*
915 (q->opacity-reconstruct_statistics[OpacityChannel].mean);
916 if (((channel & IndexChannel) != 0) &&
917 (image->colorspace == CMYKColorspace) &&
918 (reconstruct_image->colorspace == CMYKColorspace))
919 correlation_indexes[x]=(indexes[x]-
920 image_statistics[OpacityChannel].mean)*(reconstruct_indexes[x]-
921 reconstruct_statistics[OpacityChannel].mean);
926 if (SyncCacheViewAuthenticPixels(correlation_view,exception) == MagickFalse)
928 if (image->progress_monitor != (MagickProgressMonitor) NULL)
933 #if defined(MAGICKCORE_OPENMP_SUPPORT)
934 #pragma omp critical (MagickCore_GetNormalizedCrossCorrelationError)
936 proceed=SetImageProgress(image,SimilarityImageTag,progress++,
938 if (proceed == MagickFalse)
942 correlation_view=DestroyCacheView(correlation_view);
943 reconstruct_view=DestroyCacheView(reconstruct_view);
944 image_view=DestroyCacheView(image_view);
945 correlation_statistics=GetImageChannelStatistics(correlation_image,exception);
946 correlation_image=DestroyImage(correlation_image);
947 for (i=0; i < (ssize_t) AllChannels; i++)
948 distortion[i]=correlation_statistics[i].mean/(
949 image_statistics[i].standard_deviation*
950 reconstruct_statistics[i].standard_deviation);
951 distortion[AllChannels]=0.0;
952 if ((channel & RedChannel) != 0)
953 distortion[AllChannels]+=distortion[RedChannel]*distortion[RedChannel];
954 if ((channel & GreenChannel) != 0)
955 distortion[AllChannels]+=distortion[GreenChannel]*distortion[GreenChannel];
956 if ((channel & BlueChannel) != 0)
957 distortion[AllChannels]+=distortion[BlueChannel]*distortion[BlueChannel];
958 if (((channel & OpacityChannel) != 0) && (image->matte != MagickFalse))
959 distortion[AllChannels]+=distortion[OpacityChannel]*
960 distortion[OpacityChannel];
961 if (((channel & IndexChannel) != 0) &&
962 (image->colorspace == CMYKColorspace))
963 distortion[AllChannels]+=distortion[BlackChannel]*distortion[BlackChannel];
964 distortion[AllChannels]=sqrt(distortion[AllChannels])/
965 GetNumberChannels(image,AllChannels);
966 reconstruct_statistics=(ChannelStatistics *) RelinquishMagickMemory(
967 reconstruct_statistics);
968 image_statistics=(ChannelStatistics *) RelinquishMagickMemory(
970 correlation_statistics=(ChannelStatistics *) RelinquishMagickMemory(
971 correlation_statistics);
975 static MagickBooleanType GetPeakAbsoluteError(const Image *image,
976 const Image *reconstruct_image,const ChannelType channel,
977 double *distortion,ExceptionInfo *exception)
990 image_view=AcquireCacheView(image);
991 reconstruct_view=AcquireCacheView(reconstruct_image);
992 #if defined(MAGICKCORE_OPENMP_SUPPORT)
993 #pragma omp parallel for schedule(dynamic,4) shared(status)
995 for (y=0; y < (ssize_t) image->rows; y++)
998 channel_distortion[AllChannels+1];
1000 register const IndexPacket
1002 *restrict reconstruct_indexes;
1004 register const PixelPacket
1012 if (status == MagickFalse)
1014 p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
1015 q=GetCacheViewVirtualPixels(reconstruct_view,0,y,
1016 reconstruct_image->columns,1,exception);
1017 if ((p == (const PixelPacket *) NULL) || (q == (const PixelPacket *) NULL))
1022 indexes=GetCacheViewVirtualIndexQueue(image_view);
1023 reconstruct_indexes=GetCacheViewVirtualIndexQueue(reconstruct_view);
1024 (void) ResetMagickMemory(channel_distortion,0,sizeof(channel_distortion));
1025 for (x=0; x < (ssize_t) image->columns; x++)
1030 if ((channel & RedChannel) != 0)
1032 distance=QuantumScale*fabs(p->red-(double) q->red);
1033 if (distance > channel_distortion[RedChannel])
1034 channel_distortion[RedChannel]=distance;
1035 if (distance > channel_distortion[AllChannels])
1036 channel_distortion[AllChannels]=distance;
1038 if ((channel & GreenChannel) != 0)
1040 distance=QuantumScale*fabs(p->green-(double) q->green);
1041 if (distance > channel_distortion[GreenChannel])
1042 channel_distortion[GreenChannel]=distance;
1043 if (distance > channel_distortion[AllChannels])
1044 channel_distortion[AllChannels]=distance;
1046 if ((channel & BlueChannel) != 0)
1048 distance=QuantumScale*fabs(p->blue-(double) q->blue);
1049 if (distance > channel_distortion[BlueChannel])
1050 channel_distortion[BlueChannel]=distance;
1051 if (distance > channel_distortion[AllChannels])
1052 channel_distortion[AllChannels]=distance;
1054 if (((channel & OpacityChannel) != 0) &&
1055 (image->matte != MagickFalse))
1057 distance=QuantumScale*fabs(p->opacity-(double) q->opacity);
1058 if (distance > channel_distortion[OpacityChannel])
1059 channel_distortion[OpacityChannel]=distance;
1060 if (distance > channel_distortion[AllChannels])
1061 channel_distortion[AllChannels]=distance;
1063 if (((channel & IndexChannel) != 0) &&
1064 (image->colorspace == CMYKColorspace) &&
1065 (reconstruct_image->colorspace == CMYKColorspace))
1067 distance=QuantumScale*fabs(indexes[x]-(double)
1068 reconstruct_indexes[x]);
1069 if (distance > channel_distortion[BlackChannel])
1070 channel_distortion[BlackChannel]=distance;
1071 if (distance > channel_distortion[AllChannels])
1072 channel_distortion[AllChannels]=distance;
1077 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1078 #pragma omp critical (MagickCore_GetPeakAbsoluteError)
1080 for (i=0; i <= (ssize_t) AllChannels; i++)
1081 if (channel_distortion[i] > distortion[i])
1082 distortion[i]=channel_distortion[i];
1084 reconstruct_view=DestroyCacheView(reconstruct_view);
1085 image_view=DestroyCacheView(image_view);
1089 static MagickBooleanType GetPeakSignalToNoiseRatio(const Image *image,
1090 const Image *reconstruct_image,const ChannelType channel,
1091 double *distortion,ExceptionInfo *exception)
1096 status=GetMeanSquaredError(image,reconstruct_image,channel,distortion,
1098 if ((channel & RedChannel) != 0)
1099 distortion[RedChannel]=20.0*log10((double) 1.0/sqrt(
1100 distortion[RedChannel]));
1101 if ((channel & GreenChannel) != 0)
1102 distortion[GreenChannel]=20.0*log10((double) 1.0/sqrt(
1103 distortion[GreenChannel]));
1104 if ((channel & BlueChannel) != 0)
1105 distortion[BlueChannel]=20.0*log10((double) 1.0/sqrt(
1106 distortion[BlueChannel]));
1107 if (((channel & OpacityChannel) != 0) &&
1108 (image->matte != MagickFalse))
1109 distortion[OpacityChannel]=20.0*log10((double) 1.0/sqrt(
1110 distortion[OpacityChannel]));
1111 if (((channel & IndexChannel) != 0) &&
1112 (image->colorspace == CMYKColorspace))
1113 distortion[BlackChannel]=20.0*log10((double) 1.0/sqrt(
1114 distortion[BlackChannel]));
1115 distortion[AllChannels]=20.0*log10((double) 1.0/sqrt(
1116 distortion[AllChannels]));
1120 static MagickBooleanType GetRootMeanSquaredError(const Image *image,
1121 const Image *reconstruct_image,const ChannelType channel,
1122 double *distortion,ExceptionInfo *exception)
1127 status=GetMeanSquaredError(image,reconstruct_image,channel,distortion,
1129 if ((channel & RedChannel) != 0)
1130 distortion[RedChannel]=sqrt(distortion[RedChannel]);
1131 if ((channel & GreenChannel) != 0)
1132 distortion[GreenChannel]=sqrt(distortion[GreenChannel]);
1133 if ((channel & BlueChannel) != 0)
1134 distortion[BlueChannel]=sqrt(distortion[BlueChannel]);
1135 if (((channel & OpacityChannel) != 0) &&
1136 (image->matte != MagickFalse))
1137 distortion[OpacityChannel]=sqrt(distortion[OpacityChannel]);
1138 if (((channel & IndexChannel) != 0) &&
1139 (image->colorspace == CMYKColorspace))
1140 distortion[BlackChannel]=sqrt(distortion[BlackChannel]);
1141 distortion[AllChannels]=sqrt(distortion[AllChannels]);
1145 MagickExport MagickBooleanType GetImageChannelDistortion(Image *image,
1146 const Image *reconstruct_image,const ChannelType channel,
1147 const MetricType metric,double *distortion,ExceptionInfo *exception)
1150 *channel_distortion;
1158 assert(image != (Image *) NULL);
1159 assert(image->signature == MagickSignature);
1160 if (image->debug != MagickFalse)
1161 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1162 assert(reconstruct_image != (const Image *) NULL);
1163 assert(reconstruct_image->signature == MagickSignature);
1164 assert(distortion != (double *) NULL);
1166 if (image->debug != MagickFalse)
1167 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1168 if ((reconstruct_image->columns != image->columns) ||
1169 (reconstruct_image->rows != image->rows))
1170 ThrowBinaryException(ImageError,"ImageSizeDiffers",image->filename);
1172 Get image distortion.
1174 length=AllChannels+1UL;
1175 channel_distortion=(double *) AcquireQuantumMemory(length,
1176 sizeof(*channel_distortion));
1177 if (channel_distortion == (double *) NULL)
1178 ThrowFatalException(ResourceLimitFatalError,"MemoryAllocationFailed");
1179 (void) ResetMagickMemory(channel_distortion,0,length*
1180 sizeof(*channel_distortion));
1183 case AbsoluteErrorMetric:
1185 status=GetAbsoluteError(image,reconstruct_image,channel,
1186 channel_distortion,exception);
1189 case MeanAbsoluteErrorMetric:
1191 status=GetMeanAbsoluteError(image,reconstruct_image,channel,
1192 channel_distortion,exception);
1195 case MeanErrorPerPixelMetric:
1197 status=GetMeanErrorPerPixel(image,reconstruct_image,channel,
1198 channel_distortion,exception);
1201 case MeanSquaredErrorMetric:
1203 status=GetMeanSquaredError(image,reconstruct_image,channel,
1204 channel_distortion,exception);
1207 case NormalizedCrossCorrelationErrorMetric:
1209 status=GetNormalizedCrossCorrelationError(image,reconstruct_image,channel,
1210 channel_distortion,exception);
1213 case PeakAbsoluteErrorMetric:
1216 status=GetPeakAbsoluteError(image,reconstruct_image,channel,
1217 channel_distortion,exception);
1220 case PeakSignalToNoiseRatioMetric:
1222 status=GetPeakSignalToNoiseRatio(image,reconstruct_image,channel,
1223 channel_distortion,exception);
1226 case RootMeanSquaredErrorMetric:
1228 status=GetRootMeanSquaredError(image,reconstruct_image,channel,
1229 channel_distortion,exception);
1233 *distortion=channel_distortion[AllChannels];
1234 channel_distortion=(double *) RelinquishMagickMemory(channel_distortion);
1239 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1243 % 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 %
1247 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1249 % GetImageChannelDistrortion() compares the image channels of an image to a
1250 % reconstructed image and returns the specified distortion metric for each
1253 % The format of the CompareImageChannels method is:
1255 % double *GetImageChannelDistortions(const Image *image,
1256 % const Image *reconstruct_image,const MetricType metric,
1257 % ExceptionInfo *exception)
1259 % A description of each parameter follows:
1261 % o image: the image.
1263 % o reconstruct_image: the reconstruct image.
1265 % o metric: the metric.
1267 % o exception: return any errors or warnings in this structure.
1270 MagickExport double *GetImageChannelDistortions(Image *image,
1271 const Image *reconstruct_image,const MetricType metric,
1272 ExceptionInfo *exception)
1275 *channel_distortion;
1283 assert(image != (Image *) NULL);
1284 assert(image->signature == MagickSignature);
1285 if (image->debug != MagickFalse)
1286 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1287 assert(reconstruct_image != (const Image *) NULL);
1288 assert(reconstruct_image->signature == MagickSignature);
1289 if (image->debug != MagickFalse)
1290 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1291 if ((reconstruct_image->columns != image->columns) ||
1292 (reconstruct_image->rows != image->rows))
1294 (void) ThrowMagickException(&image->exception,GetMagickModule(),
1295 ImageError,"ImageSizeDiffers","`%s'",image->filename);
1296 return((double *) NULL);
1299 Get image distortion.
1301 length=AllChannels+1UL;
1302 channel_distortion=(double *) AcquireQuantumMemory(length,
1303 sizeof(*channel_distortion));
1304 if (channel_distortion == (double *) NULL)
1305 ThrowFatalException(ResourceLimitFatalError,"MemoryAllocationFailed");
1306 (void) ResetMagickMemory(channel_distortion,0,length*
1307 sizeof(*channel_distortion));
1310 case AbsoluteErrorMetric:
1312 status=GetAbsoluteError(image,reconstruct_image,AllChannels,
1313 channel_distortion,exception);
1316 case MeanAbsoluteErrorMetric:
1318 status=GetMeanAbsoluteError(image,reconstruct_image,AllChannels,
1319 channel_distortion,exception);
1322 case MeanErrorPerPixelMetric:
1324 status=GetMeanErrorPerPixel(image,reconstruct_image,AllChannels,
1325 channel_distortion,exception);
1328 case MeanSquaredErrorMetric:
1330 status=GetMeanSquaredError(image,reconstruct_image,AllChannels,
1331 channel_distortion,exception);
1334 case NormalizedCrossCorrelationErrorMetric:
1336 status=GetNormalizedCrossCorrelationError(image,reconstruct_image,
1337 AllChannels,channel_distortion,exception);
1340 case PeakAbsoluteErrorMetric:
1343 status=GetPeakAbsoluteError(image,reconstruct_image,AllChannels,
1344 channel_distortion,exception);
1347 case PeakSignalToNoiseRatioMetric:
1349 status=GetPeakSignalToNoiseRatio(image,reconstruct_image,AllChannels,
1350 channel_distortion,exception);
1353 case RootMeanSquaredErrorMetric:
1355 status=GetRootMeanSquaredError(image,reconstruct_image,AllChannels,
1356 channel_distortion,exception);
1360 return(channel_distortion);
1364 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1368 % I s I m a g e s E q u a l %
1372 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1374 % IsImagesEqual() measures the difference between colors at each pixel
1375 % location of two images. A value other than 0 means the colors match
1376 % exactly. Otherwise an error measure is computed by summing over all
1377 % pixels in an image the distance squared in RGB space between each image
1378 % pixel and its corresponding pixel in the reconstruct image. The error
1379 % measure is assigned to these image members:
1381 % o mean_error_per_pixel: The mean error for any single pixel in
1384 % o normalized_mean_error: The normalized mean quantization error for
1385 % any single pixel in the image. This distance measure is normalized to
1386 % a range between 0 and 1. It is independent of the range of red, green,
1387 % and blue values in the image.
1389 % o normalized_maximum_error: The normalized maximum quantization
1390 % error for any single pixel in the image. This distance measure is
1391 % normalized to a range between 0 and 1. It is independent of the range
1392 % of red, green, and blue values in your image.
1394 % A small normalized mean square error, accessed as
1395 % image->normalized_mean_error, suggests the images are very similar in
1396 % spatial layout and color.
1398 % The format of the IsImagesEqual method is:
1400 % MagickBooleanType IsImagesEqual(Image *image,
1401 % const Image *reconstruct_image)
1403 % A description of each parameter follows.
1405 % o image: the image.
1407 % o reconstruct_image: the reconstruct image.
1410 MagickExport MagickBooleanType IsImagesEqual(Image *image,
1411 const Image *reconstruct_image)
1430 mean_error_per_pixel;
1432 assert(image != (Image *) NULL);
1433 assert(image->signature == MagickSignature);
1434 assert(reconstruct_image != (const Image *) NULL);
1435 assert(reconstruct_image->signature == MagickSignature);
1436 if ((reconstruct_image->columns != image->columns) ||
1437 (reconstruct_image->rows != image->rows))
1438 ThrowBinaryException(ImageError,"ImageSizeDiffers",image->filename);
1441 mean_error_per_pixel=0.0;
1443 exception=(&image->exception);
1444 image_view=AcquireCacheView(image);
1445 reconstruct_view=AcquireCacheView(reconstruct_image);
1446 for (y=0; y < (ssize_t) image->rows; y++)
1448 register const IndexPacket
1450 *restrict reconstruct_indexes;
1452 register const PixelPacket
1459 p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
1460 q=GetCacheViewVirtualPixels(reconstruct_view,0,y,reconstruct_image->columns,
1462 if ((p == (const PixelPacket *) NULL) || (q == (const PixelPacket *) NULL))
1464 indexes=GetCacheViewVirtualIndexQueue(image_view);
1465 reconstruct_indexes=GetCacheViewVirtualIndexQueue(reconstruct_view);
1466 for (x=0; x < (ssize_t) image->columns; x++)
1471 distance=fabs(p->red-(double) q->red);
1472 mean_error_per_pixel+=distance;
1473 mean_error+=distance*distance;
1474 if (distance > maximum_error)
1475 maximum_error=distance;
1477 distance=fabs(p->green-(double) q->green);
1478 mean_error_per_pixel+=distance;
1479 mean_error+=distance*distance;
1480 if (distance > maximum_error)
1481 maximum_error=distance;
1483 distance=fabs(p->blue-(double) q->blue);
1484 mean_error_per_pixel+=distance;
1485 mean_error+=distance*distance;
1486 if (distance > maximum_error)
1487 maximum_error=distance;
1489 if (image->matte != MagickFalse)
1491 distance=fabs(p->opacity-(double) q->opacity);
1492 mean_error_per_pixel+=distance;
1493 mean_error+=distance*distance;
1494 if (distance > maximum_error)
1495 maximum_error=distance;
1498 if ((image->colorspace == CMYKColorspace) &&
1499 (reconstruct_image->colorspace == CMYKColorspace))
1501 distance=fabs(indexes[x]-(double) reconstruct_indexes[x]);
1502 mean_error_per_pixel+=distance;
1503 mean_error+=distance*distance;
1504 if (distance > maximum_error)
1505 maximum_error=distance;
1512 reconstruct_view=DestroyCacheView(reconstruct_view);
1513 image_view=DestroyCacheView(image_view);
1514 image->error.mean_error_per_pixel=(double) (mean_error_per_pixel/area);
1515 image->error.normalized_mean_error=(double) (QuantumScale*QuantumScale*
1517 image->error.normalized_maximum_error=(double) (QuantumScale*maximum_error);
1518 status=image->error.mean_error_per_pixel == 0.0 ? MagickTrue : MagickFalse;
1523 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1527 % S i m i l a r i t y I m a g e %
1531 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1533 % SimilarityImage() compares the reference image of the image and returns the
1534 % best match offset. In addition, it returns a similarity image such that an
1535 % exact match location is completely white and if none of the pixels match,
1536 % black, otherwise some gray level in-between.
1538 % The format of the SimilarityImageImage method is:
1540 % Image *SimilarityImage(const Image *image,const Image *reference,
1541 % RectangleInfo *offset,double *similarity,ExceptionInfo *exception)
1543 % A description of each parameter follows:
1545 % o image: the image.
1547 % o reference: find an area of the image that closely resembles this image.
1549 % o the best match offset of the reference image within the image.
1551 % o similarity: the computed similarity between the images.
1553 % o exception: return any errors or warnings in this structure.
1557 static double GetSimilarityMetric(const Image *image,const Image *reference,
1558 const ssize_t x_offset,const ssize_t y_offset,ExceptionInfo *exception)
1574 Compute the similarity in pixels between two images.
1578 image_view=AcquireCacheView(image);
1579 reference_view=AcquireCacheView(reference);
1580 for (y=0; y < (ssize_t) reference->rows; y++)
1582 register const IndexPacket
1584 *restrict reference_indexes;
1586 register const PixelPacket
1593 if (status == MagickFalse)
1595 p=GetCacheViewVirtualPixels(image_view,x_offset,y_offset+y,
1596 reference->columns,1,exception);
1597 q=GetCacheViewVirtualPixels(reference_view,0,y,reference->columns,1,
1599 if ((p == (const PixelPacket *) NULL) || (q == (const PixelPacket *) NULL))
1604 indexes=GetCacheViewVirtualIndexQueue(image_view);
1605 reference_indexes=GetCacheViewVirtualIndexQueue(reference_view);
1606 for (x=0; x < (ssize_t) reference->columns; x++)
1614 thread_similarity=0.0;
1615 distance=QuantumScale*(p->red-(MagickRealType) q->red);
1616 thread_similarity+=distance*distance;
1617 distance=QuantumScale*(p->green-(MagickRealType) q->green);
1618 thread_similarity+=distance*distance;
1619 distance=QuantumScale*(p->blue-(MagickRealType) q->blue);
1620 thread_similarity+=distance*distance;
1621 if ((image->matte != MagickFalse) && (reference->matte != MagickFalse))
1623 distance=QuantumScale*(p->opacity-(MagickRealType) q->opacity);
1624 thread_similarity+=distance*distance;
1626 if ((image->colorspace == CMYKColorspace) &&
1627 (reference->colorspace == CMYKColorspace))
1629 distance=QuantumScale*(indexes[x]-(MagickRealType)
1630 reference_indexes[x]);
1631 thread_similarity+=distance*distance;
1633 similarity+=thread_similarity;
1638 reference_view=DestroyCacheView(reference_view);
1639 image_view=DestroyCacheView(image_view);
1640 if (status == MagickFalse)
1642 similarity/=((double) reference->columns*reference->rows);
1643 similarity/=(double) GetNumberChannels(reference,AllChannels);
1644 return(sqrt(similarity));
1647 MagickExport Image *SimilarityImage(Image *image,const Image *reference,
1648 RectangleInfo *offset,double *similarity_metric,ExceptionInfo *exception)
1650 #define SimilarityImageTag "Similarity/Image"
1667 assert(image != (const Image *) NULL);
1668 assert(image->signature == MagickSignature);
1669 if (image->debug != MagickFalse)
1670 (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1671 assert(exception != (ExceptionInfo *) NULL);
1672 assert(exception->signature == MagickSignature);
1673 assert(offset != (RectangleInfo *) NULL);
1674 SetGeometry(reference,offset);
1675 *similarity_metric=1.0;
1676 if ((reference->columns > image->columns) || (reference->rows > image->rows))
1677 ThrowImageException(ImageError,"ImageSizeDiffers");
1678 similarity_image=CloneImage(image,image->columns-reference->columns+1,
1679 image->rows-reference->rows+1,MagickTrue,exception);
1680 if (similarity_image == (Image *) NULL)
1681 return((Image *) NULL);
1682 if (SetImageStorageClass(similarity_image,DirectClass) == MagickFalse)
1684 InheritException(exception,&similarity_image->exception);
1685 similarity_image=DestroyImage(similarity_image);
1686 return((Image *) NULL);
1689 Measure similarity of reference image against image.
1693 similarity_view=AcquireCacheView(similarity_image);
1694 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1695 #pragma omp parallel for schedule(dynamic,4) shared(progress,status)
1697 for (y=0; y < (ssize_t) (image->rows-reference->rows+1); y++)
1705 register PixelPacket
1708 if (status == MagickFalse)
1710 q=GetCacheViewAuthenticPixels(similarity_view,0,y,
1711 similarity_image->columns,1,exception);
1712 if (q == (const PixelPacket *) NULL)
1717 for (x=0; x < (ssize_t) (image->columns-reference->columns+1); x++)
1719 similarity=GetSimilarityMetric(image,reference,x,y,exception);
1720 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1721 #pragma omp critical (MagickCore_SimilarityImage)
1723 if (similarity < *similarity_metric)
1725 *similarity_metric=similarity;
1729 q->red=ClampToQuantum(QuantumRange-QuantumRange*similarity);
1734 if (SyncCacheViewAuthenticPixels(similarity_view,exception) == MagickFalse)
1736 if (image->progress_monitor != (MagickProgressMonitor) NULL)
1741 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1742 #pragma omp critical (MagickCore_SimilarityImage)
1744 proceed=SetImageProgress(image,SimilarityImageTag,progress++,
1746 if (proceed == MagickFalse)
1750 similarity_view=DestroyCacheView(similarity_view);
1751 return(similarity_image);