Skip to content

Problem with NDVI difference with subset image

4 messages · Sarah Goslee, ASANTOS

#
Hi,

What happens?

Do you get an error? Where?

What is your sessionInfo()?

As it says in the lssub() help, this function was written for a
particular purpose, is only known to work on linux, and may not be
widely applicable.

It's a "use at your own risk" kind of function, and you may well be
better off using one of the many other methods available for
subsetting raster data, which would save you the whole "export as
GeoTIFF, subset, reimport" sequence.

Sarah

On Sun, Jun 23, 2013 at 11:36 AM, Alexandre Santos
<alexandresantosbr at yahoo.com.br> wrote:

  
    
#
Hi Dra. Goslee,

             I work in OS linux and I like very much your lssub 
function, despites your notice. I think that the problems are in linux 
system, because:

 > landc<- stack(c("stackIm1.sample.tif","tackIm2.sample.tif"))
Error in .local(.Object, ...) :
   `/home/asantos/Documentos/Sensoriamento remoto percevejo 
bronzeador/Contornos/tackIm2.sample.tif' does not exist in the file system,
and is not recognised as a supported dataset name.

Error in .rasterObjectFromFile(x, band = band, objecttype = 
"RasterLayer",  :
   Cannot create a RasterLayer object from this file. (file does not exist)


              When I try to look the geotiff created 
(stackIm1.sample.tif) there are something wrong, because the pictures 
doesn't open (complete example below).

Thanks for your attention,

Alexandre

require(raster)
require(sp)
require(rgdal)
require(landsat)
#
#Create raster
r <- raster(nc=1000, nr=1000)
set.seed(20130622)
stackIm1 <- stack(lapply(1, function(x) setValues(r, round(runif(ncell(r))* 255))))## Simulation red band
stackIm2 <- stack(lapply(1, function(x) setValues(r, round(runif(ncell(r))* 255))))## Simulation nir

# define projection system
r.geo <- CRS("+proj=utm +zone=23 +south +datum=WGS84 +units=m +no_defs")    # geographical datum WGS84
proj4string(stackIm1) <- r.geo
proj4string(stackIm2) <- r.geo
#

#Create geotiff
writeRaster(stackIm1, filename="stackIm1.tif", format="GTiff",overwrite=TRUE)
writeRaster(stackIm2, filename="stackIm2.tif", format="GTiff",overwrite=TRUE)
#

#Subset a geotiff image 50 x 50 pixels
stackIm1.sample<-lssub("stackIm1.tif", "stackIm1.sample.tif", centerx = 0, centery = 0, widthx = 50, widthy = 50)
stackIm2.sample<-lssub("stackIm1.tif", "sstackIm2.sample.tif", centerx = 0, centery = 0, widthx = 50, widthy = 50)
#
landc<- stack(c("stackIm1.sample.tif","tackIm2.sample.tif"))





Em 23/06/2013 12:59, Sarah Goslee escreveu:

  
    
#
Hi,

On my system, it's a writeRASTER() options issue. I don't use raster
very often, but I played with the options enough to come up with a
working solution.



## this doesn't work in the right way - something about the format
## used by writeRaster isn't compatible with lssub()
#Create geotiff
writeRaster(stackIm1, filename="stackIm1.tif", format="GTiff",overwrite=TRUE)
## also, lssub() creates a file but does not return anything, so no
## need to assign it to an object
stackIm1.sample<-lssub("stackIm1.tif", "stackIm1.sample.tif", centerx
= 0, centery = 0, widthx = 50, widthy = 50)

## this works, confirming that it's a writeRaster issue
stackIm1a <- readGDAL("stackIm1.tif")
writeGDAL(stackIm1a, "stackIm1a.tif")
lssub("stackIm1a.tif", "stackIm1a.sample.tif", centerx = 0, centery =
0, widthx = 50, widthy = 50)

## This also works, but note the name of the output file
writeRaster(stackIm1, filename="stackIm1b.tif",
format="GTiff",overwrite=TRUE, bylayer=TRUE)
lssub("stackIm1b_1.tif", "stackIm1b.sample.tif", centerx = 0, centery
= 0, widthx = 50, widthy = 50)

Sarah
On Sun, Jun 23, 2013 at 4:08 PM, ASANTOS <alexandresantosbr at yahoo.com.br> wrote: