fgeo_topography() fails with edgecorrect = TRUE but passes with FALSE

[maurolepore]

From: Jian Zhang Sent: Monday, September 24, 2018 3:40:25 AM (UTC-05:00) Eastern Time (US & Canada)

I have the trouble to transform my elevation data into the required format and to use the function "fgeo_topography". The elevation data I used is attached.

> tst <- read.table("TianTong_Elev.txt", header=TRUE)
> tst <- tst[order(tst$x, tst$y),]
> elev_df <- fgeo_elevation(tst)
> fgeo_topography(elev_df, gridsize = 20, xdim = 500, ydim = 400)
Error in convex[quad.index] <- midelev - meanelev[quad.index] :
    replacement has length zero

[maurolepore]

Thanks Jian,

The problem seems to be in the chunk of code that corrects for edge effect. Notice that your code (a minimalistic version) fails with edgecorrect = TRUE but passes with edgecorrect = FALSE. I'll work on that.

library(fgeo.tool)
#> 
#> Attaching package: 'fgeo.tool'
#> The following object is masked from 'package:stats':
#> 
#>     filter

tian_tong_url <- "https://gist.githubusercontent.com/maurolepore/915fd7905e07ecc67dac5585f2c142c4/raw/db9bc25078f07072b50c561f2fe589cf394fa867/TianTong_Elev.txt"
elev <- read.table(tian_tong_url, header = TRUE)

# Fails
fgeo_topography(elev, gridsize = 20, xdim = 500, ydim = 400)
#> Error in convex[quad.index] <- midelev - meanelev[quad.index]: replacement has length zero

# Passes
fgeo_topography(elev, gridsize = 20, xdim = 500, ydim = 400, edgecorrect = FALSE)
#> # A tibble: 500 x 5
#>       gx    gy meanelev  convex slope
#>  * <dbl> <dbl>    <dbl>   <dbl> <dbl>
#>  1     0     0     389.  26.8    33.6
#>  2     0    20     378.  14.6    39.4
#>  3     0    40     366.   4.55   36.5
#>  4     0    60     355. -11.6    31.2
#>  5     0    80     366.  -9.67   60.8
#>  6     0   100     378.  -5.56   69.5
#>  7     0   120     377.  -7.95   59.0
#>  8     0   140     375.  -0.484  33.0
#>  9     0   160     368.  -2.93   24.1
#> 10     0   180     361.  -9.98   19.5
#> # ... with 490 more rows

Created on 2018-09-24 by the reprex package (v0.2.1)

[maurolepore]

I would like to test fgeo_topography() with the elevation data that the database delivers. The funciton passes with other elevation datasets but fails with elevation data from Tian Tong that I got from Jian Zang. I wonder if Jian Zang's data is different from that delivered by the database.

library(fgeo.tool)
#> 
#> Attaching package: 'fgeo.tool'
#> The following object is masked from 'package:stats':
#> 
#>     filter

# Works with data from Luquillo
fgeo_topography(fgeo.data::luquillo_elevation, gridsize = 20)
#> # A tibble: 400 x 5
#>       gx    gy meanelev convex slope
#>  * <dbl> <dbl>    <dbl>  <dbl> <dbl>
#>  1     0     0     365. -0.390  8.87
#>  2     0    20     365. -0.285 13.7 
#>  3     0    40     367.  0.210 14.4 
#>  4     0    60     367.  0.34  15.9 
#>  5     0    80     362. -0.255 20.8 
#>  6     0   100     358. -0.510 12.4 
#>  7     0   120     359. -0.28  17.1 
#>  8     0   140     360. -0.265 22.1 
#>  9     0   160     362.  0.055 24.3 
#> 10     0   180     364.  0.415 15.2 
#> # ... with 390 more rows

# Works with data from BCI
fgeo_topography(bciex::bci_elevation, gridsize = 20, xdim = 1000, ydim = 500)
#> # A tibble: 1,250 x 5
#>       gx    gy meanelev  convex slope
#>  * <dbl> <dbl>    <dbl>   <dbl> <dbl>
#>  1     0     0     123. -0.315  13.3 
#>  2     0    20     125.  1.21   17.9 
#>  3     0    40     129.  0.0975  9.65
#>  4     0    60     130. -0.132  10.9 
#>  5     0    80     132.  0.0975 11.9 
#>  6     0   100     133.  0.103  11.0 
#>  7     0   120     136. -0.06   11.1 
#>  8     0   140     138.  0.167   6.20
#>  9     0   160     140.  0.125   3.41
#> 10     0   180     141.  0.0350  2.84
#> # ... with 1,240 more rows

# Fails with data from Tian Tong
tian_tong_elevation_url <- "http://bit.ly/2DA4XMW"
tian_tong_elevation <- read.table(tian_tong_elevation_url, header = TRUE)
fgeo_topography(tian_tong_elevation, gridsize = 20, xdim = 500, ydim = 400)
#> Error in convex[quad.index] <- midelev - meanelev[quad.index]: replacement has length zero

Created on 2018-09-27 by the reprex package (v0.2.1)

[maurolepore]

I think I know what the problem is: Your elevation data is too coarse to use edgecorrect. edgecorrect works only if your elevation data is at a resolution of at least gridsize / 2.

# Fails with data from Tian Tong
tian_tong_elevation_url <- "http://bit.ly/2DA4XMW"
tian_tong_elevation <- read.table(tian_tong_elevation_url, header = TRUE)
lapply(tian_tong_elevation, function(x) head(unique(x)))
#> $x
#> [1]   0  20  40  60  80 100
#> 
#> $y
#> [1]   0  20  40  60  80 100
#> 
#> $elev
#> [1] 397.020 390.921 376.643 363.410 354.427 352.114

Created on 2018-09-27 by the reprex package (v0.2.1)

[maurolepore]

Now fgeo_topography() warns if no elevation data might be too coarse. Closes the issue.

library(fgeo.tool)
#> 
#> Attaching package: 'fgeo.tool'
#> The following object is masked from 'package:stats':
#> 
#>     filter

tian_tong_elevation_url <- "http://bit.ly/2DA4XMW"
tian_tong_elevation <- read.table(tian_tong_elevation_url, header = TRUE)

# Errs but warns that you should use `edgecorrect = FALSE`
fgeo_topography(tian_tong_elevation, gridsize = 20, xdim = 500, ydim = 400)
#> Warning: No elevation data found at `gridsize / 2`. 
#> * Is your elevation data too coarse?
#> * Do you need to use `edgecorrect = FALSE`?
#> Error in convex[quad_idx] <- midelev - meanelev[quad_idx]: replacement has length zero

# After reading the warning, `edgecorrect = TRUE`
fgeo_topography(
  tian_tong_elevation, gridsize = 20, xdim = 500, ydim = 400, edgecorrect = FALSE
)
#> # A tibble: 500 x 5
#>       gx    gy meanelev  convex slope
#>  * <dbl> <dbl>    <dbl>   <dbl> <dbl>
#>  1     0     0     389.  26.8    33.6
#>  2     0    20     378.  14.6    39.4
#>  3     0    40     366.   4.55   36.5
#>  4     0    60     355. -11.6    31.2
#>  5     0    80     366.  -9.67   60.8
#>  6     0   100     378.  -5.56   69.5
#>  7     0   120     377.  -7.95   59.0
#>  8     0   140     375.  -0.484  33.0
#>  9     0   160     368.  -2.93   24.1
#> 10     0   180     361.  -9.98   19.5
#> # ... with 490 more rows

Created on 2018-09-27 by the reprex package (v0.2.1)