Skip to content

Latest commit

 

History

History
272 lines (216 loc) · 6.91 KB

FlightGeometry.md

File metadata and controls

272 lines (216 loc) · 6.91 KB

BioSCape Flight Geometry Calculation

library(sf)
library(units)
library(tidyverse)
library(cowplot)
library(gridExtra)
library(ggforce)
library(ggmap)
library(data.table)
library(piggyback)
library(geosphere)
library(smoothr)
source("R/make_flight_line.R")

Define working projection

wproj="+proj=utm +zone=34 +south +ellps=WGS84 +datum=WGS84 +units=m +no_defs" # working projection

Aircraft Details

g3_altitude=set_units(c(28000,41000,45000),feet)
g5_altitude=set_units(c(28000,41000,51000),feet)
h_feet=g3_altitude

h=set_units(h_feet, km)  #altitude in km

Instrument Details

AVIRIS-NG

Details from https://avirisng.jpl.nasa.gov/specifications.html

AVIRISNG=data.frame(
  instrument="AVIRIS-NG",
  row=1, #just used for plotting
  tfov=set_units(36,degrees),
  ifov=set_units(1,mrad),
  spec_n=480,
  spec_res=set_units(5,nm),
  spec_min=set_units(380,nm),
  spec_max=set_units(2510,nm)
)

HyTES

Details from https://directory.eoportal.org/web/eoportal/airborne-sensors/hytes

HyTES=data.frame(
  instrument="HyTES", 
  row=2, #just used for plotting
  tfov=set_units(50,degrees),
  ifov=set_units(1.44,mrad),
  spec_n=256,
  spec_res=set_units(17.6,nm),
  spec_min=set_units(set_units(7.5,um),nm),
  spec_max=set_units(set_units(12,um),nm)
)

PRISM

Details from https://prism.jpl.nasa.gov/spectrometer_char.html

PRISM=data.frame(
  instrument="PRISM",
  row=3, #just used for plotting
  tfov=set_units(30.7,degrees),
  ifov=set_units(0.97,mrad),
  spec_n=608,
  spec_res=set_units(2.83,nm),
  spec_min=set_units(349.9,nm),
  spec_max=set_units(1053.5,nm)
)

LVIS

LVIS=data.frame(
  instrument="LVIS",
  row=4, #just used for plotting
  tfov=set_units(200,mrad) %>% set_units(degrees), #https://lvis.gsfc.nasa.gov/Data/Maps/ABoVE2019Map.html
  ifov=set_units(0.75,mrad),
  spec_n=1,
  spec_res=set_units(1,nm),
  spec_min=set_units(1064,nm),
  spec_max=set_units(1064,nm)
)

Define angular functions

get_ground_fov<- function(h,fov){
  2*set_units(h,km)*tan(set_units(fov/2,radians)) 
}

Instrument Summary

sensors=bind_rows(AVIRISNG,HyTES,PRISM,LVIS) %>% 
  rowwise() %>% 
  mutate(
    # get swath widths
    swath_width_low=get_ground_fov(min(h),tfov) %>% set_units(km),
    swath_width_mid=get_ground_fov(median(h),tfov) %>% set_units(km),  
    swath_width_high=get_ground_fov(max(h), tfov) %>% set_units(km),  
    # calculate GIFOV - Ground-projected instantaneous field of view
    gifov_low=get_ground_fov(min(h), ifov)  %>% set_units(m),  
    gifov_mid=get_ground_fov(median(h),ifov) %>% set_units(m),
    gifov_high=get_ground_fov(max(h), ifov)%>% set_units(m)
    ) %>% 
  mutate(
    swath_width_km=paste0("[",round(swath_width_low,2),",",
                          round(swath_width_high,2),"]"),
    gifov_m=paste0("[",round(gifov_low,2),",",
                                round(gifov_high),"]"))

sensors %>% 
  select(instrument,
         spectral_bands=spec_n,spectral_resolution=spec_res,
         swath_width_km,spatial_resolution_m=gifov_m) %>% 
DT::datatable()

Compare with aviris-ng data

Read data from here

avirisng=googlesheets4::read_sheet("https://docs.google.com/spreadsheets/d/1UNgRz5EZNLcKtp9cykbi-0lG_D0LWq2lRffwFX3Lpkg/edit#gid=1639195767") %>% 
  filter(grepl("Alt",.$Comments)) # keep only rows that report the altitude

avirisng$alt=set_units(as.numeric(gsub("kft.*$","",gsub("Alt = ","",avirisng$Comments)))*1000,"ft")
# subtract mean scene elevation
avirisng$alt_over_ground=avirisng$alt-set_units(avirisng$`Mean Scene Elevation`,"m")
avirisng$calculated_gifov=get_ground_fov(avirisng$alt_over_ground,set_units(1,mrad)) %>% set_units(m)

avirisng$calculated_gifov_dif=avirisng$calculated_gifov-set_units(avirisng$`Pixel Size`,m)

summary(avirisng$calculated_gifov_dif)
##      Min.   1st Qu.    Median      Mean   3rd Qu.      Max.      NA's 
##   -4.0140   -0.0034    0.0910    0.6369    0.1717 2892.7071       131

avirisng %>% 
  filter(alt<set_units(70000,ft)) %>% 
  ggplot(aes(x=`Pixel Size`,y=calculated_gifov))+
  geom_hex(bins=100)+
  scale_fill_viridis_c()+
  geom_abline(slope=1,intercept=0,col="red")+
  coord_equal()+
  xlab("Pixel size (m) reported by AVIRIS team")

AVIRIS-ng flight altitudes

avirisng %>% 
  filter(alt<set_units(60000,ft)) %>% 
  ggplot(aes(x=alt))+
  geom_histogram()
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

Typical flights

G3 https://jsc-aircraft-ops.jsc.nasa.gov/gulfstream-giii.html 6 hours / 2,600 nautical miles = 4815.2 km

d=4815*(5/6) # 6 hour flight coverage in km (minus one hour for takeoff/landing)
set_units(d,"km")
## 4012.5 [km]

Example Path

path=st_linestring(x=as.matrix(rbind(
  c(18.46929,-34.38451),
  c(18.37620,-34.08127),
  c(18.46486,-33.99350),
  c(18.77962,-33.94030),
  c(18.93478,-33.84011))),dim="XY") %>% 
  st_sfc() %>% st_sf(path=1,geom=.,crs=4326) %>% 
  smooth(method = "ksmooth")

Calculate instrument polygons

flights <- sensors %>% 
  rowwise() %>% 
  mutate(geometry=make_flight_line(path=path,
                                   swath_width = swath_width_mid)) %>% 
  st_as_sf() %>%
  mutate(instrument=factor(instrument,ordered=T,levels = c("HyTES","AVIRIS-NG","PRISM","LVIS")))


bbox=flights %>% #st_transform(flights,4326) %>% 
  st_buffer(dist = 0.2) %>% 
  st_bbox();names(bbox)=c("left","bottom","right","top")

hdf <- get_map(location=bbox,
              maptype="toner", zoom=12,
              source="stamen")
#hdf <- get_map(location=bbox,  #different background
#               maptype="satellite", zoom=12)

Visualize flight path

ggmap(hdf, extent = "normal")+
  geom_sf(data=st_transform(flights,4326),mapping=aes(fill=instrument),
          alpha=.4,color=NA,inherit.aes = F)+
  geom_sf(data=st_transform(path,4326),inherit.aes = F)+
  ylab("Latitude")+
  xlab("Longitude")
## Coordinate system already present. Adding new coordinate system, which will replace the existing one.

Write KML for Google Earth Visualization

flights %>% 
  st_segmentize(50) %>% 
#  mutate(OGR_STYLE = "BRUSH(fc:#0000FF80); PEN(c:#FF0000,w:1px)") %>%  #fiddle with kml aesthetics
  st_write("data/transect.kml",append=F,altitudeMode="clampToGround")
## Writing layer `transect' to data source `data/transect.kml' using driver `KML'
## Writing 4 features with 16 fields and geometry type Polygon.

if(F){  # if desired, push file to github release
  tag="v0.0.1"
  repo="BioSCape-io/BioSCape_flights"
  pb_new_release(repo, tag)
  pb_upload("data/transect.kml", 
          repo = repo, 
          tag = tag)
}