0

After some extensive googling, I wasn't able to find my answer (first time I couldn't surmount the issue by looking at others questions/answers). I am new to asking questions, so forgive any missteps.

I am attempting to perform what ArcGIS or QGIS performs with the slope tool, just within R. To do so I have been importing a raster that I exported from ArcGIS in GRID format with the following characteristics:

class      : RasterLayer 
dimensions : 821, 581, 477001  (nrow, ncol, ncell)
resolution : 4.996121, 4.996121  (x, y)
extent     : 2832147, 2835049, 14234048, 14238150  (xmin, xmax, ymin, ymax)
crs        : +proj=tmerc +lat_0=34.75 +lon_0=-118.583333333333 +k=0.9999 +x_0=800000.000000001 +y_0=3999999.99999999 +datum=NAD83 +units=us-ft +no_defs 
source     : rr_2020_shell 
names      : rr_2020_shell 
values     : 5623.253, 6401.356  (min, max)

It is already projected in the correct coordinate system (EPSG: 3423) but when I go to find the slope using the following code:

RR_2020_Slope = terrain(RR_2020_St1_Raster,'slope', units = 'degrees', neighbors = 8, filename = 'RR_2020_Slope.grd', overwrite = T)

The result is a slope raster that ranges from 0 to 1.28°, which is very different from what I have calculated in ArcGIS using the slope tool. Using the same DEM raster in the same projection in ArcGIS I used the slope tool with an input of 'Degree' for the output measurement, 'Planar' for the method, and 1 for Z factor and my resulting slope raster ranges from 0.001 to 73.396°.

Overall I am wondering where my mistake in R originates from, is it an elevation resolution problem? Are there issues with my projection? Forgive me, I can't necessarily include the data as they are sensitive materials but perhaps there is a clear and obvious mistake in my approach or assumptions about the functions I have used?

1 Answers1

2

The only red flag I see is that you say "it is already projected in the correct coordinate system". Projecting raster data degrades the quality. As cell values get smoothed, the slopes will get smaller. This may be particularly pronounced if the relief is at the scale of the cell size (e.g. sand dunes vs mountain chains). Have you compared with what you get with the original data?

Another source of error could be that the units of the values are different from the units of the coordinate reference system. But it would appear that in your case both are in feet.

Can you also try this with terra::terrain()?

Robert Hijmans
  • 40,301
  • 4
  • 55
  • 63
  • Hi there, thanks for your answer, I'm not sure this would provide enough of an error as the cell size is 5x5ft and the relief is generally 30ft or more. Though, I'd love to find a way to check, though I'm not quite sure how to remove the raster projection. I tried terra::terrain() but received a similar result. – William Taylor Jan 15 '22 at 18:27
  • It sounded like you changed the projection. If that is not the case, then you can ignore that comment. Without data one can only speculate. Perhaps you can provide a small part of your data? 10 x 10 cells would be plenty. And you could add some noise to it. – Robert Hijmans Jan 16 '22 at 20:09
  • Thanks for the tip, let me work on cutting a piece of the data out and I'll edit my original post. – William Taylor Jan 18 '22 at 16:32