GRASS Programmer's Manual  6.4.1(2011)
clean_nodes.c
Go to the documentation of this file.
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;   /* value not important, just for debug */
00062             int clean = 1;
00063 
00064             nlines = Vect_get_node_n_lines(Map, node);
00065             G_debug(3, "nlines = %d", nlines);
00066 
00067             for (i = 0; i < nlines; i++) {
00068                 P_LINE *Line;
00069                 int line2;
00070                 float angle2;
00071 
00072                 line2 = Vect_get_node_line(Map, node, i);
00073                 Line = Map->plus.Line[abs(line2)];
00074                 if (!Line)
00075                     continue;
00076                 G_debug(4, "  type = %d", Line->type);
00077                 if (!(Line->type & (otype & GV_LINES)))
00078                     continue;
00079 
00080                 angle2 = Vect_get_node_line_angle(Map, node, i);
00081                 if (angle2 == -9.0)
00082                     continue;   /* Degenerated line */
00083 
00084                 G_debug(4, "  line1 = %d angle1 = %e line2 = %d angle2 = %e",
00085                         line1, angle1, line2, angle2);
00086 
00087                 if (angle2 == angle1) {
00088                     int j;
00089                     double length1, length2;
00090                     int short_line;     /* line with shorter end segment */
00091                     int long_line;      /* line with longer end segment */
00092                     int new_short_line = 0;     /* line number of short line after rewrite */
00093                     int short_type, long_type, type;
00094                     double x, y, z, nx, ny, nz;
00095 
00096                     G_debug(4, "  identical angles -> clean");
00097 
00098                     /* Length of end segments for both lines */
00099                     Vect_read_line(Map, Points, NULL, abs(line1));
00100                     if (line1 > 0) {
00101                         length1 =
00102                             Vect_points_distance(Points->x[0], Points->y[0],
00103                                                  0.0, Points->x[1],
00104                                                  Points->y[1], 0.0, 0);
00105                     }
00106                     else {
00107                         int np;
00108 
00109                         np = Points->n_points;
00110                         length1 =
00111                             Vect_points_distance(Points->x[np - 1],
00112                                                  Points->y[np - 1], 0.0,
00113                                                  Points->x[np - 2],
00114                                                  Points->y[np - 2], 0.0, 0);
00115                     }
00116 
00117                     Vect_read_line(Map, Points, NULL, abs(line2));
00118                     if (line2 > 0) {
00119                         length2 =
00120                             Vect_points_distance(Points->x[0], Points->y[0],
00121                                                  0.0, Points->x[1],
00122                                                  Points->y[1], 0.0, 0);
00123                     }
00124                     else {
00125                         int np;
00126 
00127                         np = Points->n_points;
00128                         length2 =
00129                             Vect_points_distance(Points->x[np - 1],
00130                                                  Points->y[np - 1], 0.0,
00131                                                  Points->x[np - 2],
00132                                                  Points->y[np - 2], 0.0, 0);
00133                     }
00134 
00135                     G_debug(4, "  length1 = %f length2 = %f", length1,
00136                             length2);
00137 
00138                     if (length1 < length2) {
00139                         short_line = line1;
00140                         long_line = line2;
00141                     }
00142                     else {
00143                         short_line = line2;
00144                         long_line = line1;
00145                     }
00146 
00147                     /* Remove end segment from short_line */
00148                     short_type =
00149                         Vect_read_line(Map, Points, SCats, abs(short_line));
00150 
00151                     if (short_line > 0) {
00152                         x = Points->x[1];
00153                         y = Points->y[1];
00154                         z = Points->z[1];
00155                         Vect_line_delete_point(Points, 0);      /* first */
00156                     }
00157                     else {
00158                         x = Points->x[Points->n_points - 2];
00159                         y = Points->y[Points->n_points - 2];
00160                         z = Points->z[Points->n_points - 2];
00161                         Vect_line_delete_point(Points, Points->n_points - 1);   /* last */
00162                     }
00163 
00164                     /* It may happen that it is one line: node could be deleted,
00165                      * in that case we have to read the node coords first */
00166                     Vect_get_node_coor(Map, node, &nx, &ny, &nz);
00167 
00168                     if (Points->n_points > 1) {
00169                         new_short_line =
00170                             Vect_rewrite_line(Map, abs(short_line),
00171                                               short_type, Points, SCats);
00172                     }
00173                     else {
00174                         Vect_delete_line(Map, abs(short_line));
00175                     }
00176 
00177                     /* It may happen that it is one line, in that case we have to take the new
00178                      * short line as long line, orientation is not changed */
00179                     if (abs(line1) == abs(line2)) {
00180                         if (long_line > 0)
00181                             long_line = new_short_line;
00182                         else
00183                             long_line = -new_short_line;
00184                     }
00185 
00186                     /* Add new line (must be before rewrite of long_line otherwise node could be deleted) */
00187                     long_type =
00188                         Vect_read_line(Map, NULL, LCats, abs(long_line));
00189 
00190                     Vect_reset_cats(OCats);
00191                     for (j = 0; j < SCats->n_cats; j++) {
00192                         Vect_cat_set(OCats, SCats->field[j], SCats->cat[j]);
00193                     }
00194                     for (j = 0; j < LCats->n_cats; j++) {
00195                         Vect_cat_set(OCats, LCats->field[j], LCats->cat[j]);
00196                     }
00197 
00198                     if (long_type == GV_BOUNDARY || short_type == GV_BOUNDARY) {
00199                         type = GV_BOUNDARY;
00200                     }
00201                     else {
00202                         type = GV_LINE;
00203                     }
00204 
00205                     Vect_reset_line(Points);
00206                     Vect_append_point(Points, nx, ny, nz);
00207                     Vect_append_point(Points, x, y, z);
00208                     Vect_write_line(Map, type, Points, OCats);
00209 
00210                     if (Err) {
00211                         Vect_write_line(Err, type, Points, OCats);
00212                     }
00213 
00214                     /* Snap long_line to the new short_line end */
00215                     long_type =
00216                         Vect_read_line(Map, Points, LCats, abs(long_line));
00217                     if (long_line > 0) {
00218                         Points->x[0] = x;
00219                         Points->y[0] = y;
00220                         Points->z[0] = z;
00221                     }
00222                     else {
00223                         Points->x[Points->n_points - 1] = x;
00224                         Points->y[Points->n_points - 1] = y;
00225                         Points->z[Points->n_points - 1] = z;
00226                     }
00227                     Vect_line_prune(Points);
00228                     if (Points->n_points > 1) {
00229                         Vect_rewrite_line(Map, abs(long_line), long_type,
00230                                           Points, LCats);
00231                     }
00232                     else {
00233                         Vect_delete_line(Map, abs(long_line));
00234                     }
00235 
00236                     nmodif += 3;
00237                     clean = 0;
00238 
00239                     break;
00240                 }
00241 
00242                 line1 = line2;
00243                 angle1 = angle2;
00244             }
00245 
00246             if (clean || !Vect_node_alive(Map, node))
00247                 break;
00248         }
00249     }
00250 
00251     return (nmodif);
00252 }
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Defines