Kind of 3D with D3

Fig 1 – NASA Near Earth Observation - Cloud Optical Thickness using D3js.org Orthographic projection

Charts and Graphs

Growing up in the days before graphing calculators (actually before calculators at all), I’ve had a lingering antipathy toward graph and chart art. Painful evenings plotting interminable polynomials one dot at a time across blue lined graph paper left its inevitable scars. In my experience, maps are a quite different thing altogether. Growing up in Ohio at the very edge of Appalachia, I have quite pleasant memories of perusing road atlases, planning escapes to a dimly grasped world somewhere/anywhere else.

The “I’m outa here” syndrome of a mid-century youth, landed me in Colorado where I quickly found the compelling connection between USGS Topo maps and pleasant weekends climbing in Colorado’s extensive National Forests. I still fondly recall the USGS map desk at Denver’s Lakewood Federal Center where rows of flat metal cabinet drawers slid out to expose deep piles of 1:24000 scale topography maps with their rich linear contour features.

I bore you with this personal recollection to set the stage for my amazement at discovering d3.js. My enduring prejudice of charts and graphs had already begun to crack a bit after seeing Hans Rosling’s entertaining Ted lectures from 2006. (D3 version of Wealth and Health of Nations). Perhaps my pentient for the concrete posed a barrier to more abstract demands of statistical modeling. At any rate Hans Rosling’s engaging enthusiasm was able to make a dent in my prejudice, and I found myself wondering how his graphs were coded. Apart from the explicit optimism conveyed by Rosling’s energy, the visual effect of primary colored balloons rising in celebratory fashion is quite hopefully contradictory of the harsh Hobbesian dictum.

“the life of man, solitary, poor, nasty, brutish, and short.”

The revelation that my childhood experience of the 50s is recapitulated by the modern economic experience of children in Trinidad and Tobago was something of an eye opener. Although no maps are in sight, the dynamic visualization of normally depressing socio-economic statistics lent visual energy to a changing landscape of global poverty.

Fig 2 – Hans Rosling’s Wealth & Health of Nations implemented in D3

The application of technology to the human condition seems so inspiring that it’s worth the reminder that technology has its flaws. Can it really be that the computational power in hand on the streets of Lagos Nigeria now exceeds that of Redmond or Palo Alto? Will ready access to MIT courseware on Probability and Statistics actually mitigate the desperation of Nairobi’s slums? A dash of Realpolitik is in order, but the graphics are none the less beautiful and inspiring.

Fig 3 – Mobile phones transform Lagos

D3.Geo

d3js.org is an open source javascript library for data visualization. In d3, data enters, performs, and exits with a novel select pattern geared to dynamic rendering. The houselights dim, the stage is lit, data enters stage left, the data performance is exquisite, data exits stage right. This is a time oriented framework with data on the move. The dynamism of data is the key to compelling statistics and d3 examples illustrate this time after time.

With d3.js Economics may remain the Dismal Science, but its charts and graphs can now be a thing of beauty. D3 charts are visually interesting, but more than that, chart smarts are leaking into the spatial geography domain. Mike Bostock of NYT fame and Jason Davies are prolific contributors to d3js.org.

D3 version 3.0, released just a couple of weeks ago, added a wide range of projective geometry to d3.geo. Paths, Projections, and Streams combine to realize a rich html5 rendering library capable of animation effects rarely seen in GIS or on the web.

Visually these are not your grandfathers’ charts and maps.

TopoJSON

There’s more. Sean Gillies recently remarked on the advent of TopoJSON, apparently a minor side project of mbostock, fully supported by D3.Geo.Paths. In the GIS world ESRI has long held the high ground on topological data (the Arc of GIS), and now some blokes from NYT charts and graphs have squashed it into an easy to use javascript library. Wow!

TopoJSON is still more or less a compressed transport mechanism between OGC Simple Features in a server DB and client side SVG, but I imagine there will shortly be server side conversions for PostGIS and SQL Server to take advantage of major low latency possibilities. Client side simplification using the Visvalingam’s algorithm are virtually instantaneous, so zoom reduction can occur client side quite nicely.

I think you get the idea. D3.js is powerful stuff and lots of fun.

Some Experiments

It’s the New Year’s holiday with spare time to fire up a favorite data source, NASA Neo, and try out a few things. The orthographic projection is an azimuthal projection with the advantage of reproducing the visual perspective of earth from a distant vantage point.

Simple as:

var projection = d3.geo.orthographic()
    .scale(245)
    .clipAngle(90);

Fig 4 – Natural earth world data in d3js.org orthographic projection

In addition to a graticule and circle edge paths, this takes advantage of TopoJSON formatted natural earth world data published on GitHUB by, you guessed it, mbostock.

        d3.json("topojson/world-110m.json", function (error, world) {
            nasa.svg.append("path")
                .datum(topojson.object(world, world.objects.land))
                .attr("class", "boundary")
                .attr("d", nasa.path);
        });

Add some mouse event handlers to whirl the globe:

mousedown: function () {
        nasa.p0 = [d3.event.pageX, d3.event.pageY];
        nasa.context.clearRect(0, 0, nasa.width, nasa.height);
        d3.event.preventDefault();
    },
    mousemove: function () {
        if (nasa.p0) {
            nasa.p0 = [d3.event.pageX, d3.event.pageY];
            nasa.projection.rotate([nasa.λ(nasa.p0[0]), nasa.φ(nasa.p0[1])]);
            nasa.svg.selectAll("path").attr("d", nasa.path);
        }
    },
    mouseup: function () {
        if (nasa.p0) {
            nasa.mousemove();
            nasa.context.clearRect(0, 0, nasa.width, nasa.height);
            nasa.getOverlay(nasa.overlay);
            nasa.p0 = null;
        }
    },

Interactive orthographic from mbostock has a nice set of handlers adapted for this experiment. Note there are still a few quirks regarding transposing mousedown event locations to be worked out in my experimental adaptation. (My holiday free time is running out so some things have to make do.)

With mouse action, d3’s orthographic projection becomes a globe. It responds much more smoothly in Chrome than IE, apparently due to a javascript performance advantage in Chrome.

Fig 5 – javascript performance testing

This ortho globe feels 3D all because D3 affords a fast refresh through a projection on the vector continent and graticule paths.

nasa.svg.selectAll("path").attr("d", nasa.path); 

Vectors are fast, but for deeper information content I turn to NASA’s Near Earth Observation data, available as imagery from this public WMS service. The beauty of this imagery is still compelling after all these years.

Fig 6 – NASA Land night temperature overlay

However, geodetic imagery needs to be transformed to orthographic as well. D3 has all the necessary plumbing. All I added was the NASA WMS proxy with a .net WCF service.

getOverlay: function (overlay) {
        if (overlay != null && overlay.length>0) {

            nasa.Showloading();
            nasa.image = new Image;
            nasa.image.onload = nasa.onload;
            nasa.image.src = Constants.ServiceUrl + "/GetMap?url=http://neowms.sci.gsfc.nasa.gov/wms/wms?Service=WMS%26version=1.1.1%26Request=GetMap%26Layers=" + overlay + "%26SRS=EPSG:4326%26BBOX=-180.0,-90,180,90%26width="+nasa.width+"%26height="+nasa.height+"%26format=image/jpeg%26Exceptions=text/xml";
        }
    },

    onload: function () {

        var dx = nasa.image.width,
            dy = nasa.image.height;

        nasa.context.drawImage(nasa.image, 0, 0, dx, dy);

        var sourceData = nasa.context.getImageData(0, 0, dx, dy).data,
            target = nasa.context.createImageData(nasa.width, nasa.height),
            targetData = target.data;

        for (var y = 0, i = -1; y < nasa.height; ++y) {
            for (var x = 0; x < nasa.width; ++x) {
                var p = nasa.projection.invert([x, y]), λ = p[0], φ = p[1];
                if (λ > 180 || λ < -180 || φ > 90 || φ < -90) { i += 4; continue; }
                var q = ((90 - φ) / 180 * dy | 0) * dx + ((180 + λ) / 360 * dx | 0) << 2;
                targetData[++i] = sourceData[q];
                targetData[++i] = sourceData[++q];
                targetData[++i] = sourceData[++q];
                targetData[++i] = 255;
            }
        }
        nasa.context.clearRect(0, 0, nasa.width, nasa.height);
        nasa.context.putImageData(target, 0, 0);
        nasa.Hideloading();
    },

Looping through 1,182,720 pixels in javascript is not the fastest, but just to be able to do it at all with only a dozen lines of javascript is very satisfying. There may be some server side options to improve performance and PostGIS Raster is worthy of investigation. However, html5 browsers with enhanced GPU access should eventually supply higher performance raster transforms.

Fig 7 – NASA Land night temperatures transformed with D3 Orthographic projection.invert

Selection JsTree

For this experiment I also made use of JsTree for the layer selections out of the NASA WMS GetCapabilities. Layer choices are extensive and even with a tree expander approach the options overflow the available div space of IE’s jQuery UI draggable accordion panel. Scroll bars work poorly with draggable divs. A future enhancement could be manually allocated multiple expanders following NASA’s Ocean, Atmosphere, Energy, Land, and Life categories. Unfortunately this would invalidate the nice GetCapabilities layer extraction feature of the proxy service. NASA’s WMS also provides LegendURL elements for better understanding of the color ranges, which should be added to the selection panel.

MouseWheel
Since d3.geo offers projection.scale(), mousewheel events are a nice to have option that is easily bound to a browser window with jquery-mousewheel plugin.

$(window).mousewheel(function (event, delta, deltaX, deltaY) {
            var s = nasa.projection.scale();
            if (delta > 0) {
                nasa.projection.scale(s * 1.1);
            }
            else {
                nasa.projection.scale(s * 0.9);
            }
            nasa.svg.selectAll("path").attr("d", nasa.path);
            nasa.context.clearRect(0, 0, nasa.width, nasa.height);
            nasa.getOverlay(nasa.overlay);
        });

Tilted Perspective
Even though NEO data resolution doesn’t really warrant it, using tilted perspective d3.geo.satellite projection affords a space station view point. Incorporating a steerable camera is a more complex project than time allows.

var projection = d3.geo.satellite()
    .distance(1.1)
    .scale(5500)
    .rotate([76.00, -34.50, 32.12])
    .center([-2, 5])
    .tilt(25)
    .clipAngle(25);

Fig 8 – Tilted perspective d3.geo.satellite projection of NASA BlueMarbleNG-TB

Possible Extensions:
Time is the next target, since NASA Neo exposes Time extension features.
<Extent default=”2012-12-11″ name=”time”>

Example:
http://neowms.sci.gsfc.nasa.gov/wms/wms?Service=WMS&version=1.1.1&Request=GetMap&Layers=MOD_LSTAD_D&time=2000-03-06/2000-04-26/P1D&SRS=EPSG:4326&BBOX=-180.0,-85,180,85&width=1000&height=500&format=image/jpeg&Exceptions=text/xml

One approach would be a time series download as an array of image objects from the proxy service, with a step through display bound to a slider on the UI. Once downloaded and projected the image stack could be serially displayed as an animation sequence. NASA WMS Time domain is somewhat abbreviated with only about 12-14 steps available. For a longer animation affect some type of caching worker role might store images as a continuous set of PostGIS raster data types with an ongoing monthly update. PostGIS has the additional advantage of pre-projected imagery fetching for more efficient time lapse performance.

Warning:
IIS requires adding a .json application/json mime type to make use of TopoJSON formatted data.

Summary:

Mike Bostock and Jason Davies have done us all a great favor making the power of d3 freely available on GitHUB under BSD License. For sheer prodigality http://bl.ocks.org/mbostock and http://bl.ocks.org/jasondavies examples can hardly be matched. Anyone need a Guyou Hemisphere projection or perhaps a Peirce Quincuncial tessellation? d3.geo has it ready to serve up in web friendly fashion.

James Fee sums up d3 quite succinctly in his recent post “Why D3 Will Change How We Publish Maps”. I agree. We’re all joining the blokes over in charts and graphs.

links:
http://www.jasondavies.com/maps
Interrupted Bonne
World Tour
Rotating EquiRect

Looking at JPEG 2000

ER Map Viewer JPEG 2000
Fig 1 – NAIP TX JP2 in ER Viewer 7.2

I ran into a problem at the end of last week. I was attempting to use a set of NAIP county mosaics for Texas. These are very nice 1m GSD images available from the USDA Data Gateway as 4 band JPEG 2000 image files. The 4th band falls into the infrared region and allows analysts to see more vegetation type details. I don’t need it, but it helps out the Dept of Agriculture who is, after all, funding this very useful public imagery. I just wish they published NAIP in an OGC service, WMS or WCS.

Here is a metadata description. These are useful high resolution images and I was providing them for a cable routing contractor to use in his engineering ‘as built’ documents. I’ve done this in the past without much problem using NAIP images from 2006 – 2007 furnished as MrSID compressed files.

The process is to create a corridor in AutoCAD along the proposed route. Next, draw a series of page rectangles along this route making sure the entire corridor is covered, stepping through the route. This is fairly straight forward AutoCAD drawing. After exporting as DXF, these tiles can then be processed through a small Java program to create a batch file for clipping imagery.

This batch file makes use of gdal_translate to clip out a subset of the source image and create a new GeoTiff result. The resulting batch file is just a text file with a lot of lines like this:
gdal_translate -of GTiff -co “TILED=YES” -srcwin 27537 45649 24951 5046.003420858178 C ortho_1-1_1m_j_tx375_2008_2.jp2 potter.tiff

The image files are large so compression is mandatory. Especially since the images are furnished for ftp download from USDA servers. Wavelet compression is very efficient and has become the compression of choice with very large images like these. JP2 is also interesting for web applications due to its progressive transmission by pixel and resolution accuracy. Now that JPEG 2000 implementation libraries have been out for awhile, more and more jp2 files and tools are available. I didn’t anticipate any problem with the 2008 switch from MrSID to JP2.

However, I was very wrong.

A few hours later I had a list of viewers that did not work:

and a couple that did: (but didn’t export to tiff gtiff)

TatukGIS Viewer JPEG 2000
Fig2 – NAIP TX JP2 in TatukGIS Viewer
OpenEV Viewer JPEG 2000
Fig 3 – NAIP TX JP2 in OpenEV Viewer
Kakadu kdu_show JPEG 2000
Fig 4 – NAIP TX JP2 in Kakadu kdu_show
ArcView JPEG 2000
Fig 5 – NAIP TX JP2 in ArcView
Global Mapper JPEG 2000
Fig 6 – NAIP TX JP2 in Global Mapper

Hmm that is strange. It has a cimicine odor to it, but JPEG2000 has been in the works for years and has a good range of well used tools? NAIP is very popular imagery, however, a quick Google search didn’t turn up any very specific issue. This 11/19/2008 PostMortem ppt lists some JPEG2000 problems but nothing as catastrophic as my problem:

“JPEG2000 Compression format required for 4-band CCM under the NAIP contract. APFO is experiencing the following problems/issues with JPEG 2000 compressions:

  • Blurring images.
  • Difficulty with program settings
  • Limited customer support from open source software libraries
  • Image loss with advancing zoom
  • Rendering issues due to DVD exceeding ArcGIS 9.1 file capacities
  • Rendering issues on other viewers (consumer use) have been resolved.”

This link is an interesting summary from the 11/19/2008 NAIP coordination meeting. Note that the seamless NED resolution has ramifications for processing these higher resolution national imagery programs. Sounds like another reason to move to a national high resolution LiDAR survey of at least the continental US. I suppose LiDAR collection could happen simultaneously with imagery collection and save a bundle of flying costs.

USDA was prompt when I asked about my issue, but since it was viewable in their ArcInfo all was right with the world … and their data. I couldn’t really quibble with that, but I don’t have access to an ArcInfo License so in my world … all is still not right.

Next I took a look at the header records as output from gdalinfo

C:\Program Files\FWTools2.4.2\bin\gdalinfo ortho_1-1_1m_j_tx375_2008_2.jp2
Driver: JP2ECW/ERMapper JPEG2000
Files: ortho_1-1_1m_j_tx375_2008_2.jp2
       ortho_1-1_1m_j_tx375_2008_2.j2w
Size is 59292, 57976
Coordinate System is `'
Origin = (208978.447721050060000,3947556.455350011600000)
Pixel Size = (1.000000000000000,-1.000000000000000)
Corner Coordinates:
Upper Left  (  208978.448, 3947556.455)
Lower Left  (  208978.448, 3889580.455)
Upper Right (  268270.448, 3947556.455)
Lower Right (  268270.448, 3889580.455)
Center      (  238624.448, 3918568.455)
Band 1 Block=59292x1 Type=Byte, ColorInterp=Undefined
  Overviews: arbitrary
Band 2 Block=59292x1 Type=Byte, ColorInterp=Undefined
  Overviews: arbitrary
Band 3 Block=59292x1 Type=Byte, ColorInterp=Undefined
  Overviews: arbitrary
Band 4 Block=59292x1 Type=Byte, ColorInterp=Undefined
  Overviews: arbitrary

There is a small anomaly in the Coordinate System, but all seems fine otherwise. It is interesting to see the Driver listed for gdal is “Driver: JP2ECW/ERMapper JPEG2000″. Erdas ER Mapper makes ER Viewer with the same problem. I don’t know which sdk is licensed by Tatuk and AutoDesk, but possibly they also use the ER Mapper’s JPEG2000 sdk.

Next I fired up Eclipse and took a look at some JAI. First I tried rewriting the jp2 as tif:

package com.mmc.image;

import java.awt.image.BufferedImage;
import java.io.File;
import java.io.IOException;
import javax.imageio.ImageIO;·

public class J2KTest {
public static void main(String[] args) {
try {
      final BufferedImage bi = ImageIO.read(new File("C:/temp/JP2_test.jp2"));
      File dest = new File("C:/temp/test.tif");
      ImageIO.write(bi, "tiff", dest);
      System.out.println("Complete");·
    } catch (IOException e) {
        e.printStackTrace();
    }
  }
}

“no bing!” Just another clue pointing toward the MetaData though:

Exception in thread "main" java.lang.NullPointerException
	at com.sun.media.imageioimpl.plugins.jpeg2000.J2KMetadata.replace(J2KMetadata.java:962)
	at com.sun.media.imageioimpl.plugins.jpeg2000.J2KMetadata.addNode(J2KMetadata.java:631)
	at com.sun.media.imageioimpl.plugins.jpeg2000.J2KRenderedImageCodecLib.readImageMetadata(J2KRenderedImageCodecLib.java:1006)
	at com.sun.media.imageioimpl.plugins.jpeg2000.J2KRenderedImageCodecLib.createOriginalSampleModel(J2KRenderedImageCodecLib.java:673)
	at com.sun.media.imageioimpl.plugins.jpeg2000.J2KRenderedImageCodecLib.(J2KRenderedImageCodecLib.java:261)
	at com.sun.media.imageioimpl.plugins.jpeg2000.J2KImageReaderCodecLib.read(J2KImageReaderCodecLib.java:364)
	at javax.imageio.ImageIO.read(ImageIO.java:1422)
	at javax.imageio.ImageIO.read(ImageIO.java:1282)
	at com.mmc.image.J2KTest.main(J2KTest.java:19)

Time to take a look at that Coordinate System record in the header using the debugger to look at the IIOMetadata object.

package com.mmc.image;

import java.io.File;
import java.io.IOException;
import javax.imageio.ImageIO;
import javax.imageio.ImageReader;
import javax.imageio.metadata.IIOMetadata;
import javax.imageio.stream.ImageInputStream;
import com.sun.media.imageioimpl.plugins.jpeg2000.J2KImageReaderSpi;

public class JP2_metadata {
	public static void main(String[] args) {
    try {
	    //final File file = new File("C:/temp/JP2_test.jp2");
			final File file = new File("C:/temp/JP2_test.jp2");

			final ImageInputStream iis = ImageIO.createImageInputStream(file);
			J2KImageReaderSpi readerSPI = new J2KImageReaderSpi();
			final ImageReader reader = readerSPI.createReaderInstance();
			reader.setInput(iis);

			IIOMetadata iioMetadata = reader.getImageMetadata(0);

		} catch (IOException e) {
			e.printStackTrace();
		}
	}
}

“no Java!” Just the same problem reading the metadata?

Exception in thread "main" java.lang.NullPointerException
	at com.sun.media.imageioimpl.plugins.jpeg2000.J2KMetadata.replace(J2KMetadata.java:962)
	at com.sun.media.imageioimpl.plugins.jpeg2000.J2KMetadata.addNode(J2KMetadata.java:631)
	at jj2000.j2k.fileformat.reader.FileFormatReader.readFileFormat(FileFormatReader.java:279)
	at com.sun.media.imageioimpl.plugins.jpeg2000.J2KMetadata.(J2KMetadata.java:135)
	at com.sun.media.imageioimpl.plugins.jpeg2000.J2KImageReader.getImageMetadata(J2KImageReader.java:419)
	at com.mmc.image.JP2_metadata.main(JP2_metadata.java:22)

Too bad I don’t have source code for IIOMetadata. I’d like to step into IIOMetadata and look at header boxes. Perhaps I could actually see where the problem is, maybe even fix it. However, now I’ve put way too much time into this. Also I’m not sure where to find the JP2 sources for JAI ImageIO. Next stop a purchase of Global Mapper and finish the job before I get into trouble. Global Mapper is adequate and it meets my predisposition to Jurassic solutions i.e. cheap.

Summary:

I really tried to find a solution in the Open Source world. Sometimes it just doesn’t work out. Fortunately Global Mapper was able to read the 4band images as a composite and export the images as GTiff. Global Mapper would have taken days for export of full 800Mb jp2 images to GTiff. Fortunately an added bonus let me overlay my clipping page rectangle vectors and then export clipped subsets, which are much smaller and faster. The only drawback is the manual tedium of selecting one clip export at a time. Perhaps, if I knew much about Global Mapper, there is a handy scripting approach. In the meantime, I was able to finish the job on time and I was only a little disappointed about not finding an Open Source solution.

Postscript

USDA added another clue to my puzzle, noting that apparently there is some difference in the interpretation of the standards for JPEG 2000. I imagine that could be the underlying issue, but in the meantime 2008 NAIP 4 band imagery has some anomalous problems with viewer access out here in the wider world.