SharpGIS - Spatial References#GIS from a .NET developer's perspective
http://webservices.iter.dk/
http://www.rssboard.org/rss-specificationBlogEngine.NET 3.1.0.1en-UShttp://webservices.iter.dk/opml.axdhttp://www.dotnetblogengine.net/syndication.axdMorten NielsenSharpGIS0.0000000.000000Why EPSG:4326 is usually the wrong “projection”<p>
…or why spherical coordinate systems are not flat!
</p>
<p>
One of the most common ways the round world is displayed on a map is using the simplest projection we have:
</p>
<div>
<pre style="border-style: none; margin: 0em; padding: 0px; overflow: visible; font-size: 8pt; width: 100%; color: black; line-height: 12pt; font-family: consolas,'Courier New',courier,monospace; background-color: #f4f4f4">
x = longitude
y = latitude
</pre>
</div>
<p>
The name of this projection is “Plate Carree”, and is widely used because it is so simple. However we often seem to forget that we are talking about a projection. Therefore the spatial reference for this projection is very often (mis)referenced as a spherical coordinate system like the following for EPSG:4326:
</p>
<p>
GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],
<br />
PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]]
</p>
<p>
GEOGCS denotes “Geographical Coordinate System”, meaning a spherical coordinate system using angles for latitude and longitudes, not X’s and Y’s. However, the moment we project this map onto a flat screen or piece of paper, we projected the coordinate system (using the simple formula above), so how can the spatial reference of our map be the above one?
</p>
<p>
I found that in ArcMap you can manually define a projected coordinate system that renders the same map, but correctly identifies the coordinate system as being projected. Below is the Well-known-text representation of this projection string (note how PROJCS denotes projected coordinate system, and how EPSG:4326 is defined inside it):
</p>
<p>
PROJCS["World_Plate_Carree_Degrees",<em>GEOGCS["GCS_WGS_1984",
<br />
DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],
<br />
PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]</em>],
<br />
PROJECTION["Plate_Carree"],PARAMETER["False_Easting",0],
<br />
PARAMETER["False_Northing",0],PARAMETER["Central_Meridian",0],
<br />
UNIT["Degree",111319.49079327357264771338267056]]
</p>
<p>
As far as I know we don’t have a spatial reference ID for the above (<a href="http://spatialreference.org/ref/esri/54001/" target="_blank">EPSG:54001</a> comes close but uses a different linear unit), so we often incorrectly start using EPSG:4326 to describe a projection. This has been done for so long, that this is now common (mal)practice. Even OGC’s WMS specification allows you request flat maps in a spherical coordinate system, even though this doesn’t really make any sense.
</p>
<p>
So why is this a problem? Well until recently it hasn’t really been a problem, mostly because it was rare that the spherical coordinate were correctly handled, and most applications assumed that the world was flat. However, Microsoft SQL Server 2008 can correctly handle spherical coordinates, and suddenly this becomes a major problem.
</p>
<p>
Below is a query on all the counties that follow the red line defined by two points. The blue polygons is the result returned. The map claims that its spatial reference is EPSG:4326. Because of the curvature Earth, a straight line between the two endpoints looks like a curve on the projected map, so the results returned seems like they don’t match with the line, but they are in fact the counties on the line between the two endpoints if you think in a spherical coordinate system.
</p>
<p>
<img style="border: 0px none ; display: block; float: none; margin-left: auto; margin-right: auto" src="http://webservices.iter.dk/image.axd?picture=image.png" border="0" alt="image" title="image" width="624" height="203" />
</p>
<p>
The real problem here is not that the line is not really a curve (I drew the line like that because I want the features along that latitude). The problem is that the line didn’t know that it was projected, because the map it was displayed on was incorrectly set to EPSG:4326. However had the application known that this was not a spherical coordinate system, but a Plate Carree projection, it would have known that when converting from the flat coordinate system to the spherical coordinate system, it should ensure that my line follows the latitude. But because the input line already claims to be in spherical units, the application can not know that it needs to do anything to it (aside the fact that the application of course knows that this came from a flat screen map, but the business logic is of course separate from the UI).
</p>
<p>
This brings me to another even worse issue: I’d like to make the claim that almost <strong><em>all of the data you have that claims to be in EPSG:4326 has never been EPSG:4326</em></strong> ! (point data excluded).
</p>
<p>
Example: Most of the northern US border roughly follows the 49’th latitude. Let’s for sake of argument define the whole northern border as two points, from coast to coast, 49 degrees north. If I then think that my data is in EPSG:4326 (or any other geographic coordinate system), the border will NOT be following the 49th latitude, but go along a great circle that cuts into Canada.
</p>
<p>
So let us thank Microsoft for creating a database that handles spherical coordinates correctly, and for giving us major headaches when trying to handle these things correctly in a clean clear way :-). Of course there wouldn’t be any headaches if we never mixed spherical and flat coordinate systems in the first place.
</p>
http://webservices.iter.dk/post/2009/02/06/Why-EPSG4326-is-usually-the-wrong-e2809cprojectione2809d
http://webservices.iter.dk/post/2009/02/06/Why-EPSG4326-is-usually-the-wrong-e2809cprojectione2809d#commenthttp://webservices.iter.dk/post.aspx?id=8815ec4a-644c-4bce-9bbb-36ff836c8415Fri, 06 Feb 2009 09:00:00 +0000Spatial ReferencesSqlServer SpatialMortenhttp://webservices.iter.dk/pingback.axdhttp://webservices.iter.dk/post.aspx?id=8815ec4a-644c-4bce-9bbb-36ff836c84154http://webservices.iter.dk/trackback.axd?id=8815ec4a-644c-4bce-9bbb-36ff836c8415http://webservices.iter.dk/post/2009/02/06/Why-EPSG4326-is-usually-the-wrong-e2809cprojectione2809d#commenthttp://webservices.iter.dk/syndication.axd?post=8815ec4a-644c-4bce-9bbb-36ff836c8415Spherical/Web Mercator: EPSG code 3785<p>
I just received an update from the EPSG mailing list:
</p>
<p>
New to Version 6.15 are (among other things): <em><strong>Added spherical Mercator coordinate operation method and associated CRS as seen in popular web mapping and visualisation applications.</strong></em>
</p>
<p>
It looks like they FINALLY added the spherical Mercator / Web Mercator projection used in <a href="http://maps.live.com" target="_blank">Virtual Earth</a> and <a href="http://maps.google.com" target="_blank">Google Maps</a>.
</p>
<p>
This is a big surprise. EPSG’s earlier statement whether to include it was this:
</p>
<p>
"<em>We have reviewed the coordinate reference system used by Microsoft, Google, etc. and believe that it is technically flawed. We will not devalue the EPSG dataset by including such inappropriate geodesy and cartography.</em>”
</p>
<p>
Guess they changed their mind, or did they just devalue their dataset? Then again, judging from the remarks EPSG put in there, their arrogance still shines through. There´s absolute nothing wrong with using a sphere instead of a flattened sphere. Sure it's not as accurate as for instance WGS84, but then again WGS84 is not accurate either - no ellipsoid is. But we know the exact differences between the two, and as always you will need to take these things into account so I don´t see the real issue. Viisually the distortion is far less than what you would notice, and when doing area, distance and bearing calculations you would first of all never use the mercator units without taking the projection distortion into account, and if you do your calculationg in long/lat it's more or less just as easy to use WGS84 as a base for your calculations (since no datum transform is really needed).
</p>
<p>
Anyway, finally we get an official code for "Web Mercator": EPSG:3785
</p>
<p>
Here are the details of the entry:
</p>
<p>
<strong>Full WKT with authorities</strong> (untested!)<strong>:</strong><font face="Courier New"><br />
PROJCS["Popular Visualisation CRS / Mercator", GEOGCS["Popular Visualisation CRS", DATUM["Popular Visualisation Datum", SPHEROID["Popular Visualisation Sphere", 6378137, 0, AUTHORITY["EPSG",7059]], TOWGS84[0, 0, 0, 0, 0, 0, 0], AUTHORITY["EPSG",6055]], PRIMEM["Greenwich", 0, AUTHORITY["EPSG", "8901"]], UNIT["degree", 0.0174532925199433, AUTHORITY["EPSG", "9102"]], AXIS["E", EAST], AXIS["N", NORTH], AUTHORITY["EPSG",4055]], PROJECTION["Mercator"], PARAMETER["False_Easting", 0], PARAMETER["False_Northing", 0], PARAMETER["Central_Meridian", 0], PARAMETER["Latitude_of_origin", 0], UNIT["metre", 1, AUTHORITY["EPSG", "9001"]], AXIS["East", EAST], AXIS["North", NORTH], AUTHORITY["EPSG",3785]]</font>
</p>
<p>
<strong>Projected CRS</strong><br />
<font face="courier new,courier">COORD_REF_SYS_CODE: 3785<br />
COORD_REF_SYS_NAME: Popular Visualisation CRS / Mercator <br />
AREA_OF_USE_CODE: 3544 <br />
COORD_REF_SYS_KIND: projected <br />
COORD_SYS_CODE: 4499 <br />
DATUM_CODE: <br />
SOURCE_GEOGCRS_CODE: 4055 (see below)<br />
PROJECTION_CONV_CODE: 19847 <br />
CMPD_HORIZCRS_CODE:<br />
CMPD_VERTCRS_CODE:<br />
CRS_SCOPE: Certain Web mapping and visualisation applications. <br />
REMARKS: Uses spherical development. Relative to an ellipsoidal development errors of up to 800 metres in position and 0.7% in scale may arise. Some applications call this WGS 84. It is not a recognised geodetic system: see WGS 84 / World Mercator (CRS code 3395) <br />
INFORMATION_SOURCE: Microsoft.<br />
DATA_SOURCE: OGP<br />
REVISION_DATE: 3/14/2008<br />
CHANGE_ID: <br />
SHOW_CRS: TRUE<br />
DEPRECATED: FALSE</font>
</p>
<p>
<strong>Geographic CRS:</strong><br />
<font face="courier new,courier">COORD_REF_SYS_CODE: 4055<br />
COORD_REF_SYS_NAME: Popular Visualisation CRS<br />
AREA_OF_USE_CODE: 31262<br />
COORD_REF_SYS_KIND: geographic 2D<br />
COORD_SYS_CODE: 6422<br />
DATUM_CODE: 6055 (see below)<br />
SOURCE_GEOGCRS_CODE:<br />
PROJECTION_CONV_CODE:<br />
CMPD_HORIZCRS_CODE:<br />
CMPD_VERTCRS_CODE:<br />
CRS_SCOPE: Certain Web mapping and visualisation applications.<br />
REMARKS: Some applications erroneously call this WGS 84. It uses a sphere with a radius having the same value as the semi-major axis of the WGS 84 ellipsoid. There is no geodetic recognition of this system.<br />
INFORMATION_SOURCE: Microsoft.<br />
DATA_SOURCE: OGP<br />
REVISION_DATE: 3/13/2008<br />
CHANGE_ID: <br />
SHOW_CRS: TRUE<br />
DATUM_CODE: FALSE</font><br />
<br />
<strong>Datum:</strong><br />
<font face="courier new,courier">Code: 6055 <br />
Datum Name: Popular Visualisation<br />
Datum Type: geodetic<br />
Origin Description: Not specified in the classical sense of defining a geodetic datum.<br />
Datum Epoch: <br />
Ellipsoid Code: 7059 (see below)<br />
Prime Meridian Code: 8901<br />
Area Code: 1262<br />
Datum Scope : Used by certain popular Web mapping and visualisation applications.<br />
Remarks: Not recognised by geodetic authorities.<br />
Information Source: Microsoft.<br />
Data Source: OGP <br />
Revision Date: 13-Mar-08<br />
Change ID:<br />
Deprecated: No</font><br />
<br />
<strong>Ellipsoid:</strong><br />
<font face="courier new,courier">Code: 7059<br />
</font><font face="courier new,courier">Ellipsoid Name: Popular Visualisation Sphere <br />
Semi-major axis (a): 6378137 <br />
Axes units code: 9001 <br />
Inverse flattening (1/f): <br />
Semi-minor axis (b): 6378137 <br />
Ellipsoid?: No <br />
Remarks: Sphere with radius equal to the semi-major axis of the GRS80 and WGS 84 ellipsoids. Used only for Web approximate mapping and visualisation. Not recognised by geodetic authorities. <br />
Information Source: Microsoft.<br />
Data Source: OGP<br />
Revision Date: 14-Mar-08<br />
Change ID: <br />
Deprecated?: No </font>
</p>
http://webservices.iter.dk/post/2008/05/15/SphericalWeb-Mercator-EPSG-code-3785
http://webservices.iter.dk/post/2008/05/15/SphericalWeb-Mercator-EPSG-code-3785#commenthttp://webservices.iter.dk/post.aspx?id=01f8b508-9ccf-4337-b716-30143af1a42bThu, 15 May 2008 17:11:00 +0000GISSpatial ReferencesepsgmercatorprojectionsMortenhttp://webservices.iter.dk/pingback.axdhttp://webservices.iter.dk/post.aspx?id=01f8b508-9ccf-4337-b716-30143af1a42b9http://webservices.iter.dk/trackback.axd?id=01f8b508-9ccf-4337-b716-30143af1a42bhttp://webservices.iter.dk/post/2008/05/15/SphericalWeb-Mercator-EPSG-code-3785#commenthttp://webservices.iter.dk/syndication.axd?post=01f8b508-9ccf-4337-b716-30143af1a42bStraight lines on a sphere<p>
There was <a href="http://forums.microsoft.com/MSDN/ShowPost.aspx?PostID=2668748&SiteID=1">a question</a> in the <a href="http://forums.microsoft.com/MSDN/ShowForum.aspx?ForumID=1629&SiteID=1">Sql Server Spatial/Katmai forums</a> about why a straight line on a sphere didn't "look straight".
</p>
<p>
Katmai both have planar and spherical geometry types. The spherical type makes it easier to handle things like straight lines, distances and areas over large distances, and depending on which type you choose, you can get very different results.
</p>
<p>
Consider traveling on a straight line from Vancouver to Frankfurt. They are both located roughly at latitude 50°N. So the line would follow latitude 50°N right? Well that's true for a <em>planar</em> coordinate system:
</p>
<p align="center">
<img style="width: 291px; height: 246px" src="http://webservices.iter.dk/image.axd?picture=StraightLines_Planar.jpg" border="0" alt="" width="291" height="246" />
</p>
<p>
So according to <a href="http://www.theflatearthsociety.org/">The Flat Earth Society</a> (who actually believes the Earth is flat) the green line above shows the shortest distance between the two cities. But what we are actually seeing is an artifact of working with a planar coordinate system depicting a world that is not planar.
</p>
<p>
Using a <em>spheric</em> coordinate system, we get a different result. When looking at this line in a projected coordinate system, like the Mercator projection I've used here, the line doesn't look straight (left image below). However if we look at this in a 3D view (right image), it becomes apparent that the "curved line" is actually the shortest line between the two cities.
</p>
<p align="center">
<img style="width: 549px; height: 240px" src="http://webservices.iter.dk/image.axd?picture=StraightLines_Spheric.jpg" border="0" alt="" width="549" height="240" />
</p>
<p>
The difference between using the two geometry types gets really apparent when you try to find the shortest distance from New York to the line in the two coordinate systems. In the planar case, Quebec is approximately the point where the line is the closest to New York. In the spherical case, the line never gets closer to New York than from Vancouver where it started out.
</p>
<p>
Another way that I like to think of why following a latitude is not the shortest route, is imagining two guys standing a few feet from the North Pole with the pole between them. Which way is the shortest path between them: Walking straight over the pole or following the latitude at 89.999°N around it? Now imagine how this line would look on a map.
</p>
<p align="center">
<img style="width: 473px; height: 251px" src="http://webservices.iter.dk/image.axd?picture=StraightLines_NorthPole.png" border="0" alt="" width="473" height="251" />
</p>
<p align="left">
We can extend this further and have the two men go further and further away from each other. If one man gets a little further south than the other, the line wouldn't pass the pole, but go by just next to the it - it still wouldn't follow the latitude. Since the North Pole is just a point we humans "made up", we can abstract it and think of a pole in the ground anywhere on Earth. The bottom line is that the North Pole and the latitudes are all imaginary, and doesn't really exist, and definitely doesn't help describing the shortest distance.
</p>
<p align="left">
Only at Equator would it be a good idea to follow the latitude, which is because this is the only latitude that is a <a href="http://en.wikipedia.org/wiki/Great_circle">great circle</a>, and it so happens to be that any straight line between two points will describe parts of a great circle. Longitudes however are all great circles, and the line between the two men describe parts of a longitude. The reason we think that the horizontal line between Vancouver and Germany is the shortest path, is simply because we have to <em>project</em> the round Earth onto a flat piece of paper, and this causes all sorts of things to get distorted.
</p>
http://webservices.iter.dk/post/2008/01/12/Straight-lines-on-a-sphere
http://webservices.iter.dk/post/2008/01/12/Straight-lines-on-a-sphere#commenthttp://webservices.iter.dk/post.aspx?id=b21d3512-cd55-417c-96fd-c62d8a34a83eSat, 12 Jan 2008 08:39:00 +0000Spatial ReferencesSqlServer Spatialsql serverprojectionsMortenhttp://webservices.iter.dk/pingback.axdhttp://webservices.iter.dk/post.aspx?id=b21d3512-cd55-417c-96fd-c62d8a34a83e9http://webservices.iter.dk/trackback.axd?id=b21d3512-cd55-417c-96fd-c62d8a34a83ehttp://webservices.iter.dk/post/2008/01/12/Straight-lines-on-a-sphere#commenthttp://webservices.iter.dk/syndication.axd?post=b21d3512-cd55-417c-96fd-c62d8a34a83eAccurate distance calculations in UTM projections<p>
Recently a friend of mine asked me what projection would be the best to create point-buffer circles. All map projections adds distortion so you rarely end up with a circle. If you’re lucky, you’ll have an ellipse but often they get even more complex than this. In many cases a “normal” circle is sufficient for accuracy, but still it got me thinking how I would do this 100% accurate. There are many approaches to this, fx. creating a new projection for each point, do the math in spherical coordinate systems etc, but none of them are completely trivial.
</p>
<p>
The UTM projection is one of the most commonly used projections. It’s fairly accurate to measure from point to point within small distances and close to the meridian, so this post is based on that. Most of the math here can also be applies to Mercator, where the scale reductions are applied to Y instead of X.
</p>
<p>
<img style="float: right" src="http://webservices.iter.dk/image.axd?picture=utmproj_01.png" alt="" width="175" height="220" />The <em>Mercator</em> projection is created by projecting the globe onto a cylinder whose center follows the Earth’s rotation axis. The <em>Transverse</em> Mercator is basically a "tilted" cylinder parallel to equator. Along the line that the cylinder "touches" the surface the distortion is 0, and increasing away from it in the east/west direction. There is no distortion in the North/South direction. The <em>Universal</em> Transverse Mercator is slightly smaller than the globe, so it will intersect the surface two places. This also ensures that the distortion between and around these two lines are minimal. Along the meridian (the center of the two lines) the distortion for the UTM projection is 0.9996, meaning that this is the amount you have to reduce an infinitely small line in the east/west direction (For Mercator this is normally 1 along Equator). As you move away from the meridian, the scale factor increases up toward infinity. Normally UTM projections are used close to the meridian, and every time you get too far east/west, you would switch to use a new UTM projection. That’s why The Earth is divided into 60 UTM ‘zones’, one for each 6 degrees. Often for practical reasons you might want to use a UTM zone even though it’s outside its defined usage area. Here you have to be specifically careful when measuring distances. For example, Denmark is covered by zone 32 and 33, but it’s common to only use zone 32 which means that distance calculations at the east end can be very inaccurate without applying a scale reduction.
</p>
<p>
<strong><em>Warning</em></strong>: The following contains a lot of math. You can skip to the bottom and download a C# class that does all the calculations for you.
</p>
<p>
The scale factor SR(x) for any given point at distance 'x' from the meridian is given by:
</p>
<p>
<img src="http://webservices.iter.dk/image.axd?picture=utmproj_02.png" alt="" width="149" height="45" />
</p>
<p>
Since this is the scale reduction at a point, we need to sum up all of these values along a line to find the correct distance, thus the <em>east/west</em> distance from e<sub>1</sub> to e<sub>0</sub> (expressed in distance from the meridian) becomes an integral expression:
</p>
<p>
<img src="http://webservices.iter.dk/image.axd?picture=utmproj_04.png" alt="" width="339" height="64" />
</p>
<p>
To this you need to add the north/south distance using the well-known Euclidean distance formula.
</p>
<p>
For atypical UTM using WGS84 ellipsoid and a 500000m false easting it becomes:
</p>
<p>
<img src="http://webservices.iter.dk/image.axd?picture=utmproj_06.png" alt="" width="363" height="69" />
</p>
<p>
If you haven’t refreshed your integral math lately, or just want to write this in pure code, one can also write the expression as:
</p>
<p>
<img src="http://webservices.iter.dk/image.axd?picture=utmproj_08.png" alt="" width="646" height="34" />
</p>
<p>
Where m is the scalefactor at the meridian, r is the earth radius and e<sub>1</sub> and e<sub>0</sub> have the false easting subtracted.
</p>
<p>
Often we can do with an approximate value by using the scale reduction between the start and end point. This should only be used on distances <1km in easting because it gets very inaccurate for greater distances.
</p>
<p>
<img src="http://webservices.iter.dk/image.axd?picture=utmproj_10.png" alt="" width="274" height="71" />
</p>
<strong>Simpson’s Rule</strong>
<p>
One can also use <a href="http://en.wikipedia.org/wiki/Simpsons_rule" target="_blank">Simpson’s rule</a> for calculating the integral solution. My tests showed it’s a very accurate approximation for this type of integral usually giving millimeter accuracy. It’s defined as:
</p>
<p>
<img src="http://webservices.iter.dk/image.axd?picture=utmproj_12.png" alt="" width="423" height="67" />
</p>
<p>
I performed a few performance tests for 10 million calculations using all four methods for finding a distance in a UTM projection. Below is the processing time for each method:
</p>
<table border="1" cellspacing="0" cellpadding="4">
<tbody>
<tr>
<td>Euclidean</td>
<td>4.38 s</td>
</tr>
<tr>
<td>Average / Linear</td>
<td>4.49 s</td>
</tr>
<tr>
<td>Simpsons</td>
<td>6.94 s</td>
</tr>
<tr>
<td>Integral</td>
<td>9.26 s</td>
</tr>
</tbody>
</table>
<p>
Bottom line of this: If you want performance and distance is small use the averaging method, if you want accuracy, go with the Integral method, or if you want both, use the Simpsons approach.
</p>
<strong>Distance calculation - example</strong>
<p>
start point : { E=400,000 ; N=6,100,000 }
</p>
<p>
end point : { E=600,000 ; N=6,250,000 }<br />
<img src="http://webservices.iter.dk/image.axd?picture=utmproj_14.png" alt="" width="394" height="120" />
</p>
<ul>
<li>Distance using Simpsons rule: 249,942.5569m (error: 0.0003m) </li>
<li>Distance with approximate scale reduction: 249936.0046m (error: 6.54m) </li>
<li>Distance without scale reduction: 250,000m (error=57.46m)</li>
</ul>
<strong>Where is scale reduction=1?</strong>
<p>
Solving SR(x) =1 gives us the distance from the meridian where the scale reduction is exactly 1:
</p>
<p>
<img src="http://webservices.iter.dk/image.axd?picture=utmproj_18.png" alt="" width="555" height="47" />
</p>
<strong>Area calculation</strong> <span>To calculate the area of a polygon, apply the scale reduction to the line segments and then calculate the area as you normally would.</span><span><strong>Creating perfect circles</strong>
<p>
To create a circle with a fixed radius taking scale reduction into account will give you an irregular shape in the projection. You can use the approximate method (using center as scale factor) which will give you an ellipse since the scale factor for x is constant in this case. The formula for approximate circle becomes r<sup>2</sup>=(SR*x)<sup>2</sup>+y<sup>2</sup> where SR is the scale reduction at the center of the circle.
</p>
<p>
For "accurate" circles SR becomes dependent on X. To find this solution you need to isolate e<sub>1</sub> from the above equations to find the scale corrected distance. My math is too rusty to figure that one out yet, but if you know, feel free to post it in the comments.
</p>
<p>
You can download a C# class that performs all four types of distance calculations <a href="http://webservices.iter.dk/content/binary/utmdistance.zip">here</a>.
</p>
</span>
http://webservices.iter.dk/post/2007/10/14/Accurate-distance-calculations-in-UTM-projections
http://webservices.iter.dk/post/2007/10/14/Accurate-distance-calculations-in-UTM-projections#commenthttp://webservices.iter.dk/post.aspx?id=07723618-a7cd-4292-968a-9ac2eb0472abSun, 14 Oct 2007 00:54:00 +0000GISSpatial ReferencesprojectionsMortenhttp://webservices.iter.dk/pingback.axdhttp://webservices.iter.dk/post.aspx?id=07723618-a7cd-4292-968a-9ac2eb0472ab0http://webservices.iter.dk/trackback.axd?id=07723618-a7cd-4292-968a-9ac2eb0472abhttp://webservices.iter.dk/post/2007/10/14/Accurate-distance-calculations-in-UTM-projections#commenthttp://webservices.iter.dk/syndication.axd?post=07723618-a7cd-4292-968a-9ac2eb0472abThe Microsoft Live Maps and Google Maps projection<p>
I have lately seen several blogposts confused about which datum and projection <a href="http://maps.live.com/">Microsoft’s Live Maps (Virtual Earth)</a> and <a href="http://maps.google.com/">Google Maps</a> use. As most people already know by now, they render the round earth onto a flat screen using a Mercator projection.
</p>
<p>
I think the confusion comes from that they actually use two spatial reference systems at the same time:
</p>
<ol>
<li>Geographic Longitude/Latitude coordinatesystem based on the standard WGS84 datum. </li>
<li>Mercator projection using a datum based on WGS84, BUT modified to be spheric.</li>
</ol>
<p>
So when is what used?
</p>
<p>
The Javascript API’s use (1) as input when you want to add points, lines and polygons. That is, they expect you to input any geometry in geographical coordinates, and click events etc. will also return geometry in this spatial reference. This is the coordinate system most javascript developers will use. The API will automatically project it to the spheric mercator projection.
</p>
<p>
If you want to create image overlays, or roll your own tile server on top of the map, you will need to project your images into a spheric mercator projection. The JavaScript APIs are not able to do this for you.
</p>
<p>
Here’s a bit of facts about the two projections:
</p>
<p>
The valid range of (1) is: <font face="Courier New">[-180,-85.05112877980659]</font> to <font face="Courier New">[180, 85.05112877980659]</font>.
</p>
<p>
The valid range of (2) is: <font face="Courier New">[-20037508.3427892, -20037508.3427892]</font> to <font face="Courier New">[20037508.3427892, 20037508.3427892]</font>
</p>
<p>
<strong>Well-known Text for (1):</strong><br />
<font face="Courier New">GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.0174532925199433]]</font>
</p>
<p>
<strong>Well-known Text for (2):</strong><br />
<font face="Courier New">PROJCS["Mercator Spheric", GEOGCS["WGS84basedSpheric_GCS", DATUM["WGS84basedSpheric_Datum", SPHEROID["WGS84based_Sphere", 6378137, 0], TOWGS84[0, 0, 0, 0, 0, 0, 0]], PRIMEM["Greenwich", 0, AUTHORITY["EPSG", "8901"]], UNIT["degree", 0.0174532925199433, AUTHORITY["EPSG", "9102"]], AXIS["E", EAST], AXIS["N", NORTH]], PROJECTION["Mercator"], PARAMETER["False_Easting", 0], PARAMETER["False_Northing", 0], PARAMETER["Central_Meridian", 0], PARAMETER["Latitude_of_origin", 0], UNIT["metre", 1, AUTHORITY["EPSG", "9001"]], AXIS["East", EAST], AXIS["North", NORTH]]</font>
</p>
<p>
<strong>Proj.4 definition for (1):<br />
</strong><font face="Courier New">+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs</font>
</p>
<p>
<strong>Proj.4 definition for (2)</strong> (<a href="http://proj.maptools.org/faq.html">see here</a> for an explanation of the weird ’nadgrids’ parameter):<br />
<font face="Courier New">+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +no_defs</font>
</p>
<p>
So why this weird extent for the latitude? First of all, the poles in the Mercator projection extends towards infinity, so at some point they have to cut them off (and who cares about ice anyway - apart from that its melting). If you look at (2) these latitude/longitude values projected into the spheric mercator results in a perfect square, fitting perfectly with squared image tiles, that are simple to sub-divide over and over again, as you zoom in. I expect the reason for the spheric datum is for simplicity and perfomance when reprojecting points from longitude/latitude to screen coordinates. <a href="http://cfis.savagexi.com/">Charlie Savage</a> also has a <a href="http://cfis.savagexi.com/articles/2006/05/03/google-maps-deconstructed">more mathematical approach</a> to deriving these values.
</p>
<p>
For a quick introduction to projections, coordinate systems and datums <a href="http://webservices.iter.dk/2007/05/05/SpatialReferencesCoordinateSystemsProjectionsDatumsEllipsoidsConfusing.aspx">see here</a>.
</p>
<p>
If you want to know more about how these mapping api's work, keep an eye on <a href="http://deconstructingmappingapis.blogspot.com/">Jayant's blog</a>.
</p>
<p>
<strong>Update</strong>: We now have an offical EPSG code for the projection. <a href="http://webservices.iter.dk/post/2008/05/SphericalWeb-Mercator-EPSG-code-3785.aspx">See details here</a>.
</p>
http://webservices.iter.dk/post/2007/07/27/The-Microsoft-Live-Maps-and-Google-Maps-projection
http://webservices.iter.dk/post/2007/07/27/The-Microsoft-Live-Maps-and-Google-Maps-projection#commenthttp://webservices.iter.dk/post.aspx?id=686498e1-72b6-48f9-b4b4-06137c427cbbFri, 27 Jul 2007 05:06:00 +0000GISGoogle MapsSpatial ReferencesVirtual EarthmercatorprojectionsMortenhttp://webservices.iter.dk/pingback.axdhttp://webservices.iter.dk/post.aspx?id=686498e1-72b6-48f9-b4b4-06137c427cbb1http://webservices.iter.dk/trackback.axd?id=686498e1-72b6-48f9-b4b4-06137c427cbbhttp://webservices.iter.dk/post/2007/07/27/The-Microsoft-Live-Maps-and-Google-Maps-projection#commenthttp://webservices.iter.dk/syndication.axd?post=686498e1-72b6-48f9-b4b4-06137c427cbbSpatial references, coordinate systems, projections, datums, ellipsoids – confusing?<p>
People are often mixing the above as if they were one and the same, so here’s a recap of them. One of the things you often find people saying is that “my data is in the WGS84 coordinate system”. This doesn’t really make sense, but I will get back to this later.
</p>
<p>
This is a very confusing subject, and I might have gotten a few things wrong myself, so please add a comment and I’ll update it ASAP.
</p>
<p>
<strong>Coordinate systems</strong>
</p>
<p>
A coordinate system is simply put a way of describing a spatial property relative to a center. There is more than one way of doing this:
</p>
<ul>
<li>The Geocentric coordinate system is based on a normal (X,Y,Z) coordinate system with the origin at the center of Earth. This is the system that GPS uses internally for doing it calculations, but since this is very unpractical to work with as a human being ( due to the lack of well-known concepts of east, north, up, down) it is rarely displayed to the user but converted to another coordinate system.<br />
</li>
<li>The Spherical or <a href="http://webhelp.esri.com/arcgisdesktop/9.2/index.cfm?TopicName=About_geographic_coordinate_systems">Geographic coordinate system</a> is probably the most well-known. It is based on angles relative to a prime meridian and Equator usually as Longitude and Latitude. Heights are usually given relative to either the mean sea level or the datum (I’ll get back to the datum later).<br />
</li>
<li>The Cartesian coordinate system is defined as a “flat” coordinate system placed on the surface of Earth. In some projections it’s not flat in the sense that it follows the earth’s curvature in one direction and has a known scale-error in the other direction relative to the distance of the origin. The most well-known coordinate system is the Universal Transverse Mercator (UTM), but surveyors define their own little local flat coordinate systems all the time. It is very easy to work with, fairly accurate over small distances making measurements such as length, angle and area very straightforward. Cartesian coordinate systems are strongly connected to projections that I will cover later.</li>
</ul>
<p>
<em>Sidenote: The geocentric coordinate system is strictly speaking a cartesian coordinate system too, but this is the general terms I've seen used the most when talking about world coordinate systems.</em>
</p>
<p align="center">
<img style="width: 201px; height: 227px" src="http://webservices.iter.dk/image.axd?picture=Geocentric.png" border="0" alt="Geocentric.png" title="Geocentric.png" width="201" height="227" /> <img style="width: 192px; height: 215px" src="http://webservices.iter.dk/image.axd?picture=Spherical.png" border="0" alt="Spherical.png" title="Spherical.png" width="192" height="215" /> <img style="width: 203px; height: 216px" src="http://webservices.iter.dk/image.axd?picture=Cartesian.png" border="0" alt="Cartesian.png" title="Cartesian.png" width="203" height="216" />
</p>
<p>
<strong>Datums and ellipsoids</strong>
</p>
<p>
Some of the common properties of the above coordinate systems are that they are all relative to the center of Earth and except the Geocentric coordinate system, uses a height system relative to the surface of the earth.
</p>
<p>
This poses two immediate problems:
</p>
<ul>
<li>Where is the center of the earth </li>
<li>What is the shape of the earth?</li>
</ul>
<p>
<img style="float: right; width: 350px; height: 180px" src="http://webservices.iter.dk/image.axd?picture=Geoids_sm.jpg" border="0" alt="" width="350" height="180" />By now most people should know that that the earth isn’t flat (although <a href="http://en.wikipedia.org/wiki/Flat_Earth_Society">there are still some who doubts it</a>). If we define the surface of Earth as being at the mean sea level (often referred to as the Geoid), we don’t get a spheroid or even an ellipsoid. Because of gravitational changes often caused by large masses such as mountain ranges etc, Earth is actually very irregular with variations of +/- 100 meters. Since this is not very practical to work with as a model of earth, we usually use an ellipsoid for approximation. The ellipsoid is defined by its semi-major axis, and either the flattening of the ellipoid or the semi-minor axis.
</p>
<p>
The center and orientation of the ellipsoid is what we call the datum. So the datum defines an ellipsoid and through the use of a set of points on the ground that we relate to points on the ellipsoid, we define the center of the Earth. This poses another problem, because continental drift moves the points used to define the points around all the time. This is why the name of a datum usually have a year in it, often referring to the position of those points January 1st of that year (although that may vary).
</p>
<p>
There are a vast amount of datums, some used for measurements all over the world, and other local datums defined so they fit very well with a local area. Some common ones are: <a href="http://en.wikipedia.org/wiki/WGS84">World Geodetic Datum 1984</a> (WGS84), <a href="http://en.wikipedia.org/wiki/ED50">European Datum 1950</a> (ED50) and <a href="http://en.wikipedia.org/wiki/NAD83">North American Datum 1983</a> (NAD83).
</p>
<p>
The most well-known is WGS84 used by the GPS systems today. It is a good approximation of the entire world and with fix-points defined almost all over the world. When it was defined they forgot to include points in Europe though, so the Europeans now have their own ETRS89, which is usually referred to as the “realization of WGS84 in Europe”. The problem here was solely because of continental drift, so they defined some points relative to WGS84 in 1989, and keeps track of the changes. In most use-cases it is of no real importance and you can use one or the other.
</p>
<p>
I mentioned earlier that people often refer to having their data in WGS84, and you see now why this doesn’t make sense. All you know from that is that the data is defined using the WGS84 datum, but you don’t know which coordinate system it uses.
</p>
<p>
Read more on <a href="http://webhelp.esri.com/arcgisdesktop/9.2/index.cfm?TopicName=Datums">Datums</a> and <a href="http://webhelp.esri.com/arcgisdesktop/9.2/index.cfm?TopicName=Spheroids_and_spheres">Spheroids</a>.
</p>
<p>
<strong>Projections</strong>
</p>
<p>
The earth isn’t flat, and there is no simple way of putting it down on a flat paper map (or these days onto a computer screen), so people have come up with all sorts of ingenious solutions each with their pros and cons. Some preserves area, so all objects have a relative size to each other, others preserve angles (conformal) like the <a href="http://webhelp.esri.com/arcgisdesktop/9.2/index.cfm?TopicName=Mercator">Mercator projection</a>, some try to find a good intermediate mix with only little distortion on several parameters etc. Common to them all is that they transform the world onto a flat Cartesian coordinate system, and which one to choose depends on what you are trying to show.
</p>
<p>
A common statement that I hear in GIS is the following “My map doesn’t have a projection”, but this is simply not possible (unless you have a good old rotating globe). Often people are referring to data that is in longitude/latitude and displayed on a map without having specified any projection. What happens is that the system applies the simplest projection it can: Mapping Longitude directly to X and Latitude to Y. This results in an equirectangular projection, also called the <a href="http://webhelp.esri.com/arcgisdesktop/9.2/index.cfm?TopicName=Plate_Carree">“Plate Carree” projection</a>. It results in very heavy distortion making areas look squashed close to the poles. You can almost say that the “opposite” of the Plate Carree is the <a href="http://webhelp.esri.com/arcgisdesktop/9.2/index.cfm?TopicName=Mercator">Mercator projection</a> which stretches areas close to the poles in the opposite direction, making them look very big. Mercator is the type of projection you see used on <a href="http://maps.live.com">Live maps</a> and <a href="http://maps.google.com">Google maps</a>, but as many often mistakenly thinks, they do NOT use WGS84 for the projected map, although WGS84 is used when you directly input longitude/latitude values using their API (<a href="http://webservices.iter.dk/2007/07/27/TheMicrosoftLiveMapsAndGoogleMapsProjection.aspx">read more on this here</a>).
</p>
<p>
<a href="http://webhelp.esri.com/arcgisdesktop/9.2/index.cfm?id=87&pid=86&topicname=About_projected_coordinate_systems">More on projected coordinate systems</a>
</p>
<p>
<strong>Spatial reference</strong>
</p>
<p>
The spatial reference is a combination of all the above. It defines an ellipsoid, a datum using that ellipsoid, and either a geocentric, geographic or projection coordinate system. The projection also always has a geographic coordinate system associated with it. The <a href="http://www.epsg.org/">European Petroleum Survey Group</a> (EPSG) has a huge set of predefined spatial references, each given a unique ID. These ID’s are used throughout the industry and you can download an Access database with all them from their website, as well as some very good documents on projection (or see <a href="http://spatialreference.org/">the Spatial References website</a>).
</p>
<p>
So when you hear someone saying they have their data in WGS84, you can often assume that they have longitude/latitude data in WGS84 projected using Plate Carree. The spatial reference ID of this is EPSG:4326.
</p>
<p>
Spatial References are often defined in a Well-known format defining all these parameters. The Spatial Reference EPSG:4326 can therefore also be written as:
</p>
<p>
<font face="Courier New">GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.0174532925199433]]</font>
</p>
<p>
As mentioned Live/Google maps use a Mercator projection, but although their datum is based on WGS84, they use a sphere instead of an ellipsoid. This means that they use the same center and orientation as WGS84, but without applying any flattening. The spatial reference string for their projection therefore becomes:
</p>
<p>
<font face="Courier New">PROJCS["Mercator Spheric", GEOGCS["WGS84based_GCS", DATUM["WGS84based_Datum", SPHEROID["WGS84based_Sphere", 6378137, 0], TOWGS84[0, 0, 0, 0, 0, 0, 0]], PRIMEM["Greenwich", 0, AUTHORITY["EPSG", "8901"]], UNIT["degree", 0.0174532925199433, AUTHORITY["EPSG", "9102"]], AXIS["E", EAST], AXIS["N", NORTH]], PROJECTION["Mercator"], PARAMETER["False_Easting", 0], PARAMETER["False_Northing", 0], PARAMETER["Central_Meridian", 0], PARAMETER["Latitude_of_origin", 0], UNIT["metre", 1, AUTHORITY["EPSG", "9001"]], AXIS["East", EAST], AXIS["North", NORTH]]</font>
</p>
http://webservices.iter.dk/post/2007/05/05/Spatial-references2c-coordinate-systems2c-projections2c-datums2c-ellipsoids-e28093-confusing
http://webservices.iter.dk/post/2007/05/05/Spatial-references2c-coordinate-systems2c-projections2c-datums2c-ellipsoids-e28093-confusing#commenthttp://webservices.iter.dk/post.aspx?id=b25454ec-1dd8-46c3-8786-4e759e17f3c9Sat, 05 May 2007 21:44:00 +0000GISSpatial ReferencesprojectionsMortenhttp://webservices.iter.dk/pingback.axdhttp://webservices.iter.dk/post.aspx?id=b25454ec-1dd8-46c3-8786-4e759e17f3c94http://webservices.iter.dk/trackback.axd?id=b25454ec-1dd8-46c3-8786-4e759e17f3c9http://webservices.iter.dk/post/2007/05/05/Spatial-references2c-coordinate-systems2c-projections2c-datums2c-ellipsoids-e28093-confusing#commenthttp://webservices.iter.dk/syndication.axd?post=b25454ec-1dd8-46c3-8786-4e759e17f3c9Applying on-the-fly transformation in SharpMap<p>
I have received a lot of questions on how to transform data from one coordinatesystem to another on the fly in <a href="http://sharpmap.iter.dk/">SharpMap</a>. 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 <a href="http://ogr.maptools.org/">OGR </a>for this), it is easy to setup in <a href="http://sharpmap.iter.dk/">SharpMap</a>. Below are some examples on how to accomplish this.
</p>
<p>
<a href="http://sharpmap.iter.dk/">SharpMap</a> gives you the full power to specify all the parameters in a projection. The following method demonstrates how to setup a UTM projection:
</p>
.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; }
<pre class="csharpcode">
<span class="rem">/// <summary></span>
<span class="rem">/// Creates a UTM projection for the northern
/// hemisphere based on the WGS84 datum</span>
<span class="rem">/// </summary></span>
<span class="rem">/// <param name="utmZone">Utm Zone</param></span>
<span class="rem">/// <returns>Projection</returns></span>
<span class="kwrd">private</span> IProjectedCoordinateSystem CreateUtmProjection(<span class="kwrd">int</span> utmZone)
{
CoordinateSystemFactory cFac =
<span class="kwrd">new</span> SharpMap.CoordinateSystems.CoordinateSystemFactory();
<span class="rem">//Create geographic coordinate system based on the WGS84 datum</span>
IEllipsoid ellipsoid = cFac.CreateFlattenedSphere(<span class="str">"WGS 84"</span>,
6378137, 298.257223563, LinearUnit.Metre);
IHorizontalDatum datum = cFac.CreateHorizontalDatum(<span class="str">"WGS_1984"</span>,
DatumType.HD_Geocentric, ellipsoid, <span class="kwrd">null</span>);
IGeographicCoordinateSystem gcs = cFac.CreateGeographicCoordinateSystem(
<span class="str">"WGS 84"</span>, AngularUnit.Degrees, datum,
PrimeMeridian.Greenwich,
<span class="kwrd">new</span> AxisInfo(<span class="str">"Lon"</span>, AxisOrientationEnum.East),
<span class="kwrd">new</span> AxisInfo(<span class="str">"Lat"</span>, AxisOrientationEnum.North));
<span class="rem">//Create UTM projection</span>
List<ProjectionParameter> parameters = <span class="kwrd">new</span> List<ProjectionParameter>(5);
parameters.Add(<span class="kwrd">new</span> ProjectionParameter(<span class="str">"latitude_of_origin"</span>, 0));
parameters.Add(<span class="kwrd">new</span> ProjectionParameter(<span class="str">"central_meridian"</span>, -183+6*utmZone));
parameters.Add(<span class="kwrd">new</span> ProjectionParameter(<span class="str">"scale_factor"</span>, 0.9996));
parameters.Add(<span class="kwrd">new</span> ProjectionParameter(<span class="str">"false_easting"</span>, 500000));
parameters.Add(<span class="kwrd">new</span> ProjectionParameter(<span class="str">"false_northing"</span>, 0.0));
IProjection projection = cFac.CreateProjection(
<span class="str">"Transverse Mercator"</span>, <span class="str">"Transverse_Mercator"</span>, parameters);
<span class="kwrd">return</span> cFac.CreateProjectedCoordinateSystem(
<span class="str">"WGS 84 / UTM zone "</span>+utmZone.ToString() +<span class="str">"N"</span>, gcs,
projection, LinearUnit.Metre,
<span class="kwrd">new</span> AxisInfo(<span class="str">"East"</span>, AxisOrientationEnum.East),
<span class="kwrd">new</span> AxisInfo(<span class="str">"North"</span>, AxisOrientationEnum.North));
}
</pre>
<p>
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:
</p>
<p>
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"]]
</p>
<p>
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)
</p>
<!-- code formatted by http://manoli.net/csharpformat/ -->
<pre class="csharpcode">
<span class="rem">/// <summary></span>
<span class="rem">/// Create coordinatesystem based on a Well-Known text</span>
<span class="rem">/// </summary></span>
<span class="rem">/// <param name="wkt"></param></span>
<span class="rem">/// <returns></returns></span>
<span class="kwrd">private</span> ICoordinateSystem CreateCoordinateSystemFromWKT(<span class="kwrd">string</span> wkt)
{
CoordinateSystemFactory cFac = <span class="kwrd">new</span> CoordinateSystemFactory();
<span class="kwrd">return</span> cFac.CreateFromWkt(strProj);
}
</pre>
<p>
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:
</p>
<pre class="csharpcode">
((myMap.Layers[0] <span class="kwrd">as</span> VectorLayer).DataSource <span class="kwrd">as</span> ShapeFile).CoordinateSystem
</pre>
<p>
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:
</p>
<ul>
<li>Mercator 1-standard parallel (Mercator_1SP)
</li>
<li>Mercator 1-standard parallels (Mercator_2SP)
</li>
<li>Transverse mercator (Transverse_Mercator)
</li>
<li>Lambert Conic Conformal 2-standard parallel (Lambert Conic Conformal (2SP))
</li>
<li>Albers</li>
</ul>
<p>
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...)
</p>
<p>
The following shows how to create a transformation and apply it to a vectorlayer (only vector- and label-layers supports on-the-fly transformations):
</p>
<pre class="csharpcode">
<span class="rem">//Create zone UTM 32N projection</span>
IProjectedCoordinateSystem utmProj = CreateUtmProjection(32);
<span class="rem">//Create geographic coordinate system (lets just reuse the CS from the projection)</span>
IGeographicCoordinateSystem geoCS = utmProj.GeographicCoordinateSystem;
<span class="rem">//Create transformation</span>
CoordinateTransformationFactory ctFac = <span class="kwrd">new</span> CoordinateTransformationFactory();
ICoordinateTransformation transform =
ctFac.CreateFromCoordinateSystems(source, target);
<span class="rem">//Apply transformation to a vectorlayer</span>
(myMap.Layers[0] <span class="kwrd">as</span> VectorLayer).CoordinateTransformation = transform;
</pre>
<p>
Happy transforming!
</p>
http://webservices.iter.dk/post/2006/07/07/Applying-on-the-fly-transformation-in-SharpMap
http://webservices.iter.dk/post/2006/07/07/Applying-on-the-fly-transformation-in-SharpMap#commenthttp://webservices.iter.dk/post.aspx?id=24b94c5c-f29e-49ef-ab8a-2de3ebc64c62Fri, 07 Jul 2006 09:06:00 +0000SharpMapSpatial ReferencesprojectionsMortenhttp://webservices.iter.dk/pingback.axdhttp://webservices.iter.dk/post.aspx?id=24b94c5c-f29e-49ef-ab8a-2de3ebc64c623http://webservices.iter.dk/trackback.axd?id=24b94c5c-f29e-49ef-ab8a-2de3ebc64c62http://webservices.iter.dk/post/2006/07/07/Applying-on-the-fly-transformation-in-SharpMap#commenthttp://webservices.iter.dk/syndication.axd?post=24b94c5c-f29e-49ef-ab8a-2de3ebc64c62A basic survival guide to coordinate systems<p>
<a href="http://cfis.savagexi.com/">Charlie Savage</a> has <a href="http://cfis.savagexi.com/articles/2006/04/20/on-coordinate-systems">a great post</a> introducing coordinate systems in <a href="http://cfis.savagexi.com/">his blog</a>. I too have had many questions on coordinate systems, projections and datums, and Charlie have made a great introduction to the world of "spatial reference systems". He promises to post some more on the topic, so keep track of his blog if you find these topics confusing (I know I did when I first was told about datums and coordinate systems).
</p>
http://webservices.iter.dk/post/2006/04/20/A-basic-survival-guide-to-coordinate-systems
http://webservices.iter.dk/post/2006/04/20/A-basic-survival-guide-to-coordinate-systems#commenthttp://webservices.iter.dk/post.aspx?id=4149b36d-ceb9-4971-9cbe-cf9d23ac711bThu, 20 Apr 2006 08:39:00 +0000GISSpatial ReferencesprojectionsMortenhttp://webservices.iter.dk/pingback.axdhttp://webservices.iter.dk/post.aspx?id=4149b36d-ceb9-4971-9cbe-cf9d23ac711b0http://webservices.iter.dk/trackback.axd?id=4149b36d-ceb9-4971-9cbe-cf9d23ac711bhttp://webservices.iter.dk/post/2006/04/20/A-basic-survival-guide-to-coordinate-systems#commenthttp://webservices.iter.dk/syndication.axd?post=4149b36d-ceb9-4971-9cbe-cf9d23ac711b