How to create an online map with a non-Mercator Projection: Part 1

21 Sep 2009 | Articles · Geospatial · How To

Tags: · · ·

Previous articles in this series, discussed different coordinate systems and map projections available, and why different systems suit different applications (Part 1, 2). This was followed by a look at geostatistical and thematic maps where data areas and data densities are important (Part 1, 2, 3). Virtually all online map applications use the Mercator projection. This projection was designed for maritime navigation, and greatly distorts the relative areas of the polar regions. This distortion makes it an extremely poor choice for geostatistical and thematic maps, but it is currently the most popular map projection for such applications.  This is because most online map providers (eg. Google Maps and Bing Maps) only provide the Mercator projection.

I finished with a review of a number of different equal area map projections which would be suitable for thematic global maps. Unfortunately the review was theoretical. It is all very good knowing that a particular map projection should be used, but how do we use it?  How do we create a working online map that uses a map projection other than Mercator?

This two-part article shows you how to do exactly that. This first part will show you how to create a global basemap using a number of different map projections. This basemap will be rendered with UMN MapServer and delivered using the WMS protocol. The second part of the article will show you how to use OpenLayers to display overlay data from a source such as KML or GeoRSS with a matching map projection.

We will be using UMN MapServer to implement the WMS server. This is an open source map serving application that can read data from a variety of formats, and serve maps as images or as WMS tiles. We will be using it in WMS tile mode.

 

Creating a WMS MapScript file with the Mollweide Projection

This sample will use the Mollweide projection. This is a pseudo-cylindrical equal area projection. It is probably the most aesthetically pleasing equal area map projection with widespread software support.

MapServer is controlled using a MapScript file. We only produce one WMS layer for the basemap, so we only need one MapScript file. Despite supplying WMS using a relatively unusual map projection, the MapScript file is actually quite conventional.

In order to set the MapScript file to supply map tiles using the Mollweide projection, we need to know the following information:

  • Projection Type: Mollweide (“moll”) using a prime meridian of 0 degrees
  • Geoid: WGS-84  (standard for GPS and most online web maps)
  • EPSG Code:  54009 (Mollweide, Prime Meridian 0, WGS-84)
  • Map Extents: -18040096, -2267229, 18040096, 2267229 (metres W,S,N,E)

 

MapServer uses the popular Proj.4 open source library for map projection support. You should look at the Proj.4 documentation in order to determine the projection parameters. Proj.4 uses the name “moll” to indicate the Mollweide projection. The prime meridian is simply the central meridian. This could be changed, although experience has shown that MapServer can have wrap-around problems if vector data wraps around between the far west and east extremities.

The EPSG code has to be looked up. You may need to create your own dummy code if your specific set of projection parameters does not have an EPSG code. If you do this, you will need to modify the EPSG tables used by both Proj.4 and any other referencing software (eg. Proj4JS if you are using OpenLayers).

The projection results in coordinates that use metres. Ideally, the map extents are manually calculated using the projection formulae. If the formulae are not available, you will have to derive them by trial and error.

When you have determined all of the above information, you can write the MapScript header. Here is the header for our Mollweide example:

# Mollweide Projection WMS Base Map
# Based on the Mollweide MapScript file for Equal-Area-Maps.com
NAME "outline"
UNITS meters
EXTENT -18040096 -2267229 18040096 2267229
SIZE 640 480
IMAGECOLOR 255 255 255
IMAGETYPE png
SHAPEPATH "/path to your shape files/"
FONTSET "/font path/fontset.txt"
CONFIG "PROJ_LIB" "/path to your Proj4 proj directory/"
DEBUG on
STATUS ON

##################################### 
# Web object
#
WEB
     IMAGEPATH "/your_web_dir/map/tmp/"
     IMAGEURL "/map/tmp/"
     METADATA
       "wms_title"     "WMS Equal Area Map Server"
       "wms_onlineresource" "http://www.yourdomain.com/mapserv.cgi?map=/path_to_this_file/proj_mollweide.map"
       "wms_srs"       "EPSG:54009"
     END
     LOG "/log_directory/mapserv.log"
END

#####################################
# Projection Definition
#
PROJECTION
  "proj=moll"
  "lon_0=0"
  "units=m"
END

To use this example, you will need to replace the pathnames with their correct values. Many of them are dummy values required by MapScript, but not required for our configuration because we are creating a WMS server without any text output.

Notice that at the top of the file, the units are set to “meters” and the extents are given in Mollweide-projected metres.

The WEB object is a standard WMS web object definition, although the wms_srs setting is set to EPSG:54009 – our chosen Mollweide projection system. This tells MapServer to report the map tiles as conforming to EPSG:54009. This is important if our client is going to interpret them correctly.

The PROJECTION object might be unfamiliar. This simply provides the Proj.4 parameters for our chosen projection. Due to the ubiquitous nature of the open source Proj.4 library, Proj.4 projection definitions are a de facto standard, and are fully supported by MapServer. This particular definition is a simple one and tells MapServer to use the “moll” (Mollweide) projection method with a central meridian of zero degrees (ie. centred on the Prime Meridian), and to use metres as our units.

That is all that is required to produce a WMS server using a non-conventional map projection.

 

Adding the Map Data

To make a useful basemap, we need to actually plot some map data. Our sample basemap includes coastlines, national boundaries, and a graticule. These are all plotted from ESRI ‘shp’ shapefiles.

The national boundaries were created from the Mapping Hacks “World Borders” dataset. Most coastlines are also national boundaries, so simply drawing shapes for each national boundary creates most of our coastlines. However, this data misses a number of inland seas and large lakes. Coastline data for the Caspian Sea was downloaded from the Caspian Environment Programme, and the Great Lakes were downloaded from NOAA’s Coastal Geospatial Data Project. An attempt was made to find public domain data for other large seas and lakes (eg. Lake Baikal, Great Bear & Great Slave Lakes), but I could not find any.

The ShapeLib tools were used to combine the shape files as required. The data was generally of a much higher resolution than required for a global map, so the online MapShaper tool was used to reduce the resolution of all shape files. The Great Lakes and Caspian Sea datasets were especially high resolution. Although this process reduces disk space, the real gain is the significant improvement in rendering speed. Two versions of all the boundary and coastline data were created: A low resolution version for zoomed-out views, and a higher resolution version for zoomed-in views. This two tier approach improves rendering speeds when zoomed out, without sacrificing detail when zoomed in.

The graticule is implemented in the form of a polygon shapefile. Our map display has two kinds of background. First there is the background to the globe – we want this to be white. Secondly, there is the background to the land areas (ie. the sea) – we want this to be light blue. Therefore we set the MapScript background to white, and cover the entire globe with a series of polygon tiles to produce the blue sea. We are already dividing the Earth into regular shapes with the graticule, so we use one shapefile to solve both problems. This shapefile consists of a series of ‘squares’ that tile the entire Earth. To produce the sea, the tiles are drawn as blue fill. To produce the graticule, we draw the tiles with only a thin black outline.

So how do we produce this graticule? The regular nature of the graticule means we can create it using a script. ShapeLib includes utilities that can create shape files and add shapes to them. These can be automated with a Python script which uses shpcreate and dbfcreate to create the shape file and corresponding dbf file. The individual graticule polygons can then be added with shpadd and dbfadd . The DBF file is a dummy file that simply gives each shape a name. MapServer expects the DBF file to exist, even though it is not actually used in this application.

We have now created our data files. Here is the completed MapScript file:

# Mollweide Projection WMS Base Map
# Based on the Mollweide MapScript file for Equal-Area-Maps.com
NAME "outline"
UNITS meters
EXTENT -18040096 -2267229 18040096 2267229
SIZE 640 480
IMAGECOLOR 255 255 255
IMAGETYPE png
SHAPEPATH "/path to your shape files/"
FONTSET "/font path/fontset.txt"
CONFIG "PROJ_LIB" "/path to your Proj4 proj directory/"
DEBUG on
STATUS ON

##################################### 
# Web object
#
WEB
     IMAGEPATH "/your_web_dir/map/tmp/"
     IMAGEURL "/map/tmp/"
     METADATA
       "wms_title"     "WMS Equal Area Map Server"
       "wms_onlineresource" "http://www.yourdomain.com/mapserv.cgi?map=/path_to_this_file/proj_mollweide.map"
       "wms_srs"       "EPSG:54009"
     END
     LOG "/log_directory/mapserv.log"
END

#####################################
# Projection Definition
#
PROJECTION
  "proj=moll"
  "lon_0=0"
  "units=m"
END


##################################### 
# Ocean solid fill
#
LAYER
     NAME "outline"
     METADATA
       "wms_title"   "outline"
     END
     DATA "grid30"
     PROJECTION
        "init=epsg:4326"
     END
     STATUS default
     TYPE polygon
     CLASS
          STYLE
               COLOR 192 192 255
          END
     END
END


##################################### 
# Coloured Fill: Land, Zoomed out
#
LAYER
     NAME "outline"
     METADATA
       "wms_title"   "outline"
     END
     DATA "world_borders_simple"
     PROJECTION
        "init=epsg:4326"
     END
     MINSCALE 20000001
     STATUS default
     TYPE polygon
     CLASS
          STYLE
               COLOR 241 238 232
          END
     END
END

##################################### 
# Coloured Fill: Land, Zoomed in
#
LAYER
     NAME "outline"
     METADATA
       "wms_title"   "outline"
     END
     DATA "world_borders"
     PROJECTION
        "init=epsg:4326"
     END
     MAXSCALE 20000000
     STATUS default
     TYPE polygon
     CLASS
          STYLE
               COLOR 241 238 232
          END
     END
END

##################################### 
# Coloured Fill: Lakes and Inland Seas
#
LAYER
     NAME "outline"
     METADATA
       "wms_title"   "outline"
     END
     DATA "lakes_simple"
     PROJECTION
        "init=epsg:4326"
     END
     STATUS default
     MINSCALE 20000001
     TYPE polygon
     CLASS
          STYLE
               COLOR 192 192 255
          END
     END
END
LAYER
     NAME "outline"
     METADATA
       "wms_title"   "outline"
     END
     DATA "lakes"
     PROJECTION
        "init=epsg:4326"
     END
     STATUS default
     MAXSCALE 20000000
     TYPE polygon
     CLASS
          STYLE
               COLOR 192 192 255
          END
     END
END


##################################### 
# Shape boundaries to be drawn (black)
#
LAYER
     NAME "outline"
     METADATA
       "wms_title"   "outline"
     END
     DATA "world_borders_simple"
     PROJECTION
        "init=epsg:4326"
     END
     STATUS default
     MINSCALE 20000001
     TYPE line
     CLASS
          STYLE
               SIZE 1
               COLOR 0 0 0
          END
     END
END
LAYER
     NAME "outline"
     METADATA
       "wms_title"   "outline"
     END
     DATA "world_borders"
     PROJECTION
        "init=epsg:4326"
     END
     STATUS default
     MAXSCALE 20000000
     TYPE line
     CLASS
          STYLE
               SIZE 1
               COLOR 0 0 0
          END
     END
END


##################################### 
# Graticule
#
LAYER
     NAME "outline"
     METADATA
       "wms_title"   "outline"
     END
     DATA "grid30"
     PROJECTION
        "init=epsg:4326"
     END
     STATUS default
     TYPE line
     CLASS
          STYLE
               COLOR 95 95 95
          END
     END
END


END # mapfile

You will see that there are a number of data layers. Some shapefiles appear more than once, according to what exactly is being plotted. The first layer is drawn first (ie. the bottom), and the last layer is drawn last (ie. on top of everything else). Here is the drawing order:

  • Background: White (IMAGECOLOR in the header)
  • Ocean solid fill: Graticule (grid30.shp) as a blue solid fill
  • Land solid fill (Zoomed Out): Simplified world borders (world_borders_simple.shp) as off-white fill
  • Land solid fill (Zoomed In): World borders (world_borders.shp) as off-white solid fill
  • Lakes solid fill (Zoomed Out): Simplified Lakes (lakes_simple.shp) as blue solid fill
  • Lakes solid fill (Zoomed In): Lakes (lakes.shp) as blue solid fill
  • Country outlines (Zoomed Out): Simplified world borders (world_borders_simple.shp) as thin black outline
  • Country outlines (Zoomed In): World borders (world_borders.shp) as thin black outline
  • Graticule Lines: Graticule (grid30.shp) as thin dark grey outlines

 

And there we have it. In the second part of this article we will implement the OpenLayers client used to view this basemap. We will also use OpenLayers to add an overlay layer that is displayed using the correct, matching projection.

Related Posts