[gdal-dev] Problem with gdalsrsinfo on some NOAA bsb/kap files

[Jaak Laineste]

Hi Chris,

On 15.11.2012, at 1:15, Chris Barker wrote:

>> I just created index bounding boxes for
>> all the NOAA BSB charts using GDAL
> 
> Exacltly how did you use GDAL for this?

I use GDAL as library. So I have a code which scans all the files and requests bounding boxes.

>> for my Android solution, see
>> http://imgur.com/MpADB - as you can see there are quite many wrong
>> locations, Xmin-coordinate is 0.
> 
> you mean the huge pile of boxes on the prime meridian?

Yes, exactly. It is not really 0, but something close it as the x is not adjusted to real central meridian.

>> which is interesting also. I do not have the number, but there
>> seems to be significant number (1/4?) of the NOAA files positioned wrong.
> 
> you might look at a few and see if you can see a pattern.
> 
> GDAL asside, there is a lot of info in the *.kap text header -- I
> think the first few GRPs are the bounds of the "actual" map (i.e.
> without the "collar) for instance. You may be abel to get what youwant
> by directly parsing teh header.
> 
> Not that GDAL shouldn't report the correct values.

Right, the KAP has also metadata. My code is actually public (http://goo.gl/PKkmr - see function corner() ).

I have two files from same chart, which behave differently. They are both rotated, maybe this is important somehow:
gdalinfo 11315_1.KAP -> reports OK
gdalinfo 11315_2.KAP -> reports GCP projection, but not general one. Looking to bsbdataset.cpp code I do not see how this can possibly happen. Probably this is root for next problems. The Y coordinates will be correct, but Y is not adjusted to the prime meridian, as if Yorigin in the transformation matrix would be 0. 

	Coordinate System is `'
	Origin = (0.0,0.0)
	Pixel Size = (1.0,1.0)
	GCP Projection = PROJCS["Global Mercator",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 
	...


What is interesting that after I created VRT files for conversion:
gdalwarp -t_srs "EPSG:900913" -of VRT -dstnodata 0 11315_2.KAP 11315_2.vrt

And then query this:
gdalinfo 11315_2.vrt

Then I got correct and functional projection info. VRT files themselves for 11315_1 and 11315_2 are quite different: first one is based on coordinate system and transformation, second one has also copy of all GCP-s. Both files are rotated, the one which works is clockwise (SK=35.7611111), the other is rotated counter-clockwise (SK=294.1922222 in KAP terms).

The first test file KAP is also shown fine in QGIS (and my client), and second KAP is in wrong location. I have not checked QGIS source, but I suspect it uses same approach for reading as my app: creates automatically virtual VRT, so it has same issues with same files. It has also transparency issue : no data is black instead of transparent. When I open VRT file of the problematic data in QGIS, then it is also shown in correct location, so gdalwarp "fixes" it somehow.

So - what gdalwarp could do better than plain gdalinfo and AutoCreateWarpedVRT() ?

Jaak