|
| 1 | +#' Calculate the centroid of a polygon |
| 2 | +#' |
| 3 | +#' This function calculates the centroid of a polygon and returns the latitude, longitude and uncertainty of the centroid. |
| 4 | +#' |
| 5 | +#' @param sf_df A sf object with polygons |
| 6 | +#' @param id A character string with the name of the column containing the unique identifier |
| 7 | +#' |
| 8 | +#' @return A data frame with the unique identifier, latitude, longitude and uncertainty of the centroid |
| 9 | +#' |
| 10 | +#' @examples |
| 11 | +#' \dontrun{ |
| 12 | +#' # Example of how to use the calculate_polygon_centroid function |
| 13 | +#' # Load the necessary data |
| 14 | +#' boswachterijen <- boswachterijen$boswachterijen_2024 |
| 15 | +#' |
| 16 | +#' # add a unique identifier to the sf object |
| 17 | +#' boswachterijen <- boswachterijen %>% |
| 18 | +#' dplyr::mutate(UUID = as.character(row_number())) |
| 19 | +#' |
| 20 | +#' # Calculate the centroid of the polygons |
| 21 | +#' centroids_data_final <- calculate_polygon_centroid(sf_df = boswachterijen, id = "UUID") |
| 22 | +#' |
| 23 | +#' # Plot the polygons and the centroids |
| 24 | +#' library(leaflet) |
| 25 | +#' |
| 26 | +#' # Sample 1 polygon and 1 centroid to plot using id |
| 27 | +#' sample_id <- sample(centroids_data_final$UUID, 1) |
| 28 | +#' |
| 29 | +#' leaflet() %>% |
| 30 | +#' addProviderTiles("CartoDB.Positron") %>% |
| 31 | +#' addPolygons(data = boswachterijen %>% dplyr::filter(UUID == sample_id), |
| 32 | +#' weight = 1, color = "black", fillOpacity = 0.5) %>% |
| 33 | +#' addCircles(data = centroids_data_final %>% dplyr::filter(UUID == sample_id), |
| 34 | +#' lat = ~centroidLatitude, lng = ~centroidLongitude, radius = 5, |
| 35 | +#' color = "black") %>% |
| 36 | +#' addCircles(data = centroids_data_final %>% dplyr::filter(UUID == sample_id), |
| 37 | +#' lat = ~centroidLatitude, lng = ~centroidLongitude, radius = ~centroidUncertainty, |
| 38 | +#' color = "red", weight = 1) |
| 39 | +#' } |
| 40 | +#' |
| 41 | +#' @export |
| 42 | +#' @author Sander Devisscher |
| 43 | + |
| 44 | +calculate_polygon_centroid <- function(sf_df, id){ |
| 45 | + # Checks #### |
| 46 | + ## Check if the input is an sf object #### |
| 47 | + if(!inherits(sf_df, "sf")){ |
| 48 | + stop("The input should be an sf object") |
| 49 | + } |
| 50 | + |
| 51 | + ## Check if the id is a character string #### |
| 52 | + if(!is.character(id)){ |
| 53 | + id <- as.character(id) |
| 54 | + } |
| 55 | + |
| 56 | + ## Check if the id is in the sf object #### |
| 57 | + if(!(id %in% names(sf_df))){ |
| 58 | + stop("The id is not in the sf object") |
| 59 | + } |
| 60 | + |
| 61 | + ## Check if the id is unique #### |
| 62 | + if(length(unique(sf_df[[id]])) != nrow(sf_df)){ |
| 63 | + warning("The id is not unique >> the function will continue but the output will be incorrect >> try to add a unique identifier to the sf object") |
| 64 | + } |
| 65 | + |
| 66 | + # prepare data #### |
| 67 | + ## Rename the id column to "id" #### |
| 68 | + id_col <- id |
| 69 | + names(sf_df)[names(sf_df) == id_col] <- "id" |
| 70 | + |
| 71 | + ## Extract the CRS #### |
| 72 | + crs_wgs <- CRS_extracter("wgs") |
| 73 | + |
| 74 | + ## Transform the data to the correct CRS #### |
| 75 | + sf_df <- sf_df %>% |
| 76 | + sf::st_transform(crs_wgs) |
| 77 | + |
| 78 | + ## Calculate the number of vertices #### |
| 79 | + sf_df <- sf_df %>% |
| 80 | + sf::st_make_valid() %>% |
| 81 | + dplyr::mutate(NbrVertex = mapview::npts(sf_df, by_feature = TRUE)) |
| 82 | + |
| 83 | + # Caculate Centroids #### |
| 84 | + ## Calculate centroids from sp_df #### |
| 85 | + centroids <- sf_df %>% |
| 86 | + sf::st_centroid() |
| 87 | + |
| 88 | + ## Create output #### |
| 89 | + centroids_data_final <- data.frame() |
| 90 | + |
| 91 | + UUIDS <- unique(sf_df$id) |
| 92 | + |
| 93 | + ## Create a progress bar #### |
| 94 | + pb <- progress::progress_bar$new(format = " [:bar] :percent ETA: :eta", |
| 95 | + total = nrow(sf_df), |
| 96 | + clear = FALSE, |
| 97 | + width = 60) |
| 98 | + |
| 99 | + ## Calculate the distance between the centroid and the polygon #### |
| 100 | + for(u in UUIDS){ |
| 101 | + ### Update the progress bar #### |
| 102 | + pb$tick() |
| 103 | + ### Filter the sf data #### |
| 104 | + sf_df_sub <- sf_df %>% |
| 105 | + dplyr::filter(id == u) |
| 106 | + ### Check if the polygon is valid #### |
| 107 | + if(nrow(sf_df_sub) == 0){ |
| 108 | + next |
| 109 | + warning(paste0("no fortified shape for ", u)) |
| 110 | + }else{ |
| 111 | + ### split the polygons into vertrex points #### |
| 112 | + sf_df_sub <- sf_df_sub %>% |
| 113 | + st_cast("MULTIPOINT") %>% |
| 114 | + st_cast("POINT", do_split = TRUE) |
| 115 | + |
| 116 | + ### Check if the number of points is equal to the number of vertices #### |
| 117 | + if(nrow(sf_df_sub) != unique(sf_df_sub$NbrVertex)){ |
| 118 | + warning(paste0("The number of points is not equal to the number of vertices for ", u)) |
| 119 | + } |
| 120 | + } |
| 121 | + |
| 122 | + ### Filter the centroid data #### |
| 123 | + centroids_sub <- centroids %>% |
| 124 | + dplyr::filter(id == u) |
| 125 | + |
| 126 | + ### Check if the centroid is valid #### |
| 127 | + if(nrow(centroids_sub)==0){ |
| 128 | + next |
| 129 | + warning(paste0("no centroid for ", u)) |
| 130 | + } |
| 131 | + ### Calculate the distance #### |
| 132 | + distance <- st_distance(sf_df_sub, centroids_sub) %>% |
| 133 | + units::drop_units() |
| 134 | + |
| 135 | + ### Calculate the maximum distance #### |
| 136 | + maxDistance <- round(max(distance, na.rm = TRUE)) |
| 137 | + |
| 138 | + ### Set the maximum distance to 4 if it is smaller than 4 #### |
| 139 | + if(maxDistance < 4){ |
| 140 | + warning(paste0("The maximum distance is smaller than 4 for ", u, " >> setting the maximum distance to 4")) |
| 141 | + maxDistance <- 4 # reasonable accuracy of handheld GPS devices |
| 142 | + } |
| 143 | + |
| 144 | + ### Add the maximum distance to the centroid data #### |
| 145 | + centroids_sub$centroidUncertainty <- maxDistance |
| 146 | + centroids_data_final <- rbind(centroids_data_final, centroids_sub) |
| 147 | + } |
| 148 | + |
| 149 | + ## Transform the data to a data frame #### |
| 150 | + centroids_data_final <- centroids_data_final %>% |
| 151 | + dplyr::mutate(centroidLatitude = sf::st_coordinates(geometry)[, 2], |
| 152 | + centroidLongitude = sf::st_coordinates(geometry)[, 1]) %>% |
| 153 | + dplyr::select(id, |
| 154 | + centroidLatitude, |
| 155 | + centroidLongitude, |
| 156 | + centroidUncertainty) %>% |
| 157 | + sf::st_drop_geometry() |
| 158 | + |
| 159 | + ## Rename the id column to the original name #### |
| 160 | + names(centroids_data_final)[names(centroids_data_final) == "id"] <- id_col |
| 161 | + |
| 162 | + ## Return the data #### |
| 163 | + return(centroids_data_final) |
| 164 | +} |
0 commit comments