[GRASS-stats] Raster maps from GRASS to R and back to GRASS?

[Roger Bivand]

On Tue, 21 Apr 2009, Nikos Alexandris wrote:

> Roger, Dylan and grass-stats-listers,
>
>  I can't "solve" the random-NULL production of writeRAST6. I tried
> again with only 3 bands which are identical with respect to their
> extent/resolution/etc. I am not sending more data (e.g. as I promised to
> send 6 modis bands) since the "problem" occurs with simpler dataset(s).

The NAs you show in the fake2 image are in band 7 (b7) and propagate 
correctly, b7 shown after import to R; the grass_nas image is from within 
GRASS showing the vector returned by complete.cases().

They are getting there because the NODATA= argument is not given, and is 
being set to 999. I've committed a fix to sourceforge CVS, but setting 
NODATA= manually is generally a good idea.

Roger

>
>  The grass-location with the specific MODIS bands is the one I have
> already uploaded in gregis [1].
>
> (
> I also repeated the test with landsat7 images from the spearfish
> location but all works fine there!
> )
>
>  I am trying to present the problem as simple as possible (commands
> follow after step-by-step description):
>
> Step-by-step:
> 1. launch R from within a GRASS-location
> 2. load library spgrass6 and projection information (gmeta6)
> 3. load _three_ modis bands in one object
> 4. create an NA-mask using the function complete.cases()
> 5. extract "values" only form the initial 3-modis-bands object
> 6. perform pca on "values" based on svd [that is with function prcomp()]
> 7. create new columns in the initial 3-modis-bands object
> 8. fill-in with the PC's produced by prcomp()
> 9. write back to grass with writeRAST6
>
> Result:
> Viewing the produced raster maps there are many "random" NULLs where all
> original raster maps are filled with some integer value.
>
> I guess those NULLs are not "random" at all - but I have no idea what
> the reason for this behavior is.
>
> (
> Check the "Value Tool", bottom left on the screen-shots [2][3][4]. Note
> that the mouse-pointer is a bit "shifted" -- don't know why, but it
> should be over the "black pixels"!
> )
>
> # # # # # # # # # # #
> # using: grass65, R2.8.1, latest spgrass6 installed today, GDAL 1.6
> # on: Ubuntu Intrepid 64-bit
> # launch grass in location/mapset "peloponnese_hgrs87/PERMAMENT"
> # code in R
>
> # load spgrass6 and projection information
> library(spgrass6) ; G <- gmeta6()
>
> # load modis bands
> mod_pre_b267.raw <- readRAST6 (c('mod_b2', 'mod_b6', 'mod_b7'))
> # str(mod_pre_b267.raw)
> # summary(mod_pre_b267.raw)
>
> # create NA mask using complete.cases()
> mod_pre_b267.na_mask <- complete.cases(mod_pre_b267.raw at data)
> # str(mod_pre_b267.na_mask)
> # summary(mod_pre_b267.na_mask)
>
> # get values based on na_mask
> mod_pre_b267 <- mod_pre_b267.raw at data[mod_pre_b267.na_mask, ]
>
> # perform pca (svd)
> prcomp.mod_pre_b267 <- prcomp(mod_pre_b267, scale=TRUE)
> summary(prcomp.mod_pre_b267)
> prcomp.mod_pre_b267$rotation
>
> # create new columns in initial object
> mod_pre_b267.raw at data$pc1 <- NA
> mod_pre_b267.raw at data$pc2 <- NA
> mod_pre_b267.raw at data$pc3 <- NA
>
> # fill-in with PCs
> mod_pre_b267.raw at data$pc1[mod_pre_b267.na_mask] <- prcomp.mod_pre_b267
> $x[,"PC1"]
> mod_pre_b267.raw at data$pc2[mod_pre_b267.na_mask] <- prcomp.mod_pre_b267
> $x[,"PC2"]
> mod_pre_b267.raw at data$pc3[mod_pre_b267.na_mask] <- prcomp.mod_pre_b267
> $x[,"PC3"]
>
> # write back to grass
> writeRAST6(mod_pre_b267.raw, zcol=4, vname="fake1", overwrite=TRUE)
> writeRAST6(mod_pre_b267.raw, zcol=5, vname="fake2", overwrite=TRUE)
> writeRAST6(mod_pre_b267.raw, zcol=6, vname="fake3", overwrite=TRUE)
>
> # check the "random" NULL's, e.g. with qgis
> # # # # # # # # # # #
>
>
> Hopefully you will have the time to have a look at it. Kindest regards,
> Nikos
> ---
>
> References
> [1] grass location containint three MODIS surface reflectance bands:
> https://git.berlios.de/cgi-bin/gitweb.cgi?p=gregis;a=blob_plain;f=peloponnesos/modis_peloponnese_postfire07.zip;hb=peloponnese
>
> [2] http://imageshack.gr/view.php?file=owj2xu833u6678c3zvch.png
> [3] http://imageshack.gr/view.php?file=n65vfvoo07t3njywfjjo.png
> [4] http://imageshack.gr/view.php?file=sb9l315wazoeue0rgrj9.png
>
>

-- 
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no
-------------- next part --------------
A non-text attachment was scrubbed...
Name: grass_nas.png
Type: image/png
Size: 2730 bytes
Desc: 
Url : http://lists.osgeo.org/pipermail/grass-stats/attachments/20090421/20222dd8/grass_nas.png
-------------- next part --------------
A non-text attachment was scrubbed...
Name: b7nas.png
Type: image/png
Size: 6432 bytes
Desc: 
Url : http://lists.osgeo.org/pipermail/grass-stats/attachments/20090421/20222dd8/b7nas.png