0

Folks,

I'm in need of some expert opinions. I've spent the majority of my day fighting this problem and I am out of options.

I am trying to create a map that depicts time since last fire across the landscape

Obviously some spots will have never had a fire, those are represented as 0 here, but it's somewhat arbitrary

My code is as you see below:

Data can be found HERE

#example for SO
library(rgdal)
library(raster)
library(plotKML)

##load Data
Fire <- readOGR(dsn="***YOURPATHHERE******/H_Fire_Poly_Jan3_2105", "Fire_poly_KootClip")


##build raster1
r<-raster(res=300,extent(Fire))

pop<-rasterize(Fire,r,"FIRE_YEAR",fun="max",na.rm=TRUE,background=2015)
pop<-projectRaster(pop,res=300, crs="+proj=utm +zone=11 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs")

plot(pop)


##build raster 2
background<-raster(extent(pop),res=300,vals=2015)
projection(background)<-CRS(proj4string(pop))

##make sure all is good
plot(background)
plot(pop,add=T)

##algebra
final<-background-pop

##plot
plot(final)

##export
plotKML(final)

If you don't want to run the code, this is what it looks like on google earth x

This is the main issue, where:

1) it almost looks like the polygon border information is being retained, as the tan lines are the exact borders,

2)the "max" is not always being reported here, especially around the edges it appears?

3) what my goal is, is to use this map in conjunction with other rasters in species distribution modelling (so use raster stack, model and predict function to make a map). When I try to make a map with this layer, I end up with a pile of holes due to NA's...where there should be no NA's.

If I run the models without this time since fire layer, all is good. When I include it, the NA's in the resulting map are generally where you see the tan "polygon outlines" in the first pic.

How do I get this raster layer to have crisp edges and not be such a mess as it currently is???

this is what my final map species map looks like, not the missing data mirrors the issues in the fire map

c

Glorfindel
  • 21,988
  • 13
  • 81
  • 109
ctlamb
  • 131
  • 1
  • 4
  • 14

1 Answers1

0

I think your confusion comes from the visualization in Google Earth, not from the computations. This looks all good when you plot it in R. You should not use projectRaster, however.

library(rgdal)
library(raster)
Fire <- shapefile("Fire_poly_KootClip")
r <- raster(res=300,extent(Fire))
pop <- rasterize(Fire, r, "FIRE_YEAR", fun="max", na.rm=TRUE, background=2015)
x <- pop-2015
plot(x)

Simpler would be to do

Fire$since <- Fire$FIRE_YEAR - 2015
x <- rasterize(Fire, r, "FIRE_YEAR", fun="max", na.rm=TRUE, background=-100)
plot(x)

not the -100 as no fire would be similar to 'a long time ago', not to 'just now'.

Robert Hijmans
  • 40,301
  • 4
  • 55
  • 63