[Gdal-dev] warp and translate problems

[Frank Warmerdam]

Chris Hodgson wrote:
> I'm working with TIF orthophoto files, warping them from UTM to BC 
> Albers projection, and converting them to JPEG-2000 format.
> I have a few questions:
> 
> 1) it seems that it is not possible to do this in one step using 
> gdalwarp and specifying the output format -  I'm guessing the same is 
> true for gdal_translate and specifying the  output projection?  I 
> suppose this is just a limitation of the more complicated output format 
> drivers?

Chris,

Generally speaking programs like gdalwarp can't output to formats which
only support making a copy of an existing dataset.  That is, many drivers
only have a CreateCopy() method for making new files, and you need to
pass them an existing dataset to clone.  The JP2KAK driver falls into this
category.

You can actually use gdal_translate to subset to something like jpeg2000,
but what happens inside gdal_translate is that it makes up a "virtual dataset"
with the desired characteristics that fetches image data from the real
source file as needed.

Recently the virtual file format was also extended to support warping on
the fly, but gdalwarp does not (yet) utilize this transparently to write
to formats like JPEG2000.  If it is important to avoid making the intermediate
file you could use gdalwarp to write a VRT file and then gdal_translate that
to JPEG2000 format.  The .vrt file will just be an XML description of the
transformation to apply.

> 2) I'm starting with .tif and .tfw files, then gdalwarp gives me single 
> geotiff file. then when I use gdal_translate to convert to jp2 format, I 
> don't have any georeferencing (no wld file)... do I need to compile with 
> jasper as well so that I can output geojp2 files, or is there some way 
> to tell gdal_translate to output a .wld file for the output file?

The JP2KAK driver writes out GeoJP2 (tm) files as long as the input is
georeferenced.  If you really need .wld files for some reason you could
write a little python script to do it.  I don't believe the JP2KAK driver
supports generating .wld files directly.

> 3) My intent is to build a tileindex of the output files and serve them 
> up using mapserver. I'm worried about the black "collar" that I'm left 
> with around the images after I've reprojected them - will mapserver be 
> able to deal with the fact that these all overlap with data from the 
> neighboring tile? Can mapserver deal with tileindexes that have 
> overlaps? How about tiles that aren't square in the given coordinate 
> system - if I use the grid from the original orthophoto files, will it 
> use those boundaries to select which tile to get the data from? That 
> would be the optimal solution, if it works that way. Tell me if this is 
> a question better posted to the MapServer list but I know it is gdal 
> doing most of this work anyway.

Hmm, this will indeed be a problem.  MapServer doesn't mind overlapping
images, but it is just going to draw them in some arbitrary order and the
black stuff will get written over valid imagery in some cases.  Normally
I would suggest using the MapServer OFFSITE keyword to indicate that the
black should be ignored (treated as transparent), but with a lossy
compression format like JPEG2000 you will find the black pixels vary a
bit from 0/0/0 and thus you won't be able to use OFFSET effectively.

My suggestion would be to generate output tiles that are not related 1:1
to your input files.  Instead, pick a tile size and then fill it completely
by gdalwarping multiple images into it.

My understanding is you have 2TB of imagery.  If you broke your region
down into 100000 x 100000 chunks (roughly 30GB each assuming RGB uncompressed)
then you could just compress each of those chunks.  To warp the input
images into the you would need to use gdalwarp to exactly create the
correct output tiles with the -tr and -ts switches to gdalwarp.  Then just
warp each file that intersects the region into the working file (or all
the files if it is to much trouble to work out which ones intersect).
This should give you (big) resulting files with no black areas with the
possible exception of areas that have no data.

Another possibility would be to work with me to add alpha band support
to gdalwarp.  Then you could warp each file to an RGBA file with a=0
for pixels that are not set, and a=255 for pixels that are set.
MapServer knows how to utilize per-pixel alpha information when available.

Best regards,
-- 
---------------------------------------+--------------------------------------
I set the clouds in motion - turn up   | Frank Warmerdam, warmerdam at pobox.com
light and sound - activate the windows | http://pobox.com/~warmerdam
and watch the world go round - Rush    | Geospatial Programmer for Rent