Exploring the NASA-GIBS Service with Leaflet and Shiny

NASA publishes some excellent imagery and gridded data in near real time through the Global Imagery Browse Service (GIBS) service.  I want to explore this service, and in the process develop my proficiency with Leaflet and Shiny.

RStudio has released a package to bind Leaflet to Shiny page, built on their htmlwidgets .  Using this package, it is trivial to write a shiny app that contains a leaflet map.  To begin, we need to install the leaflet package from github:

	devtools::install_github('rstudio/leaflet')
	# example shiny ui.R with leaflet.js
	#
	library(shiny)
	library(leaflet)
	shinyUI(fluidPage(leafletOutput('myMap')))
	# Basic shiny server.R function to show a leaflet control using OSM tiles
	#
	library(shiny)
	shinyServer(function(input, output) {
	output$myMap = renderLeaflet(leaflet() %>% addTiles())
	})

view raw install_and _test_leaflet.R hosted with ❤ by GitHub

Shown above is a basic shiny program to display a slippy map using open street map tiles, the live example can be explored here.  So how hard is it to get those GIBS layers in there?  Looking at the NASA Earthdata map viewer we can note some issues:

• most layers have a date component to the URI
• not all layers have the same available date range
• not all layers are available in all projections
• some layers are intended to be overlays

Sticking with web mercator (EPSG:3857), the first task is to sort out all the available layers. Here is some code to read the XML capabilities file and parse out the useful bits:

	# definitions for GIBS WMTS tile layers
	library(XML)
	library(data.table)
	# read in the list of available layers for the GIBS resource
	# returns a data.table with useful parameters
	# hardwire the location of the capabilities file
	gibsRes<-"http://map1.vis.earthdata.nasa.gov/wmts-webmerc/1.0.0/WMTSCapabilities.xml"
	# we are only interested in the <Layer> sections
	doc<-xmlParse(gibsRes)
	nodes<-getNodeSet(doc,"/r:Capabilities/r:Contents/r:Layer",c(r="http://www.opengis.net/wmts/1.0"))
	gibsDF<-xmlToDataFrame(doc,nodes=nodes)
	# the resulting data.frame needs some cleanup <- better XFile code above would be cool
	# in the meantime just parse out the strings
	gibsDF$Format<-gsub("image/","",gibsDF$Format)
	# switch to data.table for better syntax
	gibsDT<-as.data.table(gibsDF)
	gibsDT[,zoomLevel:=as.numeric(substring(TileMatrixSetLink,nchar(TileMatrixSetLink)))]
	gibsDT[,startDate:=as.POSIXct(strptime(substring(Dimension,27,36),"%Y-%m-%d"))]
	gibsDT[,endDate:=as.POSIXct(strptime(substring(Dimension,38,47),"%Y-%m-%d"))]
	gibsDT[,dateFreq:=substring(Dimension,49)]
	gibsDT[,dateRange:=gsub("/"," ",substring(Dimension,27))]
	# there are a few other fields that might need to be decoded, and the <Metadata> and <ResourceURL> tags could be picked up
	# but the munging done so far is enough to access the data with leaflet
	gibsLayers <- gibsDT[Format=='jpeg',Identifier]
	glZoomlevel <- gibsDT[Format=='jpeg',zoomLevel]
	glTileset <- gibsDT[Format=='jpeg',TileMatrixSetLink]
	glImgtype <- gibsDT[Format=='jpeg',Format]
	# I dont see any explicit metadata for if a layer is intended to be an overlay
	# other than the file type, so we will use all png files as overlays
	gibsOverlays <- gibsDT[Format=='png',Identifier]
	oglZoomlevel<- gibsDT[Format=='png',zoomLevel]
	oglTileset <- gibsDT[Format=='png',TileMatrixSetLink]
	oglImgtype <- gibsDT[Format=='png',Format]

view raw defs.R hosted with ❤ by GitHub

It appears that the file type (png or jpeg) is the key to determining if the layer is a base layer or an overlay. The other bits that need to be pulled out are the maximum zoom factor and available date ranges. Using these fields, we can specify the layer(s) we want to show with the following shiny code segments:

	# Test rig for leaflet map using web mercator GIBS layers
	#
	library(shiny)
	library(leaflet)
	source('./defs.R')
	shinyUI(fluidPage(
	# Application title
	titlePanel("GIBS Leaflet Test Page"),
	# Sidebar with a slider input for number of bins
	sidebarLayout(
	sidebarPanel(
	selectInput("selLayer","Available Image Layers",gibsLayers,selected=gibsLayers[1],width="100%"),
	selectInput("selOverlay","Available Overlays",gibsOverlays,selected=gibsOverlays[17],width="100%",multiple=TRUE,selectize=FALSE),
	dateInput("selDate",
	"Viewing Date",
	min = "1979-01-01",
	max = format(Sys.Date(), "%Y-%m-%d"),
	value = format(Sys.Date(), "%Y-%m-%d")),
	h4("Base Layer Details"),
	textOutput("baseTxt")
	),
	# Show leaflet map with a text div reporting the selected date and extents
	mainPanel(
	h4(verbatimTextOutput("mapTxt")),
	leafletOutput("mapPlot",height=600),
	h4("Overlay Details"),
	tableOutput("olTxt")
	)
	)
	))

view raw ui.R hosted with ❤ by GitHub

Some hard-coding done here to select the coastlines overlay, and the calendar range was set by inspection of the current dataset.

	# test rig for leaflet app
	# install it using:
	# devtools::install_github("rstudio/leaflet")
	#
	library(shiny)
	library(leaflet)
	source('./defs.R')
	shinyServer(function(input, output, session) {
	map = leaflet()
	output$baseTxt <-renderText({
	sel<-gibsDT[Identifier==input$selLayer]
	return(paste(sel$Title,"Date Range: ",sel$dateRange))
	})
	output$olTxt <-renderTable({
	return(gibsDT[Identifier %in% input$selOverlay,list(Title,WGS84BoundingBox,zoomLevel,dateRange)])
	})
	output$mapPlot <- renderLeaflet({
	lastBounds <- isolate( input$mapPlot_bounds )
	zoomLvl <- isolate(glZoomlevel[which(gibsLayers == input$selLayer)])
	tileSet <- isolate(glTileset[which(gibsLayers == input$selLayer)])
	iType <- isolate(glImgtype[which(gibsLayers == input$selLayer)])
	# overlay selection data - needs to allow multiselect
	olzoomLvl <- isolate(oglZoomlevel[which(gibsOverlays %in% input$selOverlay)])
	oltileSet <- isolate(oglTileset[which(gibsOverlays %in% input$selOverlay)])
	oliType <- isolate(oglImgtype[which(gibsOverlays %in% input$selOverlay)])
	# add the base layer
	m <- map %>% addTiles(paste0('http://map1{s}.vis.earthdata.nasa.gov/wmts-webmerc/',input$selLayer,'/default/',input$selDate,'/',tileSet,'/{z}/{y}/{x}.',iType),
	attribution = paste(
	'<a href="https://earthdata.nasa.gov/gibs">NASA EOSDIS GIBS</a>'
	),
	options = list(
	maxZoom = zoomLvl,
	minZoom = 1,
	tileSize = 256,
	subdomains = "abc",
	noWrap = "true",
	crs = "L.CRS.EPSG3857",
	continuousWorld = "true",
	# Prevent Leaflet from retrieving non-existent tiles on the borders.
	bounds = list(list(-85.0511287776, -179.999999975),list(85.0511287776, 179.999999975))
	))
	# add all the selected overlay layers
	for(i in 1:length(input$selOverlay) ){
	m <- m %>% addTiles(paste0('http://map1.vis.earthdata.nasa.gov/wmts-webmerc/',input$selOverlay[i],'/default/',input$selDate,'/',oltileSet[i],'/{z}/{y}/{x}.',oliType[i]),
	attribution = paste('<a href="https://earthdata.nasa.gov/gibs">NASA EOSDIS GIBS</a>'),
	options = list(
	maxZoom = 9,
	maxNativeZoom = olzoomLvl[i],
	minZoom = 1,
	tileSize = 256,
	subdomains = "abc",
	noWrap = "true",
	crs = "L.CRS.EPSG3857",
	continuousWorld = "true",
	opacity = 0.3,
	# Prevent Leaflet from retrieving non-existent tiles on the borders.
	bounds = list(list(-85.0511287776, -179.999999975),list(85.0511287776, 179.999999975))
	))
	}
	# set the extents - this could use some work
	if(!is.null(lastBounds$north)) m <- m %>% fitBounds(lastBounds$east,lastBounds$south,lastBounds$west,lastBounds$north)
	m
	})
	output$mapTxt <- renderText({return(paste(input$selDate,input$mapPlot_zoom,input$mapPlot_center,input$mapPlot_bounds$north,input$mapPlot_bounds$east,input$mapPlot_bounds$south,input$mapPlot_bounds$west))})
	})

view raw server.R hosted with ❤ by GitHub

Here is a live version of the application.

The source code for this application is available at gihub:

https://github.com/bungalow-bill/GIBS-leaflet.git

Some successful flying

This weekend we got some encouraging flight results with the hexacopter.  After making some modifications to the camera mount and spending some time building confidence with the flight simulator, we were hoping to get some usable photos and maybe even work up the courage to try a flight with the autopilot in control.

Site:  Carstairs farm
Conditions:  Sunny, light breeze from the west.
Camera:  Canon Powershot A2300, shutter speed set to 1/2000, ISO 100
GPS:  uBlox LEA-6H

The first flight went well, with the exception of the landing, which was done with only five propellers.  Surprisingly, there was little damage to the copter (one arm needed to be replaced, no props broke).

The photos are not terribly interesting, unless you like looking at targets on the lawn.  Here is a close up of the resolution target from an altitude of roughly 15m.:

On the second flight we tried the autopilot, with a modest program of four waypoints.  The results were not that good, likely because I entered a radius of 2m for each waypoint which is too tight for the system to lock on to.

The photos have some more interesting subjects, but we had to curtail the flight since the autopilot was not cooperating.  Looking at the flight logs, it's clear that the problem is that the waypoint radius is too small.

The third flight didn't really have any goals, other than to do some more flying.  I suppose if we had been more ambitious we could have reprogrammed the waypoints and tried the auto mode again.   But instead we just let the autopilot fly in loiter mode, into the trees.

oops

The photos are more of the same.

Some concluding thoughts:

• I think the camera mount is providing reasonable results given the quality of the camera.  It would be nice to improve the CHDK script so that the shutter speed selection would be handled automatically
• The next flight with the autopilot needs to have a better program, we need to spend some time understanding the options and their impact on the system, specifically the waypoint radius and delay values.
• The flight simulator is a great tool for improving pilot confidence, but more practice is definitely required, particularly for landings.

FTDI Board? We don't need no stinking FTDI board!

Connecting to GPS receivers to test, configure and log data seems to be a daily thing this month.  The way we have been doing it is to connect from the computer USB to the UART on the GPS using the handy Sparkfun FTDI breakout board.  

For example, we can connect to the uBlox LEA-6H using the uCenter software to set the configuration, following these instructions.

But the goal is to get the GPS working with the APM, which has plenty of UARTs.  So, how about we just talk directly to the GPS using the APM?  Easy:
1)  Install the patched version of the Arduino development environment, following the instructions here.
2)  Connect to the APM via USB.
3)  Open the Arduino IDE and enter the following code:

/*

  Mega multple serial test

  Receives from the main serial port, sends to the others. 

 Receives from serial port 1, sends to the main serial (Serial 0).

  This example works only on the Arduino Mega

  The circuit: 

 * Any serial device attached to Serial port 1

 * Serial monitor open on Serial port 0:

  created 30 Dec. 2008

 modified 20 May 2012

 by Tom Igoe & Jed Roach

  This example code is in the public domain.

  */

void setup() {

  // initialize both serial ports:

  Serial.begin(38400);

  Serial1.begin(38400);

}

void loop() {

  // read from port 1, send to port 0:

  if (Serial1.available()) {

    int inByte = Serial1.read();

    Serial.write(inByte); 

  }

  // read from port 0, send to port 1:

  if (Serial.available()) {

    int inByte = Serial.read();

    Serial1.write(inByte); 

  }

}

Note the baud rates are set to 38.4K, which is the standard for the uBlox connection to the APM.
4)  Click the 'Upload' button to compile the code and run it on the APM
5)  Now open uCenter and connect to the appropriate COM port at 38.4K baud.
6)  Do your GPS configurating.

Note that you will need to reinstall the firmware on the APM using Mission Planner when you are done.

MODIS Image Processing Stream

A MODIS Image Processing Stream

Inputs:

• MODIS mosaics (eg Arctic, Antarctic)



Data Flow?

• Image granules from Rapidfire

Products:

• Time lapse movies
• Custom images/mosaics
• Cloud cover masks
• Motion vectors

Components

• Data retriever, per collection (may require metadata collector)
• MetaDatabase & Metadata producer
• Cloud mask tool
• Movie composer
• Optical flow detector/filter

Currently:

• Data retriever is a batch file/php for the arctic mosaic, and Watchdog/WatchFTP for granules. Good enough for now, but we need to have a tool to decide if new data should be downloaded and how.
• Metadata producer/Database - nothing
• Cloudmask tool - I have an R script that will process arctic mosaic tiles and create a fair result.  For granules, the MOD35-L2 product is available.
• Movie composer, can be done with batch scripts with varying levels of sophistication using Imagemagick, Mapserver, gdal, and ffmpeg (also sed).  It would be useful to have a web front end for specifying and managing requests and products.
• Optical flow detector/filter - combination of OpenCV and R at this point, could be batched but needs work.

Canon A2300 Mount

Here is a simple mount for the Canon A2300 camera, constructed from an aluminum bracket salvaged from an old scanner plus a piece of angle bracket.  Weight is 186g, including the bracket, fixing bolt, camera and release cable.

I used some foam weatherstrip to pad the part that touches the camera, but I doubt it will reduce any vibrations between the bracket and the camera. 

Still to be addressed:  the connection to the airframe.  It would be interesting to use some sort of shock absorbing link there, not sure what would work...

Canon Camera Shutter Release Cable

Using a Canon handheld camera with the APM 2.5 looks to be an easy low cost method of collecting some aerial imagery.  We have a couple of suitable cameras to try out:

Canon Powershot A2300:  16MP, $100 at Walmart

Canon Powershot SD800IS:  7.1MB, $25 off Craigslist

Both cameras are supported by the Canon Hack Development Kit (CHDK) which adds two valuable features to these low cost devices:  scripting and a remote shutter release via the USB port.

The remote shutter release simply monitors the presence of a voltage on the USB power pin.  Following the instructions here, I made up a test cable using a battery and a momentary contact switch.  Testing with various voltages, it seems to need at least 4.1V to reliably trigger the shutter.

CHDK Remote Cable Test Rig W/4.8V Battery

To connect to the APM, we'll need to use one of the digital I/O pins (3.3V) as a trigger, with voltage supplied by the 5V supply pin.  Testing with the battery indicates that the current drain is negligible, so here is the circuit we came up with, based on a 4049 hex inverting buffer:

CHDK Camera Circuit

Here's the completed circuit, needs some strain reliefs, heatshrink, 

A small ceramic capacitor is bridged across 5V-GND as a buffer.  On initial testing, it works!  Next steps are to test with the APM and fine tune the scripting on the camera and autopilot.

Many thanks to my friend paulb for the expert advice, parts, soldering, and correction of errors.  I did, however, supply the beer.

Here is a discussion and some code for scripting the camera.