GRASS Programmer's Manual
6.4.2(2012)
|
00001 00002 /*- 00003 * Written by H. Mitasova, I. Kosinovsky, D. Gerdes Summer 1992 00004 * Copyright 1992, H. Mitasova 00005 * I. Kosinovsky, and D.Gerdes 00006 * 00007 * modified by McCauley in August 1995 00008 * modified by Mitasova in August 1995 00009 * modified by Mitasova in August 1999 (fix for elev color) 00010 * modified by Brown in September 1999 (fix for Timestamps) 00011 * Modified by Mitasova in Nov. 1999 (write given tension into hist) 00012 * Last modified: 2006-12-13 00013 * 00014 */ 00015 00016 #include <stdio.h> 00017 #include <math.h> 00018 #include <grass/gis.h> 00019 #include <grass/bitmap.h> 00020 #include <grass/linkm.h> 00021 00022 #include <grass/glocale.h> 00023 #include <grass/interpf.h> 00024 00025 00026 #define MULT 100000 00027 00028 00029 int IL_output_2d(struct interp_params *params, struct Cell_head *cellhd, /* current region */ 00030 double zmin, double zmax, /* min,max input z-values */ 00031 double zminac, double zmaxac, double c1min, double c1max, /* min,max interpolated values */ 00032 double c2min, double c2max, double gmin, double gmax, double ertot, /* total interplating func. error */ 00033 char *input, /* input file name */ 00034 double dnorm, int dtens, int vect, int n_points) 00035 00036 /* 00037 * Creates output files as well as history files and color tables for 00038 * them. 00039 */ 00040 { 00041 FCELL *cell1; 00042 int cf1 = 0, cf2 = 0, cf3 = 0, cf4 = 0, cf5 = 0, cf6 = 0; 00043 int nrows, ncols; 00044 int i, ii; 00045 double zstep; 00046 FCELL data1, data2; 00047 struct Colors colors; 00048 struct History hist, hist1, hist2, hist3, hist4, hist5; 00049 char *type; 00050 char *mapset = NULL; 00051 int cond1, cond2; 00052 FCELL dat1, dat2; 00053 00054 cond2 = ((params->pcurv != NULL) || (params->tcurv != NULL) 00055 || (params->mcurv != NULL)); 00056 cond1 = ((params->slope != NULL) || (params->aspect != NULL) || cond2); 00057 00058 cell1 = G_allocate_f_raster_buf(); 00059 00060 /* 00061 * G_set_embedded_null_value_mode(1); 00062 */ 00063 if (params->elev != NULL) { 00064 cf1 = G_open_fp_cell_new(params->elev); 00065 if (cf1 < 0) { 00066 fprintf(stderr, "unable to create raster map %s\n", params->elev); 00067 return -1; 00068 } 00069 } 00070 00071 if (params->slope != NULL) { 00072 cf2 = G_open_fp_cell_new(params->slope); 00073 if (cf2 < 0) { 00074 fprintf(stderr, "unable to create raster map %s\n", 00075 params->slope); 00076 return -1; 00077 } 00078 } 00079 00080 if (params->aspect != NULL) { 00081 cf3 = G_open_fp_cell_new(params->aspect); 00082 if (cf3 < 0) { 00083 fprintf(stderr, "unable to create raster map %s\n", 00084 params->aspect); 00085 return -1; 00086 } 00087 } 00088 00089 if (params->pcurv != NULL) { 00090 cf4 = G_open_fp_cell_new(params->pcurv); 00091 if (cf4 < 0) { 00092 fprintf(stderr, "unable to create raster map %s\n", 00093 params->pcurv); 00094 return -1; 00095 } 00096 } 00097 00098 if (params->tcurv != NULL) { 00099 cf5 = G_open_fp_cell_new(params->tcurv); 00100 if (cf5 < 0) { 00101 fprintf(stderr, "unable to create raster map %s\n", 00102 params->tcurv); 00103 return -1; 00104 } 00105 } 00106 00107 if (params->mcurv != NULL) { 00108 cf6 = G_open_fp_cell_new(params->mcurv); 00109 if (cf6 < 0) { 00110 fprintf(stderr, "unable to create raster map %s\n", 00111 params->mcurv); 00112 return -1; 00113 } 00114 } 00115 00116 nrows = cellhd->rows; 00117 if (nrows != params->nsizr) { 00118 fprintf(stderr, "first change your rows number to nsizr! %d %d\n", 00119 nrows, params->nsizr); 00120 return -1; 00121 } 00122 00123 ncols = cellhd->cols; 00124 if (ncols != params->nsizc) { 00125 G_warning(_("First change your cols number to nsizc %d %d"), 00126 ncols, params->nsizc); 00127 return -1; 00128 } 00129 00130 if (G_set_window(cellhd) < 0) 00131 return -1; 00132 00133 if (nrows != G_window_rows()) { 00134 fprintf(stderr, "OOPS: rows changed from %d to %d\n", nrows, 00135 G_window_rows()); 00136 return -1; 00137 } 00138 00139 if (ncols != G_window_cols()) { 00140 fprintf(stderr, "OOPS: cols changed from %d to %d\n", ncols, 00141 G_window_cols()); 00142 return -1; 00143 } 00144 00145 if (params->elev != NULL) { 00146 fseek(params->Tmp_fd_z, 0L, 0); /* seek to the beginning */ 00147 for (i = 0; i < params->nsizr; i++) { 00148 /* seek to the right row */ 00149 if (fseek(params->Tmp_fd_z, (long) 00150 ((params->nsizr - 1 - 00151 i) * params->nsizc * sizeof(FCELL)), 0) 00152 == -1) { 00153 fprintf(stderr, "cannot fseek to the right spot\n"); 00154 return -1; 00155 } 00156 ii = fread(cell1, sizeof(FCELL), params->nsizc, params->Tmp_fd_z); 00157 /* 00158 * for(j=0;j<params->nsizc;j++) fprintf(stderr,"%f ",cell1[j]); 00159 * fprintf(stderr,"\n"); 00160 */ 00161 G_put_f_raster_row(cf1, cell1); 00162 00163 } 00164 } 00165 00166 if (params->slope != NULL) { 00167 fseek(params->Tmp_fd_dx, 0L, 0); /* seek to the beginning */ 00168 for (i = 0; i < params->nsizr; i++) { 00169 /* seek to the right row */ 00170 if (fseek(params->Tmp_fd_dx, (long) 00171 ((params->nsizr - 1 - 00172 i) * params->nsizc * sizeof(FCELL)), 0) 00173 == -1) { 00174 fprintf(stderr, "cannot fseek to the right spot\n"); 00175 return -1; 00176 } 00177 fread(cell1, sizeof(FCELL), params->nsizc, params->Tmp_fd_dx); 00178 G_put_f_raster_row(cf2, cell1); 00179 } 00180 } 00181 00182 if (params->aspect != NULL) { 00183 fseek(params->Tmp_fd_dy, 0L, 0); /* seek to the beginning */ 00184 for (i = 0; i < params->nsizr; i++) { 00185 /* seek to the right row */ 00186 if (fseek(params->Tmp_fd_dy, (long) 00187 ((params->nsizr - 1 - 00188 i) * params->nsizc * sizeof(FCELL)), 0) 00189 == -1) { 00190 fprintf(stderr, "cannot fseek to the right spot\n"); 00191 return -1; 00192 } 00193 fread(cell1, sizeof(FCELL), params->nsizc, params->Tmp_fd_dy); 00194 G_put_f_raster_row(cf3, cell1); 00195 } 00196 } 00197 00198 if (params->pcurv != NULL) { 00199 fseek(params->Tmp_fd_xx, 0L, 0); /* seek to the beginning */ 00200 for (i = 0; i < params->nsizr; i++) { 00201 /* seek to the right row */ 00202 if (fseek(params->Tmp_fd_xx, (long) 00203 ((params->nsizr - 1 - 00204 i) * params->nsizc * sizeof(FCELL)), 0) 00205 == -1) { 00206 fprintf(stderr, "cannot fseek to the right spot\n"); 00207 return -1; 00208 } 00209 fread(cell1, sizeof(FCELL), params->nsizc, params->Tmp_fd_xx); 00210 G_put_f_raster_row(cf4, cell1); 00211 } 00212 } 00213 00214 if (params->tcurv != NULL) { 00215 fseek(params->Tmp_fd_yy, 0L, 0); /* seek to the beginning */ 00216 for (i = 0; i < params->nsizr; i++) { 00217 /* seek to the right row */ 00218 if (fseek(params->Tmp_fd_yy, (long) 00219 ((params->nsizr - 1 - 00220 i) * params->nsizc * sizeof(FCELL)), 0) 00221 == -1) { 00222 fprintf(stderr, "cannot fseek to the right spot\n"); 00223 return -1; 00224 } 00225 fread(cell1, sizeof(FCELL), params->nsizc, params->Tmp_fd_yy); 00226 G_put_f_raster_row(cf5, cell1); 00227 } 00228 } 00229 00230 if (params->mcurv != NULL) { 00231 fseek(params->Tmp_fd_xy, 0L, 0); /* seek to the beginning */ 00232 for (i = 0; i < params->nsizr; i++) { 00233 /* seek to the right row */ 00234 if (fseek(params->Tmp_fd_xy, (long) 00235 ((params->nsizr - 1 - 00236 i) * params->nsizc * sizeof(FCELL)), 0) 00237 == -1) { 00238 fprintf(stderr, "cannot fseek to the right spot\n"); 00239 return -1; 00240 } 00241 fread(cell1, sizeof(FCELL), params->nsizc, params->Tmp_fd_xy); 00242 G_put_f_raster_row(cf6, cell1); 00243 } 00244 } 00245 00246 if (cf1) 00247 G_close_cell(cf1); 00248 if (cf2) 00249 G_close_cell(cf2); 00250 if (cf3) 00251 G_close_cell(cf3); 00252 if (cf4) 00253 G_close_cell(cf4); 00254 if (cf5) 00255 G_close_cell(cf5); 00256 if (cf6) 00257 G_close_cell(cf6); 00258 00259 00260 /* colortable for elevations */ 00261 G_init_colors(&colors); 00262 zstep = (FCELL) (zmaxac - zminac) / 5.; 00263 for (i = 1; i <= 5; i++) { 00264 data1 = (FCELL) (zminac + (i - 1) * zstep); 00265 data2 = (FCELL) (zminac + i * zstep); 00266 switch (i) { 00267 case 1: 00268 G_add_f_raster_color_rule(&data1, 0, 191, 191, 00269 &data2, 0, 255, 0, &colors); 00270 break; 00271 case 2: 00272 G_add_f_raster_color_rule(&data1, 0, 255, 0, 00273 &data2, 255, 255, 0, &colors); 00274 break; 00275 case 3: 00276 G_add_f_raster_color_rule(&data1, 255, 255, 0, 00277 &data2, 255, 127, 0, &colors); 00278 break; 00279 case 4: 00280 G_add_f_raster_color_rule(&data1, 255, 127, 0, 00281 &data2, 191, 127, 63, &colors); 00282 break; 00283 case 5: 00284 G_add_f_raster_color_rule(&data1, 191, 127, 63, 00285 &data2, 20, 20, 20, &colors); 00286 break; 00287 } 00288 } 00289 00290 if (params->elev != NULL) { 00291 mapset = G_find_file("cell", params->elev, ""); 00292 if (mapset == NULL) { 00293 fprintf(stderr, "file [%s] not found\n", params->elev); 00294 return -1; 00295 } 00296 G_write_colors(params->elev, mapset, &colors); 00297 G_quantize_fp_map_range(params->elev, mapset, 00298 (DCELL) zminac - 0.5, (DCELL) zmaxac + 0.5, 00299 (CELL) (zminac - 0.5), (CELL) (zmaxac + 0.5)); 00300 } 00301 00302 /* colortable for slopes */ 00303 if (cond1) { 00304 if (!params->deriv) { 00305 /* 00306 * smin = (CELL) ((int)(gmin*scig)); smax = (CELL) gmax; fprintf 00307 * (stderr, "min %d max %d \n", smin,smax); G_make_rainbow_colors 00308 * (&colors,smin,smax); 00309 */ 00310 G_init_colors(&colors); 00311 G_add_color_rule(0, 255, 255, 255, 2, 255, 255, 0, &colors); 00312 G_add_color_rule(2, 255, 255, 0, 5, 0, 255, 0, &colors); 00313 G_add_color_rule(5, 0, 255, 0, 10, 0, 255, 255, &colors); 00314 G_add_color_rule(10, 0, 255, 255, 15, 0, 0, 255, &colors); 00315 G_add_color_rule(15, 0, 0, 255, 30, 255, 0, 255, &colors); 00316 G_add_color_rule(30, 255, 0, 255, 50, 255, 0, 0, &colors); 00317 G_add_color_rule(50, 255, 0, 0, 90, 0, 0, 0, &colors); 00318 } 00319 else { 00320 G_init_colors(&colors); 00321 dat1 = (FCELL) - 5.0; /* replace by min dx, amin1 (c1min, 00322 * c2min); */ 00323 dat2 = (FCELL) - 0.1; 00324 G_add_f_raster_color_rule(&dat1, 127, 0, 255, 00325 &dat2, 0, 0, 255, &colors); 00326 dat1 = dat2; 00327 dat2 = (FCELL) - 0.01; 00328 G_add_f_raster_color_rule(&dat1, 0, 0, 255, 00329 &dat2, 0, 127, 255, &colors); 00330 dat1 = dat2; 00331 dat2 = (FCELL) - 0.001; 00332 G_add_f_raster_color_rule(&dat1, 0, 127, 255, 00333 &dat2, 0, 255, 255, &colors); 00334 dat1 = dat2; 00335 dat2 = (FCELL) 0.0; 00336 G_add_f_raster_color_rule(&dat1, 0, 255, 255, 00337 &dat2, 200, 255, 200, &colors); 00338 dat1 = dat2; 00339 dat2 = (FCELL) 0.001; 00340 G_add_f_raster_color_rule(&dat1, 200, 255, 200, 00341 &dat2, 255, 255, 0, &colors); 00342 dat1 = dat2; 00343 dat2 = (FCELL) 0.01; 00344 G_add_f_raster_color_rule(&dat1, 255, 255, 0, 00345 &dat2, 255, 127, 0, &colors); 00346 dat1 = dat2; 00347 dat2 = (FCELL) 0.1; 00348 G_add_f_raster_color_rule(&dat1, 255, 127, 0, 00349 &dat2, 255, 0, 0, &colors); 00350 dat1 = dat2; 00351 dat2 = (FCELL) 5.0; /* replace by max dx, amax1 (c1max, 00352 * c2max); */ 00353 G_add_f_raster_color_rule(&dat1, 255, 0, 0, 00354 &dat2, 255, 0, 200, &colors); 00355 } 00356 00357 if (params->slope != NULL) { 00358 mapset = G_find_file("cell", params->slope, ""); 00359 if (mapset == NULL) { 00360 fprintf(stderr, "file [%s] not found\n", params->slope); 00361 return -1; 00362 } 00363 G_write_colors(params->slope, mapset, &colors); 00364 G_quantize_fp_map_range(params->slope, mapset, 0., 90., 0, 90); 00365 00366 type = "raster"; 00367 G_short_history(params->slope, type, &hist1); 00368 if (params->elev != NULL) 00369 sprintf(hist1.edhist[0], "The elevation map is %s", 00370 params->elev); 00371 if (vect) 00372 sprintf(hist1.datsrc_1, "vector map %s", input); 00373 else 00374 sprintf(hist1.datsrc_1, "site file %s", input); 00375 hist1.edlinecnt = 1; 00376 00377 G_command_history(&hist1); 00378 G_write_history(params->slope, &hist1); 00379 if (params->ts) 00380 G_write_raster_timestamp(params->slope, params->ts); 00381 00382 } 00383 00384 /* colortable for aspect */ 00385 if (!params->deriv) { 00386 G_init_colors(&colors); 00387 G_add_color_rule(0, 255, 255, 255, 0, 255, 255, 255, &colors); 00388 G_add_color_rule(1, 255, 255, 0, 90, 0, 255, 0, &colors); 00389 G_add_color_rule(90, 0, 255, 0, 180, 0, 255, 255, &colors); 00390 G_add_color_rule(180, 0, 255, 255, 270, 255, 0, 0, &colors); 00391 G_add_color_rule(270, 255, 0, 0, 360, 255, 255, 0, &colors); 00392 } 00393 else { 00394 G_init_colors(&colors); 00395 dat1 = (FCELL) - 5.0; /* replace by min dy, amin1 (c1min, 00396 * c2min); */ 00397 dat2 = (FCELL) - 0.1; 00398 G_add_f_raster_color_rule(&dat1, 127, 0, 255, 00399 &dat2, 0, 0, 255, &colors); 00400 dat1 = dat2; 00401 dat2 = (FCELL) - 0.01; 00402 G_add_f_raster_color_rule(&dat1, 0, 0, 255, 00403 &dat2, 0, 127, 255, &colors); 00404 dat1 = dat2; 00405 dat2 = (FCELL) - 0.001; 00406 G_add_f_raster_color_rule(&dat1, 0, 127, 255, 00407 &dat2, 0, 255, 255, &colors); 00408 dat1 = dat2; 00409 dat2 = (FCELL) 0.0; 00410 G_add_f_raster_color_rule(&dat1, 0, 255, 255, 00411 &dat2, 200, 255, 200, &colors); 00412 dat1 = dat2; 00413 dat2 = (FCELL) 0.001; 00414 G_add_f_raster_color_rule(&dat1, 200, 255, 200, 00415 &dat2, 255, 255, 0, &colors); 00416 dat1 = dat2; 00417 dat2 = (FCELL) 0.01; 00418 G_add_f_raster_color_rule(&dat1, 255, 255, 0, 00419 &dat2, 255, 127, 0, &colors); 00420 dat1 = dat2; 00421 dat2 = (FCELL) 0.1; 00422 G_add_f_raster_color_rule(&dat1, 255, 127, 0, 00423 &dat2, 255, 0, 0, &colors); 00424 dat1 = dat2; 00425 dat2 = (FCELL) 5.0; /* replace by max dy, amax1 (c1max, 00426 * c2max); */ 00427 G_add_f_raster_color_rule(&dat1, 255, 0, 0, 00428 &dat2, 255, 0, 200, &colors); 00429 } 00430 00431 if (params->aspect != NULL) { 00432 mapset = G_find_file("cell", params->aspect, ""); 00433 if (mapset == NULL) { 00434 fprintf(stderr, "file [%s] not found\n", params->aspect); 00435 return -1; 00436 } 00437 G_write_colors(params->aspect, mapset, &colors); 00438 G_quantize_fp_map_range(params->aspect, mapset, 0., 360., 0, 360); 00439 00440 type = "raster"; 00441 G_short_history(params->aspect, type, &hist2); 00442 if (params->elev != NULL) 00443 sprintf(hist2.edhist[0], "The elevation map is %s", 00444 params->elev); 00445 if (vect) 00446 sprintf(hist2.datsrc_1, "vector map %s", input); 00447 else 00448 sprintf(hist2.datsrc_1, "site file %s", input); 00449 hist2.edlinecnt = 1; 00450 00451 G_command_history(&hist2); 00452 G_write_history(params->aspect, &hist2); 00453 if (params->ts) 00454 G_write_raster_timestamp(params->aspect, params->ts); 00455 } 00456 00457 /* colortable for curvatures */ 00458 if (cond2) { 00459 G_init_colors(&colors); 00460 dat1 = (FCELL) amin1(c1min, c2min); /* for derivatives use min 00461 * dxx,dyy,dxy */ 00462 dat2 = (FCELL) - 0.01; 00463 G_add_f_raster_color_rule(&dat1, 127, 0, 255, 00464 &dat2, 0, 0, 255, &colors); 00465 dat1 = dat2; 00466 dat2 = (FCELL) - 0.001; 00467 G_add_f_raster_color_rule(&dat1, 0, 0, 255, 00468 &dat2, 0, 127, 255, &colors); 00469 dat1 = dat2; 00470 dat2 = (FCELL) - 0.00001; 00471 G_add_f_raster_color_rule(&dat1, 0, 127, 255, 00472 &dat2, 0, 255, 255, &colors); 00473 dat1 = dat2; 00474 dat2 = (FCELL) 0.0; 00475 G_add_f_raster_color_rule(&dat1, 0, 255, 255, 00476 &dat2, 200, 255, 200, &colors); 00477 dat1 = dat2; 00478 dat2 = (FCELL) 0.00001; 00479 G_add_f_raster_color_rule(&dat1, 200, 255, 200, 00480 &dat2, 255, 255, 0, &colors); 00481 dat1 = dat2; 00482 dat2 = (FCELL) 0.001; 00483 G_add_f_raster_color_rule(&dat1, 255, 255, 0, 00484 &dat2, 255, 127, 0, &colors); 00485 dat1 = dat2; 00486 dat2 = (FCELL) 0.01; 00487 G_add_f_raster_color_rule(&dat1, 255, 127, 0, 00488 &dat2, 255, 0, 0, &colors); 00489 dat1 = dat2; 00490 dat2 = (FCELL) amax1(c1max, c2max); /* for derivatives use max 00491 * dxx,dyy,dxy */ 00492 G_add_f_raster_color_rule(&dat1, 255, 0, 0, 00493 &dat2, 255, 0, 200, &colors); 00494 00495 if (params->pcurv != NULL) { 00496 mapset = G_find_file("cell", params->pcurv, ""); 00497 if (mapset == NULL) { 00498 fprintf(stderr, "file [%s] not found\n", params->pcurv); 00499 return -1; 00500 } 00501 G_write_colors(params->pcurv, mapset, &colors); 00502 G_quantize_fp_map_range(params->pcurv, mapset, dat1, dat2, 00503 (CELL) (dat1 * MULT), 00504 (CELL) (dat2 * MULT)); 00505 00506 type = "raster"; 00507 G_short_history(params->pcurv, type, &hist3); 00508 if (params->elev != NULL) 00509 sprintf(hist3.edhist[0], "The elevation map is %s", 00510 params->elev); 00511 if (vect) 00512 sprintf(hist3.datsrc_1, "vector map %s", input); 00513 else 00514 sprintf(hist3.datsrc_1, "site file %s", input); 00515 hist3.edlinecnt = 1; 00516 00517 G_command_history(&hist3); 00518 G_write_history(params->pcurv, &hist3); 00519 if (params->ts) 00520 G_write_raster_timestamp(params->pcurv, params->ts); 00521 } 00522 00523 if (params->tcurv != NULL) { 00524 mapset = G_find_file("cell", params->tcurv, ""); 00525 if (mapset == NULL) { 00526 fprintf(stderr, "file [%s] not found\n", params->tcurv); 00527 return -1; 00528 } 00529 G_write_colors(params->tcurv, mapset, &colors); 00530 G_quantize_fp_map_range(params->tcurv, mapset, dat1, dat2, 00531 (CELL) (dat1 * MULT), 00532 (CELL) (dat2 * MULT)); 00533 00534 type = "raster"; 00535 G_short_history(params->tcurv, type, &hist4); 00536 if (params->elev != NULL) 00537 sprintf(hist4.edhist[0], "The elevation map is %s", 00538 params->elev); 00539 if (vect) 00540 sprintf(hist4.datsrc_1, "vector map %s", input); 00541 else 00542 sprintf(hist4.datsrc_1, "site file %s", input); 00543 hist4.edlinecnt = 1; 00544 00545 G_command_history(&hist4); 00546 G_write_history(params->tcurv, &hist4); 00547 if (params->ts) 00548 G_write_raster_timestamp(params->tcurv, params->ts); 00549 } 00550 00551 if (params->mcurv != NULL) { 00552 mapset = G_find_file("cell", params->mcurv, ""); 00553 if (mapset == NULL) { 00554 fprintf(stderr, "file [%s] not found\n", params->mcurv); 00555 return -1; 00556 } 00557 G_write_colors(params->mcurv, mapset, &colors); 00558 G_quantize_fp_map_range(params->mcurv, mapset, dat1, dat2, 00559 (CELL) (dat1 * MULT), 00560 (CELL) (dat2 * MULT)); 00561 00562 type = "raster"; 00563 G_short_history(params->mcurv, type, &hist5); 00564 if (params->elev != NULL) 00565 sprintf(hist5.edhist[0], "The elevation map is %s", 00566 params->elev); 00567 if (vect) 00568 sprintf(hist5.datsrc_1, "vector map %s", input); 00569 else 00570 sprintf(hist5.datsrc_1, "site file %s", input); 00571 hist5.edlinecnt = 1; 00572 00573 G_command_history(&hist5); 00574 G_write_history(params->mcurv, &hist5); 00575 if (params->ts) 00576 G_write_raster_timestamp(params->mcurv, params->ts); 00577 } 00578 } 00579 } 00580 00581 if (params->elev != NULL) { 00582 mapset = G_find_file("cell", params->elev, ""); 00583 if (mapset == NULL) { 00584 fprintf(stderr, "file [%s] not found\n", params->elev); 00585 return -1; 00586 } 00587 type = "raster"; 00588 G_short_history(params->elev, type, &hist); 00589 00590 params->dmin = sqrt(params->dmin); 00591 00592 /* 00593 * sprintf (hist.edhist[0], "tension=%f, smoothing=%f", params->fi * 00594 * dnorm / 1000., params->rsm); 00595 */ 00596 00597 if (dtens) { 00598 if (params->rsm == -1) 00599 sprintf(hist.edhist[0], "giventension=%f, smoothing att=%d", 00600 params->fi * 1000. / dnorm, params->smatt); 00601 else 00602 sprintf(hist.edhist[0], "giventension=%f, smoothing=%f", 00603 params->fi * 1000. / dnorm, params->rsm); 00604 } 00605 else { 00606 if (params->rsm == -1) 00607 sprintf(hist.edhist[0], "tension=%f, smoothing att=%d", 00608 params->fi * 1000. / dnorm, params->smatt); 00609 else 00610 sprintf(hist.edhist[0], "tension=%f, smoothing=%f", 00611 params->fi, params->rsm); 00612 } 00613 00614 sprintf(hist.edhist[1], "dnorm=%f, dmin=%f, zmult=%f", 00615 dnorm, params->dmin, params->zmult); 00616 /* 00617 * sprintf(hist.edhist[2], "segmax=%d, npmin=%d, errtotal=%f", 00618 * params->kmax,params->kmin,ertot); 00619 */ 00620 /* 00621 * sprintf (hist.edhist[2], "segmax=%d, npmin=%d, errtotal =%f", 00622 * params->kmax, params->kmin, sqrt (ertot) / n_points); 00623 */ 00624 00625 sprintf(hist.edhist[2], "segmax=%d, npmin=%d, rmsdevi=%f", 00626 params->kmax, params->kmin, sqrt(ertot / n_points)); 00627 00628 sprintf(hist.edhist[3], "zmin_data=%f, zmax_data=%f", zmin, zmax); 00629 sprintf(hist.edhist[4], "zmin_int=%f, zmax_int=%f", zminac, zmaxac); 00630 00631 if ((params->theta) && (params->scalex)) { 00632 sprintf(hist.edhist[5], "theta=%f, scalex=%f", params->theta, 00633 params->scalex); 00634 hist.edlinecnt = 6; 00635 } 00636 else 00637 hist.edlinecnt = 5; 00638 00639 if (vect) 00640 sprintf(hist.datsrc_1, "vector map %s", input); 00641 else 00642 sprintf(hist.datsrc_1, "site file %s", input); 00643 00644 G_command_history(&hist); 00645 G_write_history(params->elev, &hist); 00646 if (params->ts) 00647 G_write_raster_timestamp(params->elev, params->ts); 00648 } 00649 00650 /* 00651 * if (title) G_put_cell_title (output, title); 00652 */ 00653 return 1; 00654 }