GRASS Programmer's Manual
6.4.2(2012)
|
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 }