GRASS Programmer's Manual  6.4.2(2012)
findzc.c
Go to the documentation of this file.
00001 
00027 #include <stdio.h>
00028 #include <math.h>
00029 
00030 
00032 #define TINY    1.0e-3
00033 
00034 
00055 int
00056 G_math_findzc(double conv[], int size, double zc[], double thresh,
00057               int num_orients)
00058 {
00059     int i, j, p;
00060 
00061     /* go through entire conv image - but skip border rows and cols */
00062     for (i = 1; i < size - 1; i++) {
00063         for (p = i * size + 1, j = 1; j < size - 1; j++, p++) {
00064             int nbr[4];
00065             int ni;
00066 
00067             /* examine the 4-neighbors of position p */
00068             nbr[0] = p - 1;     /* left */
00069             nbr[1] = p + 1;     /* right */
00070             nbr[2] = p - size;  /* up */
00071             nbr[3] = p + size;  /* down */
00072 
00073             zc[p] = 0;
00074 
00075             for (ni = 0; ni < 4; ni++) {
00076                 /* condition for a zc: sign is different than a neighbor
00077                  * and the absolute value is less than that neighbor.
00078                  * Also, threshold magnitudes to eliminate noise
00079                  */
00080                 if ((((conv[p] > 0) && (conv[nbr[ni]] < 0)) ||
00081                      ((conv[p] < 0) && (conv[nbr[ni]] > 0))) &&
00082                     (fabs(conv[p]) < fabs(conv[nbr[ni]])) &&
00083                     (fabs(conv[p] - conv[nbr[ni]]) > thresh)) {
00084                     double ang;
00085                     int dir;
00086 
00087                     /* found a zc here, get angle of gradient */
00088                     if (fabs(conv[nbr[1]] - conv[nbr[0]]) < TINY) {
00089                         ang = M_PI_2;
00090 
00091                         if (conv[nbr[2]] - conv[nbr[3]] < 0)
00092                             ang = -ang;
00093                     }
00094                     else
00095                         ang = atan2(conv[nbr[2]] - conv[nbr[3]],
00096                                     conv[nbr[1]] - conv[nbr[0]]);
00097 
00098                     /* scale -PI..PI to 0..num_orients - 1 */
00099                     dir =
00100                         num_orients * ((ang + M_PI) / (M_PI * 2.0)) + 0.4999;
00101 
00102                     /* shift scale so that 0 (not 8) is straight down */
00103                     dir = (3 * num_orients / 4 + dir) % num_orients;
00104 
00105                     /* add to differentiate between no zc and an orientation */
00106                     zc[p] = 1 + dir;
00107                     break;      /* quit looking at neighbors */
00108                 }
00109             }                   /* for ni */
00110         }                       /* for p */
00111     }
00112 
00113     return 0;
00114 }
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines