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

Linux Cross Reference
Tina6/tina-libs/tina/geometry/geomImg_nonmax.c

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

  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 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 General Public License for more details.
 14  *
 15  * You should have received a copy of the GNU 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  * ANY users of TINA who require exemption from the existing licence must
 20  * negotiate a new licence with Dr. Neil.A.Thacker, the sole agent for
 21  * the University of Manchester.
 22  *
 23  **********
 24  * 
 25  * Program :    TINA
 26  * File    :  $Source: /home/tina/cvs/tina-libs/tina/geometry/geomImg_nonmax.c,v $
 27  * Date    :  $Date: 2009/03/19 18:36:26 $
 28  * Version :  $Revision: 1.3 $
 29  * CVS Id  :  $Id: geomImg_nonmax.c,v 1.3 2009/03/19 18:36:26 paul Exp $
 30  *
 31  * Author  : Legacy TINA
 32  *
 33  * Notes : nonmax.c
 34  * 
 35  * locate directional peak in gradsq image of type float using non maximal
 36  * suppresion and form appropriate edge rect
 37  * 
 38  * fit is performed to sub pixel acuity using quadratic fit across peak in
 39  * quantised direction
 40  * 
 41  * if sub pix falls outside pixel ignore it
 42  *
 43  *********
 44 */
 45 
 46 #include "geomImg_nonmax.h"
 47 
 48 #if HAVE_CONFIG_H
 49   #include <config.h>
 50 #endif
 51 
 52 #include <math.h>
 53 #include <tina/sys/sysDef.h>
 54 #include <tina/sys/sysPro.h>
 55 #include <tina/math/mathDef.h>
 56 #include <tina/math/mathPro.h>
 57 #include <tina/image/imgDef.h>
 58 #include <tina/image/imgPro.h>
 59 #include <tina/geometry/geom_EdgeDef.h>
 60 #include <tina/geometry/geom_EdgePro.h>
 61 
 62 
 63 Imrect *nonmaxsup(Imrect * gradx, Imrect * grady, Imrect * gradsq, double thres)
 64 {
 65     Imrect *edge_image;
 66     Imregion *region;
 67     int     i, j;
 68     int     lx, ux, ly, uy;
 69     int     uxm1, uym1;
 70     float   sqthres = (float)(thres * thres);
 71     float  *rowx, *rowy;
 72     float  *grad, *grad_m1, *grad_p1;
 73 
 74     if (gradx == NULL || grady == NULL || gradsq == NULL)
 75         return (NULL);
 76 
 77     region = gradsq->region;
 78 
 79     edge_image = im_alloc(gradsq->height, gradsq->width, region, ptr_v);
 80     region = edge_image->region;
 81 
 82     lx = region->lx;
 83     ux = region->ux;
 84     ly = region->ly;
 85     uy = region->uy;
 86 
 87     uxm1 = ux - 1;
 88     uym1 = uy - 1;
 89 
 90     rowx = fvector_alloc(lx, ux);
 91     rowy = fvector_alloc(lx, ux);
 92     grad = fvector_alloc(lx, ux);
 93     grad_m1 = fvector_alloc(lx, ux);
 94     grad_p1 = fvector_alloc(lx, ux);
 95 
 96     for (i = ly + 1; i < uym1; ++i)
 97     {
 98         im_get_rowf(rowy, grady, i, lx, ux);
 99         im_get_rowf(rowx, gradx, i, lx, ux);
100         im_get_rowf(grad, gradsq, i, lx, ux);
101         im_get_rowf(grad_p1, gradsq, i + 1, lx, ux);
102         im_get_rowf(grad_m1, gradsq, i - 1, lx, ux);
103 
104         for (j = lx + 1; j < uxm1; ++j)
105         {
106             float   a, b;
107             float   st;
108             float   del;
109             float   dr, dc;
110             float   gradratio, fabsgradratio;
111             Edgel  *eptr;
112             Vec2    pos = {Vec2_id};
113 
114             if (grad[j] < sqthres)
115                 continue;
116 
117             st = grad[j];
118             gradratio = rowx[j] / rowy[j];
119             fabsgradratio = (float)fabs(gradratio);
120             if (fabsgradratio < 1.5 && fabsgradratio > 0.67)    /* diagonal */
121             {
122                 if ((st < grad[j - 1] || st <= grad[j + 1]) &&
123                     (st < grad_m1[j] || st <= grad_p1[j]))
124                     continue;
125                 if (gradratio > 0)      /* right diagonal */
126                 {
127                     if (st < grad_m1[j - 1] || st <= grad_p1[j + 1])
128                         continue;
129                     st = (float)sqrt(st);
130                     a = (float)sqrt(grad_m1[j - 1]);
131                     b = (float)sqrt(grad_p1[j + 1]);
132                     del = (a - b) / ((a + b - st * 2) * 2);
133                     dr = (float)(del + 0.5);
134                     dc = (float)(del + 0.5);
135                 } else
136                 {
137                     if (st < grad_p1[j - 1] || st <= grad_m1[j + 1])
138                         continue;
139                     st = (float)sqrt(st);
140                     a = (float)sqrt(grad_p1[j - 1]);
141                     b = (float)sqrt(grad_m1[j + 1]);
142                     del = (a - b) / ((a + b - st * 2) * 2);
143                     dr = (float)(0.5 - del);
144                     dc = (float)(del + 0.5);
145                 }
146             } else if (fabs(rowx[j]) > fabs(rowy[j]))
147             {
148                 if (st < grad[j - 1] || st <= grad[j + 1])
149                     continue;
150                 st = (float)sqrt(st);
151                 a = (float)sqrt(grad[j - 1]);
152                 b = (float)sqrt(grad[j + 1]);
153                 del = (a - b) / ((a + b - st * 2) * 2);
154                 dr = (float)0.5;
155                 dc = (float)(del + 0.5);
156             } else
157             {
158                 if (st < grad_m1[j] || st <= grad_p1[j])
159                     continue;
160                 st = (float)sqrt(st);
161                 a = (float)sqrt(grad_m1[j]);
162                 b = (float)sqrt(grad_p1[j]);
163                 del = (a - b) / ((a + b - st * 2) * 2);
164                 dr = (float)(del + 0.5);
165                 dc = (float)0.5;        /* pixel centre */
166             }
167 
168             {
169                 int     ii, jj; /* actual storage location */
170 
171                 pos = vec2(j + dc, i + dr);
172                 jj = (int)floor(vec2_x(pos));
173                 ii = (int)floor(vec2_y(pos));
174                 if (ii < ly || ii >= uy || jj < lx || jj >= ux ||
175                     IM_PTR(edge_image, ii, jj) != NULL)
176                     continue;
177 
178                 eptr = edge_alloc(EDGE_RAW);
179                 eptr->pos = pos;
180                 eptr->orient = (float)atan2(rowx[j], rowy[j]);
181                 eptr->contrast = st;
182                 IM_PTR(edge_image, ii, jj) = (void *) eptr;
183             }
184         }
185     }
186     fvector_free(rowx, lx);
187     fvector_free(rowy, lx);
188     fvector_free(grad, lx);
189     fvector_free(grad_m1, lx);
190     fvector_free(grad_p1, lx);
191     return (edge_image);
192 }
193 

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