[GRASS5] Re: [Gdal-dev] ERDAS/IMG driver error for NULL values

GRASS GIS GRASS GIS developer
,
1

Hello Otto

Otto Dassau wrote:

Also there is no coordinate system specified. This should be:

> name: Transverse Mercator
> datum: potsdam
> towgs84: 590.5,69.5,411.6,-0.796,-0.052,-3.601,8.30
> proj: tmerc
> ellps: bessel
> a: 6377397.1550000003
> es: 0.0066743722
> f: 299.1528128000
> lat_0: 0.0000000000
> lon_0: 9.0000000000
> k_0: 1.0000000000
> x_0: 3500000.0000000000
> y_0: 0.0000000000

This looks like a deficiency in the r.out.gdal script. In GRASS 5.7 now
the co-ordinate system can be printed in WKT format using `g.proj -wf`.
I assume gdal_translate has a command-line option to specify the output
co-ordinate system and this could be easily added to the corresponding
line in the r.out.gdal script.

Just a suggestion,

Paul.

2

    
Paul Kelly wrote:

This looks like a deficiency in the r.out.gdal script. In GRASS 5.7 now
the co-ordinate system can be printed in WKT format using `g.proj -wf`.
I assume gdal_translate has a command-line option to specify the output
co-ordinate system and this could be easily added to the corresponding
line in the r.out.gdal script.

To follow up on this, I think it would mean changing the last line in
r.out.gdal from
gdal_translate -of $FORMAT $CREATEKEY $METAKEY $CELLHD $OUTPUT
to
gdal_translate -of $FORMAT -a_srs '`g.proj -wf`' $CREATEKEY $METAKEY
$CELLHD $OUTPUT

Paul

3

On Wed, 30 Jun 2004 13:25:39 +0100
Paul Kelly <paul-grass@stjohnspoint.co.uk> wrote:

Paul Kelly wrote:
>
>
> This looks like a deficiency in the r.out.gdal script. In GRASS 5.7 now
> the co-ordinate system can be printed in WKT format using `g.proj -wf`.
> I assume gdal_translate has a command-line option to specify the output
> co-ordinate system and this could be easily added to the corresponding
> line in the r.out.gdal script.
>

To follow up on this, I think it would mean changing the last line in
r.out.gdal from
gdal_translate -of $FORMAT $CREATEKEY $METAKEY $CELLHD $OUTPUT
to
gdal_translate -of $FORMAT -a_srs '`g.proj -wf`' $CREATEKEY $METAKEY
$CELLHD $OUTPUT

Hi Paul,

you are right:)

gdalinfo -mm baltrum.img
Driver: HFA/Erdas Imagine Images (.img)
Size is 17559, 8026
Coordinate System is:
PROJCS["Transverse Mercator",
    GEOGCS["unknown",
        DATUM["unknown",
            SPHEROID["bessel",6377397.155,299.1528128000033],
            TOWGS84[590.5,69.5,411.6,-0.796,-0.052,-3.601,8.3]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433]],
    UNIT["meters",1],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",9],
    PARAMETER["scale_factor",1],
    PARAMETER["false_easting",3500000],
    PARAMETER["false_northing",0]]
Origin = (3391622.720000,5957462.080000)
Pixel Size = (0.32000000,-0.32000000)
Corner Coordinates:
Upper Left ( 3391622.720, 5957462.080) ( 7d21'26.28"E, 53d44'19.05"N)
Lower Left ( 3391622.720, 5954893.760) ( 7d21'29.52"E, 53d42'56.00"N)
Upper Right ( 3397241.600, 5957462.080) ( 7d26'32.75"E, 53d44'23.14"N)
Lower Right ( 3397241.600, 5954893.760) ( 7d26'35.82"E, 53d43'0.09"N)
Center ( 3394432.160, 5956177.920) ( 7d24'1.09"E, 53d43'39.60"N)
Band 1 Block=64x64 Type=UInt16, ColorInterp=Palette
    Computed Min/Max=4228.000,65535.000
  Color Table (RGB with 256 entries)
    0: 0,0,0,255
    1: 7,0,0,255
    2: 15,0,0,255
    3: 23,0,0,255
    4: 31,0,0,255
    5: 41,0,0,255
    6: 49,0,0,255
    ...
    255: 255,57,0,255

thanks for your help
    Otto

4

On Wed, Jun 30, 2004 at 11:14:45AM +0100, Paul Kelly wrote:

Hello Otto

Otto Dassau wrote:
>
>
> Also there is no coordinate system specified.

Otto, are you using the original GRASS driver or the new driver
of Radim Blazek (not in CVS yet)?

This should be:
>
> > name: Transverse Mercator
> > datum: potsdam
> > towgs84: 590.5,69.5,411.6,-0.796,-0.052,-3.601,8.30
> > proj: tmerc
> > ellps: bessel
> > a: 6377397.1550000003
> > es: 0.0066743722
> > f: 299.1528128000
> > lat_0: 0.0000000000
> > lon_0: 9.0000000000
> > k_0: 1.0000000000
> > x_0: 3500000.0000000000
> > y_0: 0.0000000000

This looks like a deficiency in the r.out.gdal script.

I don't think so. Because the r.out.gdal script is simply a wrapper
around gdal_translate.

In fact gdalinfo doesn't print the projection:

gdalinfo ~/grassdata/spearfish57/PERMANENT/cellhd/elevation.dem | head
Driver: GRASS/GRASS Database Rasters
Size is 633, 466
Coordinate System is `'
Origin = (590010.000000,4928000.000000)
Pixel Size = (30.00000000,-30.00000000)
Corner Coordinates:
Upper Left ( 590010.000, 4928000.000)
Lower Left ( 590010.000, 4914020.000)
Upper Right ( 609000.000, 4928000.000)
Lower Right ( 609000.000, 4914020.000)
[...]

I am using the new driver which contains the relevant code (I'm no
expert). Somewhere this info get's lost. Or should the driver make
use of the new GRASS 5.7 projection functions (such as GPJ_grass_to_wkt
in grass57/lib/proj/convert.c )?

Once solved in GDAL also r.out.gdal will write out the projection.

Markus

5

On Wed, 30 Jun 2004 15:02:55 +0200
Markus Neteler <neteler@itc.it> wrote:

On Wed, Jun 30, 2004 at 11:14:45AM +0100, Paul Kelly wrote:
> Hello Otto
>
> Otto Dassau wrote:
> >
> >
> > Also there is no coordinate system specified.

Otto, are you using the original GRASS driver or the new driver
of Radim Blazek (not in CVS yet)?

I compiled GDAL with the patches from Radim for new raster support also in QGIS. I didn't use "old" libgrass. The patches were:

grassdataset.cpp3
configure.in.patch2

> This should be:
> >
> > > name: Transverse Mercator
> > > datum: potsdam
> > > towgs84: 590.5,69.5,411.6,-0.796,-0.052,-3.601,8.30
> > > proj: tmerc
> > > ellps: bessel
> > > a: 6377397.1550000003
> > > es: 0.0066743722
> > > f: 299.1528128000
> > > lat_0: 0.0000000000
> > > lon_0: 9.0000000000
> > > k_0: 1.0000000000
> > > x_0: 3500000.0000000000
> > > y_0: 0.0000000000
>
> This looks like a deficiency in the r.out.gdal script.

I don't think so. Because the r.out.gdal script is simply a wrapper
around gdal_translate.

In fact gdalinfo doesn't print the projection:

gdalinfo ~/grassdata/spearfish57/PERMANENT/cellhd/elevation.dem | head
Driver: GRASS/GRASS Database Rasters
Size is 633, 466
Coordinate System is `'
Origin = (590010.000000,4928000.000000)
Pixel Size = (30.00000000,-30.00000000)
Corner Coordinates:
Upper Left ( 590010.000, 4928000.000)
Lower Left ( 590010.000, 4914020.000)
Upper Right ( 609000.000, 4928000.000)
Lower Right ( 609000.000, 4914020.000)
[...]

I am using the new driver which contains the relevant code (I'm no
expert). Somewhere this info get's lost. Or should the driver make
use of the new GRASS 5.7 projection functions (such as GPJ_grass_to_wkt
in grass57/lib/proj/convert.c )?

Once solved in GDAL also r.out.gdal will write out the projection.

Markus

_______________________________________________
grass5 mailing list
grass5@grass.itc.it
http://grass.itc.it/mailman/listinfo/grass5

6

On Wed, 30 Jun 2004, Otto Dassau wrote:

On Wed, 30 Jun 2004 15:02:55 +0200
Markus Neteler <neteler@itc.it> wrote:

On Wed, Jun 30, 2004 at 11:14:45AM +0100, Paul Kelly wrote:

Hello Otto

Otto Dassau wrote:

Also there is no coordinate system specified.

Otto, are you using the original GRASS driver or the new driver
of Radim Blazek (not in CVS yet)?

I compiled GDAL with the patches from Radim for new raster support also in QGIS. I didn't use "old" libgrass. The patches were:

grassdataset.cpp3
configure.in.patch2

Oh right I had assumed you were using old libgrass.

This should be:

name: Transverse Mercator
datum: potsdam
towgs84: 590.5,69.5,411.6,-0.796,-0.052,-3.601,8.30
proj: tmerc
ellps: bessel
a: 6377397.1550000003
es: 0.0066743722
f: 299.1528128000
lat_0: 0.0000000000
lon_0: 9.0000000000
k_0: 1.0000000000
x_0: 3500000.0000000000
y_0: 0.0000000000

This looks like a deficiency in the r.out.gdal script.

I don't think so. Because the r.out.gdal script is simply a wrapper
around gdal_translate.

In fact gdalinfo doesn't print the projection:

gdalinfo ~/grassdata/spearfish57/PERMANENT/cellhd/elevation.dem | head
Driver: GRASS/GRASS Database Rasters
Size is 633, 466
Coordinate System is `'
Origin = (590010.000000,4928000.000000)
Pixel Size = (30.00000000,-30.00000000)
Corner Coordinates:
Upper Left ( 590010.000, 4928000.000)
Lower Left ( 590010.000, 4914020.000)
Upper Right ( 609000.000, 4928000.000)
Lower Right ( 609000.000, 4914020.000)
[...]

I am using the new driver which contains the relevant code (I'm no
expert). Somewhere this info get's lost. Or should the driver make
use of the new GRASS 5.7 projection functions (such as GPJ_grass_to_wkt
in grass57/lib/proj/convert.c )?

Yes I see what you mean now, that Radim's new patches call GPJ_grass_to_wkt() and the co-ordinate system should be included automatically. I'll look at it if I have some time; mightn't be for a few days though.

The change to r.out.gdal should just be considered a temporary workaround then.

Also I noticed that the datum name didn't get printed in Otto's example (Deutsches_Haupt... for potsdam) either. That looks like a bug which I will have a look at.

Paul

7

On Wednesday 30 June 2004 17:11, Paul Kelly wrote:

> grassdataset.cpp3
> configure.in.patch2

>> I am using the new driver which contains the relevant code (I'm no
>> expert). Somewhere this info get's lost. Or should the driver make
>> use of the new GRASS 5.7 projection functions (such as GPJ_grass_to_wkt
>> in grass57/lib/proj/convert.c )?

Yes I see what you mean now, that Radim's new patches call
GPJ_grass_to_wkt() and the co-ordinate system should be included
automatically. I'll look at it if I have some time; mightn't be for a few
days though.

GPJ_grass_to_wkt is disabled in grassdataset.cpp3 because of circular references,
projection is not set.

Radim

8

I compiled GDAL with the patches from Radim for new raster support
also in QGIS.

I'm just wondering if the brand new 1.2.1 version of GDAL includes
what's needed for displaying GRASS rasters in QGIS? What about PostGIS
support?

see http://article.gmane.org/gmane.comp.gis.grass.devel/3999

Hamish

9

On Thursday 01 July 2004 05:56, Hamish wrote:

> I compiled GDAL with the patches from Radim for new raster support
> also in QGIS.

I'm just wondering if the brand new 1.2.1 version of GDAL includes
what's needed for displaying GRASS rasters in QGIS? What about PostGIS
support?

The patches for 5.7 are not yet in GDAL CVS, it is also possible
to use libgrass.

Radim