PageRenderTime 31ms CodeModel.GetById 21ms RepoModel.GetById 1ms app.codeStats 0ms

/Proj4/nad_cvt.c

http://github.com/route-me/route-me
C | 71 lines | 56 code | 5 blank | 10 comment | 19 complexity | 17c0f73ecf4e786e8ae72fede8442748 MD5 | raw file
  1. #ifndef lint
  2. static const char SCCSID[]="@(#)nad_cvt.c 4.3 95/09/23 GIE REL";
  3. #endif
  4. #define PJ_LIB__
  5. #include "projects.h"
  6. #define MAX_TRY 9
  7. #define TOL 1e-12
  8. LP
  9. nad_cvt(LP in, int inverse, struct CTABLE *ct) {
  10. LP t, tb;
  11. if (in.lam == HUGE_VAL)
  12. return in;
  13. /* normalize input to ll origin */
  14. tb = in;
  15. tb.lam -= ct->ll.lam;
  16. tb.phi -= ct->ll.phi;
  17. tb.lam = adjlon(tb.lam - PI) + PI;
  18. t = nad_intr(tb, ct);
  19. if (inverse) {
  20. LP del, dif;
  21. int i = MAX_TRY;
  22. if (t.lam == HUGE_VAL) return t;
  23. t.lam = tb.lam + t.lam;
  24. t.phi = tb.phi - t.phi;
  25. do {
  26. del = nad_intr(t, ct);
  27. /* This case used to return failure, but I have
  28. changed it to return the first order approximation
  29. of the inverse shift. This avoids cases where the
  30. grid shift *into* this grid came from another grid.
  31. While we aren't returning optimally correct results
  32. I feel a close result in this case is better than
  33. no result. NFW
  34. To demonstrate use -112.5839956 49.4914451 against
  35. the NTv2 grid shift file from Canada. */
  36. if (del.lam == HUGE_VAL)
  37. {
  38. if( getenv( "PROJ_DEBUG" ) != NULL )
  39. fprintf( stderr,
  40. "Inverse grid shift iteration failed, presumably at grid edge.\n"
  41. "Using first approximation.\n" );
  42. /* return del */;
  43. break;
  44. }
  45. t.lam -= dif.lam = t.lam - del.lam - tb.lam;
  46. t.phi -= dif.phi = t.phi + del.phi - tb.phi;
  47. } while (i-- && fabs(dif.lam) > TOL && fabs(dif.phi) > TOL);
  48. if (i < 0) {
  49. if( getenv( "PROJ_DEBUG" ) != NULL )
  50. fprintf( stderr,
  51. "Inverse grid shift iterator failed to converge.\n" );
  52. t.lam = t.phi = HUGE_VAL;
  53. return t;
  54. }
  55. in.lam = adjlon(t.lam + ct->ll.lam);
  56. in.phi = t.phi + ct->ll.phi;
  57. } else {
  58. if (t.lam == HUGE_VAL)
  59. in = t;
  60. else {
  61. in.lam -= t.lam;
  62. in.phi += t.phi;
  63. }
  64. }
  65. return in;
  66. }