main.c 4.6 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186
  1. /****************************************************************************
  2. *
  3. * MODULE: r.circle
  4. *
  5. * AUTHOR(S): Bill Brown - CERL (Jan, 1993)
  6. * Markus Neteler
  7. *
  8. * PURPOSE: Creates a raster map containing concentric rings
  9. * around a given point.
  10. *
  11. * COPYRIGHT: (C) 2006-2008 by the GRASS Development Team
  12. *
  13. * This program is free software under the GNU General Public
  14. * License (>=v2). Read the file COPYING that comes with GRASS
  15. * for details.
  16. *
  17. ***************************************************************************/
  18. #include <stdlib.h>
  19. #include <strings.h>
  20. #include <math.h>
  21. #include <grass/gis.h>
  22. #include <grass/raster.h>
  23. #include <grass/glocale.h>
  24. static double distance(double *, double *, double, double, int);
  25. #ifndef HUGE_VAL
  26. #define HUGE_VAL 1.7976931348623157e+308
  27. #endif
  28. int main(int argc, char *argv[])
  29. {
  30. struct GModule *module;
  31. struct Option *coord, *out_file, *min, *max, *mult;
  32. struct Flag *flag;
  33. int *int_buf;
  34. struct Cell_head w;
  35. struct History history;
  36. int cellfile;
  37. double east, north, pt[2], cur[2], row, col, fmult;
  38. double fmin, fmax;
  39. int binary;
  40. G_gisinit(argv[0]);
  41. module = G_define_module();
  42. G_add_keyword(_("raster"));
  43. G_add_keyword(_("buffer"));
  44. G_add_keyword(_("geometry"));
  45. G_add_keyword(_("circle"));
  46. module->description =
  47. _("Creates a raster map containing concentric "
  48. "rings around a given point.");
  49. out_file = G_define_standard_option(G_OPT_R_OUTPUT);
  50. coord = G_define_standard_option(G_OPT_M_COORDS);
  51. coord->required = YES;
  52. coord->description = _("The coordinate of the center (east,north)");
  53. min = G_define_option();
  54. min->key = "min";
  55. min->type = TYPE_DOUBLE;
  56. min->required = NO;
  57. min->description = _("Minimum radius for ring/circle map (in meters)");
  58. max = G_define_option();
  59. max->key = "max";
  60. max->type = TYPE_DOUBLE;
  61. max->required = NO;
  62. max->description = _("Maximum radius for ring/circle map (in meters)");
  63. mult = G_define_option();
  64. mult->key = "multiplier";
  65. mult->type = TYPE_DOUBLE;
  66. mult->required = NO;
  67. mult->description = _("Data value multiplier");
  68. flag = G_define_flag();
  69. flag->key = 'b';
  70. flag->description = _("Generate binary raster map");
  71. if (G_parser(argc, argv))
  72. exit(EXIT_FAILURE);
  73. G_scan_easting(coord->answers[0], &east, G_projection());
  74. G_scan_northing(coord->answers[1], &north, G_projection());
  75. pt[0] = east;
  76. pt[1] = north;
  77. fmult = 1.0;
  78. if (min->answer)
  79. sscanf(min->answer, "%lf", &fmin);
  80. else
  81. fmin = 0;
  82. if (max->answer)
  83. sscanf(max->answer, "%lf", &fmax);
  84. else
  85. fmax = HUGE_VAL;
  86. if (fmin > fmax)
  87. G_fatal_error(_("Please specify a radius in which min < max"));
  88. if (mult->answer)
  89. if (1 != sscanf(mult->answer, "%lf", &fmult))
  90. fmult = 1.0;
  91. /* nonsense test */
  92. if (flag->answer && (!min->answer && !max->answer))
  93. G_fatal_error(_("Please specify min and/or max radius when "
  94. "using the binary flag"));
  95. if (flag->answer)
  96. binary = 1; /* generate binary pattern only, useful for MASK */
  97. else
  98. binary = 0;
  99. G_get_set_window(&w);
  100. cellfile = Rast_open_c_new(out_file->answer);
  101. int_buf = (int *)G_malloc(w.cols * sizeof(int));
  102. {
  103. int c;
  104. for (row = 0; row < w.rows; row++) {
  105. G_percent(row, w.rows, 2);
  106. cur[1] = Rast_row_to_northing(row + 0.5, &w);
  107. for (col = 0; col < w.cols; col++) {
  108. c = col;
  109. cur[0] = Rast_col_to_easting(col + 0.5, &w);
  110. int_buf[c] =
  111. (int)(distance(pt, cur, fmin, fmax, binary) * fmult);
  112. if (int_buf[c] == 0)
  113. Rast_set_null_value(&int_buf[c], 1, CELL_TYPE);
  114. }
  115. Rast_put_row(cellfile, int_buf, CELL_TYPE);
  116. }
  117. }
  118. G_free(int_buf);
  119. Rast_close(cellfile);
  120. Rast_short_history(out_file->answer, "raster", &history);
  121. Rast_command_history(&history);
  122. Rast_write_history(out_file->answer, &history);
  123. G_done_msg(_("Raster map <%s> created."),
  124. out_file->answer);
  125. return (EXIT_SUCCESS);
  126. }
  127. /*******************************************************************/
  128. static double distance(double from[2], double to[2], double min, double max,
  129. int binary)
  130. {
  131. static int first = 1;
  132. double dist;
  133. if (first) {
  134. first = 0;
  135. G_begin_distance_calculations();
  136. }
  137. dist = G_distance(from[0], from[1], to[0], to[1]);
  138. if ((dist >= min) && (dist <= max))
  139. if (!binary)
  140. return dist;
  141. else
  142. return 1;
  143. else
  144. return 0; /* should be NULL ? */
  145. }
  146. /**********************************************************************/