An embedded and charset-unspecified text was scrubbed... Name: not available URL: <https://stat.ethz.ch/pipermail/r-sig-geo/attachments/20130623/1de15b43/attachment.pl>
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:
Dear Members,
I'm having trouble calculating NDVI difference, first the images Geotiff can not view using Ubuntu 4.12 64-bit and therefore can not get to the NDVI, follow an example:
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)
#
#Calculate NDVI difference
multi.espc<-stack(c("stackIm1.sample.tif","stackIm2.sample.tif"))
band3<-raster(multi.espc,1)
band4<-raster(multi.espc,2)
ndvi<-(band4-band3)/(band4+band3)
#
Sarah Goslee http://www.functionaldiversity.org
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, 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:
Dear Members,
I'm having trouble calculating NDVI difference, first the images Geotiff can not view using Ubuntu 4.12 64-bit and therefore can not get to the NDVI, follow an example:
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)
#
#Calculate NDVI difference
multi.espc<-stack(c("stackIm1.sample.tif","stackIm2.sample.tif"))
band3<-raster(multi.espc,1)
band4<-raster(multi.espc,2)
ndvi<-(band4-band3)/(band4+band3)
#
======================================================================
Alexandre dos Santos
Prote??o Florestal
Coordenador do curso T?cnico em Florestas
Vice Coordenador do curso de Engenharia Florestal
IFMT - Instituto Federal de Educa??o, Ci?ncia e Tecnologia de Mato Grosso
Campus C?ceres
Caixa Postal 244
Avenida dos Ramires, s/n
Bairro: Distrito Industrial
C?ceres - MT CEP: 78.200-000
Fone: (+55) 65 8132-8112 (TIM) (+55) 65 9686-6970 (VIVO)
alexandre.santos at cas.ifmt.edu.br
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:
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, 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:
Dear Members,
I'm having trouble calculating NDVI difference, first the images Geotiff
can not view using Ubuntu 4.12 64-bit and therefore can not get to the NDVI,
follow an example:
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)
#
#Calculate NDVI difference
multi.espc<-stack(c("stackIm1.sample.tif","stackIm2.sample.tif"))
band3<-raster(multi.espc,1)
band4<-raster(multi.espc,2)
ndvi<-(band4-band3)/(band4+band3)
#
Sarah Goslee http://www.functionaldiversity.org