Serving Google Maps and Google Earth Data with CDAT

Here are some notes about using CDAT to act as a Google Maps / Google Earth server, serving pictures of climate data (esp. gridded numerical model output). This is "work in progress". Please contact JonBlower for more information.

Summary

Essentially, what is required is a Web server that can parse a URL that is generated by the Google Maps/Earth interface and return a PNG or GIF image of the data that have been requested.

I propose to use the Climate Data Analysis Tools (CDAT, a set of Python libraries) to extract the data from NetCDF files (and possibly do re-projections) and an Apache server with mod_python to serve the Web requests. I am not yet sure whether mod_python can use the CDAT libraries - I am using some of these things for the first time. If it doesn't, then I'm pretty sure I can call CDAT from a CGI script.

Installation

I am setting this up under Debian Linux (version Sarge). Actually, this is running as CoLinux under Windows XP, but this shouldn't make a difference.

Prerequisites

Before installing CDAT I installed a number of Linux packages (as root) with apt-get:

• apache2 (the Web Server)
• apache2-dev (for building mod_python for CDAT)
• build tools: make, gcc, g++ (and g77 if you want FORTRAN support, but this is not necessary)
• libx11-dev (X Windows development libraries, without which the build of Tk will fail)
• The C Shell (csh) - wasn't installed by default in my Linux

The 2.0.x series of the Apache server must be used rather than the new 2.2.x series. This is because the current version of mod_python is not compatible with Apache 2.2.x.

Note that the package names may be different in other Linux distributions. For example, in CentOS 4 the Apache Web server package names begin with "httpd-"

Installing CDAT

I downloaded CDAT version 4.0 from http://cdat.sf.net. I created a destination directory /usr/local/cdat-4.0 and gave myself (the jon user) write permissions in this directory. Then I extracted the CDAT distribution, changed to the root of the extracted files and ran (as the jon user):

./express_install /usr/local/cdat-4.0/ --cdms-only --disable-opendap

I suspect that the ioapi package will be useful (seems to do projections and includes the GDAL library) but I didn't install it (takes a long time to build). This requires Fortran (which is why I installed g77 as a pre-requisite). Note that the build process takes quite a while (especially under CoLinux) so this is a good time to have a cup of tea (or six).

Note that when I try building everything, the building of the netbpm component fails. The log files keep growing until they fill the disk with the messages:

Where is the install package you created with 'make package'?
package directory (/tmp/netpbm) ==> This does not appear to be a Netpbm install package.
A file named /tmp/netpbm/pkginfo does not exist.

and

Use of uninitialized value in string eq at ./buildtools/installnetpbm.pl line 32, <STDIN> line 9.
Use of uninitialized value in scalar chomp at ./buildtools/installnetpbm.pl line 31, <STDIN> line 9.

both of which are repeated in the log files ad infinitum. Maybe there's an issue with my version of Perl (5.8.4).

Installing mod_python

We can't use the mod_python that we would get simply by using apt-get install libapache2-mod-python because it is built for a different version of Python from that which is used in CDAT. We need to build it ourselves. This is what I did:

• Downloaded mod_python version 3.2.8 from http://www.modpython.org/
• Unzipped the tar archive and changed into the new directory
• Compiled mod_python as a shared object, using the version of Python that CDAT uses with:

./configure --with-python=/usr/local/cdat-4.0/bin/python --with-max-locks=32

• Ran make, then changed to root and ran make install.
• mod_python.so then ends up in /usr/lib/apache2/modules. In CentOS the path for mod_python.so is /etc/httpd/modules
• This also installs the mod_python libraries into CDAT (i.e. the Python installation that CDAT uses)

Configuring Apache to use mod_python

With the Debian installation of Apache2, the configuration of the web server is split up amongst lots of files. To make Apache2 load mod_python, I made sure that there was a file called /etc/apache2/mods-enabled/mod_python.load with the contents:

LoadModule python_module /usr/lib/apache2/modules/mod_python.so

In CoLinux this file is included from the main Apache configuration file, /etc/apache2/apache2/conf. In CentOS it is not necessary to have a separate mod_python.load file. Instead, the LoadModule lines are added directly to the main Apache configuration file, /etc/httpd/conf/httpd.conf. The correct line for CentOS is

LoadModule python_module modules/mod_python.so

After this I restarted the web server (/etc/init.d/apache2 restart, or service httpd restart in CentOS)

I created the directory /var/www/cdat and I edited /etc/apache2/sites-enabled/000-default (as root), adding the following lines in the <VirtualHost *> section:

<Directory /var/www/cdat>
  AddHandler mod_python .py
  PythonHandler mod_python.publisher
  PythonDebug On
</Directory>

This makes all requests for http://<myhost>/cdat/*.py get handled by the Publisher handler of mod_python. In CentOS the <VirtualHost *> sections are found in /etc/httpd/conf/httpd.conf. This file also defines the DocumentRoot path, which is /var/www/html by default. Therefore, the path for the cdat directory in httpd.conf in CentOS should be /var/www/html/cdat.

Note You have to make sure that the CDAT Python is the first Python in the PATH of the Apache user (test by running "which python")

I restarted Apache by running (as root) /etc/init.d/apache2 restart. In CentOS the command to restart Apache is service httpd restart.

Creating a test CDAT script

I created /var/www/cdat/cdattest.py with the following contents:

import cdms

def test_cdat(req, var='u'):
   req.content_type = "text/plain"

   f = cdms.open('/usr/local/cdat-4.0/sample_data/u_2000.nc')
   u = f(var)
   req.write('Units of ' + var + ': ' + u.units)
   f.close()

   return

This script loads a NetCDF data file from the CDAT samples, extracts one of its variables (specified by the var argument in the query string of the URL) and writes the units of this variable to the browser. Loaded the page http://<myhost>/cdat/cdattest.py/test_cdat in a browser. Got the string "Units of u: m/s", verifying that mod_python was working correctly with CDAT. The URL http://<myhost>/cdat/cdattest.py/test_cdat?var=u has the same result but http://<myhost>/cdat/cdattest.py/test_cdat?var=v raises an error because there is no variable called 'v' in the file.

If the above test CDAT script didn't work for you, check that mod_python is installed correctly. Maybe your mod_python is pointing to an installation of Python that doesn't contain the CDAT libraries, or maybe the mod_python module didn't get loaded into CDAT. See Installing mod_python above.

Test this by loading http://<myhost>/cdat/cdattest.py/make_pic in a browser.

Google Maps picture server

After a bit of work I managed to create a Python script that generates pictures for displaying in a Google Map. The full source is here.

Here is the top-level function, showing the main steps:

def genpic(req, x=1, y=1, zoom=15, scale_min=0, scale_max=50, opacity=100):

    # Get a Mercator grid of size 256 * 256 for the selected tile
    merc_grid = get_mercator_grid(int(x), int(y), int(zoom))
    # Get the latitude and longitude extents of the grid
    latExtent = (merc_grid.getLatitude().getValue()[0], merc_grid.getLatitude().getValue()[-1])
    lonExtent = (merc_grid.getLongitude().getValue()[0], merc_grid.getLongitude().getValue()[-1])

    # Open the dataset
    f = cdms.open('/var/www/cdat/glb.alllev.day10.000303.nc')
    # Get the temperature variable, first depth level only, and project onto the Mercator grid
    temp = f('temp', lon=lonExtent, lat=latExtent, lev=slice(0,1), squeeze=1, grid=merc_grid)
    f.close()
        
    # Calculate the scale factors to calculate the colour index from
    # the data value (index = m * value + c)
    m = 253 / (float(scale_max) - float(scale_min))
    c = 2 - m * float(scale_min)
       
    # Convert the data array into an array of colour indices. 
    # Valid pixels will have values between 0 and 255.
    temp *= m
    temp += c
        
    # Set out-of-range values to index 1
    temp = MV.choose(temp > 255.0 and temp < 2.0, (temp, 1))
        
    # Turn array into array of unsigned bytes.  This seems to turn missing
    # values into zeros automatically
    b = temp.astype(MV.UnsignedInt8)

    # Build the PNG and write it to the browser
    write_indexed_png(req, b.getValue().tolist(), get_palette(), get_alpha_channel(opacity))

Example: http://<myhost>/cdat/gmapspicserver.py/genpic?scale_min=-5&scale_max=30&x=2&y=2&opacity=66

Notes

The power of CDAT is evident here: having created a Mercator grid, we can extract the data from the source file and project it onto this new grid in a single step. This gives us a two-dimensional MV array (which is similar to a Numeric (numpy) array with extra features). In a couple of lines, we convert the values in this array into a set of colour indices between 0 and 255. Index 0 represents missing values and index 1 represents values that are out of the given scale range.

The speed of the code is at least comparable to that of the Java code that currently serves pictures to the Godiva2 site. In my test setup (running under Apache2 under Colinux on my Windows box), this code served an average of 2.12 pictures per second (tested using the Web Application Stress Tool). The Godiva2 site (which runs on different hardware and accesses different data!) serves 2.42 pictures per second, so the performance gap is not large, if it exists at all.

The PNG generation code is adapted from the minipng library by Dan Sandler: see http://dsandler.org/soft/python/minipng.py. Many thanks to Dan for this library.

To do

• The PNG pictures served by this script appear completely transparent in Internet Explorer 6. Is this because of the alpha channel?
• The MV.choose() operation (which sets the out-of-range values) is relatively slow: if this is commented out then the script can serve 3.15 pictures per second.

-- JonBlower - 11 May 2006