[gdal-dev] How do I select 3 of 5 bands for output?

[Stephen Woodbridge]

Ahh, thank you for the link to the vrttut. This has got me going in the 
right direction. I actually did try the -of VRT but didn't know it was 
supposed to be a XML file so discounted it when the file was small.

So, my next problem is that when I run gdalwarp on the VRT file I get 
assertion errors like the following on about half of the tiles I try to 
create.

$ cat a
gdalwarp -srcnodata 0 -dstnodata 0 -s_srs +init=epsg:26916 -t_srs 
EPSG:4326 -te -86.25 38.3125 -86.1875 38.375 -rb -wm 250 --config 
GDAL_ONE_BIG_READ ON -co "TILED=YES" 
/var/data/tiler-tiffs/11234-1226164462/src/tmp.vrt a.tif

$ sh a
Creating output file that is 6150P x 6150L.
Processing input file /var/data/tiler-tiffs/11234-1226164462/src/tmp.vrt.
gdalwarp: ../compressed/compressed.cpp:1724: bool 
kd_tile::read_tile_part_header(): Assertion `tpart_body_length >= 0' failed.
a: line 1: 11502 Aborted                 gdalwarp -srcnodata 0 
-dstnodata 0 -s_srs +init=epsg:26916 -t_srs EPSG:4326 -te -86.25 38.3125 
-86.1875 38.375 -rb -wm 250 --config GDAL_ONE_BIG_READ ON -co 
"TILED=YES" /var/data/tiler-tiffs/11234-1226164462/src/tmp.vrt a.tif

Is this because the tile extends outside the source image?

+-----------------+
|  src_image      |
|            +----------+
|            |  tile    |
|            |          |
|            +----------+
|                 |
+-----------------+

Or something similar to this? Any thoughts on the cause(s) and workarounds?

Thanks again for all the help.

-Steve

Lucena, Ivan wrote:
> Hi Steve,
> 
> No, VRT files are not Geotiff, see [1], and they should not take more than a bunch of Kbytes. The problem is that in *our* gdal_translate command line *we* forgot "-of vrt" and the default format is Geotiff.
> 
>    system("gdaltranslate -OF VRT $infile -b 1 -b 2 -b 3 tmp.vrt");
> 
> Good luck,
> 
> Ivan
> 
> [1]: http://www.gdal.org/gdal_vrttut.html
> 
>>  -------Original Message-------
>>  From: Stephen Woodbridge <woodbri at swoodbridge.com>
>>  Subject: Re: [gdal-dev] How do I select 3 of 5 bands for output?
>>  Sent: Nov 08 '08 03:36
>>  
>>  Stephen Woodbridge wrote:
>>  > Lucena, Ivan wrote:
>>  >> Hi Steve,
>>  >>
>>  >>>  I have a script that was processing MrSID files and need to convert it
>>  >>>  to support jp2 files also. It seems pretty straight forward except the
>>  >>>  jp2 files have bands R,G,B,IR and I need to chop these into GeoTif
>>  >>> tiles
>>  >>>  with only R,G,B bands.
>>  >>>  
>>  >>>  The script was iterating over the src image and chopping it into tiles
>>  >>>  using the following command:
>>  >>>  
>>  >>>  system("gdalwarp -srcnodata 0 -dstnodata 0 -s_srs +init=epsg:$epsg
>>  >>>  -t_srs EPSG:4326 -te $xm $ym $xx $yx -rb -wm 250 --config
>>  >>>  GDAL_ONE_BIG_READ ON -co "TILED=YES" $infile $work/$outd/$outf\n");
>>  >>>  
>>  >>>  1) How does this behave with R,G,B,IR as input and how can I get R,G,B
>>  >>>  only into the GeoTif?
>>  >>
>>  >> I would add one call to gdal_translate with "-b 1 -b 2 -b 3" first to
>>  >> build a VRT and then call gdalwarp from the VRT.
>>  >>
>>  >> {
>>  >>   system("gdaltranslate $infile -b 1 -b 2 -b 3 tmp.vrt");
>>  >>   system("gdalwarp -srcnodata 0 -dstnodata 0 -s_srs +init=epsg:$epsg
>>  >> -t_srs EPSG:4326 -te $xm $ym $xx $yx -rb -
>>  >> wm 250 --config GDAL_ONE_BIG_READ ON -co "TILED=YES" tmp.vrt
>>  >> $work/$outd/$outf\n");
>>  >> }
>>  >> system("gdalmanage delete tmp.vrt");
>>  >
>>  > OK, this might work, but it is incredibly slow. Then again probably
>>  > doing much of anything with these images will be slow because they are
>>  > HUGE. So looking for ways to speed this up. It looks like the .vrt file
>>  > is really a GeoTif would it help later if the .vrt was made with -co
>>  > "TILED=YES", can you do this with a .vrt file. I think the uncompressed
>>  > 3 band geoTif will be in the ballpark of 10GB, is this being made as a
>>  > BigTiff? or will I run into problems when I hit the 2GB limit?
>>  >
>>  > Anyway to get this to go faster? CPU is at 101%, 78MB of memory is being
>>  > used, I have LOTS more if I can get gdal to use it.
>>  
>>  OK, it took 2 hr 10 min to run:
>>  
>>  GDAL_SKIP=JPEG2000 gdal_translate
>>  /var/data/raw_data/minnesota/ortho_1-1_1m_j_mn161_2008_1.jp2  -b 1 -b 2
>>  -b 3 tmp.vrt
>>  
>>  and the output file is 5.0 GB. And on top of this I still need to
>>  gdalwarp this VRT file a bunch of times to chop it into smaller tiles. I
>>  have hundreds of files to process. Any thoughts on how to speed this up?
>>  
>>  -Steve
>>  
>>  _______________________________________________
>>  gdal-dev mailing list
>>  gdal-dev at lists.osgeo.org
>>  http://lists.osgeo.org/mailman/listinfo/gdal-dev
>>