[GRASSLIST:3469] Re: need help with a script PS

Paul Kelly paul-grass at stjohnspoint.co.uk
Sun May 23 18:06:25 EDT 2004


Hello Michael

On Sat, 22 May 2004, Michael Barton wrote:

[...]
> The script works fine if I don't enter the information. If enter the 
> information, echo the gdal command created by the script, and paste the 
> echo'ed output into the command line, it runs fine. I've tried every 
> permutation I can think of to get gdalwarp to recognize the -t_srs parameter. 
> Obviously, I haven't yet hit on the correct on.

I think it's some kind of quoting problem---there are a lot of " and I 
didn't completely follow it.
What about changing line 94 from
         tgt_srs="-t_srs "\"$tgt_srs\"
to
         tgt_srs="-t_srs \"$tgt_srs\""

BUT---in wondering why you didn't just use 'g.proj -jf' to get the target 
co-ordinate system to match that of the current GRASS location, it gave 
me a good idea for a script to automatically reproject and import an image 
into GRASS using gdalwarp and r.in.gdal, instead of having to run 
r.in.gdal on a separate location and then r.proj to re-project the image 
inside GRASS.

See attached; do you think it is useful?
The main reservation I have about it is that it removes the need for users 
to think about what they are doing when re-projecting data, e.g. the 
re-sampling algorithms used in gdalwarp are hidden from the user and they 
might not be aware of distortion caused at the edge of null areas and
coastlines and things like that. More automation means less thinking about
how exactly you are processing your data and what errors you might be
introducing.

The other possible problem I can see is due to the recently confirmed 
'feature' in PROJ whereby if no datum transformation parameters are given in 
a co-ordinate system description, +towgs84=0,0,0 is assumed and datum 
transformation performed on this basis. If the input dataset doesn't have
datum information (leading to the above assumption) but the GRASS location 
does, I suppose there is a chance of shifting errors. Not sure how likely
this is to happen in practice but it could be a problem with this idea.

Paul
-------------- next part --------------
#!/bin/sh
############################################################################
#
# MODULE:	r.in.gdalproj
# AUTHOR(S):    Michael Barton and Paul Kelly
# PURPOSE:	Import a GDAL-compatible image in a different projection 
#               from the current GRASS location, directly into the location
#               without having to use r.proj (uses gdalwarp instead)
# COPYRIGHT:	(C) 2004 by the GRASS Development Team
#
#		This program is free software under the GNU General Public
#		License (>=v2). Read the file COPYING that comes with GRASS
#		for details.
#
#############################################################################

#%Module
#%  description: re-projects and imports raster data into GRASS using gdalwarp and r.in.gdal
#%End
#%option
#% key: input
#% type: string
#% gisprompt: file,file,file
#% description: Input raster image or dataset to be imported into GRASS
#% required : yes
#%end
#%option
#% key: output
#% type: string
#% description: Name of resulting GRASS raster layer
#% required : yes
#%end


if  [ -z "$GISBASE" ]
then
	echo
	echo You must be in GRASS GIS to run this program
	echo
	exit 1
fi

if   [ "$1" != "@ARGS_PARSED@" ]
then
	exec g.parser "$0" "$@"
fi

echo
echo Re-projecting input image to current GRASS location co-ordinate system
echo

tempfile=`g.tempfile $$`.tif

gdalwarp -t_srs "`g.proj -jf`" "$GIS_OPT_input" "$tempfile"

echo 
echo Importing into GRASS
echo 

r.in.gdal input="$tempfile" output=$GIS_OPT_output

rm -f "$tempfile"


More information about the grass-user mailing list