Welcome to OGeek Q&A Community for programmer and developer-Open, Learning and Share
Welcome To Ask or Share your Answers For Others

Categories

0 votes
607 views
in Technique[技术] by (71.8m points)

r - NetLogo - misalignment with imported GIS shapefiles

I've got a NetLogo model in which each animal occupies a "territory", wherein all patches that belong to the animal are the same color as the animal. I then use the R extension in NetLogo to create a minimum convex polygon (MCP) surrounding each territory and export those polygons as a shapefile. I then import the GIS file back into NetLogo using the GIS extension. My rationale for this is related to this question I posted before (NetLogo - applying values to patches within polygons). That is, I can use gis:intersecting to give a variable to those patches that fall within the GIS-imported polygons. The process of creating the MCP, exporting, and importing is done at each time step because the territory updates each tick. All works well except that the imported shapefile does not align perfectly with the original polygons. The attached image shows this, where the blue outlines are from the imported shapefile. enter image description here. I've tried gis:set-world-envelope (list min-pxcor max-pxcor min-pycor max-pycor) but to no avail. Does anybody know if I'm doing something wrong or if it is error inherent to exporting and then importing a shapefile for which no projection exists? Any help here would be great because having them align would solve several issues, including the previous post. The whole code is pretty long so I've attached some snippets below that are most relevant. Thanks!

extensions [r gis ] 

breed [females female]

globals
[
  hr-dataset
]    

females-own 
[ 
  Name 
  X 
  Y
]

patches-own 
[ 
  is-hr?
]

to setup 
  clear-all
  r:clear 
  ...
  ask n-of 5 patches 
  [ 
    sprout-females 1
    [ 
      ...
      set X (list pxcor) 
      set Y (list pycor) 
    ] 
  ] 
  reset-ticks 
end 

to go 
  ...
  expand
  calc-homerange
  tick 
end

to expand 
  repeat 10
  [
    ask females
    [
      move-to target
      set X lput pxcor X
      set Y lput pycor Y
    ]
  ]
end

to calc-homerange
  r:eval "library(adehabitatHR)"
  r:eval "library(sp)"
  r:eval "library(rgdal)"

  ; create an empty data.frame"
  r:eval "turtles <- data.frame()"

  ; merge the Name, X- and Y-lists of all females to one big data.frame
  ask females
  [
    (r:putdataframe "turtle" "X" X "Y" Y)     
    r:eval (word "turtle <- data.frame(turtle, Name = '" Name "')")
    r:eval "turtles <- rbind(turtles, turtle)"
  ]

  ; create SpatialPointsDataFrame
  r:eval "spdf <- SpatialPointsDataFrame(turtles[1:2], turtles[3])" 
  r:eval "homerange <- mcp(spdf, percent = 100)"
  r:eval "writeOGR(homerange, '.', 'homerange-rgdal', overwrite_layer = TRUE, driver = 'ESRI Shapefile')"

mark-homeranges
end

to mark-homeranges
  clear-drawing
  ...
  set hr-dataset gis:load-dataset "C:/Program Files (x86)/NetLogo 5.0.4/homerange-rgdal.shp"
  gis:set-world-envelope (list min-pxcor max-pxcor min-pycor max-pycor)    
  gis:set-drawing-color blue
  gis:draw hr-dataset 2
  ask patches gis:intersecting hr-dataset
  [
    set is-hr? true
  ]
end
See Question&Answers more detail:os

与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
Welcome To Ask or Share your Answers For Others

1 Reply

0 votes
by (71.8m points)

I think Seth is correct that it's an off-by-one error in mapping patch coordinates to world coordinates. It may be related to the bug that was mentioned here. His proposed solution of using (list min-pxcor - 0.5 max-pxcor + 0.5 min-pycor - 0.5 max-pycor + 0.5) should work. If not, send me the model and I'll take a look if I have time.


与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
OGeek|极客中国-欢迎来到极客的世界,一个免费开放的程序员编程交流平台!开放,进步,分享!让技术改变生活,让极客改变未来! Welcome to OGeek Q&A Community for programmer and developer-Open, Learning and Share
Click Here to Ask a Question

...