SharpGIS

#GIS from a .NET developer's perspective

SharpMap v0.9 released

SharpMap v0.9 is now uploaded to the CodePlex Workspace.

Download it here.

Thanks to everyone who has supported the development of SharpMap. We have come a long way since the previous release candidate.

Unfortunately this also marks the end of my involvement with SharpMap. In a months time I'm moving from Copenhagen to California to start on a new and exciting GIS developer job. Hopefully more on that later...

Using XAML for rendering maps

I was just reading a few articles on XAML - Microsoft new UI language for rendering vector content and user interfaces - when it hit me that this probably can get you much further in web mapping interaction than SVG can. XAML is part of the "Windows Presentation Foundation" a part of .NET Framework 3.0

 

My first attempt to generate XAML was done by creating a new renderer for SharpMap, and below you can see my first results in XamlPad. The data here has been created based on a SharpMap map by my quick'n'dirty XamlRenderer for SharpMap.

 

 

 

All it takes with the new renderer is initializing the new renderer with your map, and request a new map:

SharpMap.Xaml.Renderer XamlRenderer = new SharpMap.Xaml.Renderer(map);

string xaml = XamlRenderer.GetMap();

 

If you zip the generated XAML file you even get a 20% reduction over a PNG image.

 

Next step is to add interaction like zooming, querying etc…

 

Download XAML for the above worldmap (261,36 KB)

Applying on-the-fly transformation in SharpMap

I have received a lot of questions on how to transform data from one coordinatesystem to another on the fly in SharpMap. Usually the problem is that they have data in different coordinatesystems and want to match them. Although I would recommend applying transformations once-and-for-all to increase performance (you could use OGR for this), it is easy to setup in SharpMap. Below are some examples on how to accomplish this.

SharpMap gives you the full power to specify all the parameters in a projection. The following method demonstrates how to setup a UTM projection:

.csharpcode { font-size: small; color: black; font-family: Courier New , Courier, Monospace; background-color: #ffffff; /*white-space: pre;*/ } .csharpcode pre { margin: 0em; } .csharpcode .rem { color: #008000; } .csharpcode .kwrd { color: #0000ff; } .csharpcode .str { color: #006080; } .csharpcode .op { color: #0000c0; } .csharpcode .preproc { color: #cc6633; } .csharpcode .asp { background-color: #ffff00; } .csharpcode .html { color: #800000; } .csharpcode .attr { color: #ff0000; } .csharpcode .alt { background-color: #f4f4f4; width: 100%; margin: 0em; } .csharpcode .lnum { color: #606060; }
/// <summary>
/// Creates a UTM projection for the northern
/// hemisphere based on the WGS84 datum
/// </summary>
/// <param name="utmZone">Utm Zone</param>
/// <returns>Projection</returns>
private IProjectedCoordinateSystem CreateUtmProjection(int utmZone)
{
CoordinateSystemFactory cFac = 
      new SharpMap.CoordinateSystems.CoordinateSystemFactory();
//Create geographic coordinate system based on the WGS84 datum
IEllipsoid ellipsoid = cFac.CreateFlattenedSphere("WGS 84", 
           6378137, 298.257223563, LinearUnit.Metre);
IHorizontalDatum datum = cFac.CreateHorizontalDatum("WGS_1984", 
                     DatumType.HD_Geocentric, ellipsoid, null);
IGeographicCoordinateSystem gcs = cFac.CreateGeographicCoordinateSystem(
                     "WGS 84", AngularUnit.Degrees, datum,
                     PrimeMeridian.Greenwich,
                     new AxisInfo("Lon", AxisOrientationEnum.East),
                     new AxisInfo("Lat", AxisOrientationEnum.North));
//Create UTM projection
List<ProjectionParameter> parameters = new List<ProjectionParameter>(5);
parameters.Add(new ProjectionParameter("latitude_of_origin", 0));
parameters.Add(new ProjectionParameter("central_meridian", -183+6*utmZone));
parameters.Add(new ProjectionParameter("scale_factor", 0.9996));
parameters.Add(new ProjectionParameter("false_easting", 500000));
parameters.Add(new ProjectionParameter("false_northing", 0.0));
IProjection projection = cFac.CreateProjection(
"Transverse Mercator", "Transverse_Mercator", parameters);
return cFac.CreateProjectedCoordinateSystem(
         "WGS 84 / UTM zone "+utmZone.ToString() +"N", gcs,
projection, LinearUnit.Metre,
new AxisInfo("East", AxisOrientationEnum.East),
new AxisInfo("North", AxisOrientationEnum.North));
}

If you have a well-known text-representation, you can also create a projection from this. A WKT for an UTM projection might look like this:

PROJCS["WGS 84 / UTM zone 32N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.01745329251994328,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",9],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AUTHORITY["EPSG","32632"]]

SharpMap comes with WKT parsers for parsing a WKT to a coordinate system (note: the current v0.9RC1 has a few bug in its WKT parser, but if you get problems parsing the WKT, use the current source from the repository, where these issues have been resolved)

/// <summary>
/// Create coordinatesystem based on a Well-Known text
/// </summary>
/// <param name="wkt"></param>
/// <returns></returns>
private ICoordinateSystem CreateCoordinateSystemFromWKT(string wkt)
{
CoordinateSystemFactory cFac = new CoordinateSystemFactory();
return cFac.CreateFromWkt(strProj);
}

If your data is based on shapefile data and they have a .prj file defining the coordinatesystem, you can simply retrieve the CS from the shapefile instead:

((myMap.Layers[0] as VectorLayer).DataSource as ShapeFile).CoordinateSystem

The next step is to create a transformation between two coordinate systems. SharpMap currently supports transforming between a geographic coordinate system and one of the following projections:

  • Mercator 1-standard parallel (Mercator_1SP)
  • Mercator 1-standard parallels (Mercator_2SP)
  • Transverse mercator (Transverse_Mercator)
  • Lambert Conic Conformal 2-standard parallel (Lambert Conic Conformal (2SP))
  • Albers

Unfortunately datum-shifts and transformations between two projections are still down the pipeline, but the above will be sufficient in most cases. (for those interested full transformation between all supported projections as well as datum-shifts are almost done...)

The following shows how to create a transformation and apply it to a vectorlayer (only vector- and label-layers supports on-the-fly transformations):

//Create zone UTM 32N projection
IProjectedCoordinateSystem utmProj = CreateUtmProjection(32);
//Create geographic coordinate system (lets just reuse the CS from the projection)
IGeographicCoordinateSystem geoCS = utmProj.GeographicCoordinateSystem;
//Create transformation
CoordinateTransformationFactory ctFac = new CoordinateTransformationFactory();
ICoordinateTransformation transform = 
   ctFac.CreateFromCoordinateSystems(source, target);
//Apply transformation to a vectorlayer
(myMap.Layers[0] as VectorLayer).CoordinateTransformation = transform;

Happy transforming!

Using OGR with SharpMap

Christian Graefe recently posted an update on his additions to the GDAL/OGR  C# wrapper and his OGR provider for SharpMap.

Apart from the data sources that SharpMap already supports, Christians OGR provider adds support for "MapInfo File", "TIGER", "S57", "DGN", "Memory", "CSV", "GML", "Interlis 1", "Interlis 2", "SQLite" and "ODBC". If I'm not completely wrong he is also working on a GDAL raster layer as well...

Look at http://www.volley-boy.de/home/sharpmap.aspx

You will need to install FWTools to use it.

Great work Christian !

Using NetTopologySuite in SharpMap

I've created a small demo showing how to utilize NetTopologySuite in SharpMap. The great thing about NetTopologySuite is that is has a bunch of geoprocessing methods that SharpMap doesn't have (yet anyway). The demo uses the NetTopologySuite library for the following:

  1. Apply some geoprocessing to the geometries (in this case it puts a buffer on all rivers)
  2. Uses NTS for making exact intersection queries. By clicking on an object the buffered feature(s) are highlighted and attribute data shown beneath.

The basics are:

Create a NTS datasource based on another datasource which takes a delegate to apply to all geometries:

NtsProvider ntsDataSource = new NtsProvider(myShapeFileDatasource, CreateBuffers);

Create the geometry-operation delegate method which should be applied to all features:

private void CreateBuffers(List<GisSharpBlog.NetTopologySuite.Features.Feature> features)
{
   foreach (GisSharpBlog.NetTopologySuite.Features.Feature feature in
features)
   feature.Geometry = feature.Geometry.Buffer(0.5);
}

The 'ntsDataSource' object can then be assigned as datasource to a vectorlayer.

Clicking and highlighting features are done by reacting to a click-event on an ImageButton (in ASP.NET) or a PictureBox (in Windows.Forms).

protected void imgMap_Click(object sender, ImageClickEventArgs e)
{
   //Convert click-point in image coordinates to world coordinates
   Point ClickPnt = myMap.ImageToWorld(new System.Drawing.Point(e.X, e.Y));
  
SharpMap.Data.FeatureDataSet ds = new SharpMap.Data.FeatureDataSet();
  
//Execute click-query on first layer in layers collection
  
(myMap.Layers[0] as SharpMap.Layers.VectorLayer).DataSource.ExecuteIntersectionQuery(ClickPnt, ds);
   if (ds.Tables.Count > 0) //We have a result
  
{
     
//Add clicked features to a new selection layer
     
SharpMap.Layers.VectorLayer laySelected = new SharpMap.Layers.VectorLayer("Selection");
      laySelected.DataSource =
new GeometryProvider(ds.Tables[0]);
      laySelected.Style.Fill =
new System.Drawing.SolidBrush(System.Drawing.Color.Yellow);
      myMap.Layers.Add(laySelected);
   }
   //Call method that renders a new map and displays it to the user
  
GenerateMap();
}

 

Download NetTopologySuite_Demo.zip (544,32 KB)

SharpMap v0.9 RC1 online !

SharpMap v0.9RC1 is now available for download from the SharpMap website.

Thanks to all the people who kept sending feedback and bug reports!

I've resolved a lot of bugs since beta 1, as well as added some obviously missing properties and methods in the API. Especially the Shapefile reader has been enhanced to handle a lot of special cases (in one case it was even able to render a file that ArcGIS claimed to be corrupt). On-the-fly transformation has also been completed and now supports transformation between geographic and projected coordinate systems.

Thanks to Humberto Ferreira for providing us with an Oracle Spatial provider as well. SharpMap now supports data from ESRI Shapefiles, PostGIS, Oracle, Microsoft SQL Server, Access and ECW raster.

Last but not least, I've spent a great deal of time enhancing the documentation, so go take a look at that as well. There are lots of code-samples hidden in there as well.

Make Microsoft SQL Server geospatial

I’ve always thought that on the spatial support, MSSQL was way behind many of the other database servers, in its lack of supporting storage of geometry. With the new .NET CLR you can actually add your own .NET-based object types and I’ve also tried implementing the Simple Features Specification in SQL Server. There are some limitations that made me give this up though. First of all, a user data type cannot be more than 8000 bytes. That is at most no more than 500 vertices in a geometry object, which is far too little for an accurate coastline for instance. Another problem is that SQL Server doesn’t support inheritance chains, so you can’t make a good object-oriented implementation of your datatype either.

…so yesterday I went for a completely different and much simpler approach. I decided to just store the geometry as Well-Known Binary in an image column. The reason for using an image column is that it can hold up to 2Gb of data, which should be sufficient for most geometry objects ;-). A binary field has the same 8000 byte limitation as UDT so this is no good. In addition to the geometry field, I create four real-type fields, holding the min/max values of the geometry envelope. This makes it efficient to do boundingbox based queries on the data. Any additional fields would be attributes for the geometry object.

I implemented the whole thing using SharpMap. First I created a small upload application that takes a shapefile, creates a table in the database and uploads the geometry and attributes to it. SharpMap has the necessary data readers and WKB formatters for this. The second part was to create a data provider from which SharpMap could draw from. I more or less based this on the PostGreSQL/PostGIS data provider for SharpMap, by changing the boundingbox query to use the four envelope fields. All this was not much more than an hour’s work, so it is very simple to accomplish.

I must say I was very surprised by the performance of the approach. It is just slightly faster than the shapefile data provider, which used to be the fastest data provider for SharpMap. In comparison the PostGreSQL/PostGIS is generally 4-6 times slower.

I have created a small demo web-application you can download from here. It contains two pages: one for uploading to the database, and one for rendering a layer from the database. All you need to do is to add an empty SQL Server 2005 Express database to the \App_Data\ folder and name id "GeoDatabase.mdf".

Download SharpMapSqlServer.zip (181,74 KB) (updated May 20, 2006)

Update: The MsSqlProvider is now also included in v0.9RC1, including a method for uploading from any of the SharpMap datasources to MS SQL.

Transformation problems with NetTopologySuite and GeoTools.NET

The last couple of days I've been working on doing on-the-fly transformations in SharpMap. I started out by implementing all the interfaces in the OGC "Coordinate Transformation Services" specification. These specifications are just great and makes my work of designing APIs so much eaier. Anyway when I started at the actual implementation I also looked at NetTopologySuite (which is a .NET port of Java Topology Suite and GeoTools.NET for a little help on the more complex transformation formulas.

It didn't really work out. All my results was either wrong, resulted in Overflow and DivisionByZero exceptions and NaN values. All I could then do was working through the EPSG transformation formulas, which was I was trying to avoid in the first place (oh by the way, EPSG recently released v7.10 of their coordinate system database).

Well I finally hacked it, and my conclusion is that there are several problems in the NTS code that causes very inaccurate (if not completely wrong) transformations. Problems like an incorrect constant for converting between radians and degrees, transformation formulas that are wrong (at least the Mercator projection - haven't checked the others yet), transformation-factories that doesn't respect the order of source and target coordinate systems (it will always transform from geographic to projected coordsys), and using the TransformList(array-of-points) method in many cases gives a different result than Transform(point).

I emailed my findings to Diego and Andrew who's doing a great job on NTS, so that they can benefit from this work well, but this is just to warn you people out there that NTS is probably giving your incorrect or inaccurate transformation results, until the bug is corrected. This is by no means Diego's fault; I guess he just copied it from GeoTools.NET (since the code is more or less the same), but this reminds me why we always use the national cadastras transformation libraries for transforming data...

Oh well as you can probably understand from this post, SharpMap now supports on-the-fly transforrmation! How cool is that? :-) If you can't wait for the next release, grap the source from the GotDotNet workspace. Warning though: I haven't checked the rest of the formulas thoroughly yet, but Mercator projection is checked and corrected and UTM seems to be working...

Using SharpMap as a Google Earth Network Link

Since Google Earth is so popular, I thought why not join in and have SharpMap serve maps for GE? It is quite simple to make your own GE Network Link, and here are a few pointers, and you can download a demo web application you can use as well.

The GE Network Link basically works by requesting an XML webpage (well Google calls it KML), telling GE where to get the image tiles and where they should be placed on the globe. GE also adds a 'BBOX=minx,miny,maxx,maxy' querystring to the request so you can return a custom KML response containing the most appropriate tile(s). A response could look like this:

<?xml version="1.0" encoding="UTF-8" ?>
<kml xmlns="
http://earth.google.com/kml/2.0">
<Document>
  <open>1</open>
  <GroundOverlay>
    <open>1</open>
    <Icon>
      <href>http://localhost/TileServer.ashx?BBOX=-180,-90,0,90&size=512</href>
    </Icon>
    <LatLonBox>
      <north>90</north>
      <south>-90</south>
      <east>0</east>
      <west>-180</west>
      <rotation>0</rotation>
    </LatLonBox>
  </GroundOverlay>
</Document>
</kml>

You can add as many "GroundOverlays" as you like. For instance you can cover a requested BBOX with as many tiles as you like. That way you can increase the image size/resolution or cover a larger or smaller area than initially requested while getting quicker responses. You can read more about the KML format here.

Creating the KML is pretty straight-forward in ASP.NET, so I wont cover that here (but download the example and see for yourself).

Creating the image tile is easily done using SharpMap. Add an HttpHandler to you webpage, and create the following ProcessRequest method to render and return the map:

public void ProcessRequest(HttpContext context)

{

    //Get tilesize in request url

    int tilesize = int.Parse(context.Request.QueryString["size"]);

    //Get boundingbox requested

    string[] strbox = context.Request.QueryString["BBOX"].Split(new char[] { ',' });

    SharpMap.Geometries.BoundingBox box = new SharpMap.Geometries.BoundingBox

            (double.Parse(strbox[0]), double.Parse(strbox[1]),

            double.Parse(strbox[2]), double.Parse(strbox[3]));

 

    //Call custom method that sets up the map with the requested tilesize

    SharpMap.Map myMap = MapHelper.InitializeMap(new System.Drawing.Size(tilesize, tilesize));

    //Center on requested tile and set the appropriate view width

    myMap.Center = box.GetCentroid();

    myMap.Zoom = box.Width;

 

    //Render the map

    System.Drawing.Bitmap b = (System.Drawing.Bitmap)myMap.GetMap();

    //Create a PNG image which supports transparency

    context.Response.Clear();

    context.Response.Expires = 10000000;

    context.Response.ContentType = "image/png";

    System.IO.MemoryStream MS = new System.IO.MemoryStream();

    //Set background to the transparent color which will be see-through in Google Earth

    b.MakeTransparent(myMap.BackColor);

    b.Save(MS, System.Drawing.Imaging.ImageFormat.Png);

    byte[] buffer = MS.ToArray();

    //Send image response

    context.Response.OutputStream.Write(buffer, 0, buffer.Length);

    // tidy up 

    b.Dispose();

    context.Response.End();

}

The last thing we need to do is add the network link to GE. The easiest way to do this is to do this from within GE. From the menu Select Add –> Network link. Type a name for your network link and set the location to the URL of the XML service. Under “View-Based Refresh” , set the “When” parameter to “After camera stops”, and click OK, and you should be set to go !

Unfortunately GE doesn't provide you with the same smooth image transitions that it use for its own tiles, so you will have to do with the rather crude way of showing big red squares until the tiles have been loaded. Furthermore the BBOX GE requests isn't always very accurate, especially if you are looking south at an angle you will notice the BBOX is very much off.

The small demo application you can download below shows a world map with each country’s population density shown on a colorscale from blue to red. Copy the files to a virtual directory (make sure it runs as its own web-application, at least the sub-folders are located at the root of the webapplication). Set the network-link to point to the URL of the default.aspx page.

Download GoogleEarthServer.zip (270,65 KB)

GE also supports vectordata to be served as KML. It could be interesting to see who comes up with a "SharpMap KML Vector Server" first :-)