0

I need to convert a shape file into a raster and have absolutely no clue where to even start.

If anyone could help me out I'd really appreciate!

Update: I have found out about the 'readOGR'-function, but everytime I use it I get the following message:

Warning messages:
1: In readOGR(dsn = "C:/Users/Giaco/Dropbox/Random Walk/Waterbodies",  :
  Dropping null geometries: 308, 309
2: In readOGR(dsn = "C:/Users/Giaco/Dropbox/Random Walk/Waterbodies",  :
  Z-dimension discarded

Can somebody tell me what that means?

Edit:

altdata <- raster("altitude.tif")        
plot(altdata)
Lotic <- readOGR(dsn="C:/Users/Giaco/Dropbox/Random Walk/Waterbodies",layer="Lotic")
Lentic <- readOGR(dsn="C:/Users/Giaco/Dropbox/Random Walk/Waterbodies",layer="Lentic")

How can I plot the raster "altdata", the SpatialPointsDataFrame "Lentic" and the SpatialLinesDataFrame "Lotic" all in one plot ?

Edit:

altdata
class       : RasterLayer 
dimensions  : 1286, 963, 1238418  (nrow, ncol, ncell)
resolution  : 15, 15  (x, y)
extent      : 90938.29, 105383.3, 190000, 209290  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=tmerc +lat_0=39.66666666666666 +lon_0=1 +k=1 +x_0=200000 +y_0=300000 +ellps=intl +pm=-9.131906111111112 +units=m +no_defs 
data source : C:\Users\Giaco\Dropbox\Random Walk\altitude.tif 
names       : altitude 
values      : -32768, 32767  (min, max)

> Lentic 
class       : SpatialPointsDataFrame 
features    : 182 
extent      : -108473.2, -95455.86, -107870.9, -91344.22  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=tmerc +lat_0=39.66825833333333 +lon_0=-8.133108333333334 +k=1 +x_0=0 +y_0=0 +ellps=GRS80 +units=m +no_defs 
variables   : 3
names       : Presence,     Type, Accessible 
min values  :        0, Fountain,          0 
max values  :        1,     Well,          2 
> Lotic
class       : SpatialLinesDataFrame 
features    : 317 
extent      : -108956.5, -93832.44, -108979.5, -90747.34  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=tmerc +lat_0=39.66825833333333 +lon_0=-8.133108333333334 +k=1 +x_0=0 +y_0=0 +ellps=GRS80 +units=m +no_defs 
variables   : 1
names       : Presence 
min values  :        0 
max values  :        1 
snoops
  • 37
  • 7
  • You need to provide more information and need to see your code so syntax errors can be ruled out. What exactly do you want to represents as a raster, is your data polygon or point and the values numeric or categorical, does there continuous coverage over the area? From the gdal errors is looks like you have null geometries in the data and that the shapefile has a z value. I do not believe that either of these warnings are catastrophic and you should still be getting an sp object. Because you did not provide code I cannot have you check output. – Jeffrey Evans Apr 15 '16 at 14:39
  • Thanks for your answer. It's my first time working with shape files, so please excuse if my questions seem a little stupid. I have two shape files one contains points and the other one lines. I have succeeded in loading them into R and also in plotting them. Now I would like to add them to a raster plot. I will edit my question so you can understand better. Thanks in advance. – snoops Apr 15 '16 at 14:50

1 Answers1

2

Using base plot, there are two way to overlay vector data on a raster. First plot the raster, then you can either call plot for each feature class using the add=TRUE argument. Alternately, you can use the points and lines functions which will also add to the current plot.

Create some example data

library(raster)
library(sp)

 r <- raster(nrows=180, ncols=360, xmn=571823.6, xmx=616763.6, ymn=4423540, 
             ymx=4453690, resolution=270, crs = CRS("+proj=utm +zone=12 +datum=NAD83 
             +units=m +no_defs +ellps=GRS80 +towgs84=0,0,0"))
 r[] <- runif(ncell(r))
 pts <- sampleRandom(r, 10, na.rm = TRUE, sp = TRUE)

Plot the raster then call plot again, with add=TRUE, to add points.

plot(r)
  plot(pts, pch=20, cex=2, col="red", add=TRUE)

Or, plot the raster and use points to add the point feature class.

plot(r)
  points(pts, pch=20, cex=2, col="red")

Edit: Your extents between the raster and vector feature classes do not overlap.

We can create SpatialPolygons using the extent of your objects and an example raster (with a uniform value of 1).

library(raster)
proj <- sp::CRS("+proj=tmerc +lat_0=39.66666666666666 +lon_0=1 +k=1 +x_0=200000 +y_0=300000  
                 +ellps=intl +pm=-9.131906111111112 +units=m +no_defs")
r.ext <- as(extent(90938.29, 105383.3, 190000, 209290), "SpatialPolygons")
  proj4string(r.ext) <- proj
pt.ext <- as(extent(-108473.2, -95455.86, -107870.9, -91344.22), "SpatialPolygons")
  proj4string(pt.ext) <- proj
line.ext <- as(extent(-108956.5, -93832.44, -108979.5, -90747.34),  "SpatialPolygons")
    proj4string(line.ext) <- proj
r <- raster(r.ext, resolution = c(15,15), crs = proj)
  r[] <- rep(1, ncell(r))

Here we see that if we plot the raster and then the point and line extent polygons, you cannot see them.

plot(r)
  plot(pt.ext, add=TRUE)
  plot(line.ext, add=TRUE)

However, if we plot the line and point extent polygons they overlay just fine.

plot(line.ext)
  plot(pt.ext,add=TRUE)

If we limit the raster extent to the extent of the line object, you should see the raster sub-region, but do not. And, if you try to crop the raster you will receive an "Error in .local(x, y, ...) : extents do not overlap" error.

plot(line.ext)
  plot(r, add=TRUE)
Jeffrey Evans
  • 2,325
  • 12
  • 18
  • Thanks ! Unfortunately none of is working for me. Probably this is because I am not using vector data but spatial data as I mentioned above...Have you any other idea ? – snoops Apr 16 '16 at 11:07
  • In this case vector and spatial are interchangeable. The term "vector" is merely referring to spatial objects that are represented as points, polygons or lines. Since you do not provide any information on what is going wrong I really cant help you. I can however, speculate that your data is in a different coordinate space therefore, does not overlay. Look at the proj4string of each object to see if they match. Did you work through my example, which are spatial objects (raster and points). – Jeffrey Evans Apr 16 '16 at 14:40
  • Yes I have looked at your example and it worked perfectly fine. If I trie the same however I don't get any result. The coordinate space seems to match, but I am no expert. I will edit my quuestion so you can have a look yourself (altdata is the raster to which I need to add the two other objects). Is there any way I can upload the documents that I am using so you could try yourself ? In any case thanks again for your time and answers I really appreciate it ! – snoops Apr 17 '16 at 18:56
  • @snoops, please see my edit. The extents between your raster and vector feature classes do not overlap. This is why the forum guidelines request background information and a reproducible example. Had you previously provided the information in your edit I could have saved you a few days by point out the you had an extent issue. – Jeffrey Evans Apr 17 '16 at 19:35
  • I am pretty new to the forum, I was hoping to be able to upload the data I am using so people would understand better. In addition I am pretty new to R too and sometimes it is not so easy to create a reproducible example for me. But thanks for the heads up, I will try to do better next time ! Again, thanks a lot for your answer, you really helped me out ! – snoops 11 mins ago – snoops Apr 17 '16 at 20:39
  • @snoops, I know that it can be a challenge to provide a question that provides all of the relevant information but, a well formed question yields a better response, or any for that matter. Abbreviated, vague questions tend to go unanswered and get down-voted. Your question, as it reads now, is the minimum of what we are looking for. Unfortunately, there is no forum repository so you have to use something like dropbox or Google Drive for distributing data. For issues related to spatial data I would recommend posting on GIS Stack Exchange rather than Stackoverflow. – Jeffrey Evans Apr 17 '16 at 21:22