[GRASS-user] d.out.file now does GeoTIFFs

Hamish · October 29, 2008, 8:23am

Hi,

I just added support for GeoTIFF export from d.out.file in 6.4svn.

Resolution is not as good as if you use r.out.gdal or r.out.tiff, but
it lets you export a map with multiple raster and vector overlays,
decorations, grid lines, etc. that you compose in the xmon.

e.g. for use as a backdrop in some navigation software.
(overlay bathymetry DEM of river channel over an orthophoto)

slight contributing motivation:
I am not sure of a good way to merge raster maps with different color
tables (imagery data). In the past I have tried 'r.mapcalc r# g# b#'
operators+r.patch, d.out.png+r.in.gdal+r.region, neither was much fun.
I guess g.pnmcomp could help but I fear the added complexity of passing
the data through the display drivers will introduce noise. Perhaps
gdal_merge.py? (see r.in.wms mess WRT merging paletted GIF images)

please test,
Hamish

neteler · October 29, 2008, 8:26am

On Wed, Oct 29, 2008 at 9:23 AM, Hamish <hamish_b@yahoo.com> wrote:

Hi,

I just added support for GeoTIFF export from d.out.file in 6.4svn.

This is really great, Hamish! Maybe you remember our discussion
on this - I'll test it asap.

thanks,
Markus

Nikos_Alexandris1 · October 29, 2008, 8:51am

On Wed, 2008-10-29 at 21:23 +1300, Hamish wrote:

Hi,

I just added support for GeoTIFF export from d.out.file in 6.4svn.

Resolution is not as good as if you use r.out.gdal or r.out.tiff, but
it lets you export a map with multiple raster and vector overlays,
decorations, grid lines, etc. that you compose in the xmon.

e.g. for use as a backdrop in some navigation software.
(overlay bathymetry DEM of river channel over an orthophoto)

slight contributing motivation:
I am not sure of a good way to merge raster maps with different color
tables (imagery data). In the past I have tried 'r.mapcalc r# g# b#'
operators+r.patch, d.out.png+r.in.gdal+r.region, neither was much fun.
I guess g.pnmcomp could help but I fear the added complexity of passing
the data through the display drivers will introduce noise. Perhaps
gdal_merge.py? (see r.in.wms mess WRT merging paletted GIF images)

please test,
Hamish

Can I test it without installing it? I don't want to re-configure all
possible add-onds I have currently installed.

I mean
#update
svn up

#clean previous configuration
sudo make distclean

#re-connfigure
./configure # with whatever of options

#make
make

and then run it directly from within the configured source (?)

Thanks, Nikos

Nikos_Alexandris1 · October 29, 2008, 9:09am

On Wed, 2008-10-29 at 09:51 +0100, Nikos Alexandris wrote:

On Wed, 2008-10-29 at 21:23 +1300, Hamish wrote:
> Hi,
>
> I just added support for GeoTIFF export from d.out.file in 6.4svn.
>
> Resolution is not as good as if you use r.out.gdal or r.out.tiff, but
> it lets you export a map with multiple raster and vector overlays,
> decorations, grid lines, etc. that you compose in the xmon.
>
> e.g. for use as a backdrop in some navigation software.
> (overlay bathymetry DEM of river channel over an orthophoto)
>
>
> slight contributing motivation:
> I am not sure of a good way to merge raster maps with different color
> tables (imagery data). In the past I have tried 'r.mapcalc r# g# b#'
> operators+r.patch, d.out.png+r.in.gdal+r.region, neither was much fun.
> I guess g.pnmcomp could help but I fear the added complexity of passing
> the data through the display drivers will introduce noise. Perhaps
> gdal_merge.py? (see r.in.wms mess WRT merging paletted GIF images)
>
>
> please test,
> Hamish

Can I test it without installing it? I don't want to re-configure all
possible add-onds I have currently installed.

I mean
#update
svn up

#clean previous configuration
sudo make distclean

#re-connfigure
./configure # with whatever of options

#make
make

and then run it directly from within the configured source (?)

Thanks, Nikos

Well,
so I did.

It works great. I did load 2 raster maps, 2 vector (boundaries) maps and
a north arrow. Exported as geotiff and opened in OpenEV. The image is
georeferenced (working with epsg 2100 in my case).

Thanks, Nikos

maning_sambale · October 29, 2008, 10:09am

Well,
so I did.

How? I haven't been doing svn update for the past 2 months because
its working OK.
Would love to try this out also.

It works great. I did load 2 raster maps, 2 vector (boundaries) maps and
a north arrow. Exported as geotiff and opened in OpenEV. The image is
georeferenced (working with epsg 2100 in my case).

--
|---------|----------------------------------------------------------|
| __.-._ |"Ohhh. Great warrior. Wars not make one great." -Yoda |
| '-._"7' |"Freedom is still the most radical idea of all" -N.Branden|
| /'.-c |Linux registered user #402901, http://counter.li.org/ |
| | /T |http://esambale.wikispaces.com|
| _)_/LI
|---------|----------------------------------------------------------|

Hamish · October 29, 2008, 10:17am

Nikos Alexandris:

> Well,
> so I did.

maning sambale wrote:

How? I haven't been doing svn update for the past 2
months because its working OK.
Would love to try this out also.

if you are working from a svn source tree it is just:

cd scripts/d.out.file
svn up
make

you only have to do ./configure again and rebuild everything for deep
changes. (if you haven't in a few months it is probably a good idea
thought)

and whatever you do to "make install"

but this one is just a script, you can download the latest version and
simply replace the file. (remember 'chmod a+x')

Hamish

Hamish · October 29, 2008, 10:35am

Hamish:

> I just added support for GeoTIFF export from
d.out.file in 6.4svn.

Markus:

This is really great, Hamish! Maybe you remember our discussion
on this - I'll test it asap.

I didn't, but re-reading the thread
http://thread.gmane.org/gmane.comp.gis.grass.devel/24289/
I see it has the same problem you found from the white borders of an
unmatched region/xmon aspect ratio. mmph.

way back when Hamish wrote:

oh, NSEW bounds *including* the white out-of-region bands on the side
of the display if the aspect ratios don't match.

and Markus replied:

Exactly. That's why I wanted to get the coordinates from the monitor
itself.

Unfortunately 'd.info -g' is buggy, it reports exactly the same as
'g.region -g'. But if that was fixed it would be an easy solution.

Hamish

maning_sambale · October 29, 2008, 10:37am

Hamish,

Just tried it. This makes a quick and dirty visualization really
fast! Especially for non GIS colleagues.

I can simply import the d.mon display to geotiff then run
gdal2tile.py. VoilaI I have a KML of the output map for viewing in
GoogleEarth.

Thanks!

cheers,
maning

On 10/29/08, Hamish <hamish_b@yahoo.com> wrote:

Nikos Alexandris:

> Well,
> so I did.

maning sambale wrote:

How? I haven't been doing svn update for the past 2
months because its working OK.
Would love to try this out also.

if you are working from a svn source tree it is just:

cd scripts/d.out.file
svn up
make

you only have to do ./configure again and rebuild everything for deep
changes. (if you haven't in a few months it is probably a good idea
thought)

and whatever you do to "make install"

but this one is just a script, you can download the latest version and
simply replace the file. (remember 'chmod a+x')

Hamish

--
|---------|----------------------------------------------------------|
| __.-._ |"Ohhh. Great warrior. Wars not make one great." -Yoda |
| '-._"7' |"Freedom is still the most radical idea of all" -N.Branden|
| /'.-c |Linux registered user #402901, http://counter.li.org/ |
| | /T |http://esambale.wikispaces.com|
| _)_/LI
|---------|----------------------------------------------------------|

Hamish · October 29, 2008, 10:55am

Unfortunately 'd.info -g' is buggy, it reports exactly the same as
'g.region -g'. But if that was fixed it would be an easy solution.

apparently D_d_to_u_col() and D_d_to_u_row() + R_screen_rite() etc are
needed.

one moment....

Hamish

Hamish · October 29, 2008, 12:43pm

Hamish:

> Unfortunately 'd.info -g' is buggy, it reports exactly the same as
> 'g.region -g'. But if that was fixed it would be an easy solution.

fixed in devbr6 but not trunk; I'm unsure about the changes needed.

now it reports the coords of the full window and resolution is
map_units/pixel.

apparently D_d_to_u_col() and D_d_to_u_row() + R_screen_rite() etc are
needed.

one moment....

I solved the d.out.file geotiff aspect problem by keeping the g.region
settings and removing the white strips on the edge with dimensions from
'd.info -b'. Too bad if you had decorations in the margins...

If there were demand for it, it would be possible to add a flag to leave
the margins there and use 'd.info -g' to set the GeoTIFF georef bounds to
the bigger size.

Hamish

Nikos_Alexandris1 · October 29, 2008, 12:51pm

On Wed, 2008-10-29 at 03:17 -0700, Hamish wrote:

Nikos Alexandris:
> > Well,
> > so I did.

maning sambale wrote:
> How? I haven't been doing svn update for the past 2
> months because its working OK.
> Would love to try this out also.

if you are working from a svn source tree it is just:

cd scripts/d.out.file
svn up
make

you only have to do ./configure again and rebuild everything for deep
changes. (if you haven't in a few months it is probably a good idea
thought)

and whatever you do to "make install"

but this one is just a script, you can download the latest version and
simply replace the file. (remember 'chmod a+x')

Hamish

I thought that it's just a script since it is located in the /scripts
directory. But the presence of a Makefile "told" me to re-compile. Now I
know... Anyway, even doing more than required, I tested it. I
re-compiled grass6_devel and did NOT install (sudo make install) it.
Instead I ran GRASS from within the
"/usr/local/src/grass6_devel/bin.x86_64-unknown-linux-gnu" directory.

Nikos_Alexandris1 · October 29, 2008, 1:03pm

On Wed, 2008-10-29 at 05:43 -0700, Hamish wrote:

Hamish:
> > Unfortunately 'd.info -g' is buggy, it reports exactly the same as
> > 'g.region -g'. But if that was fixed it would be an easy solution.

fixed in devbr6 but not trunk; I'm unsure about the changes needed.

now it reports the coords of the full window and resolution is
map_units/pixel.

> apparently D_d_to_u_col() and D_d_to_u_row() + R_screen_rite() etc are
> needed.
>
> one moment....

I solved the d.out.file geotiff aspect problem by keeping the g.region
settings and removing the white strips on the edge with dimensions from
'd.info -b'. Too bad if you had decorations in the margins...

If there were demand for it, it would be possible to add a flag to leave
the margins there and use 'd.info -g' to set the GeoTIFF georef bounds to
the bigger size.

Hamish

Tested again, white borders no more there Great!! But it's true,
with a grid decoration (d.grid), "numbers" falling outside of the grid
are cut-out in the exported geotiff. Not that bad though.

Kind regards, Nikos

Glynn_Clements1 · October 29, 2008, 9:03pm

Hamish wrote:

> > Unfortunately 'd.info -g' is buggy, it reports exactly the same as
> > 'g.region -g'. But if that was fixed it would be an easy solution.

fixed in devbr6 but not trunk; I'm unsure about the changes needed.

7.0 doesn't have D_get_screen_window(); D_setup() uses R_get_window()
instead (R_get_window() reports the raster clip region, which is
initially set from $GRASS_FRAME).

Also, I would suggest reporting the values returned from D_get_u_*
rather than using the current region directly (which can always be
obtained from g.region).

Although the two should match, it's the values reported by D_get_u_*
which are actually used, so using D_get_u_* would be more robust in
the face of any future extensions (e.g. allowing the region to be
enlarged to fit the window rather than shrinking the window to fit the
region).

--
Glynn Clements <glynn@gclements.plus.com>

maning_sambale · October 30, 2008, 1:33am

Hamish,

I just found out that the output from d.out.file geotiff has some distortions.

See this: http://farm4.static.flickr.com/3007/2984908579_b73cdd3b99_o.png

The lighter relief is the original raster from the GRASS dbase and the
darker one with some green vector is the d.out.file output. I'm using
qgis to visualize.

cheers,
maning

On Thu, Oct 30, 2008 at 5:03 AM, Glynn Clements
<glynn@gclements.plus.com> wrote:

Hamish wrote:

> > Unfortunately 'd.info -g' is buggy, it reports exactly the same as
> > 'g.region -g'. But if that was fixed it would be an easy solution.

fixed in devbr6 but not trunk; I'm unsure about the changes needed.

7.0 doesn't have D_get_screen_window(); D_setup() uses R_get_window()
instead (R_get_window() reports the raster clip region, which is
initially set from $GRASS_FRAME).

Also, I would suggest reporting the values returned from D_get_u_*
rather than using the current region directly (which can always be
obtained from g.region).

Although the two should match, it's the values reported by D_get_u_*
which are actually used, so using D_get_u_* would be more robust in
the face of any future extensions (e.g. allowing the region to be
enlarged to fit the window rather than shrinking the window to fit the
region).

--
Glynn Clements <glynn@gclements.plus.com>
_______________________________________________
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

--
|---------|----------------------------------------------------------|
| __.-._ |"Ohhh. Great warrior. Wars not make one great." -Yoda |
| '-._"7' |"Freedom is still the most radical idea of all" -N.Branden|
| /'.-c |Linux registered user #402901, http://counter.li.org/ |
| | /T |http://esambale.wikispaces.com|
| _)_/LI
|---------|----------------------------------------------------------|

Hamish · October 30, 2008, 1:43am

maning sambale wrote:

I just found out that the output from d.out.file geotiff
has some distortions.

See this:
http://farm4.static.flickr.com/3007/2984908579_b73cdd3b99_o.png

The lighter relief is the original raster from the GRASS dbase and the
darker one with some green vector is the d.out.file output.

Can you try again with today's SVN? You need to update the d.info
module as well. I think it should be better now.

I'm using qgis to visualize.

yes me too, "qgis known_good.tif new.tif" and setting the transparency
to 50% is a very fast way of checking.

Today I hope to enable margin removal for all d.out.file formats.
(a flag to keep them is already in place)

I'd like to get rid of the "sleep" steps though (AFAIR the PNG driver
may still be working after "d.mon stop" returns)

Is GRASS_PNG_AUTO_WRITE the solution? or RENDER_IMMEDIATE ..?

Hamish

Hamish · October 30, 2008, 12:13pm

Glynn:

7.0 doesn't have D_get_screen_window(); D_setup() uses R_get_window()
instead (R_get_window() reports the raster clip region, which is
initially set from $GRASS_FRAME).

note what we are after for the d.info -g flag is the *unclipped* geographic region.
as you say, we already have g.region -g for the clipped version.

Hamish

Glynn_Clements1 · October 30, 2008, 4:49pm

Hamish wrote:

> I just found out that the output from d.out.file geotiff
> has some distortions.
>
> See this:
> http://farm4.static.flickr.com/3007/2984908579_b73cdd3b99_o.png
>
> The lighter relief is the original raster from the GRASS dbase and the
> darker one with some green vector is the d.out.file output.

Can you try again with today's SVN? You need to update the d.info
module as well. I think it should be better now.

> I'm using qgis to visualize.

yes me too, "qgis known_good.tif new.tif" and setting the transparency
to 50% is a very fast way of checking.

Today I hope to enable margin removal for all d.out.file formats.
(a flag to keep them is already in place)

I'd like to get rid of the "sleep" steps though (AFAIR the PNG driver
may still be working after "d.mon stop" returns)

This shouldn't happen. R_kill_driver() sends GRAPH_CLOSE to the
monitor, then tries to read() one byte. As the monitor never sends
this byte, the read() will block until it sees EOF, which happens when
the monitor terminates.

Similarly, R_close_driver() won't return until the driver's Respond()
method has completed. For the PNG driver, this means writing the image
if GRASS_PNG_AUTO_WRITE=TRUE.

Is GRASS_PNG_AUTO_WRITE the solution? or RENDER_IMMEDIATE ..?

For output, you should probably be using GRASS_RENDER_IMMEDIATE rather
than the PNG driver.

In 6.3.0 and later, you can use GRASS_PNG_READ to have the driver
overlay its output onto an existing image rather than replacing it.

GRASS_PNG_READ requires that the image already exists with the correct
format and dimensions. Typical usage:

  GRASS_RENDER_IMMEDIATE=TRUE
  GRASS_PNG_READ=FALSE
  export GRASS_RENDER_IMMEDIATE GRASS_PNG_READ
  d.erase
  GRASS_PNG_READ=TRUE
  <more d.* commands>

GRASS_PNG_AUTO_WRITE only matters if you are using a PNG monitor
rather than immediate rendering, and you need the file to be written
after each d.* command. Without GRASS_PNG_AUTO_WRITE, the file will
only be written when the driver is terminated with "d.mon stop=PNG".

--
Glynn Clements <glynn@gclements.plus.com>

Glynn_Clements1 · October 30, 2008, 5:11pm

Hamish wrote:

slight contributing motivation:
I am not sure of a good way to merge raster maps with different color
tables (imagery data). In the past I have tried 'r.mapcalc r# g# b#'
operators+r.patch, d.out.png+r.in.gdal+r.region, neither was much fun.
I guess g.pnmcomp could help but I fear the added complexity of passing
the data through the display drivers will introduce noise. Perhaps
gdal_merge.py? (see r.in.wms mess WRT merging paletted GIF images)

Maybe r.composite should have an option to use a single input map
rather than separate channels? IOW, convert the input map (whether
integer or FP) to an integer map whose categories directly correspond
to the input map's colours.

It could easily be implemented as a script:

  in=$GIS_OPT_INPUT
  out=$GIS_OPT_OUTPUT
  tmp=tmp.colors.$$
  r.mapcalc <<EOF
  $tmp.r = r#$in
  $tmp.g = g#$in
  $tmp.b = b#$in
  EOF
  r.composite r=$tmp.r g=$tmp.g b=$tmp.g levels=256 output=$tmp.c
  r.mapcalc "$out = if(isnull($in),null(),$tmp.c)"
  r.colors $out raster=$tmp.c
  g.remove rast=$tmp.r,$tmp.g,$tmp.b,$tmp.c

Actually, the first r.mapcalc should probably be isolated into a
separate "uncomposite" script, which split a colourised map into R/G/B
channels:

  r.mapcalc <<EOF
  $out.r = if(isnull($in),null(),r#$in)
  $out.g = if(isnull($in),null(),g#$in)
  $out.b = if(isnull($in),null(),b#$in)
  EOF

Also, a script to "collapse" a rule-based colour table (as generated
by e.g. r.composite) to a category-based table, using r.category and
r.what.color. Left as an exercise for the reader.

--
Glynn Clements <glynn@gclements.plus.com>

Glynn_Clements1 · October 30, 2008, 5:39pm

Hamish wrote:

> 7.0 doesn't have D_get_screen_window(); D_setup() uses R_get_window()
> instead (R_get_window() reports the raster clip region, which is
> initially set from $GRASS_FRAME).

note what we are after for the d.info -g flag is the *unclipped* geographic region.
as you say, we already have g.region -g for the clipped version.

I know.

The point is that you need the geographic coordinates which correspond
to some particular set of screen coordinates.

If you're writing a script which needs to perform its own translation
from display coordinates to geographic coordinates, you need two
pieces of information: a rectangle in display coordinates, and the
matching rectangle in geographic coordinates.

In many cases, it doesn't particularly matter which display rectangle
you pick, whether it's the full screen, or the active frame, or the
portion of the active frame to which the region is mapped (i.e. the
active frame excluding the "margins"), so long as you use the same
rectangle for both display and region coordinates.

If you need the active frame, R_get_window() will provide this
information. That reports the low-level clip window; nothing will be
drawn outside of that rectangle.

OTOH, the D_get_d_* functions report the portion to which the current
region is mapped, excluding any margins. Rasters will never be drawn
outside of that region. Vectors might be drawn outside of that region,
depending upon whether display-level clipping or culling is enabled.

Also, note that mapping the display frame to geographic coordinates
may produce an invalid region, as the north/south boundaries may
exceed 90 degrees.

--
Glynn Clements <glynn@gclements.plus.com>