Hello,
missing the possibility to measure areas in the grass tcltk map display I
created a patch which adds such a function.
Could somebody test this patch and perhaps even merge it into the grass
repository?
Thank you very much,
Jonas
--
Neu: GMX FreeDSL Komplettanschluss mit DSL 6.000 Flatrate + Telefonanschluss
für nur 17,95 Euro/mtl.!*
http://dslspecial.gmx.de/freedsl-surfflat/?ac=OM.AD.PD003K11308T4569a
--- grass/grass/branches/releasebranch_6_4/gui/tcltk/gis.m/mapcanvas.tcl 2009-04-18 16:08:16.000000000 +0200
+++ mapcanvas.tcl-with-marea 2009-04-18 16:51:23.000000000 +0200
@@ -65,7 +65,12 @@
variable liney1
variable linex2
variable liney2
-
+ variable marea
+ variable close_marea
+ variable tot_marea
+ variable north_start
+ variable east_start
+
# There is a global coords # Text to display in indicator widget, indexed by mon
# Process ID for temp files
@@ -1601,6 +1606,11 @@
variable measurement_annotation_handle
variable mlength
variable totmlength
+ variable marea
+ variable close_marea
+ variable tot_marea
+ variable north_start
+ variable east_start
variable linex1
variable liney1
variable linex2
@@ -1635,6 +1645,11 @@
MapCanvas::setcursor $mon "pencil"
set mlength 0
set totmlength 0
+ set marea 0
+ set close_marea 0
+ set tot_marea 0
+ set north_start 0
+ set east_start 0
}
@@ -1679,12 +1694,17 @@
}
}
-# measure line length
+# measure line length and area
proc MapCanvas::measure { mon x y } {
variable can
variable measurement_annotation_handle
variable mlength
variable totmlength
+ variable marea
+ variable close_marea
+ variable tot_marea
+ variable north_start
+ variable east_start
variable linex1
variable liney1
variable linex2
@@ -1714,15 +1734,41 @@
set mlength [expr {sqrt(pow(($east1 - $east2), 2) + pow(($north1 - $north2), 2))}]
set totmlength [expr {$totmlength + $mlength}]
+
+ # begin area calculation
+ # formula:
+ # 2A = | Sum ( ( Yi + Yi+1 ) x ( Xi - Xi+1 ) ) |
+ # save coordinates of first point of area measurement:
+ if { ($north_start == 0) && ($east_start == 0) && ($marea == 0) && \
+ ($close_marea == 0) && ($tot_marea == 0)} {
+ set north_start $north1
+ set east_start $east1
+ }
+
+ #calculate and cumulate subareas
+ set marea [expr {($north1 + $north2)*($east1 - $east2) + $marea}]
+
+ #calculate last subarea (added by virtual polygon segment from last point to start point)
+ set close_marea [expr {($north2 + $north_start)*($east2 - $east_start) }]
+
+ #calculate result of area measurement
+ set tot_marea [expr { abs($marea + $close_marea) / 2}]
+
+ # end area calculation
+
+
# format length numbers and units in a nice way
set out_seg [ fmt_length $mlength ]
set out_tot [ fmt_length $totmlength ]
+ set out_area [ fmt_area $tot_marea ]
monitor_annotate $measurement_annotation_handle \
[G_msg " --segment length = $out_seg\n"]
monitor_annotate $measurement_annotation_handle \
[G_msg "cumulative length = $out_tot\n"]
+ monitor_annotate $measurement_annotation_handle \
+ [G_msg " area = $out_area\n"]
set linex1 $linex2
set liney1 $liney2
@@ -1779,6 +1825,68 @@
+# format area numbers and units in a nice way, as a function of area
+proc MapCanvas::fmt_area { area } {
+
+ set mapunits [MapCanvas::get_mapunits]
+
+ set outunits $mapunits
+ set divisor "1.0"
+
+ # figure out which units to use
+ if { [string equal "meters" "$mapunits"] } {
+ if { $area > 2500000 } {
+ set outunits "km2"
+ set divisor "1000000.0"
+ } elseif { $area > 25000 } {
+ set outunits "ha"
+ set divisor "10000.0"
+ } else {
+ set outunits "m2"
+ set divisor "1.0"
+ }
+ } elseif { [string first "feet" "$mapunits"] >= 0 } {
+ # nano-bug: we match any "feet", but US Survey feet is really
+ # 5279.9894 per statute mile, or 10.6' per 1000 miles. As >1000
+ # miles the tick markers are rounded to the nearest 10th of a
+ # mile (528'), the difference in foot flavours is ignored.
+ if { $area > 27878400 } {
+ set outunits "sq miles"
+ set divisor "27878400.0"
+ } elseif { $area > 43560 } {
+ set outunits "acres"
+ set divisor "43560.0"
+ } elseif { $area > 900 } {
+ set outunits "sq yards"
+ set divisor "9.0"
+ } else {
+ set outunits "sq feet"
+ set divisor "1.0"
+ }
+ } else {
+ set outunits "unit?"
+ set divisor "1.0"
+ }
+
+ # format numbers in a nice way
+ if { [expr $area/$divisor ] >= 2500 } {
+ set outfmt "%.0f"
+ } elseif { [expr $area/$divisor ] >= 1000 } {
+ set outfmt "%.1f"
+ } elseif { [expr $area/$divisor ] > 0 } {
+ set outfmt "%.[expr {int(ceil(3 - log10($area/$divisor)))}]f"
+ } else {
+ # error: no range (nan?)
+ set outfmt "%g"
+ }
+
+ set outarea [format $outfmt [expr $area/$divisor ] ]
+
+ return [concat $outarea $outunits ]
+}
+
+
+
###############################################################################
# procedures for querying
_______________________________________________
grass-dev mailing list
[email protected]
http://lists.osgeo.org/mailman/listinfo/grass-dev
