Using the ImageMosaic plugin

Introduction

This tutorial describes the process of creating a new coverage using the new ImageMosaic plugin. The ImageMosaic plugin is authored by Simone Giannecchini, and allows the creation of a mosaic from a number of georeferenced rasters. The plugin can be used with Geotiffs, as well as rasters accompanied by a world file(.pgw for png files, .jgw for jpg files, etc.).

The JAI documentation gives a good description about what a Mosaic does:

The "Mosaic" operation creates a mosaic of two or more source images. This operation could be used for example to assemble a set of overlapping geospatially rectified images into a contiguous image. It could also be used to create a montage of photographs such as a panorama.

For this example, I am using 3 GeoTIFF files in Albers projection. I am attaching them here, in case someone wants to try them too.

Part One: Preliminary work

In order to configure a new CoverageStore and a new Coverage with this plugin, some files need to be generated first. More specifically, the following files are needed:

A shapefile that contains enclosing polygons for each raster file. This shapefile needs to have a field named "location", where values are the relative path for the mosaic rasters.
A projection file (.prj) for the above-mentioned shapefile.
A properties file (.properties). This file contains properties such as cell size in x and y direction, the number of rasters for the ImageMosaic coverage, etc.

At present, you can use one of the following ways for creating the index shapefile:

Using MosaicIndexBuilder
This tool is authored by Simone Giannecchini, and can be found in the Geotools Subversion repository http://svn.geotools.org/trunk/modules/unsupported/coveragetools
Advantages: All-in-one solution--no need to manually edit and create .properties or a .prj file.
Disadvantages: Requires you to download and build the tool. However, this is easily accomplished if you have Eclipse installed, and a SVN checkout of gt-trunk.
Using GDAL
GDAL (Geospatial Data Abstraction Library) comes in handy if you already have it installed.
Advantages: No need to compile the MosaicIndexBuilder Java tool, if you already have gdal installed.
Disadvantages: Requires manual creation and editing of .properties and .prj files.

Using MosaicIndexBuilder

You can run the MosaicIndexBuilder (from the command line, Eclipse, or whatever...) with the following arguments:

\-s geoserver-1.4.x_merge_wcs/web/target/geoserver/data/coverages/modis_mosaic \-w \*.tif \-name test

where the -s switch is obligatory, and provides the path to all rasters, and the -w switch sets the wildcard pattern for all raster files for the ImageMosaic coverage.
Note that you can also change the name of the generated shapefile by using the "-name" switch.

The projection (.prj) file

If you have non-Geotiff rasters in a projection other than WGS 84, it is a good idea to have an accompanying .prj file for each raster. Otherwise, the MosaicIndexBuilder will create a .prj file that assumes WGS 84.

Here is the full help output for the MosaicIndexBuilder, you get the latest version if you execute the program without any arguments:

MosaicIndexBuilder -h -v -s <source> -w <wildcardOpt> -name <name> -p <priority> -s (--source_directory) source  path where files are located

Using GDAL

Within the directory of the rasters, run the command:
```
gdaltindex -tileindex location test *.tif
```
Here gdaltindex creates a test.shp with a field named "location" (the exact name is needed by the Imagemosaic plugin), that has the relative path for all images with .tif extension in the current directory. The same procedure will apply for rasters accompanied by a world file (png, jpg, etc).
Construct a test.properties file. The .properties file is required to have same filename as the index shapefile (test.properties) with the contents:
```
Name=test
Levels=10000.0,10000.0
LevelsNum=1
Envelope2D=-4222839.0,-1666784.0 1575818.0,3914576.0
NumFiles=3
```
Parameters, see below for details
- Name the name to be associated with the mosaic
- Levels depicts the cell size in the X an d Y directions
- LevelsNum the total number of levels
- Envelope2D contains the lower left and upper right corners' coordinates (separated by a white space)
- Numfiles is the number of files in the Coverage
  The values for "Levels" can be gathered by inspecting the output from gdalinfo (HINT: read the absolute values from the "Pixel Size" line):
```
gdalinfo passA2006260173532.tif
```
  where the line of interest reads:
```
Pixel Size = (10000.000000000000000,-10000.000000000000000)
```
  The values for "Envelope2D" can be gathered by inspecting the output from ogrinfo (HINT: read the "Extent" line all the way to the top):
```
ogrinfo  -al test.shp
```
  Note that the values for the "Levels" parameter should be positive, since our goal is to describe the cell size in SRS units (meters, feet, degrees) in the X and Y directions. If you have overviews the Levels parameter should contain a space separated list of the pixel sizes for the overviews in use. For example, if all your tiles have the same overviews generated using gdaladdo 2 4 8 16 32 64 on a coverage whose native resolution is 0.2m per pixel, the property file should contain:
```
Levels=0.2,0.2 0.4,0.4 0.8,0.8 1.6,1.6 3.2,3.2 6.4,6.4 12.8,12.8
LevelsNum=7
```

Create a test.prj file and insert the WKT text for your projection. For example:

PROJCS["Albers Equal area",GEOGCS["WGS 84",
	DATUM["World Geodetic System 1984",
		SPHEROID["WGS 84", 6378137.0, 298.257223563,
		AUTHORITY["EPSG","7030"]],
	AUTHORITY["EPSG","6326"]],
PRIMEM["Greenwich", 0.0, AUTHORITY["EPSG","8901"]],
UNIT["degree", 0.017453292519943295],
AXIS["Geodetic latitude", NORTH], AXIS["Geodetic longitude", EAST], AUTHORITY["EPSG","4326"]],
PROJECTION["Albers Equal Area", AUTHORITY["EPSG","9822"]],
PARAMETER["central_meridian", -96.0],
PARAMETER["latitude_of_origin", 37.5],
PARAMETER["standard_parallel_1", 29.833333333333336],
PARAMETER["false_easting", 0.0],
PARAMETER["false_northing", 0.0],
PARAMETER["standard_parallel_2", 45.833333333333336],
UNIT["m", 1.0], AXIS["Easting", EAST], AXIS["Northing", NORTH], AUTHORITY["EPSG","41111"]]

This is just an Albers projection with custom parameters for my rasters, you may be able to extract it from the GeoTIFFs themselves by using gdalinfo. If you need to define a custom projection, see [Custom projection definition in Geoserver].

Regardless of the method you use for creating the index shapefile, it is a good idea to test it, especially if the ImageMosaic plugin is new to you:

How to test the resulting shapefile:

Configure it as a featureType and display a map with the correct projection, maybe with another layer in the same projection
and make sure that layers line up.

Part Two: Configuring a Coverage in Geoserver

This is a process very similar to creating a FeatureType. More specifically, one has to:

Create a new CoverageStore:

Go to "Welcome | Config | Data | CoverageStores | New" via the web interface, choose "ImageMosaicking plugin" for "Coverage Store" and give it a name:
Click "New", type in the path to the shapefile, then click the Submit button:
Click on the Apply, Save, and Load buttons.

Create a new Coverage using the new ImageMosaic CoverageStore:

Go to "Welcome | Config | Data | Coverage | New" via the web interface, choose the name of the coverage you just created:
Click on the New button and you will be presented with the Coverage Editor:
Make sure there is a value for "Native CRS", then click the Submit button. Hopefully there are no errors
Click on the Apply, Save, and Load buttons.

Results

That is it!!! Now let's see the new ImageMosaic Coverage in action:

Another map, this time the generated index shapefile is superimposed over the rasters: