2

I have a satellite image as a PNG with known coordinates (bottom, left, top, right). I want to draw a shapefile on it (the limits of a county) in R.

My first thought was to convert the PNG (it's read as a three-layered matrix, for red, green and blue) to a raster. But creating a three-layered raster may be an overkill.

And simpler solution to do this?

I found the function rasterImage(), and I'm doing like that:

Map0 <- readPNG('myMap0.png')
cM <- c(xmin,ymin,xmax,ymax)
png('myMap.png',10240,8610)
par(mar=c(0,0,0,0))
plot(myShapeFile,border=NA,xlim=c(cM[1],cM[3]),ylim=c(cM[2],cM[4]))
rasterImage(Map0,cM[1],cM[2],cM[3],cM[4])
plot(myShapeFile,add=T,border='white',lwd=5)
dev.off()

It's almost fine, except for a white border around the image. Since I used the dimensions of the original image (10240x8610), this border is probably distorting the new image, because the painted area is smaller than the original. I had this problem before, but could never solve it.

I thought the problem was with png(), so I did the following:

par(mar=c(0,0,0,0),bg=NA)
plot(myShapeFile,border='white',xlim=c(cM[1],cM[3]),ylim=c(cM[2],cM[4]),lwd=5)
dev.print(width=10240,height=8610,'myBorder.png',dev=png)

So I could try to merge the two images as layers in Gimp, but I see that the border remains, causing the shapefile to be drawn in the wrong position (it aligns perfectly with myMap.png, but not with myMap0.png). In Google Earth I know that the coordinates in cM are precise, so, since I'm giving mar=0, and xlim and ylim should refer to the limits of the plot, why this margin keeps on disturbing the image created?

Following @mdsumner comment, I tried the following:

map0 <- brick('myMap0.png')
extent(map0) <- c(cM[1],cM[3],cM[2],cM[4])
png('myMap.png',10240,8610)
par(mar=c(0,0,0,0))
plotRGB(map0,maxpixels=10240*8610)
plot(myShapeFile,add=T,border='white',lwd=5)
dev.off()

Now it's almost right, except for the left and right white margins, each 27 pixels wide).

  • 1
    @Aaron Please see my edit. – Rodrigo Aug 30 '17 at 1:58
  • Try raster package. Map0 <- brick('myMap0.png'); extent(Map0) <- c(xmin, xmax, ymin, ymax); plotRGB(Map0); plot(myShapeFile, add = TRUE) – mdsumner Aug 30 '17 at 3:17
  • @mdsumner Still no success. Please see my edit. – Rodrigo Aug 30 '17 at 3:34
  • I reckon your extents are wrong, or your shapefile is in a different CRS. Fudging map data in a PNG with free floating numbers to georeference it from scratch is a very tenuous situation, and to work through possible issues will require a lot more information. At minimum, include the plot/s you see, and print summaries of the shapefile object and the raster object - the ones you plot with. – mdsumner Aug 30 '17 at 15:14
  • @mdsumner You' right. I have re-made the calculations, and there was a small error in the bottom coordinate (since it's a mosaic of Google Maps imagery, I had to remove the 50px of the logo, and forgot to update the coordinate accordingly). You may want to include your suggestion as the accepted answer. Thank you! – Rodrigo Aug 30 '17 at 15:48

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.