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
This page was automatically generated by the
LXR engine.
Visit the LXR main site for more
information.