[GRASS-SVN] r30614 - grass/trunk/raster/r.los
svn_grass at osgeo.org
svn_grass at osgeo.org
Tue Mar 18 05:34:08 EDT 2008
Author: neteler
Date: 2008-03-18 05:34:07 -0400 (Tue, 18 Mar 2008)
New Revision: 30614
Modified:
grass/trunk/raster/r.los/Makefile
grass/trunk/raster/r.los/local_proto.h
grass/trunk/raster/r.los/main.c
grass/trunk/raster/r.los/make_list.c
grass/trunk/raster/r.los/point.h
grass/trunk/raster/r.los/pts_elim.c
grass/trunk/raster/r.los/segment.c
Log:
Earth curvature support added
Modified: grass/trunk/raster/r.los/Makefile
===================================================================
--- grass/trunk/raster/r.los/Makefile 2008-03-18 07:04:22 UTC (rev 30613)
+++ grass/trunk/raster/r.los/Makefile 2008-03-18 09:34:07 UTC (rev 30614)
@@ -2,8 +2,8 @@
PGM = r.los
-LIBES = $(SEGMENTLIB) $(GISLIB)
-DEPENDENCIES = $(SEGMENTDEP) $(GISDEP)
+LIBES = $(SEGMENTLIB) $(GPROJLIB) $(GISLIB)
+DEPENDENCIES = $(SEGMENTDEP) $(GPROJDEP) $(GISDEP)
include $(MODULE_TOPDIR)/include/Make/Module.make
Modified: grass/trunk/raster/r.los/local_proto.h
===================================================================
--- grass/trunk/raster/r.los/local_proto.h 2008-03-18 07:04:22 UTC (rev 30613)
+++ grass/trunk/raster/r.los/local_proto.h 2008-03-18 09:34:07 UTC (rev 30614)
@@ -1,7 +1,7 @@
/* color_rnge.c */
double decide_color_range(double, double, double);
/* newsegment.c */
-struct point *segment(int, int, int, double, double, int, int, int, int, SEGMENT *, SEGMENT *, SEGMENT *, int, int, int);
+struct point *segment(int, int, int, double, double, int, int, int, int, SEGMENT *, SEGMENT *, SEGMENT *, int, int, int, int, double);
/* pts_elim.c */
double find_orientation(int, int, int);
-double find_inclination(int, int, int, SEGMENT *, int, int);
+double find_inclination(int, int, int, SEGMENT *, int, int, int, double);
Modified: grass/trunk/raster/r.los/main.c
===================================================================
--- grass/trunk/raster/r.los/main.c 2008-03-18 07:04:22 UTC (rev 30613)
+++ grass/trunk/raster/r.los/main.c 2008-03-18 09:34:07 UTC (rev 30614)
@@ -35,6 +35,7 @@
#include <grass/gis.h>
#define MAIN
#include <grass/segment.h>
+#include <grass/gprojects.h>
#include <grass/glocale.h>
#include "cmd_line.h"
#include "point.h"
@@ -69,8 +70,10 @@
struct point *heads[16], *SEARCH_PT;
struct GModule *module;
struct Option *opt1, *opt2, *opt3, *opt5, *opt6, *opt7;
+ struct Flag *curvature;
struct History history;
char title[128];
+ double aa, e2;
G_gisinit(argv[0]);
@@ -109,9 +112,14 @@
opt6->type = TYPE_DOUBLE;
opt6->required = NO;
opt6->answer = "10000";
- opt6->options = "0-999999";
+ opt6->options = "0-5000000"; /* observer can be in a plane, too */
opt6->description = _("Maximum distance from the viewing point (meters)");
+ /* http://mintaka.sdsu.edu/GF/explain/atmos_refr/horizon.html */
+ curvature = G_define_flag();
+ curvature->key = 'c';
+ curvature->description = _("Consider earth curvature (current ellipsoid)");
+
if (G_parser(argc, argv))
exit(EXIT_FAILURE);
@@ -227,6 +235,18 @@
close(patt_fd);
}
+ if (curvature->answer){
+ /* try to get the radius the standard GRASS way from the libs */
+ G_get_ellipsoid_parameters (&aa, &e2);
+ if (aa == 0) {
+ /* since there was a problem, take a hardcoded radius :( */
+ G_warning(_("Problem to obtain current ellipsoid parameters, using sphere (6370997.0)"));
+ aa = 6370997.00;
+ }
+ G_debug(3, "radius: %f", aa);
+ }
+ G_message(_("Using maximum distance from the viewing point (meters): %f"), max_dist);
+
/* open, initialize and segment all files */
in_fd = open(in_name, 2);
segment_init(&seg_in, in_fd, 4);
@@ -304,7 +324,7 @@
slope_1, slope_2, flip, sign_on_y,
sign_on_x, viewpt_elev, &seg_in,
&seg_out, &seg_patt, row_viewpt,
- col_viewpt, patt_flag);
+ col_viewpt, patt_flag, curvature->answer, aa);
G_percent(segment_no, 16, 5);
} /* end of for-loop over segments */
Modified: grass/trunk/raster/r.los/make_list.c
===================================================================
--- grass/trunk/raster/r.los/make_list.c 2008-03-18 07:04:22 UTC (rev 30613)
+++ grass/trunk/raster/r.los/make_list.c 2008-03-18 09:34:07 UTC (rev 30614)
@@ -18,7 +18,8 @@
struct point *make_list(struct point *head, int y, int x,
SEGMENT * seg_in_p, int viewpt_elev,
- int quadrant, int row_viewpt, int col_viewpt)
+ int quadrant, int row_viewpt, int col_viewpt,
+ int docurv, double ellps_a)
{
double del_x, del_y, dist, orientation, inclination;
static struct point *PRESENT_PT;
@@ -38,7 +39,7 @@
/* otherwise find orientation and inclination */
orientation = find_orientation(x, y, quadrant);
inclination = find_inclination(x, y, viewpt_elev, seg_in_p,
- row_viewpt, col_viewpt);
+ row_viewpt, col_viewpt, docurv, ellps_a);
if (head == NULL) { /* first point ? */
head = make_point(orientation, inclination, y, x);
Modified: grass/trunk/raster/r.los/point.h
===================================================================
--- grass/trunk/raster/r.los/point.h 2008-03-18 07:04:22 UTC (rev 30613)
+++ grass/trunk/raster/r.los/point.h 2008-03-18 09:34:07 UTC (rev 30614)
@@ -31,14 +31,14 @@
/* delete.c */
struct point *delete(struct point *, struct point *, SEGMENT *, int, int);
/* make_list.c */
-struct point *make_list(struct point *, int, int, SEGMENT *, int, int, int, int);
+struct point *make_list(struct point *, int, int, SEGMENT *, int, int, int, int, int, double);
/* mark_pts.c */
int mark_visible_points(struct point *, SEGMENT *, int, int, double, double);
/* pts_elim.c */
-struct point *hidden_point_elimination(struct point *, int, SEGMENT *, SEGMENT *, SEGMENT *, int, int, int, int, int, int);
+struct point *hidden_point_elimination(struct point *, int, SEGMENT *, SEGMENT *, SEGMENT *, int, int, int, int, int, int, int, double);
/* segment.c */
struct point *segment(int, int, int, double, double,
- int, int, int, int, SEGMENT *, SEGMENT *, SEGMENT *, int, int, int);
+ int, int, int, int, SEGMENT *, SEGMENT *, SEGMENT *, int, int, int, int, double);
#endif
/****************************************************************/
Modified: grass/trunk/raster/r.los/pts_elim.c
===================================================================
--- grass/trunk/raster/r.los/pts_elim.c 2008-03-18 07:04:22 UTC (rev 30613)
+++ grass/trunk/raster/r.los/pts_elim.c 2008-03-18 09:34:07 UTC (rev 30614)
@@ -11,6 +11,7 @@
#include <math.h>
#include <grass/gis.h>
#include <grass/segment.h>
+#include <grass/glocale.h>
#include "point.h"
#include "radians.h"
#include "local_proto.h"
@@ -42,7 +43,8 @@
struct point *head, int viewpt_elev,
SEGMENT *seg_in_p, SEGMENT *seg_out_p, SEGMENT *seg_patt_p,
int quadrant, int sign_on_y, int sign_on_x,
- int row_viewpt, int col_viewpt, int patt_flag)
+ int row_viewpt, int col_viewpt, int patt_flag,
+ int docurv, double ellps_a)
{
struct point *CHECKED_PT, *BLOCKING_PT;
@@ -92,11 +94,11 @@
inclination_neighbor_1 =
find_inclination(neighbor_1_x,neighbor_1_y,viewpt_elev,
- seg_in_p,row_viewpt,col_viewpt);
+ seg_in_p,row_viewpt,col_viewpt, docurv, ellps_a);
inclination_neighbor_2 =
find_inclination(neighbor_2_x,neighbor_2_y,viewpt_elev,
- seg_in_p,row_viewpt,col_viewpt);
+ seg_in_p,row_viewpt,col_viewpt, docurv, ellps_a);
/* check all points behind the blocking point */
@@ -278,11 +280,12 @@
/****************************************************************/
double
-find_inclination (int x, int y, int viewpt_elev, SEGMENT *seg_in_p, int row_viewpt, int col_viewpt)
+find_inclination (int x, int y, int viewpt_elev, SEGMENT *seg_in_p, int row_viewpt, int col_viewpt,
+ int docurv, double ellps_a)
{
- double del_x, del_y,dist,atan(),sqrt();
+ double del_x, del_y,dist;
int abs();
FCELL picked_pt_elev;
extern struct Cell_head window;
@@ -293,6 +296,10 @@
dist=sqrt(del_x * del_x + del_y * del_y)*window.ns_res;
segment_get(seg_in_p,&picked_pt_elev,row_viewpt-y,x+col_viewpt);
+
+ if (docurv) /* decrease height of target point */
+ picked_pt_elev = picked_pt_elev - ((dist*dist) / (2 * ellps_a));
+
return(atan((picked_pt_elev-viewpt_elev)/dist));
}
Modified: grass/trunk/raster/r.los/segment.c
===================================================================
--- grass/trunk/raster/r.los/segment.c 2008-03-18 07:04:22 UTC (rev 30613)
+++ grass/trunk/raster/r.los/segment.c 2008-03-18 09:34:07 UTC (rev 30614)
@@ -19,7 +19,7 @@
int sign_on_y, int sign_on_x, int viewpt_elev,
SEGMENT * seg_in_p, SEGMENT * seg_out_p,
SEGMENT * seg_patt_p, int row_viewpt, int col_viewpt,
- int patt_flag)
+ int patt_flag, int docurv, double ellps_a)
{
int lower_limit_y, upper_limit_y, less, x, y,
x_actual, y_actual, x_flip, y_flip;
@@ -60,7 +60,7 @@
/* add chosen point to the point list */
head = make_list(head, y_actual, x_actual, seg_in_p,
- viewpt_elev, quadrant, row_viewpt, col_viewpt);
+ viewpt_elev, quadrant, row_viewpt, col_viewpt, docurv, ellps_a);
}
} /* end of outer loop */
@@ -79,7 +79,7 @@
head = hidden_point_elimination(head, viewpt_elev,
seg_in_p, seg_out_p, seg_patt_p,
quadrant, sign_on_y, sign_on_x,
- row_viewpt, col_viewpt, patt_flag);
+ row_viewpt, col_viewpt, patt_flag, docurv, ellps_a);
}
return (head);
More information about the grass-commit
mailing list