~ [ source navigation ] ~ [ diff markup ] ~ [ identifier search ] ~ [ freetext search ] ~ [ file search ] ~

Linux Cross Reference
Tina5/tina-libs/tina/image/imgPrc_compare.c

Version: ~
Architecture: ~ [ i386 ] ~ [ alpha ] ~ [ m68k ] ~ [ mips ] ~ [ ppc ] ~ [ sparc ] ~ [ sparc64 ] ~

  1 /**********
  2  *
  3  * Copyright (c) 2003, Division of Imaging Science and Biomedical Engineering,
  4  * University of Manchester, UK.  All rights reserved.
  5  *
  6  * Redistribution and use in source and binary forms, with or without modification,
  7  * are permitted provided that the following conditions are met:
  8  *
  9  *   . Redistributions of source code must retain the above copyright notice,
 10  *     this list of conditions and the following disclaimer.
 11  *
 12  *   . Redistributions in binary form must reproduce the above copyright notice,
 13  *     this list of conditions and the following disclaimer in the documentation
 14  *     and/or other materials provided with the distribution.
 15  *
 16  *   . Neither the name of the University of Manchester nor the names of its
 17  *     contributors may be used to endorse or promote products derived from this
 18  *     software without specific prior written permission.
 19  *
 20  *
 21  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
 22  * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
 23  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
 24  * ARE DISCLAIMED.  IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE
 25  * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
 26  * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
 27  * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
 28  * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
 29  * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
 30  * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
 31  * POSSIBILITY OF SUCH DAMAGE.
 32  *
 33  **********
 34  * 
 35  * Program :    TINA
 36  * File    :  $Source: /home/tina/cvs/tina-libs/tina/image/imgPrc_compare.c,v $
 37  * Date    :  $Date: 2007/02/15 01:52:29 $
 38  * Version :  $Revision: 1.5 $
 39  * CVS Id  :  $Id: imgPrc_compare.c,v 1.5 2007/02/15 01:52:29 paul Exp $
 40  *
 41  * Author  : NAT
 42  *
 43 */
 44 /**
 45  * @file  rountines to support regression testing im_compare() , im_transpose(). 
 46  * @brief Rows and columns transposed for char, integer, float and complex image types. 
 47  *        Compare returns number of voxels failing specified tolerance, defaults to zero
 48  *        if image types not supported for comparison.          
 49  *
 50  *
 51  *********
 52 */
 53 
 54 #include "imgPrc_compare.h"
 55 
 56 #if HAVE_CONFIG_H
 57   #include <config.h>
 58 #endif
 59 
 60 #include <math.h>
 61 #include <tina/sys/sysDef.h>
 62 #include <tina/sys/sysPro.h>
 63 #include <tina/math/mathDef.h>
 64 #include <tina/math/mathPro.h>
 65 #include <tina/image/img_GenDef.h>
 66 #include <tina/image/img_GenPro.h>
 67 #include <tina/image/imgPrc_suptype.h>
 68 #include <tina/image/img_PrcPro.h>
 69 
 70 int imi_compare(int tol, Imrect * im1, Imrect * im2)
 71 {
 72     Imregion *roi;
 73     int    *row1, *row2;
 74     int     diff;
 75     int     lx, ux, ly, uy;
 76     int     i, j, sum = 0;
 77 
 78     if (im1 == NULL || im2 == NULL)
 79         return (0);
 80 
 81     roi = roi_inter(im1->region, im2->region);
 82     if (roi == NULL)
 83         return (0);
 84     lx = roi->lx;
 85     ux = roi->ux;
 86     ly = roi->ly;
 87     uy = roi->uy;
 88 
 89     row1 = ivector_alloc(lx, ux);
 90     row2 = ivector_alloc(lx, ux);
 91 
 92     for (i = ly; i < uy; ++i)
 93     {
 94         im_get_row(row1, im1, i, lx, ux);
 95         im_get_row(row2, im2, i, lx, ux);
 96         for (j = lx; j < ux; ++j)
 97         {
 98             diff = row1[j] - row2[j];
 99             if (fabs(diff) > tol) sum++;
100         }
101     }
102 
103     ivector_free(row1, lx);
104     ivector_free(row2, lx);
105     rfree((void *) roi);
106     return (sum);
107 }
108 
109 int imf_compare(double tol, Imrect * im1, Imrect * im2)
110 {
111     Imregion *roi;
112     float  *row1, *row2;
113     float   diff;
114     int     lx, ux, ly, uy;
115     int     i, j, sum = 0;
116 
117     if (im1 == NULL || im2 == NULL)
118         return (0);
119 
120     roi = roi_inter(im1->region, im2->region);
121     if (roi == NULL)
122         return (0);
123     lx = roi->lx;
124     ux = roi->ux;
125     ly = roi->ly;
126     uy = roi->uy;
127 
128     row1 = fvector_alloc(lx, ux);
129     row2 = fvector_alloc(lx, ux);
130 
131     for (i = ly; i < uy; ++i)
132     {
133         im_get_rowf(row1, im1, i, lx, ux);
134         im_get_rowf(row2, im2, i, lx, ux);
135         for (j = lx; j < ux; ++j)
136         {
137             diff = row1[j] - row2[j];
138             if (fabs(diff) > tol) sum++;
139         }
140     }
141 
142     fvector_free(row1, lx);
143     fvector_free(row2, lx);
144     rfree((void *) roi);
145     return (sum);
146 }
147 
148 int imz_compare(double tol, Imrect * im1, Imrect * im2)
149 {
150     Imregion *roi;
151     Complex *row1;
152     Complex *row2;
153     Complex diff;
154     int     lx, ux, ly, uy;
155     int     i, j, sum = 0;
156 
157     if (im1 == NULL || im2 == NULL)
158         return (0);
159 
160     roi = roi_inter(im1->region, im2->region);
161     if (roi == NULL)
162         return (0);
163     lx = roi->lx;
164     ux = roi->ux;
165     ly = roi->ly;
166     uy = roi->uy;
167 
168     row1 = zvector_alloc(lx, ux);
169     row2 = zvector_alloc(lx, ux);
170 
171     for (i = ly; i < uy; ++i)
172     {
173         im_get_rowz(row1, im1, i, lx, ux);
174         im_get_rowz(row2, im2, i, lx, ux);
175         for (j = lx; j < ux; ++j)
176         {
177             diff = cmplx_diff(row1[j], row2[j]);
178             if (cmplx_mod(diff)>tol) sum++;
179         }
180     }
181 
182     zvector_free(row1, lx);
183     zvector_free(row2, lx);
184     rfree((void *) roi);
185     return (sum);
186 }
187 
188 int im_compare(double tol, Imrect * im1, Imrect * im2)
189 {
190     if (im1 == NULL || im2 == NULL)
191         return (0);
192     switch (im_sup_vtype(im1->vtype, im2->vtype))
193     {
194     case uchar_v:
195     case char_v:
196     case ushort_v:
197     case short_v:
198     case int_v:
199         return (imi_compare((int)tol, im1, im2));
200     case float_v:
201         return (imf_compare(tol, im1, im2));
202     case complex_v:
203         return (imz_compare(tol, im1, im2));
204     default:
205         return (0);
206     }
207 }
208 
209 Imrect         *imc_transpose(Imrect * im)
210 {
211         Imrect         *im2;
212         Imregion       *roi, roi2;
213         int            *row1, *row2;
214         int             lx, ux, ly, uy;
215         int             width, height;
216         int             i, j;
217 
218         if (im == NULL)
219                 return (NULL);
220 
221         roi = im->region;
222         if (roi == NULL)
223                 return (NULL);
224         lx = roi->lx;
225         ux = roi->ux;
226         ly = roi->ly;
227         uy = roi->uy;
228         height = im->width;
229         width = im->height;
230 
231         roi2.lx = ly;
232         roi2.ly = lx;
233         roi2.ux = uy;
234         roi2.uy = ux;
235 
236         if (im->vtype == uchar_v)
237            im2 = im_alloc(height, width, &roi2, uchar_v);
238         else
239            im2 = im_alloc(height, width, &roi2, char_v);
240 
241         row1 = ivector_alloc(lx, ux);
242         row2 = ivector_alloc(lx, ux);
243 
244         for (i = ly; i < uy; ++i)
245         {
246                 im_get_row(row1, im, i, lx, ux);
247                 for (j = lx; j < ux; ++j)
248                 {
249                          row2[j] = row1[j];
250                 }
251                 im_put_col(row2, im2, i, lx, ux);
252         }
253 
254         ivector_free(row1, lx);
255         ivector_free(row2, lx);
256         return (im2);
257 }
258 
259 Imrect         *imi_transpose(Imrect * im)
260 {
261         Imrect         *im2;
262         Imregion       *roi, roi2;
263         int            *row1, *row2;
264         int             lx, ux, ly, uy;
265         int             width, height;
266         int             i, j;
267 
268         if (im == NULL)
269                 return (NULL);
270 
271         roi = im->region;
272         if (roi == NULL)
273                 return (NULL);
274         lx = roi->lx;
275         ux = roi->ux;
276         ly = roi->ly;
277         uy = roi->uy;
278         height = im->width;
279         width = im->height;
280 
281         roi2.lx = ly;
282         roi2.ly = lx;
283         roi2.ux = uy;
284         roi2.uy = ux;
285 
286         im2 = im_alloc(height, width, &roi2, int_v);
287         row1 = ivector_alloc(lx, ux);
288         row2 = ivector_alloc(lx, ux);
289 
290         for (i = ly; i < uy; ++i)
291         {
292                 im_get_row(row1, im, i, lx, ux);
293                 for (j = lx; j < ux; ++j) 
294                         row2[j] = row1[j];
295                 im_put_col(row2, im2, i, lx, ux);
296         }
297 
298         ivector_free(row1, lx);
299         ivector_free(row2, lx);
300         return (im2);
301 }
302 
303 Imrect         *imf_transpose(Imrect * im)
304 {
305         Imrect         *im2;
306         Imregion       *roi, roi2;
307         float          *row1, *row2;
308         int             lx, ux, ly, uy;
309         int             width, height;
310         int             i, j;
311 
312         if (im == NULL)
313                 return (NULL);
314 
315         roi = im->region;
316         if (roi == NULL)
317                 return (NULL);
318         lx = roi->lx;
319         ux = roi->ux;
320         ly = roi->ly;
321         uy = roi->uy;
322 
323         height = im->width;
324         width = im->height;
325 
326         roi2.lx = ly;
327         roi2.ly = lx;
328         roi2.ux = uy;
329         roi2.uy = ux;
330 
331         im2 = im_alloc(height, width, &roi2, float_v);
332         row1 = fvector_alloc(lx, ux);
333         row2 = fvector_alloc(lx, ux);
334 
335         for (i = ly; i < uy; ++i)
336         {
337                 im_get_rowf(row1, im, i, lx, ux);
338                 for (j = lx; j < ux; ++j) 
339                         row2[j] = row1[j];
340                 im_put_colf(row2, im2, i, lx, ux);
341         }
342 
343         fvector_free(row1, lx);
344         fvector_free(row2, lx);
345         return (im2);
346 }
347 
348 Imrect         *imz_transpose(Imrect * im)
349 {
350         Imrect         *im2;
351         Imregion       *roi, roi2;
352         Complex        *row1, *row2;
353         int             lx, ux, ly, uy;
354         int             width, height;
355         int             i, j;
356 
357         if (im == NULL)
358                 return (NULL);
359 
360         roi = im->region;
361         if (roi == NULL)
362                 return (NULL);
363         lx = roi->lx;
364         ux = roi->ux;
365         ly = roi->ly;
366         uy = roi->uy;
367 
368         height = im->width;
369         width = im->height;
370 
371         roi2.lx = ly;
372         roi2.ly = lx;
373         roi2.ux = uy;
374         roi2.uy = ux;
375 
376         im2 = im_alloc(height, width, &roi2, complex_v);
377 
378         row1 = zvector_alloc(lx, ux);
379         row2 = zvector_alloc(lx, ux);
380 
381         for (i = ly; i < uy; ++i)
382         {
383                 im_get_rowz(row1, im, i, lx, ux);
384                 for (j = lx; j < ux; ++j)
385                 {
386                         row2[j].x = row1[j].x;
387                         row2[j].y = row1[j].y;
388                 }
389                 im_put_colz(row2, im2, i, lx, ux);
390         }
391 
392         zvector_free(row1, lx);
393         zvector_free(row2, lx);
394         return (im2);
395 }
396 
397 Imrect         *im_transpose(Imrect *im)
398 {
399         if (im == NULL)
400                 return (NULL);
401         switch (im->vtype)
402         {
403         case uchar_v:
404         case char_v:
405                 return (imc_transpose(im));
406         case short_v:
407         case ushort_v:
408         case uint_v:
409         case int_v:
410                 return (imi_transpose(im));
411         case float_v:
412                 return (imf_transpose(im));
413         case complex_v: 
414                 return (imz_transpose(im));
415         default:
416                 return (NULL);
417         }
418 }
419 
420 /**
421 * @ brief Function to test image processing library build.
422 * @param Image types of uchar_v, char_v, ushort_v, short_v, uint_v, int_v and float_v are supported.
423 *
424 * Algebraic, trigonametric and fft functions properly tested with data being promoted
425 * to higher cast levels where appropriate. Complex image data processing is implicitly tested
426 * but complex representation is generally unsuitable for image prcessing.
427 * More complicated routines tested by less stringent transform tests which test if the library
428 * is broken rather than any functional definition.
429 */
430 int im_test(Vartype vtype)
431 {
432      Imrect *imtest1=NULL, *imtest2=NULL, *imsqr1=NULL, *imsqr2=NULL, *immod2=NULL;
433      Imrect *imlog1=NULL, *imlog2=NULL, *imfft1=NULL, *imdx1=NULL, *imtrans1=NULL;
434      Imrect *im2=NULL, *im3=NULL, *im4=NULL, *im5=NULL;
435      int itest, tflag=1;
436 
437      if (vtype == complex_v) 
438      {
439          format("Inappropriate representation for image processing routines. \n");
440          return(tflag);
441      }
442 /* random noise centred at zero and width 64 */
443      imtest1 = imf_norm_noise(256, 256, 1, 1, 0.0, 64.0);
444      imtest2 = imf_norm_noise(256, 256, 1, 1, 0.0, 64.0);
445 
446      itest = im_compare(0.01, imtest1,imtest2);
447      if (itest > 0) format ("random image generator passes test \n");
448      else 
449      {
450           format("random image generator failed cannot test further \n");
451           return(0);
452      }
453  
454 /*  set up requested data type */
455      if (vtype == uchar_v || vtype == ushort_v || vtype == uint_v)
456      {
457          im2 = im_mod(imtest1);
458          im_free(imtest1);
459          imtest1 = im2;
460          im2 = im_mod(imtest2);
461          im_free(imtest2);
462          imtest2 = im2;
463      }
464      im2 = im_cast(imtest1, vtype);
465      im_free(imtest1);
466      imtest1 = im2;
467 
468      im2 = im_cast(imtest2, vtype);
469      im_free(imtest2);
470      imtest2 = im2;
471 
472 /* begin testing in earnest, pay particular attention to first faliure */
473      im2 = im_copy(imtest1);
474      itest = im_compare(0.01, imtest1,im2);
475      im_free(im2);
476      if (itest == 0) format ("copy passes test \n");
477      else
478      {
479           format("copy failed cannot test further \n");
480           return(0);
481      }
482 
483      imsqr1 = im_prod(imtest1,imtest1);
484      im2 = im_sqr(imtest1);
485      itest = im_compare(0.01, imsqr1,im2);
486      im_free(im2);
487      if (itest == 0) format ("sqr and prod pass algebraic test \n");
488      else
489      {
490           format("sqr and prod failed algebraic test %d times \n");
491           tflag = 0;
492      }
493 
494      imsqr2 = im_sqr(imtest2);
495      im2 = im_sqrt(imsqr2);
496      immod2 = im_mod(imtest2);
497      itest = im_compare(0.01, im2, immod2);
498      im_free(im2);
499      if (itest == 0) format ("complex sqrt and mod pass algebraic test \n");
500      else
501      {
502           format("float sqr and sqrt failed algebraic test %d times \n", itest);
503           tflag = 0;
504      }
505 
506      im2 = imz_sqrt(imtest1);
507      im3 = im_prod(im2,im2);
508      itest = im_compare(0.01, im3, imtest1);
509      im_free(im2);
510      im_free(im3);
511      if (itest == 0) format ("complex sqr and sqrt pass algebraic test \n");
512      else
513      {
514           format("compex sqr and sqrt failed algebrainc test %d times \n", itest);
515           tflag = 0;
516      }
517 
518      imlog1 = im_log(imsqr1);
519      imlog2 = im_log(imsqr2);
520      im2 = im_sum(imlog1, imlog2);
521      im3 = im_prod(imsqr1,imsqr2);
522      im4 = im_log(im3); 
523 
524      itest = im_compare(0.01, im2, im4);
525      im_free(im2);
526      im_free(im3);
527      im_free(im4);
528      if (itest == 0) format ("log passes algebraic test \n");
529      else
530      {
531           format("log failed algebraic test %d . < 20 O.K. for integers \n", itest);
532           if (itest > 20 ) tflag = 0;
533      }
534 
535      im2 = im_times(0.5,imlog2);
536      im3 = im_exp(im2);
537      itest = im_compare(0.01, im3, immod2);
538      im_free(im2);
539      im_free(im3);
540      if (itest == 0) format ("times and exp pass algebraic test \n");
541      else
542      {
543           format("exp failed algebraic test %d times \n", itest);
544           tflag = 0;
545      }
546 
547      im_free(immod2);
548      im_free(imlog1);
549      im_free(imlog2);
550      im_free(imsqr1);
551      im_free(imsqr2);
552 
553      imfft1 = im_fft(imtest1, NULL);
554      im2 = im_fft_inverse(imfft1, NULL);
555 
556      itest = im_compare(0.01, imtest1, im2);
557      im_free(im2);
558      if (itest == 0) format ("fft and inverse passes algebraic test \n");
559      else
560      {
561           format("fft failed test %d \n", itest);
562           tflag = 0;
563      }
564 
565      imtrans1 = im_transpose(imtest1);
566      im2 = im_transpose(im2);
567      itest = im_compare(0.01, imtest1, im2);
568      im_free(im2);
569      if (itest == 0) format ("transpose passes test \n");
570      else
571      {
572           format("transpose failed test %d \n", itest);
573           tflag = 0;
574      }
575 
576      imdx1 = imf_diffx(imtest1);
577      im2 = im_copy(imtest1);
578      im_shift(im2, 0, 1);
579      im3 = im_copy(im2);
580      im_shift(im2, 0, -2);
581      im4 = imf_wsum(0.5, -0.5, im2, im3);
582      itest = im_compare(0.01, im4, imdx1);
583      im_free(im2);
584      im_free(im3);
585      im_free(im4);
586      im_free(im5);
587      if (itest == 0) format ("x derivative and wsum pass algebraic test \n");
588      else
589      {
590           format("x derivative failed algebraic test %d times \n", itest);
591           tflag = 0;
592      }
593 
594      im2 = im_transpose(imdx1);
595      im3 = imf_diffy(imtrans1);
596      itest = im_compare(0.01, im2, im3);
597      im_free(im2);
598      im_free(im3);
599      if (itest == 0) format ("y derivative passes transpose test \n");
600      else
601      {
602           format("y derivative failed transpose test %d times \n", itest);
603           tflag = 0;
604      }
605 
606 /* now less specific transpose tests  which test that an algorithm runs but not if correct */
607 
608      im2 = im_median(imtest1);
609      im3 = im_transpose(im2);
610      im4 = im_median(imtrans1);
611      itest = im_compare(0.01, im3, im4);
612      im_free(im2);
613      im_free(im3);
614      im_free(im4);
615      if (itest == 0) format ("median passes transpose test \n");
616      else
617      {
618           format("median failed transpose test %d times \n", itest);
619           tflag = 0;
620      }
621 
622      im2 = im_tsmooth(imtest1);
623      im3 = im_transpose(im2);
624      im4 = im_tsmooth(imtrans1);
625      itest = im_compare(0.01, im3, im4);
626      im_free(im2);
627      im_free(im3);
628      im_free(im4);
629      if (itest == 0) format ("tsmooth passes transpose test \n");
630      else
631      {
632           format("tsmooth failed transpose test %d times \n", itest);
633           tflag = 0;
634      }
635 
636      im2 = imf_lsf_smooth(imtest1, 5.0);
637      im3 = im_transpose(im2);
638      im4 = imf_lsf_smooth(imtrans1, 5.0);
639      itest = im_compare(0.01, im3, im4);
640      im_free(im2);
641      im_free(im3);
642      im_free(im4);
643      if (itest == 0) format ("lsf smooth passes transpose test \n");
644      else
645      {
646           format("lsf smooth failed transpose test %d times \n", itest);
647           tflag = 0;
648      }
649 
650      im2 = imf_gauss(imtest1, 5.0, 0.01);
651      im3 = im_transpose(im2);
652      im4 = imf_gauss(imtrans1, 5.0, 0.01);
653      itest = im_compare(0.01, im3, im4);
654      im_free(im2);
655      im_free(im3);
656      im_free(im4);
657      if (itest == 0) format ("gaussian smooth passes transpose test \n");
658      else
659      {
660           format("gaussian smooth failed transpose test %d times \n", itest);
661           tflag = 0;
662      }
663 
664      im2 = imf_laplacian(imtest1);
665      im3 = im_transpose(im2);
666      im4 = imf_laplacian(imtrans1);
667      itest = im_compare(0.01, im3, im4);
668      im_free(im2);
669      im_free(im3);
670      im_free(im4);
671      if (itest == 0) format ("laplacian passes transpose test \n");
672      else
673      {
674           format("laplacian failed transpose test %d times \n", itest);
675           tflag = 0;
676      }
677 
678 /*   the following routines have less numerical precision */
679      im2 = imf_ddn(imtest1);
680      im3 = im_transpose(im2);
681      im4 = imf_ddn(imtrans1);
682      itest = im_compare(1.0, im3, im4);
683      im_free(im2);
684      im_free(im3);
685      im_free(im4);
686      if (itest == 0) format ("2nd derivative ddn passes transpose test \n");
687      else
688      {
689           format("2nd derivative ddn failed transpose test %d times. < 10 O.K. for floats \n", itest);
690           if (itest > 10) tflag=0;
691      }
692 
693      im2 = imf_ddt(imtest1);
694      im3 = im_transpose(im2);
695      im4 = imf_ddt(imtrans1);
696      itest = im_compare(1.0, im3, im4);
697      im_free(im2);
698      im_free(im3);
699      im_free(im4);
700      if (itest == 0) format ("2nd derivative ddt passes transpose test \n");
701      else
702      {
703           format("2nd derivative ddt failed transpose test %d times. < 10 O.K. for floats \n", itest);
704           if (itest > 10) tflag=0;
705      }
706 
707      im_free(imtrans1);
708      im_free(imdx1);
709      im_free(imfft1);
710      im_free(imtest1);
711      im_free(imtest2);
712 
713      if (tflag) 
714          format("Library has PASSED tests \n");
715      else 
716          format("Library has FAILED tests, check library build \n");
717 
718 
719      return(tflag);
720 
721 
722 }
723 

~ [ source navigation ] ~ [ diff markup ] ~ [ identifier search ] ~ [ freetext search ] ~ [ file search ] ~

This page was automatically generated by the LXR engine.
Visit the LXR main site for more information.