vandg.js 3.1 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120
  1. var adjust_lon = require('../common/adjust_lon');
  2. var HALF_PI = Math.PI/2;
  3. var EPSLN = 1.0e-10;
  4. var asinz = require('../common/asinz');
  5. /* Initialize the Van Der Grinten projection
  6. ----------------------------------------*/
  7. exports.init = function() {
  8. //this.R = 6370997; //Radius of earth
  9. this.R = this.a;
  10. };
  11. exports.forward = function(p) {
  12. var lon = p.x;
  13. var lat = p.y;
  14. /* Forward equations
  15. -----------------*/
  16. var dlon = adjust_lon(lon - this.long0);
  17. var x, y;
  18. if (Math.abs(lat) <= EPSLN) {
  19. x = this.x0 + this.R * dlon;
  20. y = this.y0;
  21. }
  22. var theta = asinz(2 * Math.abs(lat / Math.PI));
  23. if ((Math.abs(dlon) <= EPSLN) || (Math.abs(Math.abs(lat) - HALF_PI) <= EPSLN)) {
  24. x = this.x0;
  25. if (lat >= 0) {
  26. y = this.y0 + Math.PI * this.R * Math.tan(0.5 * theta);
  27. }
  28. else {
  29. y = this.y0 + Math.PI * this.R * -Math.tan(0.5 * theta);
  30. }
  31. // return(OK);
  32. }
  33. var al = 0.5 * Math.abs((Math.PI / dlon) - (dlon / Math.PI));
  34. var asq = al * al;
  35. var sinth = Math.sin(theta);
  36. var costh = Math.cos(theta);
  37. var g = costh / (sinth + costh - 1);
  38. var gsq = g * g;
  39. var m = g * (2 / sinth - 1);
  40. var msq = m * m;
  41. var con = Math.PI * this.R * (al * (g - msq) + Math.sqrt(asq * (g - msq) * (g - msq) - (msq + asq) * (gsq - msq))) / (msq + asq);
  42. if (dlon < 0) {
  43. con = -con;
  44. }
  45. x = this.x0 + con;
  46. //con = Math.abs(con / (Math.PI * this.R));
  47. var q = asq + g;
  48. con = Math.PI * this.R * (m * q - al * Math.sqrt((msq + asq) * (asq + 1) - q * q)) / (msq + asq);
  49. if (lat >= 0) {
  50. //y = this.y0 + Math.PI * this.R * Math.sqrt(1 - con * con - 2 * al * con);
  51. y = this.y0 + con;
  52. }
  53. else {
  54. //y = this.y0 - Math.PI * this.R * Math.sqrt(1 - con * con - 2 * al * con);
  55. y = this.y0 - con;
  56. }
  57. p.x = x;
  58. p.y = y;
  59. return p;
  60. };
  61. /* Van Der Grinten inverse equations--mapping x,y to lat/long
  62. ---------------------------------------------------------*/
  63. exports.inverse = function(p) {
  64. var lon, lat;
  65. var xx, yy, xys, c1, c2, c3;
  66. var a1;
  67. var m1;
  68. var con;
  69. var th1;
  70. var d;
  71. /* inverse equations
  72. -----------------*/
  73. p.x -= this.x0;
  74. p.y -= this.y0;
  75. con = Math.PI * this.R;
  76. xx = p.x / con;
  77. yy = p.y / con;
  78. xys = xx * xx + yy * yy;
  79. c1 = -Math.abs(yy) * (1 + xys);
  80. c2 = c1 - 2 * yy * yy + xx * xx;
  81. c3 = -2 * c1 + 1 + 2 * yy * yy + xys * xys;
  82. d = yy * yy / c3 + (2 * c2 * c2 * c2 / c3 / c3 / c3 - 9 * c1 * c2 / c3 / c3) / 27;
  83. a1 = (c1 - c2 * c2 / 3 / c3) / c3;
  84. m1 = 2 * Math.sqrt(-a1 / 3);
  85. con = ((3 * d) / a1) / m1;
  86. if (Math.abs(con) > 1) {
  87. if (con >= 0) {
  88. con = 1;
  89. }
  90. else {
  91. con = -1;
  92. }
  93. }
  94. th1 = Math.acos(con) / 3;
  95. if (p.y >= 0) {
  96. lat = (-m1 * Math.cos(th1 + Math.PI / 3) - c2 / 3 / c3) * Math.PI;
  97. }
  98. else {
  99. lat = -(-m1 * Math.cos(th1 + Math.PI / 3) - c2 / 3 / c3) * Math.PI;
  100. }
  101. if (Math.abs(xx) < EPSLN) {
  102. lon = this.long0;
  103. }
  104. else {
  105. lon = adjust_lon(this.long0 + Math.PI * (xys - 1 + Math.sqrt(1 + 2 * (xx * xx - yy * yy) + xys * xys)) / 2 / xx);
  106. }
  107. p.x = lon;
  108. p.y = lat;
  109. return p;
  110. };
  111. exports.names = ["Van_der_Grinten_I", "VanDerGrinten", "vandg"];