stere.js 5.2 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166
  1. var HALF_PI = Math.PI/2;
  2. var EPSLN = 1.0e-10;
  3. var sign = require('../common/sign');
  4. var msfnz = require('../common/msfnz');
  5. var tsfnz = require('../common/tsfnz');
  6. var phi2z = require('../common/phi2z');
  7. var adjust_lon = require('../common/adjust_lon');
  8. exports.ssfn_ = function(phit, sinphi, eccen) {
  9. sinphi *= eccen;
  10. return (Math.tan(0.5 * (HALF_PI + phit)) * Math.pow((1 - sinphi) / (1 + sinphi), 0.5 * eccen));
  11. };
  12. exports.init = function() {
  13. this.coslat0 = Math.cos(this.lat0);
  14. this.sinlat0 = Math.sin(this.lat0);
  15. if (this.sphere) {
  16. if (this.k0 === 1 && !isNaN(this.lat_ts) && Math.abs(this.coslat0) <= EPSLN) {
  17. this.k0 = 0.5 * (1 + sign(this.lat0) * Math.sin(this.lat_ts));
  18. }
  19. }
  20. else {
  21. if (Math.abs(this.coslat0) <= EPSLN) {
  22. if (this.lat0 > 0) {
  23. //North pole
  24. //trace('stere:north pole');
  25. this.con = 1;
  26. }
  27. else {
  28. //South pole
  29. //trace('stere:south pole');
  30. this.con = -1;
  31. }
  32. }
  33. this.cons = Math.sqrt(Math.pow(1 + this.e, 1 + this.e) * Math.pow(1 - this.e, 1 - this.e));
  34. if (this.k0 === 1 && !isNaN(this.lat_ts) && Math.abs(this.coslat0) <= EPSLN) {
  35. this.k0 = 0.5 * this.cons * msfnz(this.e, Math.sin(this.lat_ts), Math.cos(this.lat_ts)) / tsfnz(this.e, this.con * this.lat_ts, this.con * Math.sin(this.lat_ts));
  36. }
  37. this.ms1 = msfnz(this.e, this.sinlat0, this.coslat0);
  38. this.X0 = 2 * Math.atan(this.ssfn_(this.lat0, this.sinlat0, this.e)) - HALF_PI;
  39. this.cosX0 = Math.cos(this.X0);
  40. this.sinX0 = Math.sin(this.X0);
  41. }
  42. };
  43. // Stereographic forward equations--mapping lat,long to x,y
  44. exports.forward = function(p) {
  45. var lon = p.x;
  46. var lat = p.y;
  47. var sinlat = Math.sin(lat);
  48. var coslat = Math.cos(lat);
  49. var A, X, sinX, cosX, ts, rh;
  50. var dlon = adjust_lon(lon - this.long0);
  51. if (Math.abs(Math.abs(lon - this.long0) - Math.PI) <= EPSLN && Math.abs(lat + this.lat0) <= EPSLN) {
  52. //case of the origine point
  53. //trace('stere:this is the origin point');
  54. p.x = NaN;
  55. p.y = NaN;
  56. return p;
  57. }
  58. if (this.sphere) {
  59. //trace('stere:sphere case');
  60. A = 2 * this.k0 / (1 + this.sinlat0 * sinlat + this.coslat0 * coslat * Math.cos(dlon));
  61. p.x = this.a * A * coslat * Math.sin(dlon) + this.x0;
  62. p.y = this.a * A * (this.coslat0 * sinlat - this.sinlat0 * coslat * Math.cos(dlon)) + this.y0;
  63. return p;
  64. }
  65. else {
  66. X = 2 * Math.atan(this.ssfn_(lat, sinlat, this.e)) - HALF_PI;
  67. cosX = Math.cos(X);
  68. sinX = Math.sin(X);
  69. if (Math.abs(this.coslat0) <= EPSLN) {
  70. ts = tsfnz(this.e, lat * this.con, this.con * sinlat);
  71. rh = 2 * this.a * this.k0 * ts / this.cons;
  72. p.x = this.x0 + rh * Math.sin(lon - this.long0);
  73. p.y = this.y0 - this.con * rh * Math.cos(lon - this.long0);
  74. //trace(p.toString());
  75. return p;
  76. }
  77. else if (Math.abs(this.sinlat0) < EPSLN) {
  78. //Eq
  79. //trace('stere:equateur');
  80. A = 2 * this.a * this.k0 / (1 + cosX * Math.cos(dlon));
  81. p.y = A * sinX;
  82. }
  83. else {
  84. //other case
  85. //trace('stere:normal case');
  86. A = 2 * this.a * this.k0 * this.ms1 / (this.cosX0 * (1 + this.sinX0 * sinX + this.cosX0 * cosX * Math.cos(dlon)));
  87. p.y = A * (this.cosX0 * sinX - this.sinX0 * cosX * Math.cos(dlon)) + this.y0;
  88. }
  89. p.x = A * cosX * Math.sin(dlon) + this.x0;
  90. }
  91. //trace(p.toString());
  92. return p;
  93. };
  94. //* Stereographic inverse equations--mapping x,y to lat/long
  95. exports.inverse = function(p) {
  96. p.x -= this.x0;
  97. p.y -= this.y0;
  98. var lon, lat, ts, ce, Chi;
  99. var rh = Math.sqrt(p.x * p.x + p.y * p.y);
  100. if (this.sphere) {
  101. var c = 2 * Math.atan(rh / (0.5 * this.a * this.k0));
  102. lon = this.long0;
  103. lat = this.lat0;
  104. if (rh <= EPSLN) {
  105. p.x = lon;
  106. p.y = lat;
  107. return p;
  108. }
  109. lat = Math.asin(Math.cos(c) * this.sinlat0 + p.y * Math.sin(c) * this.coslat0 / rh);
  110. if (Math.abs(this.coslat0) < EPSLN) {
  111. if (this.lat0 > 0) {
  112. lon = adjust_lon(this.long0 + Math.atan2(p.x, - 1 * p.y));
  113. }
  114. else {
  115. lon = adjust_lon(this.long0 + Math.atan2(p.x, p.y));
  116. }
  117. }
  118. else {
  119. lon = adjust_lon(this.long0 + Math.atan2(p.x * Math.sin(c), rh * this.coslat0 * Math.cos(c) - p.y * this.sinlat0 * Math.sin(c)));
  120. }
  121. p.x = lon;
  122. p.y = lat;
  123. return p;
  124. }
  125. else {
  126. if (Math.abs(this.coslat0) <= EPSLN) {
  127. if (rh <= EPSLN) {
  128. lat = this.lat0;
  129. lon = this.long0;
  130. p.x = lon;
  131. p.y = lat;
  132. //trace(p.toString());
  133. return p;
  134. }
  135. p.x *= this.con;
  136. p.y *= this.con;
  137. ts = rh * this.cons / (2 * this.a * this.k0);
  138. lat = this.con * phi2z(this.e, ts);
  139. lon = this.con * adjust_lon(this.con * this.long0 + Math.atan2(p.x, - 1 * p.y));
  140. }
  141. else {
  142. ce = 2 * Math.atan(rh * this.cosX0 / (2 * this.a * this.k0 * this.ms1));
  143. lon = this.long0;
  144. if (rh <= EPSLN) {
  145. Chi = this.X0;
  146. }
  147. else {
  148. Chi = Math.asin(Math.cos(ce) * this.sinX0 + p.y * Math.sin(ce) * this.cosX0 / rh);
  149. lon = adjust_lon(this.long0 + Math.atan2(p.x * Math.sin(ce), rh * this.cosX0 * Math.cos(ce) - p.y * this.sinX0 * Math.sin(ce)));
  150. }
  151. lat = -1 * phi2z(this.e, Math.tan(0.5 * (HALF_PI + Chi)));
  152. }
  153. }
  154. p.x = lon;
  155. p.y = lat;
  156. //trace(p.toString());
  157. return p;
  158. };
  159. exports.names = ["stere"];