[GRASS5] r.in.gdal and negative values

[Markus Neteler]

On Tue, May 31, 2005 at 05:06:26PM +1200, Hamish wrote:
> re. http://grass.itc.it/pipermail/grass5/2005-January/017029.html
> 
> 
> Currently r.in.gdal imports unprojected data in the negative, i.e. 0,0
> is the top right of the image.
> 
> (unlike r.in.tiff|png in GRASS 5.4 which makes 0,0 bottom left)
> 
> This can be "fixed" with r.region by the user after import.
> 
> 
> some comments from the source code:
> 
> /* TODO: most unreferenced formats are read in with negative 
>  * coordinates - desired??
>  */
> 
> [...]
> 
>         /* use negative XY coordinates per default for unprojected data
>            but not for all formats... (MN: I don't like it at all) */

... so you can see that we are two people not being in favor of
negative xy coordinates...

>         /* for hDriver names see gdal/frmts/gdalallregister.cpp */
> 
>         if ( l1bdriver || 
>              ( strcmp(GDALGetDriverShortName(hDriver),"GTiff") || 
>                strcmp(GDALGetDriverShortName(hDriver),"JPEG") ) == 0 )
>         {
>           /* e.g. L1B - NOAA/AVHRR data must be treated differently */
>           /* define positive xy coordinate system to avoid GCPs confusion */
>           cellhd.north  = cellhd.rows;
>           cellhd.south  = 0.0;
> 
> [...]
>         else
>         {
>           /* for all other unprojected data ... */
>           /* define negative xy coordinate system to avoid GCPs confusion */
>           cellhd.north  = 0.0;
>           cellhd.south  = (-1) * cellhd.rows;
> 
> 
> 
> I don't understand the logic or intention here: who's supposed to be 
> getting set to negative and why?

As written: L1B data have a strange GCP organization. We should
probably change the test to refuse L1B data completely as
i.rectify is not suitable to recitfy AVHRR data.

The new gdalwarp/thin plate splines warper performs better (so
a new error message could suggest exactly this).

> 
> FWIW the above strcmp() part of the test will always fail.

Mhh, it worked at least on a couple of platforms.
Since I wrote it it could be certainly programmed in a better way.

> i.e. strcmp() returns 0 if the strings match. It will never match both a 
> GeoTIFF and a JPEG, so one side or the other of the "OR" will be always 
> be non-zero and the ==0 part will always be false. So it really reduces 
> to if(l1bdriver) {make positive}  else {make negative}.

This was probably intended (Frank wanted always negative coordinates which
didn't work for AVHRR). So I implemented the exception for L1B.
 
> (l1bdriver is AVHRR; http://www.gdal.org/frmt_l1b.html
>   Will this data ever be without GCPs & so touch this part of the code?!)

It's always with GCPs which are written into a POINTS file during
import.

> [And PNG should be doing the same thing as TIFF & JPEG.]

Yes.

> So what is the intent of the test?
> 
> Are GeoTIFFs and JPEGs supposed to be ending up positive or negative?
>   (positive please)

Frank wanted negative (I don't know why - maybe there is no reason).
I suggest to change to positive values for all and to refuse L1B or,
at least, issue a warning.
 
> Who's getting what confused with GCPs? Does that situation still exist?
> 
> 
> I prefer to make everything positive as it's less confusing.

Agreed.

> gdalinfo reports this on an unprojected TIFF:
> Driver: GTiff/GeoTIFF
> Size is 1309, 4920
> ...
> Corner Coordinates:
> Upper Left  (    0.0,    0.0)
> Lower Left  (    0.0, 4920.0)
> Upper Right ( 1309.0,    0.0)
> Lower Right ( 1309.0, 4920.0)
> Center      (  654.5, 2460.0)
> 
> 
> So in GDAL 0,0 is top left (image coordinates (or satellite scanlines)
> as opposed to cartesian coordinates)

Note that GDAL uses *corners* of pixels for coordinates.

> GRASS XY uses cartesian coordinates so 0,0 should be bottom left if we
> want to keep map units positive. (external support for idea: this is the 
> reason for having false eastings and northings)

Well, false eastings and northings do not appear in XY locations AFAIK.

> So is this the "GCP confusion" the code speaks of? Solution is to make 
> it neither? :/  As we are GRASS, matching GRASS XY conventions seems 
> like logical thing to do.

The story is:
- Frank wrote the module with negative values in the XY case
- I added special treatment for AVHRR to make it work. Then I
  discovered that polynomials perform bad for AVHRR. Finally
  I organized partial fundings for the thin spline warper in GDAL
  which works way better.


> I would like to use i.point's $GROUP/POINTS file as a template for 
> setting GCPs in a GeoTIFF directly with gdal_translate or geotifcp,
> but it seems like I'll have to fuddle the i.points results around
> slightly to do so? (flip y-axis) If everything is positive it is 
> easy enough to write a script to form a "gdal_translate -gcp"
> command string from i.points output.
> 
> 
> It seems it is trivial to make everything positive in the r.in.gdal 
> code; can anyone justify a reason not to? 
> 
> Aka any reason I should not just go ahead and do this in 6.1-cvs?

Please go ahead...!

Markus