[GRASS5] datum transforms...

[Eric G. Miller]

After some prodding by Markus, I began looking at the datum transform
problem again.  And, I see pj_transform() seems to already do what we
need (inverse project, datum shift, forward project) with any part
possibly being not applicable.  Does anyone know of a good reason not to
use it?

After a little hacking, I have a modified s.proj that seems to work
well (not heavily tested).

Changes required:

pj_datums.c: s/ntv1 can.dat/ntv1_can.dat/
do_proj.c:
   * LL<->LL is no longer a no-op.
   * change all variations of pj_fwd and pj_inv to one call to
     pj_transform
get_proj.c:
   pj_get_kv()
     * don't just return 1 on strncmp(proj_in,"Lat",3) == 0 ||
                            strcmp(info->proj,"ll") == 0
     * when in_proj_keys->key[i] = "proj" 
            && in_proj_keys->value[i] = "ll"
	    set buffa = "proj=latlong" 
       (libproj uses "latlong" or "longlat" vs. GRASS's "ll")

Before using pj_transform(), set the environment variable PROJ_LIB to
$GISBASE/etc/nad

#define GRIDDIR "/etc/nad"

static void set_proj_lib_env (void)
{
    const char *gisbase = G_gisbase();
    char *value = NULL;
    size_t len = strlen(gisbase) * sizeof(GRIDDIR) + 1;
    value = G_malloc (len);
    sprintf (value, "%s%s", gisbase, GRIDDIR);
    if (setenv ("PROJ_LIB", value, 1) == 0)
        fprintf(stderr, "set PROJ_LIB='%s'\n", value);
    else
        fprintf(stderr, "failed to set PROJ_LIB='%s'\n", value);

    free (value);
}

That was about all I needed to do to get s.proj2 working with
pj_transform (it took me the longest to figure out the space vs.
underscore in pj_datums.c was the cause of a failure...).

Since pj_transform operates on arrays and has height values (required),
a slightly different interface from do_proj might be appropriate (to
take advantage of transforming arrays in one "pass").

Anyway, comments?

-- 
"...the plural of anecdote is [not?] data."  - attrib. to George Stigler