Search code examples
rgisrastershapefiler-raster

georeference png to a shapefile in R

New to GIS, apologies if this is basic. I have a png of a radio signal strength map:

enter image description here

I want to extract county-level data on average signal strength in China.

I can load a shapefile of Chinese counties and can create a raster file of signal strength, but I don't know how to link the two so I can get the signal strength in each county. I have to do this to thousands of pngs. Thank you!

library(raster)
ChinaRD3<-getData('GADM',country="CHN",level=3) #county level
library(rgeos)
ChinaRD3<-gSimplify(ChinaRD3,tol=0.01, topologyPreserve=TRUE)

raster = raster(paste(datadir,"screenshot2.png",sep="")) # the png in my post
v <- extract(raster, ChinaRD3)

# simplify to display mean values
output = unlist(lapply(v, function(x) if (!is.null(x)) mean(x, na.rm=TRUE) else NA ))
Solution

  • This is how you can georeference the image

    library(terra)
r <- rast("https://i.sstatic.net/auWQt.png")
# you can drop the alpha channel 
r <- r[[1:3]]
# declare the R,G,B channels
RGB(r) <- 1:3

# remove the areas around the earth
# this could be made more precise, I am sure
x <- crop(r, c(0.00683, 0.872, 0.041, 0.984))

# set the extent
ext(x) <- c(-180,180,-90,90)

# show that it works
library(geodata)
chn1 <- geodata::gadm(country="CHN", level=1, path=".")

plot(z)
lines(chn1)

    enter image description here

    Now comes the hard part. You need to translate the colors to values. That would not have been so hard if there was just the legend with or without shading and the country borders, gridlines and coastlines. But if you zoom in you can see that the borders are blurred (perhaps these were jpegs earlier in their life?), creating very many colors. With terra you can reduce the colors to <= 255, and then perhaps select the dominant colors. But there probably are better tools (image magick?) to first fix the image so that it only has the colors you care about.

    z <- colorize(x, "col")
f <- freq(z)[,-1]
f <- f[rev(order(f[,2])), ]
head(f)