Loading...
Searching...
No Matches
get_ellipse.c
Go to the documentation of this file.
/*!
\file lib/gis/get_ellipse.c
\brief GIS Library - Getting ellipsoid parameters from the database.
This routine returns the ellipsoid parameters from the database.
If the PROJECTION_FILE exists in the PERMANENT mapset, read info
from that file, otherwise return WGS 84 values.
New 05/2000 by al: for datum shift the f parameter is needed too.
This all is not a clean design, but it keeps backward-
compatibility.
Looks up ellipsoid in ellipsoid table and returns the
a, e2 and f parameters for the ellipsoid
(C) 2001-2009 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.
\author CERL
*/
#include <unistd.h>
#include <ctype.h>
#include <string.h>
#include <stdlib.h>
#include <math.h> /* for sqrt() */
#include <grass/gis.h>
#include <grass/glocale.h>
static const char PERMANENT[] = "PERMANENT";
static struct table {
struct ellipse {
char *name;
char *descr;
double a;
double e2;
double f;
} *ellipses;
int count;
int size;
int initialized;
} table;
/* static int get_a_e2 (char *, char *, double *,double *); */
static int get_a_e2_f(const char *, const char *, double *, double *, double *);
static int compare_ellipse_names(const void *, const void *);
static int get_ellipsoid_parameters(struct Key_Value *, double *, double *);
/*!
* \brief get ellipsoid parameters
*
* This routine returns the semi-major axis <b>a</b> (in meters) and
* the eccentricity squared <b>e2</b> for the ellipsoid associated
* with the database. If there is no ellipsoid explicitly associated
* with the database, it returns the values for the WGS 84 ellipsoid.
*
* \param[out] a semi-major axis
* \param[out] e2 eccentricity squared
*
* \return 1 success
* \return 0 default values used
*/
int G_get_ellipsoid_parameters(double *a, double *e2)
{
int stat;
char ipath[GPATH_MAX];
struct Key_Value *proj_keys;
proj_keys = NULL;
G_file_name(ipath, "", PROJECTION_FILE, PERMANENT);
if (access(ipath, 0) != 0) {
*a = 6378137.0;
*e2 = .006694385;
return 0;
}
proj_keys = G_read_key_value_file(ipath);
stat = get_ellipsoid_parameters(proj_keys, a, e2);
G_free_key_value(proj_keys);
return stat;
}
/*!
* \brief Get ellipsoid parameters by name
*
* This routine returns the semi-major axis <i>a</i> (in meters) and
* eccentricity squared <i>e2</i> for the named ellipsoid.
*
* \param name ellipsoid name
* \param[out] a semi-major axis
* \param[out] e2 eccentricity squared
*
* \return 1 on success
* \return 0 if ellipsoid not found
*/
int G_get_ellipsoid_by_name(const char *name, double *a, double *e2)
{
int i;
G_read_ellipsoid_table(0);
for (i = 0; i < table.count; i++) {
if (G_strcasecmp(name, table.ellipses[i].name) == 0) {
*a = table.ellipses[i].a;
*e2 = table.ellipses[i].e2;
return 1;
}
}
return 0;
}
/*!
* \brief Get ellipsoid name
*
* This function returns a pointer to the short name for the
* <i>n</i><i>th</i> ellipsoid. If <i>n</i> is less than 0 or greater
* than the number of known ellipsoids, it returns a NULL pointer.
*
* \param n ellipsoid identificator
*
* \return ellipsoid name
* \return NULL if no ellipsoid found
*/
const char *G_ellipsoid_name(int n)
{
G_read_ellipsoid_table(0);
return n >= 0 && n < table.count ? table.ellipses[n].name : NULL;
}
/*!
* \brief Get spheroid parameters by name
*
* This function returns the semi-major axis <i>a</i> (in meters), the
* eccentricity squared <i>e2</i> and the inverse flattening <i>f</i>
* for the named ellipsoid.
*
* \param name spheroid name
* \param[out] a semi-major axis
* \param[out] e2 eccentricity squared
* \param[out] f inverse flattening
*
* \return 1 on success
* \return 0 if no spheroid found
*/
int G_get_spheroid_by_name(const char *name, double *a, double *e2, double *f)
{
int i;
G_read_ellipsoid_table(0);
for (i = 0; i < table.count; i++) {
if (G_strcasecmp(name, table.ellipses[i].name) == 0) {
*a = table.ellipses[i].a;
*e2 = table.ellipses[i].e2;
*f = table.ellipses[i].f;
return 1;
}
}
return 0;
}
/*!
* \brief Get description for nth ellipsoid
*
* This function returns a pointer to the description text for the
* <i>n</i>th ellipsoid. If <i>n</i> is less than 0 or greater
* than the number of known ellipsoids, it returns a NULL pointer.
*
* \param n ellipsoid identificator
*
* \return pointer to ellipsoid description
* \return NULL if no ellipsoid found
*/
const char *G_ellipsoid_description(int n)
{
G_read_ellipsoid_table(0);
return n >= 0 && n < table.count ? table.ellipses[n].descr : NULL;
}
static int get_a_e2_f(const char *s1, const char *s2, double *a, double *e2,
double *f)
{
double b, recipf;
if (sscanf(s1, "a=%lf", a) != 1)
return 0;
if (*a <= 0.0)
return 0;
if (sscanf(s2, "e=%lf", e2) == 1) {
*f = (double)1.0 / -sqrt(((double)1.0 - *e2)) + (double)1.0;
return (*e2 >= 0.0);
}
if (sscanf(s2, "f=1/%lf", f) == 1) {
if (*f <= 0.0)
return 0;
recipf = (double)1.0 / (*f);
*e2 = recipf + recipf - recipf * recipf;
return (*e2 >= 0.0);
}
if (sscanf(s2, "b=%lf", &b) == 1) {
if (b <= 0.0)
return 0;
if (b == *a) {
*f = 0.0;
*e2 = 0.0;
}
else {
recipf = ((*a) - b) / (*a);
*f = (double)1.0 / recipf;
*e2 = recipf + recipf - recipf * recipf;
}
return (*e2 >= 0.0);
}
return 0;
}
static int compare_ellipse_names(const void *pa, const void *pb)
{
const struct ellipse *a = pa;
const struct ellipse *b = pb;
return G_strcasecmp(a->name, b->name);
}
/*!
\brief Read ellipsoid table
\param fatal non-zero value for G_fatal_error(), otherwise
G_warning() is used
\return 1 on success
\return 0 on error
*/
int G_read_ellipsoid_table(int fatal)
{
FILE *fd;
char file[GPATH_MAX];
char buf[1024];
char badlines[256];
int line;
int err;
if (G_is_initialized(&table.initialized))
return 1;
snprintf(file, sizeof(file), "%s/etc/proj/ellipse.table", G_gisbase());
fd = fopen(file, "r");
if (fd == NULL) {
(fatal ? G_fatal_error : G_warning)(
_("Unable to open ellipsoid table file <%s>"), file);
G_initialize_done(&table.initialized);
return 0;
}
err = 0;
*badlines = 0;
for (line = 1; G_getl2(buf, sizeof buf, fd); line++) {
char name[100], descr[100], buf1[100], buf2[100];
struct ellipse *e;
G_strip(buf);
if (*buf == 0 || *buf == '#')
continue;
if (sscanf(buf, "%s \"%99[^\"]\" %s %s", name, descr, buf1, buf2) !=
4) {
err++;
snprintf(buf, sizeof(buf), " %d", line);
if (*badlines)
G_strlcat(badlines, ",", sizeof(badlines));
G_strlcat(badlines, buf, sizeof(badlines));
continue;
}
if (table.count >= table.size) {
table.size += 60;
table.ellipses =
G_realloc(table.ellipses, table.size * sizeof(struct ellipse));
}
e = &table.ellipses[table.count];
e->name = G_store(name);
e->descr = G_store(descr);
if (get_a_e2_f(buf1, buf2, &e->a, &e->e2, &e->f) ||
get_a_e2_f(buf2, buf1, &e->a, &e->e2, &e->f))
table.count++;
else {
err++;
snprintf(buf, sizeof(buf), " %d", line);
if (*badlines)
G_strlcat(badlines, ",", sizeof(badlines));
G_strlcat(badlines, buf, sizeof(badlines));
continue;
}
}
fclose(fd);
if (!err) {
/* over correct typed version */
qsort(table.ellipses, table.count, sizeof(struct ellipse),
compare_ellipse_names);
G_initialize_done(&table.initialized);
return 1;
}
(fatal ? G_fatal_error : G_warning)(
n_(("Line%s of ellipsoid table file <%s> is invalid"),
("Lines%s of ellipsoid table file <%s> are invalid"), err),
badlines, file);
G_initialize_done(&table.initialized);
return 0;
}
static int get_ellipsoid_parameters(struct Key_Value *proj_keys, double *a,
double *e2)
{
const char *str, *str1;
if (!proj_keys) {
return -1;
}
if ((str = G_find_key_value("ellps", proj_keys)) != NULL) {
if (strncmp(str, "sphere", 6) == 0) {
str = G_find_key_value("a", proj_keys);
if (str != NULL) {
if (sscanf(str, "%lf", a) != 1)
G_fatal_error(_("Invalid a: field '%s' in file %s in <%s>"),
str, PROJECTION_FILE, PERMANENT);
}
else
*a = 6370997.0;
*e2 = 0.0;
return 0;
}
else {
if (G_get_ellipsoid_by_name(str, a, e2) == 0)
G_fatal_error(_("Invalid ellipsoid '%s' in file %s in <%s>"),
str, PROJECTION_FILE, PERMANENT);
else
return 1;
}
}
else {
str = G_find_key_value("a", proj_keys);
str1 = G_find_key_value("es", proj_keys);
if ((str != NULL) && (str1 != NULL)) {
if (sscanf(str, "%lf", a) != 1)
G_fatal_error(_("Invalid a: field '%s' in file %s in <%s>"),
str, PROJECTION_FILE, PERMANENT);
if (sscanf(str1, "%lf", e2) != 1)
G_fatal_error(_("Invalid es: field '%s' in file %s in <%s>"),
str, PROJECTION_FILE, PERMANENT);
return 1;
}
else {
str = G_find_key_value("proj", proj_keys);
if ((str == NULL) || (strcmp(str, "ll") == 0)) {
*a = 6378137.0;
*e2 = .006694385;
return 0;
}
else
G_fatal_error(_("No ellipsoid info given in file %s in <%s>"),
PROJECTION_FILE, PERMANENT);
}
}
return 1;
}
NULL
#define NULL
Definition ccmath.h:32
G_initialize_done
void G_initialize_done(int *p)
Definition counter.c:77
G_is_initialized
int G_is_initialized(int *p)
Definition counter.c:60
b
double b
Definition driver/set_window.c:5
G_file_name
char * G_file_name(char *path, const char *element, const char *name, const char *mapset)
Builds full path names to GIS data files.
Definition file_name.c:61
G_read_ellipsoid_table
int G_read_ellipsoid_table(int fatal)
Read ellipsoid table.
Definition get_ellipse.c:246
G_ellipsoid_description
const char * G_ellipsoid_description(int n)
Get description for nth ellipsoid.
Definition get_ellipse.c:182
G_get_spheroid_by_name
int G_get_spheroid_by_name(const char *name, double *a, double *e2, double *f)
Get spheroid parameters by name.
Definition get_ellipse.c:153
G_ellipsoid_name
const char * G_ellipsoid_name(int n)
Get ellipsoid name.
Definition get_ellipse.c:132
G_get_ellipsoid_parameters
int G_get_ellipsoid_parameters(double *a, double *e2)
get ellipsoid parameters
Definition get_ellipse.c:66
G_get_ellipsoid_by_name
int G_get_ellipsoid_by_name(const char *name, double *a, double *e2)
Get ellipsoid parameters by name.
Definition get_ellipse.c:104
PERMANENT
#define PERMANENT
Definition get_projinfo.c:19
G_getl2
int G_getl2(char *buf, int n, FILE *fd)
Gets a line of text from a file of any pedigree.
Definition getl.c:60
G_fatal_error
void G_fatal_error(const char *msg,...)
Print a fatal error message to stderr.
Definition gis/error.c:159
G_warning
void G_warning(const char *msg,...)
Print a warning message to stderr.
Definition gis/error.c:203
G_gisbase
const char * G_gisbase(void)
Get full path name of the top level module directory.
Definition gisbase.c:39
G_free_key_value
void G_free_key_value(struct Key_Value *kv)
Free allocated Key_Value structure.
Definition key_value1.c:104
G_find_key_value
const char * G_find_key_value(const char *key, const struct Key_Value *kv)
Find given key (case sensitive).
Definition key_value1.c:85
G_read_key_value_file
struct Key_Value * G_read_key_value_file(const char *file)
Read key/values pairs from file.
Definition key_value3.c:55
file
#define file
name
const char * name
Definition named_colr.c:6
G_strcasecmp
int G_strcasecmp(const char *x, const char *y)
String compare ignoring case (upper or lower).
Definition strings.c:47
G_store
char * G_store(const char *s)
Copy string to allocated memory.
Definition strings.c:87
G_strip
void G_strip(char *buf)
Removes all leading and trailing white space from string.
Definition strings.c:300
G_strlcat
size_t G_strlcat(char *dst, const char *src, size_t dsize)
Size-bounded string concatenation.
Definition strlcat.c:59
err
SYMBOL * err(FILE *fp, SYMBOL *s, char *msg)
Definition symbol/read.c:216