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

Linux Cross Reference
Tina4/src/vision/improc/nonmax.c

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

  1 /**@(#)
  2 **/
  3 #include <math.h>
  4 #include <tina/sys.h>
  5 #include <tina/sysfuncs.h>
  6 #include <tina/math.h>
  7 #include <tina/mathfuncs.h>
  8 #include <tina/vision.h>
  9 
 10 
 11 /* nonmax.c
 12  * 
 13  * locate directional peak in gradsq image of type float using non maximal
 14  * suppresion and form appropriate edge rect
 15  * 
 16  * fit is performed to sub pixel acuity using quadratic fit across peak in
 17  * quantised direction
 18  * 
 19  * if sub pix falls outside pixel ignore it
 20  * 
 21  */
 22 
 23 Imrect *nonmaxsup(Imrect * gradx, Imrect * grady, Imrect * gradsq, double thres)
 24 {
 25     Imrect *edge_image;
 26     Imregion *region;
 27     int     i, j;
 28     int     lx, ux, ly, uy;
 29     int     uxm1, uym1;
 30     float   sqthres = (float)(thres * thres);
 31     float  *rowx, *rowy;
 32     float  *grad, *grad_m1, *grad_p1;
 33 
 34     if (gradx == NULL || grady == NULL || gradsq == NULL)
 35         return (NULL);
 36 
 37     region = gradsq->region;
 38 
 39     edge_image = im_alloc(gradsq->height, gradsq->width, region, ptr_v);
 40     region = edge_image->region;
 41 
 42     lx = region->lx;
 43     ux = region->ux;
 44     ly = region->ly;
 45     uy = region->uy;
 46 
 47     uxm1 = ux - 1;
 48     uym1 = uy - 1;
 49 
 50     rowx = fvector_alloc(lx, ux);
 51     rowy = fvector_alloc(lx, ux);
 52     grad = fvector_alloc(lx, ux);
 53     grad_m1 = fvector_alloc(lx, ux);
 54     grad_p1 = fvector_alloc(lx, ux);
 55 
 56     for (i = ly + 1; i < uym1; ++i)
 57     {
 58         im_get_rowf(rowy, grady, i, lx, ux);
 59         im_get_rowf(rowx, gradx, i, lx, ux);
 60         im_get_rowf(grad, gradsq, i, lx, ux);
 61         im_get_rowf(grad_p1, gradsq, i + 1, lx, ux);
 62         im_get_rowf(grad_m1, gradsq, i - 1, lx, ux);
 63 
 64         for (j = lx + 1; j < uxm1; ++j)
 65         {
 66             float   a, b;
 67             float   st;
 68             float   del;
 69             float   dr, dc;
 70             float   gradratio, fabsgradratio;
 71             Edgel  *eptr;
 72             Edgel  *edge_alloc();
 73             Vec2    pos = {Vec2_id};
 74 
 75             if (grad[j] < sqthres)
 76                 continue;
 77 
 78             st = grad[j];
 79             gradratio = rowx[j] / rowy[j];
 80             fabsgradratio = (float)fabs(gradratio);
 81             if (fabsgradratio < 1.5 && fabsgradratio > 0.67)    /* diagonal */
 82             {
 83                 if ((st < grad[j - 1] || st <= grad[j + 1]) &&
 84                     (st < grad_m1[j] || st <= grad_p1[j]))
 85                     continue;
 86                 if (gradratio > 0)      /* right diagonal */
 87                 {
 88                     if (st < grad_m1[j - 1] || st <= grad_p1[j + 1])
 89                         continue;
 90                     st = (float)sqrt(st);
 91                     a = (float)sqrt(grad_m1[j - 1]);
 92                     b = (float)sqrt(grad_p1[j + 1]);
 93                     del = (a - b) / ((a + b - st * 2) * 2);
 94                     dr = (float)(del + 0.5);
 95                     dc = (float)(del + 0.5);
 96                 } else
 97                 {
 98                     if (st < grad_p1[j - 1] || st <= grad_m1[j + 1])
 99                         continue;
100                     st = (float)sqrt(st);
101                     a = (float)sqrt(grad_p1[j - 1]);
102                     b = (float)sqrt(grad_m1[j + 1]);
103                     del = (a - b) / ((a + b - st * 2) * 2);
104                     dr = (float)(0.5 - del);
105                     dc = (float)(del + 0.5);
106                 }
107             } else if (fabs(rowx[j]) > fabs(rowy[j]))
108             {
109                 if (st < grad[j - 1] || st <= grad[j + 1])
110                     continue;
111                 st = (float)sqrt(st);
112                 a = (float)sqrt(grad[j - 1]);
113                 b = (float)sqrt(grad[j + 1]);
114                 del = (a - b) / ((a + b - st * 2) * 2);
115                 dr = (float)0.5;
116                 dc = (float)(del + 0.5);
117             } else
118             {
119                 if (st < grad_m1[j] || st <= grad_p1[j])
120                     continue;
121                 st = (float)sqrt(st);
122                 a = (float)sqrt(grad_m1[j]);
123                 b = (float)sqrt(grad_p1[j]);
124                 del = (a - b) / ((a + b - st * 2) * 2);
125                 dr = (float)(del + 0.5);
126                 dc = (float)0.5;        /* pixel centre */
127             }
128 
129             {
130                 int     ii, jj; /* actual storage location */
131 
132                 pos = vec2(j + dc, i + dr);
133                 jj = (int)floor(vec2_x(pos));
134                 ii = (int)floor(vec2_y(pos));
135                 if (ii < ly || ii >= uy || jj < lx || jj >= ux ||
136                     IM_PTR(edge_image, ii, jj) != NULL)
137                     continue;
138 
139                 eptr = edge_alloc(EDGE_RAW);
140                 eptr->pos = pos;
141                 eptr->orient = (float)atan2(rowx[j], rowy[j]);
142                 eptr->contrast = st;
143                 IM_PTR(edge_image, ii, jj) = (void *) eptr;
144             }
145         }
146     }
147     fvector_free((void *) rowx, lx);
148     fvector_free((void *) rowy, lx);
149     fvector_free((void *) grad, lx);
150     fvector_free((void *) grad_m1, lx);
151     fvector_free((void *) grad_p1, lx);
152     return (edge_image);
153 }
154 

~ [ 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.