[Gdal-dev] modis sinusoidal tiles: bug in gdal?

[Vincent Schut]

Frank Warmerdam wrote:
> Vincent Schut wrote:
>> The puzzle:
>> Then I tried to merge some modis tiles by using VRT, using vrt_merge.py.
>> I was very puzzeled when I got an error message telling me it wouldn't
>> merge my tiles because they did not have the same scale. 
>
> Vincent,
>
> vrt_merge.py isn't part of GDAL - could you provide a pointer to the
> version you are using?
I know that. It is however a very valuable tool for me :-)
I think I picked it up from the ml somewhere sometime long ago... Not
easy to provide a pointer. Might have been this one:
http://www.vso.cape.com/~nhv/files/gdal/gdal_vrtmerge.py, but I'm not
sure. And I've adapted it to correct for an index error in the
colortable support somewhere, if I remember correctly. And for
ngbindings instead of ogbindings. I can send you a copy if you'd like.
>
>> After some
>> investigation, it indeed appeared they didn't! Pixel size (as calculated
>> by gdal) is different per tile. At around (0,0) degrees latlong, pixels
>> are square and of size 463.312716527916677 meters. The further one goes
>> from the equator or 0 deg longitude, the more the pixel size in that
>> direction differs from this 463.312716527916677 and pixels are not
>> square anymore.
>
> I see the pixels are just very slightly not-square.  The actually bounds
> in sinusoidal projection are apparently computed from the lat/long bounds
> for the image so they are sensitive to small differences due to floating
> point error, and issues in the sinusoidal equations.  I see the
> differences
> are roughly 12 digits into the pixel size and represent a difference in
> pixel size of roughly one billionth of a meter or (I think) roughly one
> millionth of a meter in the origin/corner values.
That is true. I see now that the metadata doesn't provide projected
corner coordinates, only pixel size
(CHARACTERISTICBINSIZE=463.312716527778)... Hmm, NASA, bad boy...
However, I have read somewhere that the sinusoidal projection would be
one of the most accurate for worldwide coverage, so I still am a bit
surprised at the differences. As I was used to work in unprojected
lat/lon spaces till now, I admit I have been a bit oversensitive to
nanometer scales... ;-)
OK, now to find out whether CHARACTERISTICBINSIZE is more accurate, or
the latlon corner coordinates are more accurate...
>
> I believe the issue is an undue sensitivity to slight differences in
> pixel size in vrt_merge.py.
I think so too, now. I think I'll create my own script to mosaic these
tiles, just forcing them right into the grid... :)
>
>> the gdalinfo's:
>>
>> tile h18 v09 (ul = 0,0, note the square pixels)
>> Origin = (0.000000000000000,-0.000000000000000)
>> Pixel Size = (463.312716527916677,-463.312716527916677)
>
>>
>> tile: h17 v01
>
>> Origin = (-1111950.519667000044137,8895604.157332999631763)
>> Pixel Size = (463.312716527916677,-463.312716527499731)
>
>>
>> tile h08 v06
>> Origin = (-11119505.196667000651360,3335851.558999999891967)
>> Pixel Size = (463.312716527917246,-463.312716527916677)
>
> Best regards,
Thanks!
VS.