00001
00020 #include <stdlib.h>
00021 #include <grass/gis.h>
00022 #include <grass/Vect.h>
00023 #include <grass/glocale.h>
00024
00038 int
00039 Vect_clean_small_angles_at_nodes(struct Map_info *Map, int otype,
00040 struct Map_info *Err)
00041 {
00042 int node;
00043 int nmodif = 0;
00044 struct line_pnts *Points;
00045 struct line_cats *SCats, *LCats, *OCats;
00046
00047 Points = Vect_new_line_struct();
00048 SCats = Vect_new_cats_struct();
00049 LCats = Vect_new_cats_struct();
00050 OCats = Vect_new_cats_struct();
00051
00052 for (node = 1; node <= Vect_get_num_nodes(Map); node++) {
00053 int i, nlines;
00054
00055 G_debug(3, "node = %d", node);
00056 if (!Vect_node_alive(Map, node))
00057 continue;
00058
00059 while (1) {
00060 float angle1 = -100;
00061 int line1 = -999;
00062 int clean = 1;
00063
00064
00065 nlines = Vect_get_node_n_lines(Map, node);
00066 G_debug(3, "nlines = %d", nlines);
00067
00068 for (i = 0; i < nlines; i++) {
00069 P_LINE *Line;
00070 int line2;
00071 float angle2;
00072
00073 line2 = Vect_get_node_line(Map, node, i);
00074 Line = Map->plus.Line[abs(line2)];
00075 if (!Line)
00076 continue;
00077 G_debug(4, " type = %d", Line->type);
00078 if (!(Line->type & (otype & GV_LINES)))
00079 continue;
00080
00081 angle2 = Vect_get_node_line_angle(Map, node, i);
00082 if (angle2 == -9.0)
00083 continue;
00084
00085 G_debug(4, " line1 = %d angle1 = %e line2 = %d angle2 = %e",
00086 line1, angle1, line2, angle2);
00087
00088 if (angle2 == angle1) {
00089 int j;
00090 double length1, length2;
00091 int short_line;
00092 int long_line;
00093 int new_short_line = 0;
00094 int short_type, long_type, type;
00095 double x, y, z, nx, ny, nz;
00096
00097 G_debug(4, " identical angles -> clean");
00098
00099
00100 Vect_read_line(Map, Points, NULL, abs(line1));
00101 if (line1 > 0) {
00102 length1 =
00103 Vect_points_distance(Points->x[0], Points->y[0],
00104 0.0, Points->x[1],
00105 Points->y[1], 0.0, 0);
00106 }
00107 else {
00108 int np;
00109
00110 np = Points->n_points;
00111 length1 =
00112 Vect_points_distance(Points->x[np - 1],
00113 Points->y[np - 1], 0.0,
00114 Points->x[np - 2],
00115 Points->y[np - 2], 0.0, 0);
00116 }
00117
00118 Vect_read_line(Map, Points, NULL, abs(line2));
00119 if (line2 > 0) {
00120 length2 =
00121 Vect_points_distance(Points->x[0], Points->y[0],
00122 0.0, Points->x[1],
00123 Points->y[1], 0.0, 0);
00124 }
00125 else {
00126 int np;
00127
00128 np = Points->n_points;
00129 length2 =
00130 Vect_points_distance(Points->x[np - 1],
00131 Points->y[np - 1], 0.0,
00132 Points->x[np - 2],
00133 Points->y[np - 2], 0.0, 0);
00134 }
00135
00136 G_debug(4, " length1 = %f length2 = %f", length1,
00137 length2);
00138
00139 if (length1 < length2) {
00140 short_line = line1;
00141 long_line = line2;
00142 }
00143 else {
00144 short_line = line2;
00145 long_line = line1;
00146 }
00147
00148
00149 short_type =
00150 Vect_read_line(Map, Points, SCats, abs(short_line));
00151
00152 if (short_line > 0) {
00153 x = Points->x[1];
00154 y = Points->y[1];
00155 z = Points->z[1];
00156 Vect_line_delete_point(Points, 0);
00157 }
00158 else {
00159 x = Points->x[Points->n_points - 2];
00160 y = Points->y[Points->n_points - 2];
00161 z = Points->z[Points->n_points - 2];
00162 Vect_line_delete_point(Points, Points->n_points - 1);
00163 }
00164
00165
00166
00167 Vect_get_node_coor(Map, node, &nx, &ny, &nz);
00168
00169 if (Points->n_points > 1) {
00170 new_short_line =
00171 Vect_rewrite_line(Map, abs(short_line),
00172 short_type, Points, SCats);
00173 }
00174 else {
00175 Vect_delete_line(Map, abs(short_line));
00176 }
00177
00178
00179
00180 if (abs(line1) == abs(line2)) {
00181 if (long_line > 0)
00182 long_line = new_short_line;
00183 else
00184 long_line = -new_short_line;
00185 }
00186
00187
00188
00189 long_type =
00190 Vect_read_line(Map, NULL, LCats, abs(long_line));
00191
00192 Vect_reset_cats(OCats);
00193 for (j = 0; j < SCats->n_cats; j++) {
00194 Vect_cat_set(OCats, SCats->field[j], SCats->cat[j]);
00195 }
00196 for (j = 0; j < LCats->n_cats; j++) {
00197 Vect_cat_set(OCats, LCats->field[j], LCats->cat[j]);
00198 }
00199
00200 if (long_type == GV_BOUNDARY || short_type == GV_BOUNDARY) {
00201 type = GV_BOUNDARY;
00202 }
00203 else {
00204 type = GV_LINE;
00205 }
00206
00207 Vect_reset_line(Points);
00208 Vect_append_point(Points, nx, ny, nz);
00209 Vect_append_point(Points, x, y, z);
00210 Vect_write_line(Map, type, Points, OCats);
00211
00212 if (Err) {
00213 Vect_write_line(Err, type, Points, OCats);
00214 }
00215
00216
00217 long_type =
00218 Vect_read_line(Map, Points, LCats, abs(long_line));
00219 if (long_line > 0) {
00220 Points->x[0] = x;
00221 Points->y[0] = y;
00222 Points->z[0] = z;
00223 }
00224 else {
00225 Points->x[Points->n_points - 1] = x;
00226 Points->y[Points->n_points - 1] = y;
00227 Points->z[Points->n_points - 1] = z;
00228 }
00229 Vect_line_prune(Points);
00230 if (Points->n_points > 1) {
00231 Vect_rewrite_line(Map, abs(long_line), long_type,
00232 Points, LCats);
00233 }
00234 else {
00235 Vect_delete_line(Map, abs(long_line));
00236 }
00237
00238 nmodif += 3;
00239 clean = 0;
00240
00241 break;
00242 }
00243
00244 line1 = line2;
00245 angle1 = angle2;
00246 }
00247
00248 if (clean || !Vect_node_alive(Map, node))
00249 break;
00250 }
00251 }
00252
00253 return (nmodif);
00254 }