<!DOCTYPE html PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN">
<html>
<head>
<meta content="text/html;charset=ISO-8859-1" http-equiv="Content-Type">
<title></title>
</head>
<body bgcolor="#ffffff" text="#000000">
Hi Roger,<br>
<br>
No, I didn't use the calculate statistics trick. I just opened the
rasters and when they didn't display properly, I simply changed to
"Stretched" and then back to "Classified" and it looked fine. I should
mention though, that I added each band separately rather than adding
both bands at once. Zev<br>
<br>
Roger Bivand wrote:
<blockquote
cite="midPine.LNX.4.44.0605100959010.5583-100000@reclus.nhh.no"
type="cite">
<pre wrap="">On Tue, 9 May 2006, Zev Ross wrote:
</pre>
<blockquote type="cite">
<pre wrap="">Hi Roger and Andy,
Thanks for the help. I'm using ArcGIS 9.1. writeGDAL works as promised.
The one odd behavior was the initial display in ArcGIS once I first added
the two bands. Using the "stretched" classification, the map initially
looked all black. Under properties, I changed the symbology to
"Classified" and it looked fine -- then I changed back to "Stretched" and
again it looked fine. Quirky behavior like this from ArcGIS is nothing
new.
</pre>
</blockquote>
<pre wrap=""><!---->
Did you do the:
ArcToolbox | Data Management Tools | Raster | Calculate Statistics
trick too? I think that is what Arc (sometimes) does on the fly when
changing symbology (dreadful term!). But when NAs are present, it can get
stuck.
</pre>
<blockquote type="cite">
<pre wrap="">I tried the suggested write ascii command from maptools on my own data
set and it didn't work initially due to the fact that my grid cells
apparently are not perfect squares (perhaps this is a rounding issue). No
need to keep fiddling with that, though, as writeGDAL does the trick.
Thanks again,
</pre>
</blockquote>
<pre wrap=""><!---->
Yes, they need to be equal for the Arc ASCII grid route to work. If they
are, it works (when the decimal point in FLOAT is correct too).
Roger
</pre>
<blockquote type="cite">
<pre wrap="">Zev
Roger Bivand wrote:
On Thu, 4 May 2006, Zev Ross wrote:
Hi All,
I'd like to write the results of a kriging call to an ascii grid using
GSTAT -- but I would like to write BOTH the predictions and the variances
to grids. Is there a more elegant (and less dangerous) way to do it than
what I have below?
Dangerous for whom? Elegant, would be nice but life is short? The answer
depends on the software that is going to read the output grids, and how it
treats locales, etc., since both predictions and variances will be
floating point.
depth_uk <- krige(DEPTH~slope.asc, depth, slope, vgm_depth_r)
# write the predictions
write.asciigrid(depth_uk "c:/junk/rk/predictions.asc")
# replace predictions with variances and then write the variances
depth_uk$var1.pred<-depth_uk$var1.var
write.asciigrid(depth_uk, "c:/junk/rk/.asc")
If the software on the other side reads GeoTiff, my preference would be:
library(rgdal)
writeGDAL(depth_uk, "depth_uk.tif")
which I have seen work with ENVI, but not with ArcGIS 9.1; this preserves
coordinate reference system metadata if set.
In a forthcoming release of rgdal, you should be able to pass options= to
writeGDAL() - specifically INTERLEAVE=PIXEL, see:
<a class="moz-txt-link-freetext" href="http://grass.itc.it/grass61/manuals/html61_user/r.out.gdal.html">http://grass.itc.it/grass61/manuals/html61_user/r.out.gdal.html</a>
and I have seen this help with ArcGIS 9.1, although it wasn't predictable
(the legend scale showed correct values but the visualisation was wrong
sometimes - I tried on a Wednesday if that helps!). ArcGIS only accepted
single band GeoTiff files, it thought 3-band were coloured images. ENVI
simply read the GDAL-generated GeoTiffs (with 4 bands in the case we
tried - point pattern kernel densities at different bandwidths) correctly
without making any assumptions.
Depending on your locales, the ASCII grid route is being maintained in the
maptools package and functions in the sp package will be deprecated. So
library(maptools)
writeAsciiGrid(depth_uk, "preds.txt", attr="var1.pred", dec=<your choice>)
writeAsciiGrid(depth_uk, "vars.txt", attr="var1.var", dec=<your choice>)
should get the values into ArcGIS 9.1 through the Toolbox (it is very
sensitive to the "."/"," dec= setting). [The intention is to gather
input/output functions in maptools and rgdal, freeing the other packages
from having often older, duplicate copies of functions that do not get
maintained.]
Again, how to do it does depend on what software is going to read the
output ASCII grids, and what assumptions (often undocumented) it makes
about the files.
Please let us know how you get on,
Roger
PS this sample code comes from Tomislav Hengl's page
(<a class="moz-txt-link-freetext" href="http://spatial-analyst.net/regkriging.php">http://spatial-analyst.net/regkriging.php</a>)
Nice link!
Zev
--
Zev Ross
ZevRoss Spatial Analysis
303 Fairmount Ave
Ithaca, NY 14850
(607) 277-0004 (phone)
(866) 877-3690 (fax toll-free)
<a class="moz-txt-link-abbreviated" href="mailto:zev@zevross.com">zev@zevross.com</a>
<a class="moz-txt-link-abbreviated" href="http://www.zevross.com">www.zevross.com</a>
--
Zev Ross
ZevRoss Spatial Analysis
303 Fairmount Ave
Ithaca, NY 14850
(607) 277-0004 (phone)
(866) 877-3690 (fax toll-free)
<a class="moz-txt-link-abbreviated" href="mailto:zev@zevross.com">zev@zevross.com</a>
<a class="moz-txt-link-abbreviated" href="http://www.zevross.com">www.zevross.com</a>
</pre>
</blockquote>
<pre wrap=""><!---->
</pre>
</blockquote>
<br>
<div class="moz-signature">-- <br>
<title>Untitled Document</title>
<meta http-equiv="Content-Type" content="text/html; ">
<font face="Arial, Helvetica, sans-serif, Trebuchet MS" size="-1">Zev
Ross<br>
<strong><font color="#660000">ZevRoss Spatial Analysis</font></strong><br>
303 Fairmount Ave<br>
Ithaca, NY 14850<br>
(607) 277-0004 (phone)<br>
(866) 877-3690 (fax toll-free)<br>
<a class="moz-txt-link-abbreviated" href="mailto:zev@zevross.com">zev@zevross.com</a><br>
<a class="moz-txt-link-abbreviated" href="http://www.zevross.com">www.zevross.com</a></font><br>
</div>
</body>
</html>