1 /**********
2 *
3 * This file is part of the TINA Open Source Image Analysis Environment
4 * henceforth known as TINA
5 *
6 * TINA is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU Lesser General Public License as
8 * published by the Free Software Foundation.
9 *
10 * TINA is distributed in the hope that it will be useful,
11 * but WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 * GNU Lesser General Public License for more details.
14 *
15 * You should have received a copy of the GNU Lesser General Public License
16 * along with TINA; if not, write to the Free Software Foundation, Inc.,
17 * 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
18 *
19 **********
20 *
21 * Program : TINA
22 * File : $Source: /home/tina/cvs/tina-libs/tina/vision/visPgh_histfuncs.c,v $
23 * Date : $Date: 2003/11/12 15:44:37 $
24 * Version : $Revision: 1.4 $
25 * CVS Id : $Id: visPgh_histfuncs.c,v 1.4 2003/11/12 15:44:37 neil Exp $
26 *
27 * Author : NAT and A.Ashbrook
28 */
29 /**
30 * @file Wrapper routines for building and comparing Pairwise Geometric Histogrammes.
31 * @brief Includes normalisation and comparison routines for fixed size histogrammes.
32 */
33 #include "visPgh_histfuncs.h"
34
35 #include <stdio.h>
36 #include <math.h>
37 #include <tina/sys/sysDef.h>
38 #include <tina/sys/sysPro.h>
39 #include <tina/math/mathDef.h>
40 #include <tina/math/mathPro.h>
41 #include <tina/image/img_GenPro.h>
42 #include <tina/geometry/geom_LinePro.h>
43 #include <tina/vision/visPgh_hist.h>
44
45
46 /* ---------- Functions ---------- */
47
48 static void sqr_root_pairwise(Imrect *hist)
49 {
50 int i,j, min_x, min_y, max_x, max_y;
51 double hist_data;
52
53 min_x = hist->region->lx;
54 min_y = hist->region->ly;
55 max_x = hist->region->ux;
56 max_y = hist->region->uy;
57
58 /* Each histogram entry is replaced by its square root */
59
60 for(i=min_y;i<max_y;i++)
61 for(j=min_x;j<max_x;j++)
62 {
63 IM_PIX_GET(hist,i,j,hist_data);
64 if (hist_data < 0.0) hist_data = 0.0;
65 IM_PIX_SET(hist,i,j,(double)sqrt((double)hist_data));
66 }
67 }
68
69 /* ---------- */
70
71 void sqr_root_and_normalize_pairwise(Imrect *hist)
72 {
73 int i,j, min_x, min_y, max_x, max_y;
74 double sum=0.0;
75 double gl;
76
77 /* (1) Find the sum of the histogram enties.
78 (2) Divide each histogram entry by this value.
79 (3) Replace each histogram entry by its sqr root.
80 This is done for statistical reasons. */
81
82 min_x = hist->region->lx;
83 min_y = hist->region->ly;
84 max_x = hist->region->ux;
85 max_y = hist->region->uy;
86
87 for(i=min_y;i<max_y;i++)
88 for(j=min_x;j<max_x;j++)
89 {
90 IM_PIX_GET(hist,i,j,gl);
91 sum += gl;
92 }
93
94 if (sum!=0.0)
95 {
96 for(i=min_y;i<max_y;i++)
97 for(j=min_x;j<max_x;j++)
98 {
99 IM_PIX_GET(hist,i,j,gl)
100 IM_PIX_SET(hist,i,j,gl/(double)sum);
101 }
102
103 sqr_root_pairwise(hist);
104 }
105 }
106
107 /* ---------- Functions ---------- */
108
109 double sum_pairwise(Imrect *hist)
110 {
111 int i,j, min_x, min_y, max_x, max_y;
112 double gl;
113 double sum = 0.0;
114
115 min_x = hist->region->lx;
116 min_y = hist->region->ly;
117 max_x = hist->region->ux;
118 max_y = hist->region->uy;
119
120 for(i=min_y;i<max_y;i++)
121 for(j=min_x;j<max_x;j++)
122 {
123 IM_PIX_GET(hist,i,j,gl);
124 sum += gl;
125 }
126
127 return (sum);
128 }
129
130 /* ---------- */
131
132 Imrect *build_pairwise(Line2 *ref_line, Model_poly_header *ref_model, List *line_list,
133 double window_r)
134 {
135 Imrect *pairwise;
136 List *ptr;
137 Line2 *lptr;
138 Hist_ref *ref;
139 double l, l0;
140 int type = -1;
141 int i,j;
142
143 /* (1) Allocate space for the histogram.
144 (2) Compile the histogram from the ref line and
145 the line list.
146 (3) Sqrt the hist to get the correct statistics.
147 (4) Add a Hist_ref to the histograms props list.
148 */
149
150 pairwise = im_pairs_alloc();
151 if (pairwise==NULL)
152 {
153 format("Error: Not space for pairwise histogram.\n");
154 return NULL;
155 }
156
157 for(i = pairwise->region->ly; i < pairwise->region->uy; i++)
158 for(j = pairwise->region->lx; j < pairwise->region->ux; j++)
159 IM_PIX_SET(pairwise,i, j, 0.0);
160
161 for(ptr=line_list;ptr!=NULL;ptr=ptr->next)
162 {
163 lptr=ptr->to;
164 if (vec2_dist(ref_line->p, lptr->p) < window_r)
165 compare_lines(pairwise, ref_line, lptr, &type);
166 }
167
168 l0 = vec2_dist(ref_line->p1, ref_line->p2);
169 l = sum_pairwise(pairwise)/l0;
170
171 sqr_root_pairwise(pairwise);
172
173 ref = ralloc(sizeof(Hist_ref));
174 if (ref==NULL)
175 {
176 format("Error: No space for Hist_ref.\n");
177 im_free(pairwise);
178 return NULL;
179 }
180
181 ref->model = ref_model;
182 ref->matches = NULL;
183 ref->line_seg = ref_line;
184 ref->l0 = l0;
185 ref->l = l;
186 ref->type = type;
187 ref->scale_factor = 1.0;
188
189 pairwise->props = proplist_add(pairwise->props, ref, HIST_REF_TYPE,
190 hist_ref_free);
191 /* remove fast indexing scheme for now NAT 15/5/95
192 findex = nonzero_index(pairwise);
193 pairwise->props = proplist_add(pairwise->props, findex, FAST_INDEX1,
194 vector_free);
195 if (type==DIRECTED)
196 {
197 findex = nonzero_index2(pairwise);
198 pairwise->props = proplist_add(pairwise->props, findex, FAST_INDEX2,
199 vector_free);
200 }
201 pairwise->label = ref_line->label;
202 */
203 return (pairwise);
204
205 }
206
207 /* ---------- */
208
209 Imrect *build_normalized_pairwise(Line2 *ref_line, Model_poly_header *ref_model,
210 List *line_list, double window_r)
211 {
212 Imrect *pairwise;
213 List *ptr;
214 Line2 *lptr;
215 Hist_ref *ref;
216 double l, l0;
217 int type = -1;
218 int i,j;
219
220 /* (1) Allocate space for the histogram.
221 (2) Compile the histogram from the ref line and
222 the line list.
223 (3) Normalize the histogram.
224 (4) Sqrt the hist to get the correct statistics.
225 (5) Add a Hist_ref to the histograms props list.
226 */
227
228 pairwise = im_pairs_alloc();
229 if (pairwise==NULL)
230 {
231 printf("Error: Not space for pairwise histogram.\n");
232 return NULL;
233 }
234
235 for(i = pairwise->region->ly; i < pairwise->region->uy; i++)
236 for(j = pairwise->region->lx; j < pairwise->region->ux; j++)
237 IM_PIX_SET(pairwise,i, j, 0.0);
238
239 for(ptr=line_list;ptr!=NULL;ptr=ptr->next)
240 {
241 lptr=ptr->to;
242 if (vec2_dist(ref_line->p, lptr->p) < window_r)
243 compare_lines(pairwise, ref_line, lptr, &type);
244 }
245
246 l0 = vec2_dist(ref_line->p1, ref_line->p2);
247 l = sum_pairwise(pairwise)/l0;
248
249 sqr_root_and_normalize_pairwise(pairwise);
250
251 ref = ralloc(sizeof(Hist_ref));
252 if (ref==NULL)
253 {
254 printf("Error: No space for Hist_ref.\n");
255 im_free(pairwise);
256 return NULL;
257 }
258
259 ref->model = ref_model;
260 ref->matches = NULL;
261 ref->line_seg = ref_line;
262 ref->l0 = l0;
263 ref->l = l;
264 ref->type = type;
265 ref->scale_factor = 1.0;
266
267 pairwise->props = proplist_add(pairwise->props, ref, HIST_REF_TYPE,
268 hist_ref_free);
269 /*
270 pairwise->label = ref_line->label;
271 */
272 return (pairwise);
273
274 }
275
276 /* ---------- */
277
278 void hist_ref_free(Hist_ref *hist_ref, int type)
279 {
280
281 /* (1) Free the matches list. (Only the links,
282 not the data which they point to).
283 (2) Free the copied line seg.
284 (3) Free the hist_ref itself. */
285
286 if (type!=HIST_REF_TYPE)
287 {
288 format("Error: Ilegal use of 'hist_ref_free'.\n");
289 return;
290 }
291
292 list_rm_links(hist_ref->matches);
293 /* line2_free(hist_ref->line_seg);
294 not to be done here NAT
295 */
296 rfree(hist_ref);
297
298 }
299
300
301 /* ---------- */
302
303 Imrect *pairs_build_norm_hist_scale(Line2 *ref_line, Model_poly_header *ref_model, List *line_list, double window_r, double scale)
304 {
305 Imrect *pairwise;
306 List *ptr;
307 Line2 *lptr, *scaled_ref, *scaled_line, *clip_line;
308 Hist_ref *ref;
309 double l, l0;
310 int type = -1;
311 Transform2 trans;
312 int i,j;
313
314
315 /* (1) Allocate space for the histogram.
316 (2) Compile the histogram from the scaled ref line and
317 the scaled line list.
318 (3) Normalize the histogram.
319 (4) Sqrt the hist to get the correct statistics.
320 (5) Add a Hist_ref to the histograms props list.
321 */
322
323 trans.R = rot2_with_scale(0.0, scale);
324 trans.t = vec2_zero();
325
326 scaled_ref = line2_copy(ref_line);
327 line2_transform(scaled_ref, trans);
328
329 pairwise = im_pairs_alloc();
330 if (pairwise==NULL)
331 {
332 printf("Error: Not space for pairwise histogram.\n");
333 return NULL;
334 }
335
336 //init pairwaise
337
338 for(i = pairwise->region->ly; i < pairwise->region->uy; i++)
339 for(j = pairwise->region->lx; j < pairwise->region->ux; j++)
340 IM_PIX_SET(pairwise,i, j, 0.0);
341
342
343
344 for(ptr=line_list;ptr!=NULL;ptr=ptr->next)
345 {
346 lptr= (Line2 *) ptr->to;
347 scaled_line = line2_copy(lptr);
348 line2_transform(scaled_line, trans);
349
350 clip_line = line2_circ_inter(scaled_line, scaled_ref->p, window_r);
351
352
353 if (clip_line!=NULL)
354 {
355 compare_lines(pairwise, scaled_ref, clip_line, &type);
356 }
357
358 line2_free(scaled_line);
359 line2_free(clip_line);
360 }
361
362 l0 = vec2_dist(scaled_ref->p1, scaled_ref->p2);
363 l = sum_pairwise(pairwise)/l0;
364
365
366 sqr_root_and_normalize_pairwise(pairwise);
367
368 ref = ralloc(sizeof(Hist_ref));
369 if (ref==NULL)
370 {
371 printf("Error: No space for Hist_ref.\n");
372 im_free(pairwise);
373 return NULL;
374 }
375
376 line2_free(scaled_ref);
377
378 ref->model = ref_model;
379 ref->matches = NULL;
380 ref->line_seg = ref_line;
381 ref->l0 = l0;
382 ref->l = l;
383 ref->type = type;
384 ref->scale_factor = scale;
385
386 pairwise->props = proplist_add(pairwise->props, ref, HIST_REF_TYPE,
387 hist_ref_free);
388 /* No longer supported APA 19/6/95
389 pairwise->label = ref_line->label; */
390
391 return (pairwise);
392
393 }
394
395 /* ---------- */
396
397 Hist_ref *hist_ref_get(Imrect *hist)
398 {
399 Hist_ref *hist_ref;
400
401 hist_ref = prop_get(hist->props, HIST_REF_TYPE);
402
403 if (hist_ref==NULL)
404 {
405 printf("Error: Histogram has no hist_ref.\n");
406 return NULL;
407 }
408
409 return hist_ref;
410 }
411
412 /* ---------- */
413
414 Line2 *ref_line_from_hist(Imrect *hist)
415 {
416 Hist_ref *hist_ref;
417 Line2 *ref_line;
418
419 if ((hist_ref = hist_ref_get(hist))==NULL) return NULL;
420
421 if ((ref_line = hist_ref->line_seg)==NULL)
422 {
423 printf("Error: Histogram has no reference line.\n");
424 return NULL;
425 }
426
427 return ref_line;
428 }
429
430 /* ---------- */
431
432 List *model_geom_from_hist(Imrect *hist)
433 {
434 Hist_ref *hist_ref;
435 Model_poly_header *model_poly_head;
436 List *model_geom;
437
438 if ((hist_ref = hist_ref_get(hist))==NULL) return NULL;
439
440 if ((model_poly_head = hist_ref->model)==NULL)
441 {
442 printf("Error: Histogram has no model_poly_header.\n");
443 return NULL;
444 }
445
446 if ((model_geom = model_poly_head->model_poly_data)==NULL)
447 {
448 printf("Error: Histogram has no model geometry.\n");
449 return NULL;
450 }
451
452 return model_geom;
453 }
454
455 /* ---------- */
456
457 List *matches_list_from_hist(Imrect *hist)
458 {
459 Hist_ref *hist_ref;
460 List *matches_list;
461
462 if ((hist_ref = hist_ref_get(hist))==NULL) return NULL;
463
464 if ((matches_list = hist_ref->matches)==NULL)
465 {
466 printf("Error: Histogram has no matches.\n");
467 return NULL;
468 }
469
470 return matches_list;
471 }
472
473 /* ---------- */
474
475 Hist_ref *hist_ref_copy(Hist_ref *hist_ref)
476 {
477 Hist_ref *copy;
478
479 copy = ralloc(sizeof(Hist_ref));
480
481 copy->line_seg = hist_ref->line_seg;
482 copy->model = hist_ref->model;
483 copy->l0 = hist_ref->l0;
484 copy->l = hist_ref->l;
485 copy->scale_factor = hist_ref->scale_factor;
486 copy->type = hist_ref->type;
487
488 return copy;
489 }
490
491 /* ---------- */
492
493 double dot_product(Imrect *im1, Imrect *im2, int type1, int type2)
494 {
495 double dp = 0.0;
496 Vector *index;
497 Imregion *roi;
498 float *row1, *row2;
499 int lx, ux, ly, uy;
500 int i, j;
501
502 if (im1 == NULL || im2 == NULL)
503 return (-1.0);
504 if (type1!=type2) return(-1.0);
505
506 if ((index = prop_get(im1->props,FAST_INDEX1)))
507 {
508 /* remove next bit until finished NAT 15/5/95
509 dp = fast_dot_prod(im1,im2,type1,type2,index);
510 */
511 return(dp);
512 }
513 else if ((index = prop_get(im2->props,FAST_INDEX1)))
514 {
515 /* remove next bit until finished NAT 15/5/95
516 dp = fast_dot_prod(im1,im2,type1,type2,index);
517 */
518 return(dp);
519 }
520
521 roi = roi_inter(im1->region, im2->region);
522 if (roi == NULL)
523 return (-1.0);
524 lx = roi->lx;
525 ux = roi->ux;
526 ly = roi->ly;
527 uy = roi->uy;
528
529 for (i = ly; i < uy; ++i)
530 {
531 row1 = IM_ROW(im1,i);
532 row2 = IM_ROW(im2,i);
533 for (j = lx; j < ux; ++j)
534 dp += row1[j]*row2[j];
535 }
536
537 rfree(roi);
538 return (dp);
539 }
540
541 double dot_product2(Imrect *im1, Imrect *im2, int type1, int type2)
542 {
543 double dp = 0.0;
544 Imregion *roi;
545 Vector *index1,*index2;
546 float *row1, *row2;
547 int lx, ux, ly, uy;
548 int ishift;
549 int i, j;
550
551 if (im1 == NULL || im2 == NULL)
552 return (-1.0);
553 if (type1!=type2) return(-1.0);
554
555 if ((index1 = prop_get(im1->props,FAST_INDEX1))
556 && (index2 = prop_get(im1->props,FAST_INDEX2)))
557 {
558 /* remove until ready NAT 16/5/95
559 dp = fast_dot_prod2(im1,im2,type1,type2,index1,index2);
560 */
561 return(dp);
562 }
563 else if ((index1 = prop_get(im2->props,FAST_INDEX1))
564 && (index2 = prop_get(im2->props,FAST_INDEX2)))
565 {
566 /* dp = fast_dot_prod2(im2,im1,type1,type2,index1,index2);
567 */
568 return(dp);
569 }
570
571 roi = roi_inter(im1->region, im2->region);
572 if (roi == NULL)
573 return (-1.0);
574 lx = roi->lx;
575 ux = roi->ux;
576 ly = roi->ly;
577 uy = roi->uy;
578 ishift = (int)((ux-lx)/2.0);
579 for (i = ly; i < uy; ++i)
580 {
581 row1 = IM_ROW(im1,i);
582 row2 = IM_ROW(im2,-i-1);
583 for (j = lx; j < ishift; ++j)
584 dp += row1[j]*row2[j+ishift];
585 for (j = ishift; j < ux; ++j)
586 dp += row1[j]*row2[j-ishift];
587 }
588
589 rfree(roi);
590 return (dp);
591 }
592
593 /* ---------- */
594
This page was automatically generated by the
LXR engine.
Visit the LXR main site for more
information.