Thursday, September 15, 2005

Overlaying Mercator Projected WMS Layers on GoogleMaps

One troubling aspect of the GoogleMaps API is their use of the Mercator Projection. This projection is not so bad for looking at large-scale maps (zoomed in), but is awful for visualizing small-scale (zoomed out) data, due to the large degree of exaggeration of land masses near the poles. Making matters more complicated is the fact that most Mercator Projections deal in meters as units while GoogleMaps delivers content to and from their API in decimal degrees. This makes things difficult especially when layering global WMS layers on GoogleMaps.

The screenshot to the left illustrates what happens when layering country boundaries from a WMS using a geographic projection on top of GoogleMaps. The interesting thing is that if we were to zoom in on a particular portion of the world, the boundaries would line up. This is why many people are able to use a geographic projection with Google Maps if only using for large-scale mapping projects.

OK, so what if we wanted to utilize a Mercator Projection from a WMS service on top of GoogleMaps? Here is how I did it, using Mapserver and Proj4 for the WMS content (details at the mapserver website):

1. Add the appropriate Mercator definition in the Proj4 EPSG file:
# Mercator
<54004> +proj=merc +lat_ts=0 +lon_0=0 +k=1.000000 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m no_defs <>

2. Alter the Mapserver Map File so the Output Projection looks like:
PROJECTION
"init=epsg:54004"
END

3. Add the epsg:54004 projection definition to the WEB section of the Mapfile like:
WEB
METADATA
"wms_srs" "EPSG:4326 EPSG:54004"
.....

4. Make sure that the layers in the mapfile reflect their native projection.

5. Alter the Javascript in the GMaps_WMSSpec_0.1.js code (references to the source of this code appear in an earlier post). Primarily, I added two functions to convert decimal degrees to Mercator Meters so the GoogleMaps API can call the WMS service in a language it can understand:

// dd2MercMeters Conversion Functions
// from http://search.cpan.org/src/RRWO/GPS-Lowrance-0.31/lib/Geo/Coordinates/MercatorMeters.pm
//
var MAGIC_NUMBER=6356752.3142;
var DEG2RAD=0.0174532922519943;
var PI=3.14159267;
function dd2MercMetersLng(p_lng) {
return MAGIC_NUMBER*(p_lng*DEG2RAD);
}
function dd2MercMetersLat(p_lat) {
if (p_lat >= 85) p_lat=85;
if (p_lat <= -85) p_lat=-85;
return MAGIC_NUMBER*Math.log(Math.tan(((p_lat*DEG2RAD)+(PI/2)) /2));
}

6. Make a couple of alterations in the WMSSpec.prototype.getWMSURL function. Attached here is the whole function as I code it:

WMSSpec.prototype.getWMSURL=function(a,b,c, url, layers, format) {
if (!this.Styles) this.Styles="default";
var ts = this.tileSize;
var ul = this.getLatLng(a*ts,(b+1)*ts, c);
var lr = this.getLatLng((a+1)*ts, b*ts, c);
// user mercator for small scale
if (map.getZoomLevel() > 12) {
var bbox = dd2MercMetersLng(ul.x) + "," + dd2MercMetersLat(ul.y) + "," + dd2MercMetersLng(lr.x)+ "," + dd2MercMetersLat(lr.y);
var l_strSRS="SRS=EPSG:54004";
} else {
var bbox = ul.x + "," + ul.y + "," + lr.x + "," + lr.y;
var l_strSRS="SRS=EPSG:4326";
}
var url = url + "REQUEST=GetMap&SERVICE=WMS&VERSION=1.1.1&LAYERS=" + layers + "&STYLES=" + this.Styles + "&FORMAT=" + format + "&BGCOLOR=0xFFFFFF&TRANSPARENT=TRUE&"+l_strSRS+"&BBOX=" + bbox + "&WIDTH=" + ts + "&HEIGHT=" + ;
this.cURL=url;
return url;
}


You'll notice in step #6 above that I don't use Mercator for resolutions below 13. This is because I believe the nifty little javascript re-project is not quite as accurate as it can be (hey-- what do you expect from javascript doing reprojections!). Anyway, reverting to latlong at lower zoom levels will fix this.

Thats it! The resulting map will looksomething like the screenshot below. The lines going across the top of the Russian Federation, by the way, I think are an artifact of drawing country boundaries across the international dateline. Since I'm only using point data for my real application i didn't worry about this. Not sure if this is an issue but I think it crept in because my source data is stored in geographic units.

Tuesday, September 06, 2005

Measuring flow in a small creek using a simple weir

Recently, I constructed a simple weir for measuring flow in Owen's Creek at the point where our ranch-road crosses the creek. There is an old metal pipe, 8' in diameter which the water flows through, under the road. In front of this pipe I laid a rough-cut 2x10 at the base by digging out gravel in the stream bed. On top of this I placed another rough-cut 2x10 with a four foot rectangular notch. By measuring the head (height of water flowing over the notch) one can determine the flow rate of the water in the stream (See table below).

One improvement I want to make is placing a piece of metal in front of the board and up an inch or so from the level of the notch. This will make a thinner area for the water to flow over which gives a more accurate measurement.

There are still some leaks from the bottom of the weir, running about 25gpm (a conservative guess on my part).
Following are some estimates of streamflow using this method, with an estimated accuracy of +/- 30gpm:


  • July 30th, 2005: 400gpm

  • September 4th, 2005: 315gpm

  • September 23rd, 2005: 250gpm

  • October 1, 2005: 577gpm (1 inch rain in last 24 hours after no rain all summer)

  • Assuming irrigation rates of 96-146gpm, I will be leaving a downstream flow of 165-215gpm flow (using September 4th measurements).

    DISCHARGE FROM A RECTANGULAR WEIR WITH A FOUR-FOOT WIDE CREST GIVEN IN US GALLONS PER MINUTE






















    Gauge Reading in InchesGPM
    0.547
    0.7586
    1145
    1.25199
    1.5258
    1.75320
    2387
    2.25476
    2.5552
    2.75628
    3708
    3.25792
    3.5880
    3.75968
    41080
    4.251176
    4.51268
    4.751364
    51472