GRASS Programmer's Manual
6.4.2(2012)
|
00001 00022 #include <stdio.h> 00023 #include <stdlib.h> 00024 00025 #include <grass/gis.h> 00026 #include <grass/gstypes.h> 00027 00028 #include "rowcol.h" 00029 00033 #define TFAST_PTS 800 00034 00038 #define MFAST_LNS 400 00039 00040 static geoline *copy_line(geoline *); 00041 static geoline *thin_line(geoline *, float); 00042 00051 static geoline *copy_line(geoline * gln) 00052 { 00053 geoline *newln; 00054 int i, np; 00055 00056 newln = (geoline *) G_malloc(sizeof(geoline)); /* G_fatal_error */ 00057 if (!newln) { 00058 return (NULL); 00059 } 00060 00061 np = newln->npts = gln->npts; 00062 00063 if (2 == (newln->dims = gln->dims)) { 00064 newln->p2 = (Point2 *) G_calloc(np, sizeof(Point2)); /* G_fatal_error */ 00065 if (!newln->p2) { 00066 return (NULL); 00067 } 00068 00069 for (i = 0; i < np; i++) { 00070 newln->p2[i][X] = gln->p2[i][X]; 00071 newln->p2[i][Y] = gln->p2[i][Y]; 00072 } 00073 } 00074 else { 00075 newln->p3 = (Point3 *) G_calloc(np, sizeof(Point3)); /* G_fatal_error */ 00076 if (!newln->p3) { 00077 return (NULL); 00078 } 00079 00080 for (i = 0; i < np; i++) { 00081 newln->p3[i][X] = gln->p3[i][X]; 00082 newln->p3[i][Y] = gln->p3[i][Y]; 00083 newln->p3[i][Z] = gln->p3[i][Z]; 00084 } 00085 } 00086 00087 newln->next = NULL; 00088 00089 return (newln); 00090 } 00091 00092 00104 static geoline *thin_line(geoline * gln, float factor) 00105 { 00106 geoline *newln; 00107 int i, nextp, targp; 00108 00109 newln = (geoline *) G_malloc(sizeof(geoline)); /* G_fatal_error */ 00110 if (!newln) { 00111 return (NULL); 00112 } 00113 00114 targp = (int)(gln->npts / factor); 00115 00116 if (targp < 2) { 00117 targp = 2; 00118 } 00119 00120 newln->npts = targp; 00121 00122 if (2 == (newln->dims = gln->dims)) { 00123 newln->p2 = (Point2 *) G_calloc(targp, sizeof(Point2)); /* G_fatal_error */ 00124 if (!newln->p2) { 00125 return (NULL); 00126 } 00127 00128 for (i = 0; i < targp; i++) { 00129 if (i == targp - 1) { 00130 nextp = gln->npts - 1; /* avoid rounding error */ 00131 } 00132 else { 00133 nextp = (int)((i * (gln->npts - 1)) / (targp - 1)); 00134 } 00135 00136 newln->p2[i][X] = gln->p2[nextp][X]; 00137 newln->p2[i][Y] = gln->p2[nextp][Y]; 00138 } 00139 } 00140 else { 00141 newln->p3 = (Point3 *) G_calloc(targp, sizeof(Point3)); /* G_fatal_error */ 00142 if (!newln->p3) { 00143 return (NULL); 00144 } 00145 00146 for (i = 0; i < targp; i++) { 00147 if (i == targp - 1) { 00148 nextp = gln->npts - 1; /* avoid rounding error */ 00149 } 00150 else { 00151 nextp = (int)((i * (gln->npts - 1)) / (targp - 1)); 00152 } 00153 00154 newln->p3[i][X] = gln->p3[nextp][X]; 00155 newln->p3[i][Y] = gln->p3[nextp][Y]; 00156 newln->p3[i][Z] = gln->p3[nextp][Z]; 00157 } 00158 } 00159 00160 newln->next = NULL; 00161 00162 return (newln); 00163 } 00164 00172 float gv_line_length(geoline * gln) 00173 { 00174 int n; 00175 float length = 0.0; 00176 00177 for (n = 0; n < gln->npts - 1; n++) { 00178 if (gln->p2) { 00179 length += GS_P2distance(gln->p2[n + 1], gln->p2[n]); 00180 } 00181 else { 00182 length += GS_distance(gln->p3[n + 1], gln->p3[n]); 00183 } 00184 } 00185 00186 return (length); 00187 } 00188 00196 int gln_num_points(geoline * gln) 00197 { 00198 int np = 0; 00199 geoline *tln; 00200 00201 for (tln = gln; tln; tln = tln->next) { 00202 np += tln->npts; 00203 } 00204 00205 return (np); 00206 } 00207 00215 int gv_num_points(geovect * gv) 00216 { 00217 return (gln_num_points(gv->lines)); 00218 } 00219 00220 00221 00232 int gv_decimate_lines(geovect * gv) 00233 { 00234 int T_pts, A_ppl, N_s; 00235 float decim_factor, slength[MFAST_LNS], T_slength, A_slength; 00236 geoline *gln, *prev; 00237 00238 /* should check if already exists & free if != gv->lines */ 00239 if (TFAST_PTS > (T_pts = gv_num_points(gv))) { 00240 gv->fastlines = gv->lines; 00241 00242 return (1); 00243 } 00244 00245 N_s = 0; 00246 T_slength = 0.0; 00247 decim_factor = T_pts / TFAST_PTS; 00248 A_ppl = T_pts / gv->n_lines; /* (int) Average points per line */ 00249 00250 prev = NULL; 00251 00252 for (gln = gv->lines; gln; gln = gln->next) { 00253 if (gln->npts > A_ppl) { 00254 if (prev) { 00255 prev->next = thin_line(gln, decim_factor); 00256 prev = prev->next; 00257 } 00258 else { 00259 prev = gv->fastlines = thin_line(gln, decim_factor); 00260 } 00261 } 00262 else if (N_s < MFAST_LNS) { 00263 T_slength += slength[N_s++] = gv_line_length(gln); 00264 } 00265 } 00266 00267 A_slength = T_slength / N_s; 00268 N_s = 0; 00269 00270 for (gln = gv->lines; gln; gln = gln->next) { 00271 if (gln->npts <= A_ppl) { 00272 if (N_s < MFAST_LNS) { 00273 if (slength[N_s++] > A_slength) { 00274 if (prev) { 00275 prev->next = copy_line(gln); 00276 prev = prev->next; 00277 } 00278 else { 00279 prev = gv->fastlines = copy_line(gln); 00280 } 00281 } 00282 } 00283 } 00284 } 00285 00286 G_debug(3, "Decimated lines have %d points.", 00287 gln_num_points(gv->fastlines)); 00288 00289 return (1); 00290 }