Re: [GRASS-dev] Re: [GRASS GIS] #73: r.out.gdal tiff output does not work

[Glynn Clements]

Markus Metz wrote:

> This is now a mix of r.in.gdal and r.out.gdal. The two modules 
> complement each other, and I guess a minimum requirement would be that 
> something exported with r.out.gdal and then imported again with 
> r.in.gdal should be identical (taking into consideration the region 
> settings during export) and not have any data loss, be it due to 
> clipping data or due to nodata problems.

That's only possible if the exported data is in a format which can
hold the full range of the data. A UInt32 TIFF can't hold the full
range of CELL data (it can't hold any negative values).

> >>> The only "fool-proof" nodata value for CELL is 0x8000; all other
> >>> values may be actual data values.
> >>>   
> >> For grass rasters of type CELL, not for gdal export as Int32. gdalinfo 
> >> says NoData Value=2147483648, but cells that should be nodata are 
> >> -2147483648 and are thus not recognized as nodata.
> >> 
> >
> > It worked as of a few days ago, in the version I have checked out:
> >
> > raster/r.in.gdal/main.c:780:
> >
> >if (((GInt32 *) cell)[indx] == (GInt32) dfNoData) {
> >
> > (GInt32)2147483648 and (GInt32)-2147483648 are equal.
>
> I think the problem is that gdal allows the nodata value to be out of 
> range with regard to the raster type,

Yes.

> e.g. in gdal, a nodata value of 
> - with GDT_Byte is legal and can be used to effectively disable 
> nodata information. No idea if this conforms to GeoTIFF specs (probably 
> yes), but I'm pretty sure that this is the current policy of gdal 
> (possible with e.g. gdal_translate).
> That would also mean that r.out.gdal with type = GDT_Int32 and nodata = 
> 2147483648, NULL cells become -2147483648 but the nodata value in the 
> metadata stays 2147483648 (gdalinfo on the export output), which in turn 
> means that other software, also QGIS, does not see any nodata cells in 
> the export output ->information loss about NULL cells. That's why I 
> wouldn't use 0x8000 as default nodata value for GDT_Int32 in 
> r.out.gdal even if r.in.gdal has no problems with it.

You're misinterpretating my use of "0x8000" here. I'm talking
about the C interpretation, which is the (signed) integer value
-2147483648.

> > If it's been changed, use nodata=-0x8000. 
> 
> That's currently done in r.out.gdal in trunk and devbr6.

> > Or cast the nodata value
> > to the source data type then back to double.
> 
> For r.out.gdal, there was discussion about not to override user options 
> and instead issue an error or a warning (going for error now) that the 
> nodata value is out of range.

This isn't about overriding user options, but interpreting them
correctly.

It's debatable whether just using atof(nodataopt->answer) directly is
actually a correct interpretation when the input map is CELL, or
whether the code should use e.g. atoi(nodataopt->answer) or
(CELL)atof(nodataopt->answer) for CELL inputs. Reading as int (or
casting to it) will interpret 0x8000 as -2^31, while reading as FP
will interpret it as 2^31.

OTOH, reading as FP allows you to specify +2^31 as the nodata value
when exporting CELL as UInt32. Maybe it would be better to read as
double then cast to the destination type.

> Currently, all default nodata values are 
> within the range of the selected GDAL data type. User-defined nodata 
> values are tested whether they are within the range of the selected data 
> type, if not, r.out.gdal in trunk and devbr6 exits with an error 
> message. r.out.gdal also exits with an error message if the nodata value 
> is present in the raster to be exported in order to avoid data loss.
> This is a rather strict and conservative approach because of the 
> reported problems with nodata in r.out.gdal exports.

Right.

Even so, 0x8000 is still the ideal value for exporting CELL data
to either UInt32 or Int32. You just have to pick the correct
interpretation of it (unsigned int or signed int respectively).

-- 
Glynn Clements 
___
grass-dev mailing list
grass-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-dev

[GRASS-dev] Re: [GRASS GIS] #560: WinGRASS not deleting temp ppm files from map display

[GRASS GIS]

#560: WinGRASS not deleting temp ppm files from map display
---+
  Reporter:  isaacullah|   Owner:  grass-dev@lists.osgeo.org
  Type:  defect|  Status:  new  
  Priority:  major |   Milestone:  6.4.0
 Component:  Display   | Version:  6.4.0 RCs
Resolution:|Keywords:  v.digit ppm temp 
  Platform:  MSWindows XP  | Cpu:  x86-32   
---+
Comment (by glynn):

 Replying to [comment:1 hamish]:
 > the only thing I can think of is that the files are still in use, so can
 not be swept up by clean_temp on exit. ?

 AFAICT, clean_temp doesn't work (at all) on Windows. See the "TODO"s in
 the code.

 clean_temp only removes files if they are owned by the current user, and
 either the PID is no longer in use or the file is more than 4 days old.

 But it cannot determine the current user on Windows (no getuid()), so it
 uses -1, but MSVCRT's stat() always reports zero for the st_uid field.
 Also, the code which checks whether a given PID is still in use always
 fails (which is handled the same as if the PID *is* in use).

 If the current user check was skipped, it should at least delete files
 more than 4 days old.

 Implementing the ownership check would be non-trivial, due to the
 complexity of the NT security model (e.g. a file can be owned by a group
 rather than a user, there's a difference between owner and creator, a
 "user" consists of both a domain and a username, etc).

-- 
Ticket URL: 
GRASS GIS 
___
grass-dev mailing list
grass-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-dev

Re: [GRASS-dev] v.net.visibility memory leak? - was Re: [GRASS-user] traveling salesman problem in air

[Glynn Clements]

Markus Metz wrote:

> >> I remember a discussion in the dev ml, I think with the osgeo4w 
> >> installer, that memory allocated with malloc/calloc/realloc should be 
> >> free-d with free, not G_free. If this is the case, then there is a 
> >> problem in the vector libs, because space for a struct line_pnts * is 
> >> allocated with calloc and free-d with G_free, and relevant (or all) 
> >> calls to malloc/calloc/realloc/free in diglib should be replaced with 
> >> G_malloc etc. What to do?
> >> 
> >
> > GRASS code should always use the G_* versions unless there is a clear
> > reason not to, e.g. if the pointer will be passed to an external
> > library which will assume ownership of it.
> 
> OK, so lib/vector/diglib needs to be updated. What about dglib and 
> rtree? AFAICT they are written such that the code can be used as 
> external standalone libraries and not a single G_*() version is used, I 
> think.

Two choices: either rewrite the libraries to use the G_* versions, or
ensure that their callers use free() for pointers returned from those
libraries (if the callers only use deallocation functions provided by
the libraries, this point is moot).

The choice would be influenced by whether the library is available
elsewhere, and if so the extent to which the GRASS version has been
modified, and whether the library requires the caller to allocate or
free memory.

For an unmodified (or minimally modified) version of an external
library which performs all allocation and deallocation itself, I would
leave it using malloc() and free().

If the library isn't available elsewhere, or has been heavily
modified, or it requires the caller to allocate and/or free memory, I
would modify it to use the G_* versions.

> >> If Vect_set_release_support() does not cause a severe time penalty, I 
> >> would recommend to use it. IMHO 73,854,014 - 444,026 (previous result) 
> >> bytes = ca. 70 MB lingering around in memory is a bit unnecessary.
> >
> > That should only be done if vectors will be closed signficantly prior
> > to termination, e.g.:
> >
> > for (i = 0; i < nvects; i++) {
> > vect = Vect_open_old(...);
> > ...
> > Vect_close(vect);
> > }
> 
> Can we add this to the Programmer's manual under GRASS 6 Vector 
> Architecture where Vect_set_release_support() is mentioned or have I 
> missed something?

That would be a good idea, IMHO.

-- 
Glynn Clements 
___
grass-dev mailing list
grass-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-dev

[GRASS-dev] add area measurement

[yuhuu]

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.0 +0200
+++ mapcanvas.tcl-with-marea	2009-04-18 16:51:23.0 +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 > 250 } {
+	set outunits "km2"
+	set divisor "100.0"
+	} elseif { $area > 25000 } {
+	set outunits "ha"
+	set divisor "1.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 ]
+}
+
+
+
 ###

[GRASS-dev] Re: [GRASS GIS] #562: patch proposal: do not show delimiters for empty labels in case of multiple cats per feature

[GRASS GIS]

#562: patch proposal: do not show delimiters for empty labels in case of 
multiple
cats per feature
--+-
  Reporter:  mlennert |   Owner:  grass-dev@lists.osgeo.org
  Type:  enhancement  |  Status:  new  
  Priority:  normal   |   Milestone:  7.0.0
 Component:  default  | Version:  svn-trunk
Resolution:   |Keywords:   
  Platform:  Unspecified  | Cpu:  Unspecified  
--+-
Comment (by hamish):

 what's the command line option to do that? I'm not sure I know what that
 does.


 my gut feeling is that those empty delimiters actually do show
 information: positional. if you remove them and only one field has a
 value, which field did it come from, a, b, or c?


 Hamish

-- 
Ticket URL: 
GRASS GIS 
___
grass-dev mailing list
grass-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-dev

[GRASS-dev] [GRASS GIS] #562: patch proposal: do not show delimiters for empty labels in case of multiple cats per feature

[GRASS GIS]

#562: patch proposal: do not show delimiters for empty labels in case of 
multiple
cats per feature
-+--
 Reporter:  mlennert |   Owner:  grass-dev@lists.osgeo.org
 Type:  enhancement  |  Status:  new  
 Priority:  normal   |   Milestone:  7.0.0
Component:  default  | Version:  svn-trunk
 Keywords:   |Platform:  Unspecified  
  Cpu:  Unspecified  |  
-+--
 The attached patch (against svn7) slightly modifies d.vect in that it only
 displays delimiters ("/") for those category values where label values
 actually exist. I attache two screenshots to make clearer what the patch
 does.

 I would just submit the patch without posting, but I'm not sure if some
 people find it important to see the existance of null values for some
 categories. However, if no one objects, I'll commit.

 Moritz

-- 
Ticket URL: 
GRASS GIS 
___
grass-dev mailing list
grass-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-dev

[GRASS-dev] Re: [GRASS GIS] #554: wxGUI location wizard: Search in description with 0 hits

[GRASS GIS]

#554: wxGUI location wizard: Search in description with 0 hits
---+
  Reporter:  hamish|   Owner:  martinl
  Type:  defect|  Status:  assigned   
  Priority:  major |   Milestone:  6.4.0  
 Component:  wxGUI | Version:  6.4.0 RCs  
Resolution:|Keywords:  location wizard
  Platform:  MSWindows XP  | Cpu:  x86-32 
---+
Comment (by martinl):

 Hi,

 Replying to [comment:2 neteler]:
 >
 
http://grass.osgeo.org/grass64/binary/mswindows/native/WinGRASS-6.4.0RC4-1-Setup.exe

 could we release also weekly updated native installer? - to enable users
 to test recently fixed bugs?

 Martin

-- 
Ticket URL: 
GRASS GIS 
___
grass-dev mailing list
grass-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-dev

[GRASS-dev] Re: [GRASS GIS] #551: WinGRASS-6.4_RC3 error creating location with EPSG codes.

[GRASS GIS]

#551: WinGRASS-6.4_RC3 error creating location with EPSG codes.
---+
  Reporter:  isaacullah|   Owner:  grass-dev@lists.osgeo.org
  Type:  defect|  Status:  closed   
  Priority:  major |   Milestone:  6.4.0
 Component:  Python| Version:  6.4.0 RCs
Resolution:  fixed |Keywords:  create location epsg 
  Platform:  MSWindows XP  | Cpu:  Unspecified  
---+
Changes (by martinl):

  * status:  new => closed
  * resolution:  => fixed

Comment:

 Replying to [ticket:551 isaacullah]:
 > Trying to create a location with EPSG codes in WinGRAS-6.4_RC3 using the
 WX-GUI gives following error: "Unable to read EPSG codes: string index out
 of range"
 >
 > The path to epsg file is correct, and the file is normal. It seems
 wingrass cannot parse the epsg file.
 >
 > Work around: Change directory symbol in path to epsg file from "\" to
 "/". Enter a known epsg code (eg. search spatialrefernce.org for the
 correct code) in the EPSG box, and press "next". Grass creates new
 location.
 >
 > This work around does not let one "Browse EPSG codes", however. If one
 trys to do so, the same error "Unable to read EPSG codes: string index out
 of range" pops up.
 >
 > This bug occurs in both osgeo4W version and the native installer version
 of WinGRASS-6.4_RC3

 Fixed in RC4, see also #554 for additional fixes.

 M.

-- 
Ticket URL: 
GRASS GIS 
___
grass-dev mailing list
grass-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-dev

[GRASS-dev] Re: [GRASS GIS] #554: wxGUI location wizard: Search in description with 0 hits

[GRASS GIS]

#554: wxGUI location wizard: Search in description with 0 hits
---+
  Reporter:  hamish|   Owner:  martinl
  Type:  defect|  Status:  assigned   
  Priority:  major |   Milestone:  6.4.0  
 Component:  wxGUI | Version:  6.4.0 RCs  
Resolution:|Keywords:  location wizard
  Platform:  MSWindows XP  | Cpu:  x86-32 
---+
Changes (by martinl):

  * status:  new => assigned
  * owner:  grass-dev@lists.osgeo.org => martinl
 * cc: grass-dev@lists.osgeo.org (added)

Comment:

 Replying to [ticket:554 hamish]:
 > - Search in description: ___ something with no hits, eg 4326 or typo

 you search in the description, not in the codes - so you get empty list.
 Use "Reload EPGS codes" to reload the list.

 > At this point the list of codes is empty and the following error message
 appears on the DOS console window:
 >
 > {{{
 > Traceback (most recent call last):
 >  File
 
"C:\OSGeo4W\apps\grass\grass-6.4.0svn\etc\wxpython\gui_modules\location_wizard.py",
 line 1365, in OnSearch
 > self.epsgcode = self.epsglist.Search(index=1, str=value)[0]
 > TypeError: 'NoneType' object is unsubscriptable
 > }}}

 Fixed in r36775.

 > so I guess two things are needed:
 > - only try to grab epsg code to fill in "EPSG code: " box if >=1
 result was found, and
 > - automatically reload EPSG codes before each new search (if 0 hits then
 the 'search in desc' box is unusable until you press the Reload EPSG Codes
 button.

 Done in r36776.

 M.

-- 
Ticket URL: 
GRASS GIS 
___
grass-dev mailing list
grass-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-dev

[GRASS-dev] Re: [GRASS GIS] #559: wingrass stand-alone installer: MSys console doesn't work

[GRASS GIS]

#559: wingrass stand-alone installer: MSys console doesn't work
---+
  Reporter:  hamish|   Owner:  cnielsen 
  Type:  defect|  Status:  assigned 
  Priority:  minor |   Milestone:  6.4.0
 Component:  Installation  | Version:  6.4.0 RCs
Resolution:|Keywords:  wingrass 
  Platform:  MSWindows XP  | Cpu:  x86-32   
---+
Comment (by mlennert):

 Replying to [comment:2 cnielsen]:
 > Replying to [ticket:559 hamish]:
 > > ps- is there a way to automatically minimize the DOS box when
 launching GRASS? it just sort of sits there idle.
 >
 > Yes this can be done in the creation of the shortcut (*.lnk) files.
 >
 > {{{SW_SHOWMINIMIZED}}} will start it minimized
 >
 > {{{SH_HIDE}}} should (I haven't tested it yet) start it but run it in
 the background (ie. hide it completely). Although there would be no extra
 window, it may hide error messages like the one that says wxdigit isn't
 working.
 >
 > Which option should I implement?

 Can this be altered by the user via the properties of the link ? If yes,
 then I would plead for SH_HIDE and then in case of problems we can tell
 the user to change that. If not, then we probably have not much of a
 choice and have to leave it with SW_SHOWMINIMIZED, unless there is a
 different way of starting GRASS which would allow debugging.

 Other option: if there is a GRASS menu in the program menu, add a "Debug
 mode" entry which starts with the window apparent.

 Moritz

-- 
Ticket URL: 
GRASS GIS 
___
grass-dev mailing list
grass-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-dev

[GRASS-dev] Re: [GRASS GIS] #557: wxGUI startup gui: shade start and double click mapset

[GRASS GIS]

#557: wxGUI startup gui: shade start and double click mapset
--+-
  Reporter:  hamish   |   Owner:  martinl
  Type:  enhancement  |  Status:  assigned   
  Priority:  trivial  |   Milestone:  6.4.0  
 Component:  wxGUI| Version:  6.4.0 RCs  
Resolution:   |Keywords:  startup
  Platform:  Linux| Cpu:  Unspecified
--+-
Changes (by martinl):

  * status:  new => assigned
  * owner:  grass-dev@lists.osgeo.org => martinl
 * cc: grass-dev@lists.osgeo.org (added)

Comment:

 Done in r36771 (trunk), r36772 (devbr6) and r36773 (relbr64).

 Martin

-- 
Ticket URL: 
GRASS GIS 
___
grass-dev mailing list
grass-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-dev

[GRASS-dev] Re: [GRASS GIS] #561: meta: add a new "normal" bug priority

[GRASS GIS]

#561: meta: add a new "normal" bug priority
--+-
  Reporter:  hamish   |   Owner:  grass-dev@lists.osgeo.org
  Type:  enhancement  |  Status:  closed   
  Priority:  minor|   Milestone:  Website  
 Component:  Website  | Version:   
Resolution:  fixed|Keywords:  trac 
  Platform:  Unspecified  | Cpu:  Unspecified  
--+-
Changes (by martinl):

  * status:  new => closed
  * resolution:  => fixed

Comment:

 Done.

-- 
Ticket URL: 
GRASS GIS 
___
grass-dev mailing list
grass-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-dev