[R-br] ?Problema para mudar projeção em um raster

[ASANTOS]

Obrigado Mauro,

    Resolvi o erro na rotina e mesmo assim o spTransform() não esta 
convertendo de utm na zona 23K para latlong no hemisfério sul, segue CRM 
completo ate o recorte da área de interesse no geoTiff do topodata:

require(sp)
require(rgdal)
require(raster)
# Criar uma área de menor dimensão que a imagem inteira que abarque a 
região de interesse
xcc<-773759.1
ycc<-7841546
p.central<-cbind(xcc,ycc)

###Criando os vértices da área
coordV<-NULL
coordV 
<-rbind(coordV,cbind(p.central[,1]+c(-1500,1500,1500,-1500,-1500),p.central[,2]+c(1500,1500,-1500,-1500,1500)))
plot(coordV[,1],coordV[,2])
points(p.central[,1],p.central[,2], col="red")
#

# Cria um polígono com o contorno definido
bnds <- cbind(x=c(coordV[,1]), y=c(coordV[,2]))

# CRS UTM
SP <- SpatialPolygons(list(Polygons(list(Polygon(bnds)), "1")))
proj4string(SP) = CRS("+proj=utm +zone=23+south +datum=WGS84 +units=m 
+no_defs") ## Projeção

# CRS em lat long para o polígomno ficar na mesma projeção que o GeoTif
SPlatlong<- spTransform(SP, CRS("+proj=longlat +south +datum=WGS84"))
#
#
#Obtenção do objeto declividade do TOPODATA
url=("http://www.dsr.inpe.br/topodata/data/geotiff/19S42_SN.zip")
download.file(url, destfile = "19S42_SN.zip")
### descompacta
system("unzip 19S42_SN.zip")
demlatlog<-raster('19S42_SN.tif')### Abre o slope TOPODATA
#
#Cortar a região de interesse
demlatlog.crop <- crop(demlatlog, extent(SPlatlong), snap='out')
#

Obrigado,

Alexandre





  Em 22/07/2013 22:09, Mauro Sznelwar escreveu:
> Tentei rodar o seu CRM e deu problema aqui
> > bnds <- cbind(x=c(coordV2[21:25,1]), y=c(coordV2[21:25,2]))
> Error in coordV2[21:25, 1] : subscript out of bounds
> >
>
> Segue CRM:
>
> require(sp)
> require(rgdal)
> # Criar uma área de menor dimensão que a imagem inteira que abarque a
> região de interesse
> xcc<-773759.1
> ycc<-7841546
> p.central<-cbind(xcc,ycc)
>
> ###Criando os vértices da área
> coordV <-
> rbind(coordV,cbind(p.central[,1]+c(-1500,1500,1500,-1500,-1500),p.central[,2]+c(1500,1500,-1500,-1500,1500)))
> coordV
> coordV2=coordV[!is.na(coordV[,1]),]
> #plot(coordV2[21:25,1],coordV2[21:25,2])
> #points(p.central[,1],p.central[,2], col="red")
> #
>
> # Cria um polígono com o contorno definido
> bnds <- cbind(x=c(coordV2[21:25,1]), y=c(coordV2[21:25,2]))
>
> # CRS UTM
> SP <- SpatialPolygons(list(Polygons(list(Polygon(bnds)), "1")))
> proj4string(SP) = CRS("+proj=utm +zone=23+south +datum=WGS84 +units=m
> +no_defs") ## Projeção
>
> # CRS em lat long para o polígomno ficar na mesma projeção que o GeoTif
> SPlatlong<- spTransform(SP, CRS("+proj=longlat +south +datum=WGS84"))
> #
>
>
> Obrigado,
>
> Alexandre
>
>
>
>

-- 
======================================================================
Alexandre dos Santos
Proteção 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)
e-mails:alexandresantosbr em yahoo.com.br
         alexandre.santos em cas.ifmt.edu.br
Lattes: http://lattes.cnpq.br/1360403201088680
======================================================================

-------------- Próxima Parte ----------
Um anexo em HTML foi limpo...
URL: <http://listas.inf.ufpr.br/pipermail/r-br/attachments/20130723/b59bab9b/attachment.html>