My EPSG:54004 mystery solved!


EPSG:54004 problem fig 1
Fig 1 – DIA Looks a Lot Better!

With a helpful comment from SharpGIS I was able to finally pin down my baffling problem with EPSG:54004.

The problem is in the datums.

ESRI:54004

PROJCS["World_Mercator",
    GEOGCS["GCS_WGS_1984",
        DATUM["WGS_1984",
            SPHEROID["WGS_1984",6378137,298.257223563]],
        PRIMEM["Greenwich",0],
        UNIT["Degree",0.017453292519943295]],
    PROJECTION["Mercator_1SP"],
    PARAMETER["False_Easting",0],
    PARAMETER["False_Northing",0],
    PARAMETER["Central_Meridian",0],
    PARAMETER["Standard_Parallel_1",0],
    UNIT["Meter",1],
    AUTHORITY["EPSG","54004"]]

As Morten pointed out the 54004 datum includes a flattening, 298.257223563 :
SPHEROID["WGS_1984",6378137,298.257223563]],

So 54004 should be treated as an Ellipsoid rather than a Sphere.

There is a subtle difference in 900913. If you notice 900913 also includes a flattening:
SPHEROID["WGS_1984",6378137,298.257223563]],

EPSG:900913

PROJCS["Google Mercator",
    GEOGCS["WGS 84",
        DATUM["World Geodetic System 1984",
            SPHEROID["WGS 84",6378137.0,298.257223563,
                AUTHORITY["EPSG","7030"]],
            AUTHORITY["EPSG","6326"]],
        PRIMEM["Greenwich",0.0,
            AUTHORITY["EPSG","8901"]],
        UNIT["degree",0.017453292519943295],
        AXIS["Geodetic latitude",NORTH],
        AXIS["Geodetic longitude",EAST],
        AUTHORITY["EPSG","4326"]],
    PROJECTION["Mercator_1SP"],
    PARAMETER["semi_minor",6378137.0],
    PARAMETER["latitude_of_origin",0.0],
    PARAMETER["central_meridian",0.0],
    PARAMETER["scale_factor",1.0],
    PARAMETER["false_easting",0.0],
    PARAMETER["false_northing",0.0],
    UNIT["m",1.0],
    AXIS["Easting",EAST],
    AXIS["Northing",NORTH],
    AUTHORITY["EPSG","900913"]]

However, you might not notice in addition it includes an explicit minor axis parameter.
PARAMETER["semi_minor",6378137.0],
And this minor axis is identical to the major axis. The axis definition overrides the flattening in the Datum and is probably technically incorrect, but the idea was just to get a spherical mercator into a definition that people could use to match Google Maps. I’ve seen this definition in PostGIS, GeoServer, and OpenLayers.

I had already noticed this and played with a MercatorEllipsoid function to see if that would fix my problem. However, sadly, I made an error and miscalculated eccentricity. The real equation for e goes something like this:
double f = 1 / 298.257223563;
e = Math.Sqrt(2*f – Math.Pow(f,2));

resulting in e = 0.081819190842621486;

Once I made the correction for the proper eccentricity in MercatorEllipsoid, ESRI:54004 lines up with the EPSG:3857. DIA is back in its rightful place.

My MercatorEllipsoid function now calculates correct BBOX parameters for GetMap requests, but only for com.esri.wms.Esrimap services. Looks like ESRI is the expert and correctly produces ESRI:54004 with Datum as defined. However, not so with GeoServer.

GeoServer seems to ignore the flattening 298.257223563 or else assume it is like the 900913 where flattening is overridden by a minor axis parameter:
semi_minor axis = major.PARAMETER[semi_minor,6378137.0]

This leads to some problems. My WMS client now has to decide which service correctly interprets the DATUM on 54004. For now I just check for “com.esri.wms.Esrimap” in the WMS url and change datums accordingly. This will undoubtedly lead to problems with other services since I don’t yet know how MapServer or others treat 54004.

Summary

  1. ESRI is right again!
  2. Always check your math one more time
  3. Community responses produce answers

Thanks everyone!

EPSG Problems


EPSG:54004 problem fig 1
Fig 1 – Something is not right with DIA?

Life is full of problems. With eight kids we have lots of problems, especially as school gets started up again. Life in the chaos lane has its moments, but the problem that has been baffling me this week has to do with EPSG coordinate systems and in particular the venerable but unofficial EPSG:54004.

EPSG:54004 is really not official. You won’t find it in the official EPSG registry. It is commonly used, however, ever since it was added as a “user extension” by Google, originally for use in Google maps.

Here are the basics of EPSG:54004

PROJCS["World_Mercator",
    GEOGCS["GCS_WGS_1984",
        DATUM["WGS_1984",
            SPHEROID["WGS_1984",6378137,298.257223563]],
        PRIMEM["Greenwich",0],
        UNIT["Degree",0.017453292519943295]],
    PROJECTION["Mercator_1SP"],
    PARAMETER["False_Easting",0],
    PARAMETER["False_Northing",0],
    PARAMETER["Central_Meridian",0],
    PARAMETER["Standard_Parallel_1",0],
    UNIT["Meter",1],
    AUTHORITY["EPSG","54004"]]
It is similar to EPSG:900913

PROJCS["Google Mercator",
    GEOGCS["WGS 84",
        DATUM["World Geodetic System 1984",
            SPHEROID["WGS 84",6378137.0,298.257223563,
                AUTHORITY["EPSG","7030"]],
            AUTHORITY["EPSG","6326"]],
        PRIMEM["Greenwich",0.0,
            AUTHORITY["EPSG","8901"]],
        UNIT["degree",0.017453292519943295],
        AXIS["Geodetic latitude",NORTH],
        AXIS["Geodetic longitude",EAST],
        AUTHORITY["EPSG","4326"]],
    PROJECTION["Mercator_1SP"],
    PARAMETER["semi_minor",6378137.0],
    PARAMETER["latitude_of_origin",0.0],
    PARAMETER["central_meridian",0.0],
    PARAMETER["scale_factor",1.0],
    PARAMETER["false_easting",0.0],
    PARAMETER["false_northing",0.0],
    UNIT["m",1.0],
    AXIS["Easting",EAST],
    AXIS["Northing",NORTH],
    AUTHORITY["EPSG","900913"]]

And both are similar to the official, official EPSG:3857

They should be equivalent. All are using earth radius 6378137.0m with a spherical mercator approximation that happens to be very efficient for quad tree tile pyramids. It has been so successful that it’s used in all the online map services, Google, Bing, Yahoo, and Open Street Map. It is also pretty simple to convert latitude, longitude into spherical Mercator coordinates. Using my dog eared John Snyder’s Manual I get something like this:

private Point Mercator(double lon, double lat)
{

        /* spherical mercator for Google, Bing, Yahoo
        * epsg:900913, epsg:54004, epsg:3857 R= 6378137
        * x = longitude
        * y= R*ln(tan(pi/4 +latitude/2)
        */············

        double R = 6378137.0;
        double x = R * Math.PI / 180 * lon;
        double y = R * Math.Log(Math.Tan(Math.PI / 180 *
         (45 + lat / 2.0)));
        return new Point(x, y);
}

The Problem

Now on to the problem, apparently not all EPSG:54004s are equal?

Looking at the screen shot above you can see that the MRLC NLCD 2001 impervious Surface layer using EPSG:54004 shows DIA a few miles south of the Bing Maps DIA? Bing Maps, as noted above, uses the EPSG:3857 spherical mercator which should be the same, but apparently is not? A quick check against OSM and Yahoo base maps verifies the same discrepancy.

I’ve been known to get tangled up in coordinate systems, so just as a sanity check here is a screenshot with DRCOG’s Network layer also rendered as EPSG:54004:

EPSG:54004 problem fig 2
Fig 2 – DRCOG EPSG:54004 looks right with Bing Maps but wrong with MRLC?

The DRCOG EPSG:54004 and Bing Maps align correctly so it is now 2 to 1 against the MRLC EPSG:54004. But what is going on here?

Taking a peek at DRCOG’s GetCapabilities I can see it is served out of GeoServer. “http://gis.drcog.org:80/geoserver/wms” Here is the user epsg section of GeoServer’s epsg.properties file from my local copy of GeoServer:

54004=PROJCS["WGS84 / Simple Mercator", GEOGCS["WGS 84", DATUM["WGS_1984",
SPHEROID["WGS_1984", 6378137.0, 298.257223563]],
PRIMEM["Greenwich", 0.0], UNIT["degree", 0.017453292519943295]],
PROJECTION["Mercator_1SP_Google"],
PARAMETER["latitude_of_origin", 0.0],
PARAMETER["central_meridian", 0.0],
PARAMETER["scale_factor", 1.0],
PARAMETER["false_easting", 0.0],
PARAMETER["false_northing", 0.0], UNIT["m", 1.0],
AXIS["x", EAST], AXIS["y", NORTH], AUTHORITY["EPSG","54004"]]

It looks familiar.

The MRLC GetCapabilites shows an ESRI server,
“http://ims.cr.usgs.gov:80/wmsconnector/com.esri.wms.Esrimap.”
Perhaps this is a clue. I know of several other ESRI WMS services so going to the National Atlas, I can see, yes, it publishes an EPSG:54004 SRS: <SRS>EPSG:54004</SRS>

I plug the National Atlas into my handy WMSClient and take a look. I selected the major highways and interstates which are visible, with opacity turned up, against the Impervious Surface backdrop. You can see a high degree of correlation between the two layers, and both are no where near the Bing base map?

EPSG:54004 problem fig 3
Fig 3 – National Atlas EPSG:54004 agrees with impervious surface?

OK, now it’s 2 to 2 and the end of the ninth. I am clear that the EPSG:54004 discrepancy has something to do with a difference in interpretation between WMS servers. It appears that ESRI is consistently different from Bing Maps and GeoServer. ESRI is the expert in this area, so my thought is there may be some unusual EPSG consideration that ESRI has caught and others possibly have not.

A bit of Googling brings up an older ESRI technical article that may have some bearing. However, nothing very clear about this large of a discrepancy. Going with a hunch I tried pretending ESRI EPSG:54004 was ellipsoidal rather than spherically based. Doing a quick switch to an Ellipsoidal Mercator calculation, however, makes no discernable difference.

private Point MercatorEllipsoid(double lon, double lat)
{
        double a = 6378137.0;// meters equitorial radius            
        double e = 1 / 298.257223563;
        double ctrmerid = 0.0;

        double px = lon * Math.PI / 180;
        double py = lat * Math.PI / 180;
        double pc = ctrmerid * Math.PI / 180;

        double x = a * (px - pc);
        double y = a * (Math.Log(Math.Tan((45.0 * Math.PI / 180) + py / 2.0)
        * Math.Pow((1.0 - e * Math.Sin(py))
        / (1.0 + e * Math.Sin(py)), e / 2.0)));
        return new Point(x, y);

}

Switching to the universal EPSG:4326 lined things up well enough. Unfortunately it does nothing to solve the mystery.

EPSG:4326  fig 3
Fig 4 – Switching to EPSG:4326 brings layers into alignment at this zoom level

Summary

Mystery is not resolved. If others have some insight, I’d like to hear it, especially from experts with ESRI com.esri.wms.Esrimap experience. If anyone wants to play with these EPSG:54004 layers I have a simple WMS client here: Silverlight WMSClient.

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.

A look at Silverlight 3 Projection Properties


Silverlight Map Control
Fig 1 – Example of Silverlight 3 Projection

Silverlight 3 introduces some new features. One that I was anxious to try out is the ‘Projection’ property. Charles Petzold took some time out from reading English novels to post some interesting experiments with this new feature.

Since I have been playing with the LidarServer WMS, I decided to apply some projection property enhancement to my LidarServer WMS viewer. The folks at USGS had partnered with some other government agencies to collect a set of LiDAR for the entire Denver metro area. This was needed for planning and security arrangements surrounding the Democratic National Convention last summer. This data was all made available to the public, presumably not terrorists. Mark Eaton, the NSDI liaison for Colorado, was kind enough to load this rather large data set onto a hard drive for me and I passed it over to the LidarServer folks to publish as a sample. This is a great resource and I hope other metro governments around the country can do the same. The really nice thing about the LiDAR for Denver, at least for my experiment, is that it includes lots of highrise buildings.

This is important only because LidarServer has incorporated a very useful extension to their WMS service. They added an additional non-spec request called GetProfileView. Basically this type of request gets the point cloud perpendicular to a profile line. The direction is determined by a kind of right hand rule pointing your thumb along a profile line with fingers curled in the direction of the cloud profile. Like many standards there are edges that are usefully broken. By adding the GetProfileView request to the WMS repertoire LidarServer opens a whole new dimension to viewing LidAr, the popularly known 3rd dimension.

Of course the desktop GIS folks have enjoyed this for some time, however, we in the browser interface world find this all fairly novel. The advent of Silverlight 3 starts things rolling down the road to 3D, but as Charles Petzold points out Silverlight 3 applies 3D transforms to 2D elements. Turning 2D elements into 3D objects is still problematic. I used Petzold’s first approach, utilizing the PlaneProjection class to set up my experiment.

First I had to have a way to draw profile lines and set GetProfileView properties, so I added a new expander panel. In the process I ran across a really helpful bit of code: http://silverlight.net/forums/p/12467/40404.aspx. There are more sophisticated DragDrop controls like DragDrop Manager, but the first code I ran across did the job well and was easy to understand. Adding any Silverlight control element inside of Canvas makes it draggable:

<Canvas x:Name="SpinProfileCanvas" Grid.Column="0" Grid.Row="0" Grid.RowSpan="2"
 HorizontalAlignment="Left" VerticalAlignment="Top">
    <local:DragDropPanel x:Name="SpinProfileDragPanel">
    .
    .
    </local:DragDropPanel>
</Canvas>

I quickly made all my panels draggable using this code and changing all my menu panels from MapControl.MapLayer to plain ordinary Canvas layouts coextensive to the Map Control.

I started the draw profile code using the Polyline shape class. However, the transition from a Point collection to individual segments required by the GetProfileView request was cumbersome, especially since I wanted a rubber line hint while drawing. I finally switched to just a List collection of Line shapes. I could then iterate through the List and build all the requests for a series of profile line segments. These then appear in a ScrollViewer wrapped StackPanel inside a DragDropPanel as a horizontal series of Image elements.

Image img = new Image();
img.Source = new BitmapImage(new Uri(profileServer + lidarname + "?
  SERVICE=WMS&VERSION=1.3&REQUEST=GetProfileView
  &FORMAT=image/png&EXCEPTIONS=INIMAGE&CRS=EPSG:3785
  &LEFT_XY=" + p0.X + "%2C" + p0.Y + "&RIGHT_XY=" + p1.X + "%2C" + p1.Y +
  "&DEPTH=" + DepthText.Text + "&SHOW_DRAPELINE=" + drapeline +
  "&SHOW_GRATICULE=" + graticule + "&WIDTH=800&HEIGHT=" + profileHeight +
  "&COLORIZATION=" + colorization + "&FILTER=" + filter));
bdr.Child = img;
profileImages.Children.Add(bdr);

Because these are image sources it isn’t necessary to write a proxy service to get around the non-site of origin security restrictions. The result is a set of profiles, one for each line segment. I also added a rollover scroll affect so that as you roll over a line segment the associated profile scrolls into view.

Silverlight Map Control
Fig 2 – First Iteration of GetProfileView using a simple ScrollView

User Note: LidarServer imposes a 1000m limit on GetProfileView segment lengths. You’ll need to zoom to within a few blocks. Sorry, getting a profile of the Denver skyline isn’t possible. I already tried that. It does make sense to limit things to reasonable sizes.

This is interesting but not really Silverlight 3 Projection yet. The next enhancement added a rectangle draw tool to get four sides. These four sides could then be plugged into a set of Grids with Silverlight 3 Projection attributes:

<Grid x:Name="frontgrid"
  Grid.Column="0" Grid.Row="1"  Grid.ColumnSpan="2" >
  <Grid.Projection>
    <PlaneProjection x:Name="front"
      LocalOffsetZ="100" RotationY="0" />
  </Grid.Projection>

This Grid can then be populated with the profile image similar to the previous ScrollViewer , but the beauty here is that now the sides of my rectangle profiles can each be transformed into its own plane. Finally these planes can all be part of a Storyboard animation:

<Storyboard x:Name="RotationAnimation">
  <DoubleAnimation Storyboard.TargetName="front"
    Storyboard.TargetProperty="RotationY"
    By="360" Duration="0:0:20"
    RepeatBehavior="Forever" />

Now I can draw a rectangle around a building such as Colorado’s Capitol, and see the four elevations of LiDAR point clouds rotating in a projection cube, actually ‘almost cube’. There is no top or bottom. I suppose I could grab a WMS GetMap for the top and even set a GetMap Filter to ‘Ground’ in order to fill in the bottom of a cube for a complete isometric view arrangement. That will be another project. I am also thinking about an n sided polygon profile tool, rather than the limited rectangle profile.

The Rectangle shape has some problems for stretching an area of interest. The Xaml Rectangle shape is defined by an upper left point along with width and height int attributes. However, the width and height throw out of range errors when negative which complicates life for an arbitrary stretchy rectangle draw. You’ll see from the example that it only draws from upper left to lower right. I decided to keep it this way for the time being, rather than using a 4sided polygon. One beneficial side affect is that the resulting rectangle shape has well defined north, south, east, west sides without resorting to a clockwise algorithm and a bit more code.

In the meantime I added some buttons for pause/resume on my animation as well as a grow/shrink on the elevation size. This lets me see some of the taller buildings where the lower profiles just don’t cut it. It may be worth adding a stretch capability to the Image Panel to let these adjustments be arbitrary sizes. You notice that the Projection attributes can also nest other elements like buttons and text. Buttons still have their full event driven nature. You gotta love recursion! So I popped a full screen button onto my rotating profile image. You have to be quick to click a moving button, but the reward is a new “draggable” panel with a 2500px max dimension profile request. Lots of detail. Besides you can always pause the rotation nonsense first.

Not bad for a few days work! I liked this little exercise because it illustrates the power of a declarative XML graphics grammar like xaml. We are not talking tons of development effort to get a rather sophisticated browser interface. Silverlight is raising the bar on user expectations and it will be interesting to see how this works out in the more traditional world of GIS.

Summary:
Silverlight 3 raises the bar on expectations, but also portends change by dragging more stuff into the client. Traditional GIS was very much a desktop deal, however, you can see changes coming, and perhaps more rapidly than we suppose. We may be nearing a technology tipping point here.

This all reminds me of an earlier transition in GIS, when the data layer broke away from monolithic GIS systems, becoming an independent ‘plug & play’ component with various DB vendors taking chunks of GIS functionality for themselves. In this case we’re seeing the browser client starting to take chunks out of the traditional GIS middle tier. With Silverlight, vectors have moved out to the client. The necessary number of round trips found in older web technologies made these types of interfaces prohibitive. (Anybody still remember ArcIMS?) Silverlight shrinks the desktop performance advantage and moves some more chunks of GIS out of the server. I don’t know if this is good or bad, but it sure is fun!