This is probably getting a bit off-topic now, maybe better keep this converation going on python-gis-sig list? I'm cc'ing there....

2009/10/27 John [H2O] <washakie@gmail.com>
What I want to do is quite simple (I think). I simply would like to be able

1) Download any resolution of the hdf files available via :
http://rapidfire.sci.gsfc.nasa.gov/realtime/2009XXX/  where XXX is day of
Basically you have the MOD021K, MOD02HK and MOD02QK data granules, plus the MOD03 data granule. The former are the 1k, 500m and 250m data (hence the name) and the latter is a collection of geolocation terms (imaging geometries, etc).

You realise that these are swath data don't you? So they have no direct information on the location of each pixel, although you can work it out from the sensor height, attitude and so on. 

2) Import the hdf into python, access the data (reproject if necessary).
Yes, it's the merged RGB level 1 data I want to ultimately plot.

Importing into python is easy with GDAL. Your reprojection is a problem because the data aren't even projected. As I said above, they are swath observations that need to be gridded (Level 2 product).

2b) produce geotiffs for use elsewhere

Again, easy using GDAL (provided you don't care about using swath data, OR that you know how to grid it).

3) Plot the data in basemap, plot some other data (my own, xy tracks,
points, etc) over the data.

Again, you need to know the projection of the MODIS data for this to be of any use.

I don't have an off-hand way of gridding the Level 1 data, but someone may have. Other than the MS2GT <http://nsidc.org/data/modis/ms2gt/> (MODIS Swath to Grid toolbox), there's also the information given in <http://www.mcst.ssai.biz/mcstweb/L1B/product.html>, and I guess you could use that as a way of developing a python MS2GT ;D