[GRASS-SVN] r50005 - grass-addons/grass6/raster/r.viewshed

svn_grass at osgeo.org svn_grass at osgeo.org
Sat Dec 31 08:43:47 EST 2011


Author: mmetz
Date: 2011-12-31 05:43:47 -0800 (Sat, 31 Dec 2011)
New Revision: 50005

Modified:
   grass-addons/grass6/raster/r.viewshed/eventlist.cpp
   grass-addons/grass6/raster/r.viewshed/statusstructure.cpp
Log:
r.viewshed: add latlong support (backport from trunk)

Modified: grass-addons/grass6/raster/r.viewshed/eventlist.cpp
===================================================================
--- grass-addons/grass6/raster/r.viewshed/eventlist.cpp	2011-12-31 13:15:54 UTC (rev 50004)
+++ grass-addons/grass6/raster/r.viewshed/eventlist.cpp	2011-12-31 13:43:47 UTC (rev 50005)
@@ -522,13 +522,28 @@
 get_square_distance_from_viewpoint(const AEvent & a, const Viewpoint & vp)
 {
 
-    double eventy, eventx;
+    double eventy, eventx, dist;
 
     calculate_event_position(a, vp.row, vp.col, &eventy, &eventx);
 
-    double dist = (eventx - vp.col) * (eventx - vp.col) +
-	(eventy - vp.row) * (eventy - vp.row);
-    /*don't take sqrt, it is expensive; suffices for comparison */
+    if (G_projection() == PROJECTION_LL) {
+	struct Cell_head window;
+	
+	G_get_window(&window);
+	
+	dist = G_distance(G_col_to_easting(vp.col + 0.5, &window),
+			  G_row_to_northing(vp.row + 0.5, &window),
+			  G_col_to_easting(eventx + 0.5, &window),
+			  G_row_to_northing(eventy + 0.5, &window));
+
+	dist = dist * dist;
+    }
+    else {
+	/*don't take sqrt, it is expensive; suffices for comparison */
+	dist = (eventx - vp.col) * (eventx - vp.col) +
+	    (eventy - vp.row) * (eventy - vp.row);
+    }
+
     return dist;
 }
 
@@ -539,13 +554,28 @@
 					      const Viewpoint & vp)
 {
 
-    double eventy, eventx;
+    double eventy, eventx, dist;
 
     calculate_event_position(a, vp.row, vp.col, &eventy, &eventx);
-    double dist = (eventx - vp.col) * (eventx - vp.col) +
-	(eventy - vp.row) * (eventy - vp.row);
-    /*don't take sqrt, it is expensive; suffices for comparison */
 
+    if (G_projection() == PROJECTION_LL) {
+	struct Cell_head window;
+	
+	G_get_window(&window);
+	
+	dist = G_distance(G_col_to_easting(vp.col + 0.5, &window),
+			  G_row_to_northing(vp.row + 0.5, &window),
+			  G_col_to_easting(eventx + 0.5, &window),
+			  G_row_to_northing(eventy + 0.5, &window));
+
+	dist = dist * dist;
+    }
+    else {
+	/*don't take sqrt, it is expensive; suffices for comparison */
+	dist = (eventx - vp.col) * (eventx - vp.col) +
+	    (eventy - vp.row) * (eventy - vp.row);
+    }
+
     print_event(a, 2);
     G_debug(2, " pos= (%.3f. %.3f) sqdist=%.3f", eventx, eventy, dist);
 
@@ -578,7 +608,7 @@
 
 
 /* ------------------------------------------------------------ 
-   //note: this is expensive because distance is not storedin the event
+   //note: this is expensive because distance is not stored in the event
    //and must be computed on the fly */
 int DistanceCompare::compare(const AEvent & a, const AEvent & b)
 {
@@ -594,12 +624,44 @@
     double eventy, eventx;
 
     calculate_event_position(a, globalVP.row, globalVP.col, &eventy, &eventx);
-    da = (eventx - globalVP.col) * (eventx - globalVP.col) +
-	(eventy - globalVP.row) * (eventy - globalVP.row);
+    if (G_projection() == PROJECTION_LL) {
+	struct Cell_head window;
+	
+	G_get_window(&window);
+	
+	da = G_distance(G_col_to_easting(globalVP.col + 0.5, &window),
+			  G_row_to_northing(globalVP.row + 0.5, &window),
+			  G_col_to_easting(eventx + 0.5, &window),
+			  G_row_to_northing(eventy + 0.5, &window));
+
+	da = da * da;
+    }
+    else {
+	/*don't take sqrt, it is expensive; suffices for comparison */
+	da = (eventx - globalVP.col) * (eventx - globalVP.col) +
+	    (eventy - globalVP.row) * (eventy - globalVP.row);
+    }
+
+
     calculate_event_position(b, globalVP.row, globalVP.col, &eventy, &eventx);
-    db = (eventx - globalVP.col) * (eventx - globalVP.col) +
-	(eventy - globalVP.row) * (eventy - globalVP.row);
+    if (G_projection() == PROJECTION_LL) {
+	struct Cell_head window;
+	
+	G_get_window(&window);
+	
+	db = G_distance(G_col_to_easting(globalVP.col + 0.5, &window),
+			  G_row_to_northing(globalVP.row + 0.5, &window),
+			  G_col_to_easting(eventx + 0.5, &window),
+			  G_row_to_northing(eventy + 0.5, &window));
 
+	db = db * db;
+    }
+    else {
+	/*don't take sqrt, it is expensive; suffices for comparison */
+	db = (eventx - globalVP.col) * (eventx - globalVP.col) +
+	    (eventy - globalVP.row) * (eventy - globalVP.row);
+    }
+
     if (da > db) {
 	return 1;
     }

Modified: grass-addons/grass6/raster/r.viewshed/statusstructure.cpp
===================================================================
--- grass-addons/grass6/raster/r.viewshed/statusstructure.cpp	2011-12-31 13:15:54 UTC (rev 50004)
+++ grass-addons/grass6/raster/r.viewshed/statusstructure.cpp	2011-12-31 13:43:47 UTC (rev 50005)
@@ -116,11 +116,22 @@
        //sn->gradient = (sn->elev  - vp->elev)/(sn->dist2vp); */
        
     double diffElev = elev - vp->elev;
-    double dx = ((double)sn->col - vp->col) * hd.ew_res;
-    double dy = ((double)sn->row - vp->row) * hd.ns_res;
     
-    sn->dist2vp = (dx * dx) + (dy * dy);
+    if (G_projection() == PROJECTION_LL) {
+	double dist = G_distance(G_col_to_easting(sn->col + 0.5, &(hd.window)),
+				 G_row_to_northing(sn->row + 0.5, &(hd.window)),
+				 G_col_to_easting(vp->col + 0.5, &(hd.window)),
+				 G_row_to_northing(vp->row + 0.5, &(hd.window)));
 
+	sn->dist2vp = dist * dist;
+    }
+    else {
+	double dx = ((double)sn->col - vp->col) * hd.ew_res;
+	double dy = ((double)sn->row - vp->row) * hd.ns_res;
+	
+	sn->dist2vp = (dx * dx) + (dy * dy);
+    }
+
     if (diffElev == 0) {
 	sn->gradient[1] = 0;
 	return;
@@ -162,11 +173,23 @@
        //sn->gradient = (sn->elev  - vp->elev)/(sn->dist2vp); */
        
     double diffElev = elev - vp->elev;
-    double dx = (col - vp->col) * hd.ew_res;
-    double dy = (row - vp->row) * hd.ns_res;
-    double dist2vp = (dx * dx) + (dy * dy);
+    double dist2vp;
 
+    if (G_projection() == PROJECTION_LL) {
+	double dist = G_distance(G_col_to_easting(col + 0.5, &(hd.window)),
+				 G_row_to_northing(row + 0.5, &(hd.window)),
+				 G_col_to_easting(vp->col + 0.5, &(hd.window)),
+				 G_row_to_northing(vp->row + 0.5, &(hd.window)));
 
+	dist2vp = dist * dist;
+    }
+    else {
+	double dx = (col - vp->col) * hd.ew_res;
+	double dy = (row - vp->row) * hd.ns_res;
+	
+	dist2vp = (dx * dx) + (dy * dy);
+    }
+
     /* PI / 2 above, - PI / 2 below */
     sn->gradient[e_idx] = atan(diffElev / sqrt(dist2vp));
 



More information about the grass-commit mailing list