【问题标题】:R: Creating a covariate on a projected coordinate systemR:在投影坐标系上创建协变量
【发布时间】:2017-05-24 22:38:53
【问题描述】:

我正在尝试使用包 Spatstat 创建协变量,我有一个 shp 文件作为多边形,以便在点位于多边形之外时为协变量分配一个值。我的问题是坐标投影有问题,因为当我完成这个过程时,我有一个奇怪的数字,有人可以帮助我吗?我的代码:

#Creating covariate Z:

W<- owin(xrange=c(767768.3,768773.6),yrange=c(517335.8,518044.8))
xp<-c(768773.6) #cambia!
yp<-c(518044.8)
X<-ppp(x=xp,y=yp,W)
par(mfrow=c(2,2))
Z<-density(X,500)/(7e-06)

#reading and projecting the shp:

Prodes1<-readOGR(dsn="forestsseparate.shp", layer="forestsseparate", dropNULLGeometries=TRUE)
projection(Prodes1)<-CRS("+proj=utm +zone=18 +ellps=intl +towgs84=307,304,-318,0,0,0,0 +units=m +no_defs")

#Extracting the shp polygon coordinates:

xcoords<-rev(Prodes1@polygons[[1]]@Polygons[[1]]@coords[,1])
ycoords<-rev(Prodes1@polygons[[1]]@Polygons[[1]]@coords[,2])
pruebalu<-owin(poly=list(x=xcoords,y=ycoords))

mx<-(Z$xrange[2]-Z$xrange[1])/Z$xstep
my<-(Z$yrange[2]-Z$yrange[1])/Z$ystep

xcoorde<-matrix(NA,nrow=1,ncol=mx)
ycoorde<-matrix(NA,nrow=1,ncol=my)

#Finding each pixel center

for(i in 1:mx){
  xcoorde[i]<-(Z$xrange[1]+i*Z$xstep) -Z$xstep/2
}

for(i in 1:my){
  ycoorde[i]<-(Z$yrange[1]+i*Z$ystep) -Z$ystep/2
}

#Changing the pixel value if is outside the polygon:

for(i in 1:length(Z$xcol)){
    for(j in 1: length(Z$yrow)){

        if(inside.owin(xcoorde[i],ycoorde[j], pruebalu)==FALSE){
            whichPixel<-nearest.raster.point(xcoorde[i], ycoorde[j], Z)
            Z$v[whichPixel$col,whichPixel$row]<-0.4
        }
    }
}

这是我的结果:

【问题讨论】:

    标签: r pixel coordinate-systems spatstat


    【解决方案1】:

    对不起,我弄错了 Z$v 元素的顺序,这里是更正:

    for(i in 1:length(Z$xcol)){
        for(j in 1: length(Z$yrow)){
    
            if(inside.owin(xcoorde[i],ycoorde[j], pruebalu)==FALSE){
                whichPixel<-nearest.raster.point(xcoorde[i], ycoorde[j], Z)
                Z$v[whichPixel$row,whichPixel$col]<-0.4
            }
        }
    }
    

    现在完美了!

    【讨论】:

      猜你喜欢
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 2010-11-14
      • 1970-01-01
      • 1970-01-01
      • 2020-11-04
      相关资源
      最近更新 更多