Please, visit the README for general information about this package
Create a raster from QuadKey-identified data for a specified map area
and zoom level. A detailed explanation of what tile coordinates
(tileX
, tileY
) and pixel coordinates
(pixelX
, pixelY
) are can be found in
quadkeyr
vignette From QuadKey to Simple Features data.frame and other conversions
and in the official
documentation.
If you are working with a data.frame with a quadkey
column, as shown as INPUT in the example image, you can visualize the
data by combining the QuadKeys’ spatial information with any other
reported variable. For example, let’s read an example dataset:
head(data)
#> X lat lon quadkey date_time n_crisis percent_change
#> 1 1 -34.58574 -58.55438 2103213001233302 2020-04-15 NA 2.863940
#> 2 2 -34.51335 -58.57635 2103213001213202 2020-04-15 NA -2.600991
#> 3 3 -34.59026 -58.57086 2103213001233221 2020-04-15 NA 1.460141
#> 4 4 -34.52692 -58.54340 2103213001231110 2020-04-15 NA 2.607693
#> 5 5 -34.53145 -58.52692 2103213001320003 2020-04-15 NA 3.242862
#> 6 6 -34.52692 -58.58734 2103213001230110 2020-04-15 NA 1.172856
#> day hour
#> 1 2020-04-15 0
#> 2 2020-04-15 0
#> 3 2020-04-15 0
#> 4 2020-04-15 0
#> 5 2020-04-15 0
#> 6 2020-04-15 0
To create the raster, it is preferable to work with a regular grid. That is why we need to ensure that we have all the QuadKeys covering the bounding box of the area we want to study.
The function add_regular_polygon_grid()
retrieves all
the QuadKeys required to fill the bounding box created with the QuadKeys
provided in the quadkey
column. The additional QuadKeys
added in the process will have NA
values for all columns
except quadkey
and geometry
.
If you are unsure whether your grid is complete, you can still run
the function, and it will return the same sf
POLYGON
data.frame without any extra rows.
complete_polygon_grid <- add_regular_polygon_grid(data = data)
complete_polygon_grid
#> $data
#> Simple feature collection with 396 features and 9 fields
#> Geometry type: POLYGON
#> Dimension: XY
#> Bounding box: xmin: -58.60107 ymin: -34.60156 xmax: -58.5022 ymax: -34.50203
#> Geodetic CRS: WGS 84
#> # A tibble: 396 × 10
#> quadkey X date_time n_crisis percent_change day hour
#> * <chr> <dbl> <dttm> <dbl> <dbl> <date> <dbl>
#> 1 210321300… 150 2020-04-15 00:00:00 NA 2.29 2020-04-15 0
#> 2 210321300… 149 2020-04-15 00:00:00 212. 0.848 2020-04-15 0
#> 3 210321300… 148 2020-04-15 00:00:00 NA 28.4 2020-04-15 0
#> 4 210321300… 147 2020-04-15 00:00:00 179. 3.75 2020-04-15 0
#> 5 210321300… 146 2020-04-15 00:00:00 20.3 2.81 2020-04-15 0
#> 6 210321300… 145 2020-04-15 00:00:00 NA -12.8 2020-04-15 0
#> 7 210321300… 144 2020-04-15 00:00:00 127. -5.37 2020-04-15 0
#> 8 210321300… 143 2020-04-15 00:00:00 NA 6.60 2020-04-15 0
#> 9 210321300… 142 2020-04-15 00:00:00 NA 3.16 2020-04-15 0
#> 10 210321300… 141 2020-04-15 00:00:00 NA -2.84 2020-04-15 0
#> # ℹ 386 more rows
#> # ℹ 3 more variables: tileX <dbl>, tileY <dbl>, geometry <POLYGON [°]>
#>
#> $num_cols
#> [1] 18
#>
#> $num_rows
#> [1] 22
As you can see in the following plot, we obtain a regular polygon grid where only the QuadKeys that are reported with values for the variable of interest are colored.
The function add_regular_polygon_grid()
returns three
outputs in a list:
data
num_cols
and;num_rows
The number of columns and rows will be useful information for the raster creation step.
ggplot() +
geom_sf(data = complete_polygon_grid$data,
aes(fill = percent_change)) +
scale_fill_viridis_c(na.value = "transparent") +
theme_minimal()
After obtaining the regular polygon grid,
we can create a STARS
object using create_stars_raster()
.
stars_object <- create_stars_raster(data = complete_polygon_grid$data,
template = complete_polygon_grid$data,
var = "percent_change",
nx = complete_polygon_grid$num_cols,
ny = complete_polygon_grid$num_rows)
stars_object
#> stars object with 2 dimensions and 1 attribute
#> attribute(s):
#> Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
#> percent_change -23.94377 -2.50754 1.109488 0.6168087 3.838209 28.42059 246
#> dimension(s):
#> from to offset delta refsys point x/y
#> x 1 18 -58.6 0.005493 WGS 84 FALSE [x]
#> y 1 22 -34.5 -0.004524 WGS 84 FALSE [y]
If you want to visualize the QuadKeys without creating a raster
image, you also have the option of converting the QuadKeys into an
sf
POLYGON data.frame using the function
quadkey_df_to_polygon()
, and then plotting it.
# create polygon column
poly_data <- quadkey_df_to_polygon(data = data)
poly_data
#> Simple feature collection with 150 features and 9 fields
#> Geometry type: POLYGON
#> Dimension: XY
#> Bounding box: xmin: -58.60107 ymin: -34.60156 xmax: -58.5022 ymax: -34.50203
#> Geodetic CRS: WGS 84
#> First 10 features:
#> X lat lon quadkey date_time n_crisis percent_change
#> 1 1 -34.58574 -58.55438 2103213001233302 2020-04-15 NA 2.8639400
#> 2 2 -34.51335 -58.57635 2103213001213202 2020-04-15 NA -2.6009911
#> 3 3 -34.59026 -58.57086 2103213001233221 2020-04-15 NA 1.4601414
#> 4 4 -34.52692 -58.54340 2103213001231110 2020-04-15 NA 2.6076928
#> 5 5 -34.53145 -58.52692 2103213001320003 2020-04-15 NA 3.2428621
#> 6 6 -34.52692 -58.58734 2103213001230110 2020-04-15 NA 1.1728565
#> 7 7 -34.51787 -58.59283 2103213001212321 2020-04-15 13.07492 -0.3100404
#> 8 8 -34.58574 -58.59833 2103213001232302 2020-04-15 NA 11.4323758
#> 9 9 -34.56764 -58.52142 2103213001322012 2020-04-15 NA -18.2960744
#> 10 10 -34.55407 -58.59833 2103213001230320 2020-04-15 NA -15.7189196
#> day hour geometry
#> 1 2020-04-15 0 POLYGON ((-58.55713 -34.588...
#> 2 2020-04-15 0 POLYGON ((-58.5791 -34.5156...
#> 3 2020-04-15 0 POLYGON ((-58.57361 -34.592...
#> 4 2020-04-15 0 POLYGON ((-58.54614 -34.529...
#> 5 2020-04-15 0 POLYGON ((-58.52966 -34.533...
#> 6 2020-04-15 0 POLYGON ((-58.59009 -34.529...
#> 7 2020-04-15 0 POLYGON ((-58.59558 -34.520...
#> 8 2020-04-15 0 POLYGON ((-58.60107 -34.588...
#> 9 2020-04-15 0 POLYGON ((-58.52417 -34.569...
#> 10 2020-04-15 0 POLYGON ((-58.60107 -34.556...
The function add_regular_polygon_grid()
serves as a
wrapper for three functions: get_qk_coords()
,
regular_qk_grid()
, and grid_to_polygon()
.
We will demonstrate some of the uses of these functions in isolation:
The regular_qk_grid()
function generates a regular grid
by utilizing create_qk_grid()
. To accomplish this, it
necessitates defining the bounding box, specified by xmin
,
xmax
, ymin
, and ymax
, which
delineates the area for QuadKey grid creation.
regular_qk_grid()
supplies these arguments to
create_qk_grid()
by estimating the bounding box from the
QuadKeys provided as an argument. However, if desired, you can manually
input these arguments and use the function by itself.
For this example, we have selected:
xmin
= -59xmax
= -57ymin
= -35ymax
= -34Let’s plot them as points.
The QuadKey grid can have a zoom level between 1 (less detail) to 23 (more detail).
The function create_qk_grid()
will return three
outputs:
grid$data
a dataframe with tileX
,
tileY
and the QuadKey value for each element of the
grid
grid$num_rows
the number of rows and
grid$num_cols
the number of columns of the
grid.
grid <- create_qk_grid(
xmin = -59,
xmax = -57,
ymin = -35,
ymax = -34,
zoom = 12)
head(grid$data)
#> # A tibble: 6 × 3
#> tileX tileY quadkey
#> <dbl> <dbl> <chr>
#> 1 1376 2473 210321302002
#> 2 1377 2473 210321302003
#> 3 1378 2473 210321302012
#> 4 1379 2473 210321302013
#> 5 1380 2473 210321302102
#> 6 1381 2473 210321302103
grid$num_rows
#> [1] 15
grid$num_cols
#> [1] 24
The coordinates can be extracted from the QuadKeys using the function
get_qk_coords()
.
grid_coords <- get_qk_coord(data = grid$data)
head(grid_coords)
#> Simple feature collection with 6 features and 3 fields
#> Geometry type: POINT
#> Dimension: XY
#> Bounding box: xmin: -59.0625 ymin: -34.958 xmax: -58.62305 ymax: -34.958
#> Geodetic CRS: WGS 84
#> # A tibble: 6 × 4
#> tileX tileY quadkey geometry
#> <dbl> <dbl> <chr> <POINT [°]>
#> 1 1376 2473 210321302002 (-59.0625 -34.958)
#> 2 1377 2473 210321302003 (-58.97461 -34.958)
#> 3 1378 2473 210321302012 (-58.88672 -34.958)
#> 4 1379 2473 210321302013 (-58.79883 -34.958)
#> 5 1380 2473 210321302102 (-58.71094 -34.958)
#> 6 1381 2473 210321302103 (-58.62305 -34.958)
We can visualize the points in the map to understand better the results.
We have a grid of points representing the QuadKeys. Remember that these points represent the upper-left corner of each QuadKey, which might give the impression that they do not precisely cover the entire area defined by the initial points.
As we are creating the polygons from a sf
POINT
data.frame we can use grid_to_polygon()
.
polygrid <- grid_to_polygon(grid_coords)
polygrid
#> Simple feature collection with 360 features and 3 fields
#> Geometry type: POLYGON
#> Dimension: XY
#> Bounding box: xmin: -59.0625 ymin: -35.03 xmax: -56.95312 ymax: -33.94336
#> Geodetic CRS: WGS 84
#> # A tibble: 360 × 4
#> # Rowwise:
#> tileX tileY quadkey geometry
#> * <dbl> <dbl> <chr> <POLYGON [°]>
#> 1 1399 2459 210321132133 ((-57.04102 -34.01624, -56.95312 -34.01624, -56.953…
#> 2 1398 2459 210321132132 ((-57.12891 -34.01624, -57.04102 -34.01624, -57.041…
#> 3 1397 2459 210321132123 ((-57.2168 -34.01624, -57.12891 -34.01624, -57.1289…
#> 4 1396 2459 210321132122 ((-57.30469 -34.01624, -57.2168 -34.01624, -57.2168…
#> 5 1395 2459 210321132033 ((-57.39258 -34.01624, -57.30469 -34.01624, -57.304…
#> 6 1394 2459 210321132032 ((-57.48047 -34.01624, -57.39258 -34.01624, -57.392…
#> 7 1393 2459 210321132023 ((-57.56836 -34.01624, -57.48047 -34.01624, -57.480…
#> 8 1392 2459 210321132022 ((-57.65625 -34.01624, -57.56836 -34.01624, -57.568…
#> 9 1391 2459 210321123133 ((-57.74414 -34.01624, -57.65625 -34.01624, -57.656…
#> 10 1390 2459 210321123132 ((-57.83203 -34.01624, -57.74414 -34.01624, -57.744…
#> # ℹ 350 more rows
It worked! As you can see here, the coordinates we randomly selected
as a starting point for the bounding box are situated within the
polygons, but not at a specific position inside each polygon. This was
expected, you can read about Quadkey conversion in the vignette From a QuadKey to a Simple Features data.frame and other conversions
.
If you want to see the grid, you can also check the app: