GRASS Programmer's Manual  6.4.2(2012)
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
proj/datum.c
Go to the documentation of this file.
1 
16 #include <unistd.h>
17 #include <string.h>
18 #include <ctype.h>
19 #include <stdlib.h>
20 
21 #include <grass/gis.h>
22 #include <grass/glocale.h>
23 #include <grass/gprojects.h>
24 #include "local_proto.h"
25 
37 int GPJ_get_datum_by_name(const char *name, struct gpj_datum *dstruct)
38 {
39  struct datum_list *list, *listhead;
40 
41  list = listhead = read_datum_table();
42 
43  while (list != NULL) {
44  if (G_strcasecmp(name, list->name) == 0) {
45  dstruct->name = G_store(list->name);
46  dstruct->longname = G_store(list->longname);
47  dstruct->ellps = G_store(list->ellps);
48  dstruct->dx = list->dx;
49  dstruct->dy = list->dy;
50  dstruct->dz = list->dz;
51  free_datum_list(listhead);
52  return 1;
53  }
54  list = list->next;
55  }
56  free_datum_list(listhead);
57  return -1;
58 }
59 
85 int GPJ_get_default_datum_params_by_name(const char *name, char **params)
86 {
87  struct gpj_datum_transform_list *list, *old;
88  int count = 1;
89 
91 
92  if (list == NULL) {
93  *params = NULL;
94  return -1;
95  }
96 
97  /* Take the first parameter set in the list as the default
98  * (will normally be a 3-parameter transformation) */
99  *params = G_store(list->params);
100 
101  while (list->next != NULL) {
102  count++;
103  old = list;
104  list = list->next;
105  G_free(old);
106  }
107 
108  G_free(list);
109  return count;
110 
111 }
112 
136 int GPJ_get_datum_params(char **name, char **params)
137 {
138  int ret;
139  struct Key_Value *proj_keys = G_get_projinfo();
140 
141  ret = GPJ__get_datum_params(proj_keys, name, params);
142  G_free_key_value(proj_keys);
143 
144  return ret;
145 }
146 
174 int GPJ__get_datum_params(struct Key_Value *projinfo,
175  char **datumname, char **params)
176 {
177  int returnval = -1;
178 
179  if (NULL != G_find_key_value("datum", projinfo)) {
180  *datumname = G_store(G_find_key_value("datum", projinfo));
181  returnval = 1;
182  }
183  else
184  *datumname = NULL;
185 
186  if (G_find_key_value("datumparams", projinfo) != NULL) {
187  *params = G_store(G_find_key_value("datumparams", projinfo));
188  returnval = 2;
189  }
190  else if (G_find_key_value("nadgrids", projinfo) != NULL) {
191  const char *gisbase = G_gisbase();
192 
193  G_asprintf(params, "nadgrids=%s%s/%s", gisbase, GRIDDIR,
194  G_find_key_value("nadgrids", projinfo));
195  returnval = 2;
196  }
197  else if (G_find_key_value("towgs84", projinfo) != NULL) {
198  G_asprintf(params, "towgs84=%s",
199  G_find_key_value("towgs84", projinfo));
200  returnval = 2;
201  }
202  else if (G_find_key_value("dx", projinfo) != NULL
203  && G_find_key_value("dy", projinfo) != NULL
204  && G_find_key_value("dz", projinfo) != NULL) {
205  G_asprintf(params, "towgs84=%s,%s,%s",
206  G_find_key_value("dx", projinfo),
207  G_find_key_value("dy", projinfo),
208  G_find_key_value("dz", projinfo));
209  returnval = 2;
210  }
211  else
212  *params = NULL;
213 
214  return returnval;
215 
216 }
217 
240 int GPJ_ask_datum_params(const char *datumname, char **params)
241 {
242  char buff[1024], answer[100];
243  char *Tmp_file;
244  FILE *Tmp_fd = NULL;
245  struct gpj_datum_transform_list *list, *listhead, *old;
246  int transformcount, currenttransform;
247 
248  if (G_strcasecmp(datumname, "custom") != 0) {
249  Tmp_file = G_tempfile();
250  if (NULL == (Tmp_fd = fopen(Tmp_file, "w"))) {
251  G_warning(_("Unable to open temporary file"));
252  }
253 
254  fprintf(Tmp_fd, "Number\tDetails\t\n---\n");
255  listhead = GPJ_get_datum_transform_by_name(datumname);
256  list = listhead;
257  transformcount = 0;
258  while (list != NULL) {
259  /* Count how many sets of transformation paramters have been
260  * defined for this datum and print them to a temporary file
261  * in case the user asks for them to be displayed */
262  fprintf(Tmp_fd,
263  "%d\tUsed in %s\n\t(PROJ.4 Params %s)\n\t%s\n---\n",
264  list->count, list->where_used, list->params,
265  list->comment);
266  list = list->next;
267  transformcount++;
268  }
269  fclose(Tmp_fd);
270 
271  for (;;) {
272  do {
273  fprintf(stderr,
274  ("\nNow select Datum Transformation Parameters\n"));
275  fprintf(stderr,
276  ("Please think carefully about the area covered by your data\n"
277  "and the accuracy you require before making your selection.\n"));
278  fprintf(stderr,
279  ("\nEnter 'list' to see the list of available Parameter sets\n"));
280  fprintf(stderr,
281  ("Enter the corresponding number, or <RETURN> to cancel request\n"));
282  fprintf(stderr, ">");
283  } while (!G_gets(answer));
284  G_strip(answer);
285  if (strlen(answer) == 0) {
286  remove(Tmp_file);
287  G_free(Tmp_file);
288  return -1;
289  }
290  if (strcmp(answer, "list") == 0) {
291  char *pager;
292 
293  pager = getenv("GRASS_PAGER");
294  if (!pager || strlen(pager) == 0)
295  pager = "cat";
296 
297  /* Always print interactive output to stderr */
298  sprintf(buff, "%s \"%s\" 1>&2", pager,
299  G_convert_dirseps_to_host(Tmp_file));
300  G_system(buff);
301  }
302  else {
303  if ((sscanf(answer, "%d", &currenttransform) != 1) ||
304  currenttransform > transformcount ||
305  currenttransform < 1) {
306 
307  /* If a number was not typed, was less than 0 or greater
308  * than the number of sets of parameters, ask again */
309  fprintf(stderr, ("\ninvalid transformation number\n"));
310  }
311  else
312  break;
313  }
314 
315  }
316  remove(Tmp_file);
317  G_free(Tmp_file);
318 
319  list = listhead;
320  while (list != NULL) {
321  /* Search through the linked list to find the parameter string
322  * that corresponds to the number entered */
323  if (list->count == currenttransform)
324  G_asprintf(params, list->params);
325 
326  /* Continue to end of list even after we find it, to free all
327  * the memory used */
328  old = list;
329  list = old->next;
330  G_free(old);
331  }
332  }
333  else {
334  /* Here we ask the user to enter customised parameters */
335  for (;;) {
336  do {
337  fprintf(stderr,
338  ("\nPlease specify datum transformation parameters in PROJ.4 syntax. Examples:\n"));
339  fprintf(stderr,
340  ("\ttowgs84=dx,dy,dz\t(3-parameter transformation)\n"));
341  fprintf(stderr,
342  ("\ttowgs84=dx,dy,dz,rx,ry,rz,m\t(7-parameter transformation)\n"));
343  fprintf(stderr,
344  ("\tnadgrids=alaska\t(Tables-based grid-shifting transformation)\n"));
345  fprintf(stderr, _("Hit RETURN to cancel request\n"));
346  fprintf(stderr, ">");
347  } while (!G_gets(answer));
348  G_strip(answer);
349  if (strlen(answer) == 0)
350  return -1;
351  G_asprintf(params, answer);
352  sprintf(buff,
353  "Parameters to be used are:\n\"%s\"\nIs this correct?",
354  *params);
355  if (G_yes(buff, 1))
356  break;
357 
358  }
359 
360  }
361 
362  return 1;
363 
364 }
365 
378 struct gpj_datum_transform_list *GPJ_get_datum_transform_by_name(const char
379  *inputname)
380 {
381  FILE *fd;
382  char file[GPATH_MAX];
383  char buf[1024];
384  int line;
385  struct gpj_datum_transform_list *current = NULL, *outputlist = NULL;
386  struct gpj_datum dstruct;
387  int count = 0;
388 
389  GPJ_get_datum_by_name(inputname, &dstruct);
390  if (dstruct.dx < 99999 && dstruct.dy < 99999 && dstruct.dz < 99999) {
391  /* Include the old-style dx dy dz parameters from datum.table at the
392  * start of the list, unless these have been set to all 99999 to
393  * indicate only entries in datumtransform.table should be used */
394  if (current == NULL)
395  current = outputlist =
396  G_malloc(sizeof(struct gpj_datum_transform_list));
397  else
398  current = current->next =
399  G_malloc(sizeof(struct gpj_datum_transform_list));
400  G_asprintf(&(current->params), "towgs84=%.3f,%.3f,%.3f", dstruct.dx,
401  dstruct.dy, dstruct.dz);
402  G_asprintf(&(current->where_used), "whole %s region", inputname);
403  G_asprintf(&(current->comment),
404  "Default 3-Parameter Transformation (May not be optimum for "
405  "older datums; use this only if no more appropriate options "
406  "are available.)");
407  count++;
408  current->count = count;
409  current->next = NULL;
410  }
411  GPJ_free_datum(&dstruct);
412 
413  /* Now check for additional parameters in datumtransform.table */
414 
415  sprintf(file, "%s%s", G_gisbase(), DATUMTRANSFORMTABLE);
416 
417  fd = fopen(file, "r");
418  if (!fd) {
419  G_warning(_("Unable to open datum table file <%s>"), file);
420  return outputlist;
421  }
422 
423  for (line = 1; G_getl2(buf, sizeof(buf), fd); line++) {
424  char name[100], params[1024], where_used[1024], comment[1024];
425 
426  G_strip(buf);
427  if (*buf == '\0' || *buf == '#')
428  continue;
429 
430  if (sscanf(buf, "%99s \"%1023[^\"]\" \"%1023[^\"]\" \"%1023[^\"]\"",
431  name, params, where_used, comment) != 4) {
432  G_warning(_("Error in datum table file <%s>, line %d"), file,
433  line);
434  continue;
435  }
436 
437  if (G_strcasecmp(inputname, name) == 0) {
438  /* If the datum name in this line matches the one we are
439  * looking for, add an entry to the linked list */
440  if (current == NULL)
441  current = outputlist =
442  G_malloc(sizeof(struct gpj_datum_transform_list));
443  else
444  current = current->next =
445  G_malloc(sizeof(struct gpj_datum_transform_list));
446  current->params = G_store(params);
447  current->where_used = G_store(where_used);
448  current->comment = G_store(comment);
449  count++;
450  current->count = count;
451  current->next = NULL;
452  }
453  }
454 
455  fclose(fd);
456 
457  return outputlist;
458 
459 }
460 
472 {
473  FILE *fd;
474  char file[GPATH_MAX];
475  char buf[4096];
476  int line;
477  struct datum_list *current = NULL, *outputlist = NULL;
478  int count = 0;
479 
480  sprintf(file, "%s%s", G_gisbase(), DATUMTABLE);
481 
482  fd = fopen(file, "r");
483  if (!fd) {
484  G_warning(_("Unable to open datum table file <%s>"), file);
485  return NULL;
486  }
487 
488  for (line = 1; G_getl2(buf, sizeof(buf), fd); line++) {
489  char name[100], descr[1024], ellps[100];
490  double dx, dy, dz;
491 
492  G_strip(buf);
493  if (*buf == '\0' || *buf == '#')
494  continue;
495 
496  if (sscanf(buf, "%s \"%1023[^\"]\" %s dx=%lf dy=%lf dz=%lf",
497  name, descr, ellps, &dx, &dy, &dz) != 6) {
498  G_warning(_("Error in datum table file <%s>, line %d"), file,
499  line);
500  continue;
501  }
502 
503  if (current == NULL)
504  current = outputlist = G_malloc(sizeof(struct datum_list));
505  else
506  current = current->next = G_malloc(sizeof(struct datum_list));
507  current->name = G_store(name);
508  current->longname = G_store(descr);
509  current->ellps = G_store(ellps);
510  current->dx = dx;
511  current->dy = dy;
512  current->dz = dz;
513  current->next = NULL;
514 
515  count++;
516  }
517 
518  fclose(fd);
519 
520  return outputlist;
521 }
522 
529 void GPJ_free_datum(struct gpj_datum *dstruct)
530 {
531  G_free(dstruct->name);
532  G_free(dstruct->longname);
533  G_free(dstruct->ellps);
534  return;
535 }
536 
543 void free_datum_list(struct datum_list *dstruct)
544 {
545  struct datum_list *old;
546 
547  while (dstruct != NULL) {
548  G_free(dstruct->name);
549  G_free(dstruct->longname);
550  G_free(dstruct->ellps);
551  old = dstruct;
552  dstruct = old->next;
553  G_free(old);
554  }
555 
556  return;
557 }