Commit 12611ff1 authored by luroth's avatar luroth
Browse files

change request: automatic GCP placement based on frequencies

parent bdd91d6e
......@@ -67,11 +67,6 @@ calc_depth_of_field_far <- function(f, aperture, circle_of_confusion, focus_dist
return(ifelse(focus_distance <= hyperfocal_distance, (focus_distance / ( (f - focus_distance) / (hyperfocal_distance - f) + 1)), Inf))
}
#calc_shutter_speed <- function(exposure_value, aperture, iso=100) {
# return(as.numeric(iso)/(25.0 * 2.0^(2.0+as.numeric(exposure_value))/(as.numeric(aperture)^2.0)))
#}
calc_shutter_speed <- function(exposure_value, aperture, iso=100) {
return( (25.0 * aperture^2.0 * 2.0^(2-exposure_value) ) / iso )
......@@ -175,6 +170,67 @@ calc_pixel_freq_xy <- function(plot_size_x, plot_size_y, spacing_x, spacing_y, n
return(list(pixel_frequency_data_xy, pixel_frequency_data_x, pixel_frequency_data_y))
}
calc_gcp_positions <- function(mapping_area_x, mapping_area_y, gcp_n_in_x, gcp_n_in_y, gcp_arrangement_pattern) {
gcp_x <- rep(seq(from=0, to=mapping_area_x, length.out= gcp_n_in_x), each=gcp_n_in_y)
gcp_y <- rep(seq(from=0, to=mapping_area_y, length.out= gcp_n_in_y), gcp_n_in_x)
# Crossproduct according to pattern
if(gcp_arrangement_pattern == "quad") {
gcp <- data.frame(x=gcp_x, y=gcp_y)
} else {
gcp_y <- ifelse(rep(seq(from=1, to=gcp_n_in_x, by = 1), each=gcp_n_in_y) %% 2 != 0, gcp_y, gcp_y + (mapping_area_y/((gcp_n_in_y-1)*2)))
gcp <- data.frame(x=gcp_x, y=gcp_y)
gcp <- gcp %>% filter(y <= mapping_area_y)
}
return(gcp)
}
calc_gcp_recover_frequency <- function(mapping_area_x, mapping_area_y, gcp_n_in_x, gcp_n_in_y, gcp_arrangement_pattern, field_of_view_x, field_of_view_y, photo_positions) {
gcp <- calc_gcp_positions(mapping_area_x, mapping_area_y, gcp_n_in_x, gcp_n_in_y, gcp_arrangement_pattern)
gcp_x_min <- gcp[,1] - field_of_view_x/2
gcp_x_max <- gcp[,1] + field_of_view_x/2
gcp_y_min <- gcp[,2] - field_of_view_y/2
gcp_y_max <- gcp[,2] + field_of_view_y/2
# Get exposure stations
photos_x <- photo_positions[,1]
photos_y <- photo_positions[,2]
# Template for the vector we excpect to get returned
return_template <- rep(TRUE, length(photos_x))
# Test if exposure station is in ground field of view
counts <- vapply(gcp_x_min, FUN = function(x){x <= photos_x}, return_template) &
vapply(gcp_x_max, FUN = function(x){x >= photos_x}, return_template) &
vapply(gcp_y_min, FUN = function(x){x <= photos_y}, return_template) &
vapply(gcp_y_max, FUN = function(x){x >= photos_y}, return_template)
# Count hits
hit_per_gcp <- data.frame("x"=gcp[,1],
"y"=gcp[,2],
"Recover frequency"=colSums(counts))
row.names(hit_per_gcp) <- paste("GCP", seq(1, nrow(hit_per_gcp)))
hit_per_image <- rowSums(counts)
hist_ <- hist(hit_per_image, breaks=seq(from=0, to=max(hit_per_image), by = 1), plot=FALSE)
return(list(gcp, hit_per_gcp, hit_per_image, hist_))
}
opt_gcp_recover_frequency <- function(gcp_n_in_x, gcp_rec_number, gcp_rec_frequency, mapping_area_x, mapping_area_y, gcp_arrangement_pattern, field_of_view_x, field_of_view_y, photo_positions) {
diff_to_int <- gcp_n_in_x %% 1
gcp_n_in_x <- round(gcp_n_in_x)
gcp_n_in_y <- ceiling(mapping_area_y / (mapping_area_x / gcp_n_in_x))
# Caluclate difference to requested frequency
freq <- calc_gcp_recover_frequency(mapping_area_x, mapping_area_y, gcp_n_in_x, gcp_n_in_y, gcp_arrangement_pattern, field_of_view_x, field_of_view_y, photo_positions)
sum_dens <- do.call(sum, tail(as.list(freq[[4]]$density), -gcp_rec_number))
diff_freq <- sum_dens - (gcp_rec_frequency/100)
# Penaltize for odd numbers and too high bins
score <- abs(diff_freq) + diff_to_int/100
return(score)
}
long2UTM <- function(long) {
(floor((long + 180)/6) %% 60) + 1
......
......@@ -104,6 +104,10 @@ validate_inputs <- function(input) {
need(input$gcp_n_in_x<=20, "Number of GCP along mapping area width > 20"),
need(input$gcp_n_in_y>0, "Number of GCP along mapping area depth missing or invalid"),
need(input$gcp_n_in_y<=20, "Number of GCP along mapping area depth > 20"),
need(input$gcp_rec_frequency>0, "Frequency of images (%) not > 0"),
need(input$gcp_rec_frequency<=100, "Frequency of images (%) > 100"),
need(input$gcp_rec_number>=1, "Frequency of images with number of visible GCPs not >= 1"),
need(input$gcp_rec_number<=3, "Frequency of images with number of visible GCPs not <= 3"),
need(input$position_edge1_lat, "Field edge (Latitude) missing"),
need(input$position_edge1_long, "Field edge (Longitude) missing"),
need(input$position_edge2_lat, "Flight direction (Latitude) missing"),
......@@ -568,53 +572,66 @@ server_ <- function(input, output, session) {
coord_fixed() +
ggtheme_blank
})
# Toggle whether GCP arrangement or recover frequency
# Default:
shinyjs::toggleState("gcp_rec_number")
shinyjs::toggleState("gcp_rec_frequency")
# Dynamic:
observeEvent({input$edit_gcp_arrangement}, {
shinyjs::toggleState("gcp_rec_number")
shinyjs::toggleState("gcp_rec_frequency")
shinyjs::toggleState("gcp_n_in_x")
shinyjs::toggleState("gcp_n_in_y")
})
# GCP location
# Calculate GCP locations and recovery frequency
observe({
req(input$gcp_n_in_x>0, input$gcp_n_in_y>0, input$mapping_area_x>0, input$mapping_area_y>0, input$gcp_arrangement_pattern)
# Positions along axis
gcp_x <- rep(seq(from=0, to=input$mapping_area_x, length.out= input$gcp_n_in_x), each=input$gcp_n_in_y)
gcp_y <- rep(seq(from=0, to=input$mapping_area_y, length.out= input$gcp_n_in_y), input$gcp_n_in_x)
# Crossproduct according to pattern
if(input$gcp_arrangement_pattern == "quad") {
derived_values$gcp <- data.frame(x=gcp_x, y=gcp_y)
req(input$mapping_area_x>0, input$mapping_area_y>0, input$gcp_arrangement_pattern)
req(derived_values$field_of_view_x>0, derived_values$field_of_view_y>0, derived_values$photo_positions)
req(input$gcp_rec_frequency, input$gcp_rec_number)
if(input$edit_gcp_arrangement) {
values <- calc_gcp_recover_frequency(input$mapping_area_x, input$mapping_area_y, input$gcp_n_in_x, input$gcp_n_in_y,
input$gcp_arrangement_pattern, derived_values$field_of_view_x, derived_values$field_of_view_y, derived_values$photo_positions)
derived_values$gcp <- values[[1]]
derived_values$hit_per_gcp <- values[[2]]
derived_values$hit_per_image <- values[[3]]
} else {
gcp_y <- ifelse(rep(seq(from=1, to=input$gcp_n_in_x, by = 1), each=input$gcp_n_in_y) %% 2 != 0, gcp_y, gcp_y + (input$mapping_area_y/((input$gcp_n_in_y-1)*2)))
gcp <- data.frame(x=gcp_x, y=gcp_y)
derived_values$gcp <- gcp %>% filter(y <= input$mapping_area_y)
sr <- optimize(opt_gcp_recover_frequency, 1:20,
input$gcp_rec_number, input$gcp_rec_frequency, input$mapping_area_x, input$mapping_area_y, input$gcp_arrangement_pattern,
derived_values$field_of_view_x, derived_values$field_of_view_y, derived_values$photo_positions, tol=0.05)
gcp_n_in_x <- round(sr$minimum)
gcp_n_in_y <- ceiling(input$mapping_area_y / (input$mapping_area_x / gcp_n_in_x))
values <- calc_gcp_recover_frequency(input$mapping_area_x, input$mapping_area_y, gcp_n_in_x, gcp_n_in_y,
input$gcp_arrangement_pattern, derived_values$field_of_view_x, derived_values$field_of_view_y, derived_values$photo_positions)
if(max(values[[3]] ) > 12) {
gcp_n_in_x <- 2
gcp_n_in_y <- 2
values <- calc_gcp_recover_frequency(input$mapping_area_x, input$mapping_area_y, gcp_n_in_x, gcp_n_in_y,
input$gcp_arrangement_pattern, derived_values$field_of_view_x, derived_values$field_of_view_y, derived_values$photo_positions)
showNotification(type = "error", id="warning_gcp", "Automatic GCP placement not possible",duration = NULL)
} else {
removeNotification(id="warning_gcp")
}
derived_values$gcp <- values[[1]]
derived_values$hit_per_gcp <- values[[2]]
derived_values$hit_per_image <- values[[3]]
updateNumericInput(session, "gcp_n_in_x", value= gcp_n_in_x)
updateNumericInput(session, "gcp_n_in_y", value= gcp_n_in_y)
}
# Also save number of GCP
derived_values$number_of_gcp <- nrow(derived_values$gcp)
})
# Calculate GCP recovery frequency
observe({
req(derived_values$field_of_view_x>0, derived_values$field_of_view_y>0, derived_values$gcp, derived_values$photo_positions)
# Set ground field of view for all GCP
gcp_x_min <- derived_values$gcp[,1] - derived_values$field_of_view_x/2
gcp_x_max <- derived_values$gcp[,1] + derived_values$field_of_view_x/2
gcp_y_min <- derived_values$gcp[,2] - derived_values$field_of_view_y/2
gcp_y_max <- derived_values$gcp[,2] + derived_values$field_of_view_y/2
# Get exposure stations
photos_x <- derived_values$photo_positions[,1]
photos_y <- derived_values$photo_positions[,2]
# Template for the vector we excpect to get returned
return_template <- rep(TRUE, length(photos_x))
# Test if exposure station is in ground field of view
counts <- vapply(gcp_x_min, FUN = function(x){x < photos_x}, return_template) &
vapply(gcp_x_max, FUN = function(x){x > photos_x}, return_template) &
vapply(gcp_y_min, FUN = function(x){x < photos_y}, return_template) &
vapply(gcp_y_max, FUN = function(x){x > photos_y}, return_template)
# Count hits
hit_per_gcp <- data.frame("x"=derived_values$gcp[,1],
"y"=derived_values$gcp[,2],
"Recover frequency"=colSums(counts))
row.names(hit_per_gcp) <- paste("GCP", seq(1, nrow(hit_per_gcp)))
derived_values$hit_per_gcp <- hit_per_gcp
derived_values$hit_per_image <- rowSums(counts)
})
......
......@@ -49,6 +49,8 @@ default_ISO_max<- NULL
default_freq_max<- NULL
default_flight_max<- NULL
default_max_number_of_wp <- 99
default_gcp_rec_frequency <- 50
default_gcp_rec_number <- 1
default_position_edge1_lat <- 47.450812627526901
default_position_edge1_long <- 8.682496912397921
......@@ -179,7 +181,13 @@ ui_ <- fluidPage(
),
tabPanel("GCPs",
h4("Recover frequency"),
fluidRow(
column(6, numericInput("gcp_rec_frequency", "Frequency of images (%)", width = "100%", value = default_gcp_rec_frequency, step=1)),
column(6, sliderInput("gcp_rec_number", "... with min. number of visible GCPs", width = "100%", value = default_gcp_rec_number, step=1, min=1, max=3))
),
h4("Arrangement"),
checkboxInput("edit_gcp_arrangement", "Set manually"),
fluidRow(
column(6, numericInput("gcp_n_in_x", "Number of GCPs, width", width = "100%", value = default_no_gcp_x, step=1)),
column(6, numericInput("gcp_n_in_y", "Number of GCPs, depth", width = "100%", value = default_no_gcp_y, step=1))
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment