1 | /* |
---|
2 | proj4js.js -- Javascript reprojection library. |
---|
3 | |
---|
4 | Authors: Mike Adair madairATdmsolutions.ca |
---|
5 | Richard Greenwood richATgreenwoodmap.com |
---|
6 | Didier Richard didier.richardATign.fr |
---|
7 | Stephen Irons |
---|
8 | License: LGPL as per: http://www.gnu.org/copyleft/lesser.html |
---|
9 | Note: This program is an almost direct port of the C library |
---|
10 | Proj4. |
---|
11 | */ |
---|
12 | /* ====================================================================== |
---|
13 | proj4js.js |
---|
14 | ====================================================================== */ |
---|
15 | |
---|
16 | /* |
---|
17 | Author: Mike Adair madairATdmsolutions.ca |
---|
18 | Richard Greenwood rich@greenwoodmap.com |
---|
19 | License: LGPL as per: http://www.gnu.org/copyleft/lesser.html |
---|
20 | |
---|
21 | $Id: Proj.js 2956 2007-07-09 12:17:52Z steven $ |
---|
22 | */ |
---|
23 | |
---|
24 | /** |
---|
25 | * Namespace: Proj4js |
---|
26 | * |
---|
27 | * Proj4js is a JavaScript library to transform point coordinates from one |
---|
28 | * coordinate system to another, including datum transformations. |
---|
29 | * |
---|
30 | * This library is a port of both the Proj.4 and GCTCP C libraries to JavaScript. |
---|
31 | * Enabling these transformations in the browser allows geographic data stored |
---|
32 | * in different projections to be combined in browser-based web mapping |
---|
33 | * applications. |
---|
34 | * |
---|
35 | * Proj4js must have access to coordinate system initialization strings (which |
---|
36 | * are the same as for PROJ.4 command line). Thes can be included in your |
---|
37 | * application using a <script> tag or Proj4js can load CS initialization |
---|
38 | * strings from a local directory or a web service such as spatialreference.org. |
---|
39 | * |
---|
40 | * Similarly, Proj4js must have access to projection transform code. These can |
---|
41 | * be included individually using a <script> tag in your page, built into a |
---|
42 | * custom build of Proj4js or loaded dynamically at run-time. Using the |
---|
43 | * -combined and -compressed versions of Proj4js includes all projection class |
---|
44 | * code by default. |
---|
45 | * |
---|
46 | * Note that dynamic loading of defs and code happens ascynchrously, check the |
---|
47 | * Proj.readyToUse flag before using the Proj object. If the defs and code |
---|
48 | * required by your application are loaded through script tags, dynamic loading |
---|
49 | * is not required and the Proj object will be readyToUse on return from the |
---|
50 | * constructor. |
---|
51 | * |
---|
52 | * All coordinates are handled as points which have a .x and a .y property |
---|
53 | * which will be modified in place. |
---|
54 | * |
---|
55 | * Override Proj4js.reportError for output of alerts and warnings. |
---|
56 | * |
---|
57 | * See http://trac.osgeo.org/proj4js/wiki/UserGuide for full details. |
---|
58 | */ |
---|
59 | |
---|
60 | /** |
---|
61 | * Global namespace object for Proj4js library |
---|
62 | */ |
---|
63 | Proj4js = { |
---|
64 | |
---|
65 | /** |
---|
66 | * Property: defaultDatum |
---|
67 | * The datum to use when no others a specified |
---|
68 | */ |
---|
69 | defaultDatum: 'WGS84', //default datum |
---|
70 | |
---|
71 | /** |
---|
72 | * Method: transform(source, dest, point) |
---|
73 | * Transform a point coordinate from one map projection to another. This is |
---|
74 | * really the only public method you should need to use. |
---|
75 | * |
---|
76 | * Parameters: |
---|
77 | * source - {Proj4js.Proj} source map projection for the transformation |
---|
78 | * dest - {Proj4js.Proj} destination map projection for the transformation |
---|
79 | * point - {Object} point to transform, may be geodetic (long, lat) or |
---|
80 | * projected Cartesian (x,y), but should always have x,y properties. |
---|
81 | */ |
---|
82 | transform: function(source, dest, point) { |
---|
83 | if (!source.readyToUse) { |
---|
84 | this.reportError("Proj4js initialization for:"+source.srsCode+" not yet complete"); |
---|
85 | return point; |
---|
86 | } |
---|
87 | if (!dest.readyToUse) { |
---|
88 | this.reportError("Proj4js initialization for:"+dest.srsCode+" not yet complete"); |
---|
89 | return point; |
---|
90 | } |
---|
91 | |
---|
92 | // Workaround for Spherical Mercator |
---|
93 | if ((source.srsProjNumber =="900913" && dest.datumCode != "WGS84" && !dest.datum_params) || |
---|
94 | (dest.srsProjNumber == "900913" && source.datumCode != "WGS84" && !source.datum_params)) { |
---|
95 | var wgs84 = Proj4js.WGS84; |
---|
96 | this.transform(source, wgs84, point); |
---|
97 | source = wgs84; |
---|
98 | } |
---|
99 | |
---|
100 | // DGR, 2010/11/12 |
---|
101 | if (source.axis!="enu") { |
---|
102 | this.adjust_axis(source,false,point); |
---|
103 | } |
---|
104 | |
---|
105 | // Transform source points to long/lat, if they aren't already. |
---|
106 | if ( source.projName=="longlat") { |
---|
107 | point.x *= Proj4js.common.D2R; // convert degrees to radians |
---|
108 | point.y *= Proj4js.common.D2R; |
---|
109 | } else { |
---|
110 | if (source.to_meter) { |
---|
111 | point.x *= source.to_meter; |
---|
112 | point.y *= source.to_meter; |
---|
113 | } |
---|
114 | source.inverse(point); // Convert Cartesian to longlat |
---|
115 | } |
---|
116 | |
---|
117 | // Adjust for the prime meridian if necessary |
---|
118 | if (source.from_greenwich) { |
---|
119 | point.x += source.from_greenwich; |
---|
120 | } |
---|
121 | |
---|
122 | // Convert datums if needed, and if possible. |
---|
123 | point = this.datum_transform( source.datum, dest.datum, point ); |
---|
124 | |
---|
125 | // Adjust for the prime meridian if necessary |
---|
126 | if (dest.from_greenwich) { |
---|
127 | point.x -= dest.from_greenwich; |
---|
128 | } |
---|
129 | |
---|
130 | if( dest.projName=="longlat" ) { |
---|
131 | // convert radians to decimal degrees |
---|
132 | point.x *= Proj4js.common.R2D; |
---|
133 | point.y *= Proj4js.common.R2D; |
---|
134 | } else { // else project |
---|
135 | dest.forward(point); |
---|
136 | if (dest.to_meter) { |
---|
137 | point.x /= dest.to_meter; |
---|
138 | point.y /= dest.to_meter; |
---|
139 | } |
---|
140 | } |
---|
141 | |
---|
142 | // DGR, 2010/11/12 |
---|
143 | if (dest.axis!="enu") { |
---|
144 | this.adjust_axis(dest,true,point); |
---|
145 | } |
---|
146 | |
---|
147 | return point; |
---|
148 | }, // transform() |
---|
149 | |
---|
150 | /** datum_transform() |
---|
151 | source coordinate system definition, |
---|
152 | destination coordinate system definition, |
---|
153 | point to transform in geodetic coordinates (long, lat, height) |
---|
154 | */ |
---|
155 | datum_transform : function( source, dest, point ) { |
---|
156 | |
---|
157 | // Short cut if the datums are identical. |
---|
158 | if( source.compare_datums( dest ) ) { |
---|
159 | return point; // in this case, zero is sucess, |
---|
160 | // whereas cs_compare_datums returns 1 to indicate TRUE |
---|
161 | // confusing, should fix this |
---|
162 | } |
---|
163 | |
---|
164 | // Explicitly skip datum transform by setting 'datum=none' as parameter for either source or dest |
---|
165 | if( source.datum_type == Proj4js.common.PJD_NODATUM |
---|
166 | || dest.datum_type == Proj4js.common.PJD_NODATUM) { |
---|
167 | return point; |
---|
168 | } |
---|
169 | |
---|
170 | // If this datum requires grid shifts, then apply it to geodetic coordinates. |
---|
171 | if( source.datum_type == Proj4js.common.PJD_GRIDSHIFT ) |
---|
172 | { |
---|
173 | alert("ERROR: Grid shift transformations are not implemented yet."); |
---|
174 | /* |
---|
175 | pj_apply_gridshift( pj_param(source.params,"snadgrids").s, 0, |
---|
176 | point_count, point_offset, x, y, z ); |
---|
177 | CHECK_RETURN; |
---|
178 | |
---|
179 | src_a = SRS_WGS84_SEMIMAJOR; |
---|
180 | src_es = 0.006694379990; |
---|
181 | */ |
---|
182 | } |
---|
183 | |
---|
184 | if( dest.datum_type == Proj4js.common.PJD_GRIDSHIFT ) |
---|
185 | { |
---|
186 | alert("ERROR: Grid shift transformations are not implemented yet."); |
---|
187 | /* |
---|
188 | dst_a = ; |
---|
189 | dst_es = 0.006694379990; |
---|
190 | */ |
---|
191 | } |
---|
192 | |
---|
193 | // Do we need to go through geocentric coordinates? |
---|
194 | if( source.es != dest.es || source.a != dest.a |
---|
195 | || source.datum_type == Proj4js.common.PJD_3PARAM |
---|
196 | || source.datum_type == Proj4js.common.PJD_7PARAM |
---|
197 | || dest.datum_type == Proj4js.common.PJD_3PARAM |
---|
198 | || dest.datum_type == Proj4js.common.PJD_7PARAM) |
---|
199 | { |
---|
200 | |
---|
201 | // Convert to geocentric coordinates. |
---|
202 | source.geodetic_to_geocentric( point ); |
---|
203 | // CHECK_RETURN; |
---|
204 | |
---|
205 | // Convert between datums |
---|
206 | if( source.datum_type == Proj4js.common.PJD_3PARAM || source.datum_type == Proj4js.common.PJD_7PARAM ) { |
---|
207 | source.geocentric_to_wgs84(point); |
---|
208 | // CHECK_RETURN; |
---|
209 | } |
---|
210 | |
---|
211 | if( dest.datum_type == Proj4js.common.PJD_3PARAM || dest.datum_type == Proj4js.common.PJD_7PARAM ) { |
---|
212 | dest.geocentric_from_wgs84(point); |
---|
213 | // CHECK_RETURN; |
---|
214 | } |
---|
215 | |
---|
216 | // Convert back to geodetic coordinates |
---|
217 | dest.geocentric_to_geodetic( point ); |
---|
218 | // CHECK_RETURN; |
---|
219 | } |
---|
220 | |
---|
221 | // Apply grid shift to destination if required |
---|
222 | if( dest.datum_type == Proj4js.common.PJD_GRIDSHIFT ) |
---|
223 | { |
---|
224 | alert("ERROR: Grid shift transformations are not implemented yet."); |
---|
225 | // pj_apply_gridshift( pj_param(dest.params,"snadgrids").s, 1, point); |
---|
226 | // CHECK_RETURN; |
---|
227 | } |
---|
228 | return point; |
---|
229 | }, // cs_datum_transform |
---|
230 | |
---|
231 | /** |
---|
232 | * Function: adjust_axis |
---|
233 | * Normalize or de-normalized the x/y/z axes. The normal form is "enu" |
---|
234 | * (easting, northing, up). |
---|
235 | * Parameters: |
---|
236 | * crs {Proj4js.Proj} the coordinate reference system |
---|
237 | * denorm {Boolean} when false, normalize |
---|
238 | * point {Object} the coordinates to adjust |
---|
239 | */ |
---|
240 | adjust_axis: function(crs, denorm, point) { |
---|
241 | var xin= point.x, yin= point.y, zin= point.z || 0.0; |
---|
242 | var v, t; |
---|
243 | for (var i= 0; i<3; i++) { |
---|
244 | if (denorm && i==2 && point.z===undefined) { continue; } |
---|
245 | if (i==0) { v= xin; t= 'x'; } |
---|
246 | else if (i==1) { v= yin; t= 'y'; } |
---|
247 | else { v= zin; t= 'z'; } |
---|
248 | switch(crs.axis[i]) { |
---|
249 | case 'e': |
---|
250 | point[t]= v; |
---|
251 | break; |
---|
252 | case 'w': |
---|
253 | point[t]= -v; |
---|
254 | break; |
---|
255 | case 'n': |
---|
256 | point[t]= v; |
---|
257 | break; |
---|
258 | case 's': |
---|
259 | point[t]= -v; |
---|
260 | break; |
---|
261 | case 'u': |
---|
262 | if (point[t]!==undefined) { point.z= v; } |
---|
263 | break; |
---|
264 | case 'd': |
---|
265 | if (point[t]!==undefined) { point.z= -v; } |
---|
266 | break; |
---|
267 | default : |
---|
268 | alert("ERROR: unknow axis ("+crs.axis[i]+") - check definition of "+src.projName); |
---|
269 | return null; |
---|
270 | } |
---|
271 | } |
---|
272 | return point; |
---|
273 | }, |
---|
274 | |
---|
275 | /** |
---|
276 | * Function: reportError |
---|
277 | * An internal method to report errors back to user. |
---|
278 | * Override this in applications to report error messages or throw exceptions. |
---|
279 | */ |
---|
280 | reportError: function(msg) { |
---|
281 | //console.log(msg); |
---|
282 | }, |
---|
283 | |
---|
284 | /** |
---|
285 | * |
---|
286 | * Title: Private Methods |
---|
287 | * The following properties and methods are intended for internal use only. |
---|
288 | * |
---|
289 | * This is a minimal implementation of JavaScript inheritance methods so that |
---|
290 | * Proj4js can be used as a stand-alone library. |
---|
291 | * These are copies of the equivalent OpenLayers methods at v2.7 |
---|
292 | */ |
---|
293 | |
---|
294 | /** |
---|
295 | * Function: extend |
---|
296 | * Copy all properties of a source object to a destination object. Modifies |
---|
297 | * the passed in destination object. Any properties on the source object |
---|
298 | * that are set to undefined will not be (re)set on the destination object. |
---|
299 | * |
---|
300 | * Parameters: |
---|
301 | * destination - {Object} The object that will be modified |
---|
302 | * source - {Object} The object with properties to be set on the destination |
---|
303 | * |
---|
304 | * Returns: |
---|
305 | * {Object} The destination object. |
---|
306 | */ |
---|
307 | extend: function(destination, source) { |
---|
308 | destination = destination || {}; |
---|
309 | if(source) { |
---|
310 | for(var property in source) { |
---|
311 | var value = source[property]; |
---|
312 | if(value !== undefined) { |
---|
313 | destination[property] = value; |
---|
314 | } |
---|
315 | } |
---|
316 | } |
---|
317 | return destination; |
---|
318 | }, |
---|
319 | |
---|
320 | /** |
---|
321 | * Constructor: Class |
---|
322 | * Base class used to construct all other classes. Includes support for |
---|
323 | * multiple inheritance. |
---|
324 | * |
---|
325 | */ |
---|
326 | Class: function() { |
---|
327 | var Class = function() { |
---|
328 | this.initialize.apply(this, arguments); |
---|
329 | }; |
---|
330 | |
---|
331 | var extended = {}; |
---|
332 | var parent; |
---|
333 | for(var i=0; i<arguments.length; ++i) { |
---|
334 | if(typeof arguments[i] == "function") { |
---|
335 | // get the prototype of the superclass |
---|
336 | parent = arguments[i].prototype; |
---|
337 | } else { |
---|
338 | // in this case we're extending with the prototype |
---|
339 | parent = arguments[i]; |
---|
340 | } |
---|
341 | Proj4js.extend(extended, parent); |
---|
342 | } |
---|
343 | Class.prototype = extended; |
---|
344 | |
---|
345 | return Class; |
---|
346 | }, |
---|
347 | |
---|
348 | /** |
---|
349 | * Function: bind |
---|
350 | * Bind a function to an object. Method to easily create closures with |
---|
351 | * 'this' altered. |
---|
352 | * |
---|
353 | * Parameters: |
---|
354 | * func - {Function} Input function. |
---|
355 | * object - {Object} The object to bind to the input function (as this). |
---|
356 | * |
---|
357 | * Returns: |
---|
358 | * {Function} A closure with 'this' set to the passed in object. |
---|
359 | */ |
---|
360 | bind: function(func, object) { |
---|
361 | // create a reference to all arguments past the second one |
---|
362 | var args = Array.prototype.slice.apply(arguments, [2]); |
---|
363 | return function() { |
---|
364 | // Push on any additional arguments from the actual function call. |
---|
365 | // These will come after those sent to the bind call. |
---|
366 | var newArgs = args.concat( |
---|
367 | Array.prototype.slice.apply(arguments, [0]) |
---|
368 | ); |
---|
369 | return func.apply(object, newArgs); |
---|
370 | }; |
---|
371 | }, |
---|
372 | |
---|
373 | /** |
---|
374 | * The following properties and methods handle dynamic loading of JSON objects. |
---|
375 | */ |
---|
376 | |
---|
377 | /** |
---|
378 | * Property: scriptName |
---|
379 | * {String} The filename of this script without any path. |
---|
380 | */ |
---|
381 | scriptName: "proj4js-combined.js", |
---|
382 | |
---|
383 | /** |
---|
384 | * Property: defsLookupService |
---|
385 | * AJAX service to retreive projection definition parameters from |
---|
386 | */ |
---|
387 | defsLookupService: 'http://spatialreference.org/ref', |
---|
388 | |
---|
389 | /** |
---|
390 | * Property: libPath |
---|
391 | * internal: http server path to library code. |
---|
392 | */ |
---|
393 | libPath: null, |
---|
394 | |
---|
395 | /** |
---|
396 | * Function: getScriptLocation |
---|
397 | * Return the path to this script. |
---|
398 | * |
---|
399 | * Returns: |
---|
400 | * Path to this script |
---|
401 | */ |
---|
402 | getScriptLocation: function () { |
---|
403 | if (this.libPath) return this.libPath; |
---|
404 | var scriptName = this.scriptName; |
---|
405 | var scriptNameLen = scriptName.length; |
---|
406 | |
---|
407 | var scripts = document.getElementsByTagName('script'); |
---|
408 | for (var i = 0; i < scripts.length; i++) { |
---|
409 | var src = scripts[i].getAttribute('src'); |
---|
410 | if (src) { |
---|
411 | var index = src.lastIndexOf(scriptName); |
---|
412 | // is it found, at the end of the URL? |
---|
413 | if ((index > -1) && (index + scriptNameLen == src.length)) { |
---|
414 | this.libPath = src.slice(0, -scriptNameLen); |
---|
415 | break; |
---|
416 | } |
---|
417 | } |
---|
418 | } |
---|
419 | return this.libPath||""; |
---|
420 | }, |
---|
421 | |
---|
422 | /** |
---|
423 | * Function: loadScript |
---|
424 | * Load a JS file from a URL into a <script> tag in the page. |
---|
425 | * |
---|
426 | * Parameters: |
---|
427 | * url - {String} The URL containing the script to load |
---|
428 | * onload - {Function} A method to be executed when the script loads successfully |
---|
429 | * onfail - {Function} A method to be executed when there is an error loading the script |
---|
430 | * loadCheck - {Function} A boolean method that checks to see if the script |
---|
431 | * has loaded. Typically this just checks for the existance of |
---|
432 | * an object in the file just loaded. |
---|
433 | */ |
---|
434 | loadScript: function(url, onload, onfail, loadCheck) { |
---|
435 | var script = document.createElement('script'); |
---|
436 | script.defer = false; |
---|
437 | script.type = "text/javascript"; |
---|
438 | script.id = url; |
---|
439 | script.src = url; |
---|
440 | script.onload = onload; |
---|
441 | script.onerror = onfail; |
---|
442 | script.loadCheck = loadCheck; |
---|
443 | if (/MSIE/.test(navigator.userAgent)) { |
---|
444 | script.onreadystatechange = this.checkReadyState; |
---|
445 | } |
---|
446 | document.getElementsByTagName('head')[0].appendChild(script); |
---|
447 | }, |
---|
448 | |
---|
449 | /** |
---|
450 | * Function: checkReadyState |
---|
451 | * IE workaround since there is no onerror handler. Calls the user defined |
---|
452 | * loadCheck method to determine if the script is loaded. |
---|
453 | * |
---|
454 | */ |
---|
455 | checkReadyState: function() { |
---|
456 | if (this.readyState == 'loaded') { |
---|
457 | if (!this.loadCheck()) { |
---|
458 | this.onerror(); |
---|
459 | } else { |
---|
460 | this.onload(); |
---|
461 | } |
---|
462 | } |
---|
463 | } |
---|
464 | }; |
---|
465 | |
---|
466 | /** |
---|
467 | * Class: Proj4js.Proj |
---|
468 | * |
---|
469 | * Proj objects provide transformation methods for point coordinates |
---|
470 | * between geodetic latitude/longitude and a projected coordinate system. |
---|
471 | * once they have been initialized with a projection code. |
---|
472 | * |
---|
473 | * Initialization of Proj objects is with a projection code, usually EPSG codes, |
---|
474 | * which is the key that will be used with the Proj4js.defs array. |
---|
475 | * |
---|
476 | * The code passed in will be stripped of colons and converted to uppercase |
---|
477 | * to locate projection definition files. |
---|
478 | * |
---|
479 | * A projection object has properties for units and title strings. |
---|
480 | */ |
---|
481 | Proj4js.Proj = Proj4js.Class({ |
---|
482 | |
---|
483 | /** |
---|
484 | * Property: readyToUse |
---|
485 | * Flag to indicate if initialization is complete for this Proj object |
---|
486 | */ |
---|
487 | readyToUse: false, |
---|
488 | |
---|
489 | /** |
---|
490 | * Property: title |
---|
491 | * The title to describe the projection |
---|
492 | */ |
---|
493 | title: null, |
---|
494 | |
---|
495 | /** |
---|
496 | * Property: projName |
---|
497 | * The projection class for this projection, e.g. lcc (lambert conformal conic, |
---|
498 | * or merc for mercator). These are exactly equivalent to their Proj4 |
---|
499 | * counterparts. |
---|
500 | */ |
---|
501 | projName: null, |
---|
502 | /** |
---|
503 | * Property: units |
---|
504 | * The units of the projection. Values include 'm' and 'degrees' |
---|
505 | */ |
---|
506 | units: null, |
---|
507 | /** |
---|
508 | * Property: datum |
---|
509 | * The datum specified for the projection |
---|
510 | */ |
---|
511 | datum: null, |
---|
512 | /** |
---|
513 | * Property: x0 |
---|
514 | * The x coordinate origin |
---|
515 | */ |
---|
516 | x0: 0, |
---|
517 | /** |
---|
518 | * Property: y0 |
---|
519 | * The y coordinate origin |
---|
520 | */ |
---|
521 | y0: 0, |
---|
522 | /** |
---|
523 | * Property: localCS |
---|
524 | * Flag to indicate if the projection is a local one in which no transforms |
---|
525 | * are required. |
---|
526 | */ |
---|
527 | localCS: false, |
---|
528 | |
---|
529 | /** |
---|
530 | * Property: queue |
---|
531 | * Buffer (FIFO) to hold callbacks waiting to be called when projection loaded. |
---|
532 | */ |
---|
533 | queue: null, |
---|
534 | |
---|
535 | /** |
---|
536 | * Constructor: initialize |
---|
537 | * Constructor for Proj4js.Proj objects |
---|
538 | * |
---|
539 | * Parameters: |
---|
540 | * srsCode - a code for map projection definition parameters. These are usually |
---|
541 | * (but not always) EPSG codes. |
---|
542 | */ |
---|
543 | initialize: function(srsCode, callback) { |
---|
544 | this.srsCodeInput = srsCode; |
---|
545 | |
---|
546 | //Register callbacks prior to attempting to process definition |
---|
547 | this.queue = []; |
---|
548 | if( callback ){ |
---|
549 | this.queue.push( callback ); |
---|
550 | } |
---|
551 | |
---|
552 | //check to see if this is a WKT string |
---|
553 | if ((srsCode.indexOf('GEOGCS') >= 0) || |
---|
554 | (srsCode.indexOf('GEOCCS') >= 0) || |
---|
555 | (srsCode.indexOf('PROJCS') >= 0) || |
---|
556 | (srsCode.indexOf('LOCAL_CS') >= 0)) { |
---|
557 | this.parseWKT(srsCode); |
---|
558 | this.deriveConstants(); |
---|
559 | this.loadProjCode(this.projName); |
---|
560 | return; |
---|
561 | } |
---|
562 | |
---|
563 | // DGR 2008-08-03 : support urn and url |
---|
564 | if (srsCode.indexOf('urn:') == 0) { |
---|
565 | //urn:ORIGINATOR:def:crs:CODESPACE:VERSION:ID |
---|
566 | var urn = srsCode.split(':'); |
---|
567 | if ((urn[1] == 'ogc' || urn[1] =='x-ogc') && |
---|
568 | (urn[2] =='def') && |
---|
569 | (urn[3] =='crs')) { |
---|
570 | srsCode = urn[4]+':'+urn[urn.length-1]; |
---|
571 | } |
---|
572 | } else if (srsCode.indexOf('http://') == 0) { |
---|
573 | //url#ID |
---|
574 | var url = srsCode.split('#'); |
---|
575 | if (url[0].match(/epsg.org/)) { |
---|
576 | // http://www.epsg.org/# |
---|
577 | srsCode = 'EPSG:'+url[1]; |
---|
578 | } else if (url[0].match(/RIG.xml/)) { |
---|
579 | //http://librairies.ign.fr/geoportail/resources/RIG.xml# |
---|
580 | //http://interop.ign.fr/registers/ign/RIG.xml# |
---|
581 | srsCode = 'IGNF:'+url[1]; |
---|
582 | } |
---|
583 | } |
---|
584 | this.srsCode = srsCode.toUpperCase(); |
---|
585 | if (this.srsCode.indexOf("EPSG") == 0) { |
---|
586 | this.srsCode = this.srsCode; |
---|
587 | this.srsAuth = 'epsg'; |
---|
588 | this.srsProjNumber = this.srsCode.substring(5); |
---|
589 | // DGR 2007-11-20 : authority IGNF |
---|
590 | } else if (this.srsCode.indexOf("IGNF") == 0) { |
---|
591 | this.srsCode = this.srsCode; |
---|
592 | this.srsAuth = 'IGNF'; |
---|
593 | this.srsProjNumber = this.srsCode.substring(5); |
---|
594 | // DGR 2008-06-19 : pseudo-authority CRS for WMS |
---|
595 | } else if (this.srsCode.indexOf("CRS") == 0) { |
---|
596 | this.srsCode = this.srsCode; |
---|
597 | this.srsAuth = 'CRS'; |
---|
598 | this.srsProjNumber = this.srsCode.substring(4); |
---|
599 | } else { |
---|
600 | this.srsAuth = ''; |
---|
601 | this.srsProjNumber = this.srsCode; |
---|
602 | } |
---|
603 | |
---|
604 | this.loadProjDefinition(); |
---|
605 | }, |
---|
606 | |
---|
607 | /** |
---|
608 | * Function: loadProjDefinition |
---|
609 | * Loads the coordinate system initialization string if required. |
---|
610 | * Note that dynamic loading happens asynchronously so an application must |
---|
611 | * wait for the readyToUse property is set to true. |
---|
612 | * To prevent dynamic loading, include the defs through a script tag in |
---|
613 | * your application. |
---|
614 | * |
---|
615 | */ |
---|
616 | loadProjDefinition: function() { |
---|
617 | //check in memory |
---|
618 | if (Proj4js.defs[this.srsCode]) { |
---|
619 | this.defsLoaded(); |
---|
620 | return; |
---|
621 | } |
---|
622 | |
---|
623 | //else check for def on the server |
---|
624 | var url = Proj4js.getScriptLocation() + 'defs/' + this.srsAuth.toUpperCase() + this.srsProjNumber + '.js'; |
---|
625 | Proj4js.loadScript(url, |
---|
626 | Proj4js.bind(this.defsLoaded, this), |
---|
627 | Proj4js.bind(this.loadFromService, this), |
---|
628 | Proj4js.bind(this.checkDefsLoaded, this) ); |
---|
629 | }, |
---|
630 | |
---|
631 | /** |
---|
632 | * Function: loadFromService |
---|
633 | * Creates the REST URL for loading the definition from a web service and |
---|
634 | * loads it. |
---|
635 | * |
---|
636 | */ |
---|
637 | loadFromService: function() { |
---|
638 | //else load from web service |
---|
639 | var url = Proj4js.defsLookupService +'/' + this.srsAuth +'/'+ this.srsProjNumber + '/proj4js/'; |
---|
640 | Proj4js.loadScript(url, |
---|
641 | Proj4js.bind(this.defsLoaded, this), |
---|
642 | Proj4js.bind(this.defsFailed, this), |
---|
643 | Proj4js.bind(this.checkDefsLoaded, this) ); |
---|
644 | }, |
---|
645 | |
---|
646 | /** |
---|
647 | * Function: defsLoaded |
---|
648 | * Continues the Proj object initilization once the def file is loaded |
---|
649 | * |
---|
650 | */ |
---|
651 | defsLoaded: function() { |
---|
652 | this.parseDefs(); |
---|
653 | this.loadProjCode(this.projName); |
---|
654 | }, |
---|
655 | |
---|
656 | /** |
---|
657 | * Function: checkDefsLoaded |
---|
658 | * This is the loadCheck method to see if the def object exists |
---|
659 | * |
---|
660 | */ |
---|
661 | checkDefsLoaded: function() { |
---|
662 | if (Proj4js.defs[this.srsCode]) { |
---|
663 | return true; |
---|
664 | } else { |
---|
665 | return false; |
---|
666 | } |
---|
667 | }, |
---|
668 | |
---|
669 | /** |
---|
670 | * Function: defsFailed |
---|
671 | * Report an error in loading the defs file, but continue on using WGS84 |
---|
672 | * |
---|
673 | */ |
---|
674 | defsFailed: function() { |
---|
675 | Proj4js.reportError('failed to load projection definition for: '+this.srsCode); |
---|
676 | Proj4js.defs[this.srsCode] = Proj4js.defs['WGS84']; //set it to something so it can at least continue |
---|
677 | this.defsLoaded(); |
---|
678 | }, |
---|
679 | |
---|
680 | /** |
---|
681 | * Function: loadProjCode |
---|
682 | * Loads projection class code dynamically if required. |
---|
683 | * Projection code may be included either through a script tag or in |
---|
684 | * a built version of proj4js |
---|
685 | * |
---|
686 | */ |
---|
687 | loadProjCode: function(projName) { |
---|
688 | if (Proj4js.Proj[projName]) { |
---|
689 | this.initTransforms(); |
---|
690 | return; |
---|
691 | } |
---|
692 | |
---|
693 | //the URL for the projection code |
---|
694 | var url = Proj4js.getScriptLocation() + 'projCode/' + projName + '.js'; |
---|
695 | Proj4js.loadScript(url, |
---|
696 | Proj4js.bind(this.loadProjCodeSuccess, this, projName), |
---|
697 | Proj4js.bind(this.loadProjCodeFailure, this, projName), |
---|
698 | Proj4js.bind(this.checkCodeLoaded, this, projName) ); |
---|
699 | }, |
---|
700 | |
---|
701 | /** |
---|
702 | * Function: loadProjCodeSuccess |
---|
703 | * Loads any proj dependencies or continue on to final initialization. |
---|
704 | * |
---|
705 | */ |
---|
706 | loadProjCodeSuccess: function(projName) { |
---|
707 | if (Proj4js.Proj[projName].dependsOn){ |
---|
708 | this.loadProjCode(Proj4js.Proj[projName].dependsOn); |
---|
709 | } else { |
---|
710 | this.initTransforms(); |
---|
711 | } |
---|
712 | }, |
---|
713 | |
---|
714 | /** |
---|
715 | * Function: defsFailed |
---|
716 | * Report an error in loading the proj file. Initialization of the Proj |
---|
717 | * object has failed and the readyToUse flag will never be set. |
---|
718 | * |
---|
719 | */ |
---|
720 | loadProjCodeFailure: function(projName) { |
---|
721 | Proj4js.reportError("failed to find projection file for: " + projName); |
---|
722 | //TBD initialize with identity transforms so proj will still work? |
---|
723 | }, |
---|
724 | |
---|
725 | /** |
---|
726 | * Function: checkCodeLoaded |
---|
727 | * This is the loadCheck method to see if the projection code is loaded |
---|
728 | * |
---|
729 | */ |
---|
730 | checkCodeLoaded: function(projName) { |
---|
731 | if (Proj4js.Proj[projName]) { |
---|
732 | return true; |
---|
733 | } else { |
---|
734 | return false; |
---|
735 | } |
---|
736 | }, |
---|
737 | |
---|
738 | /** |
---|
739 | * Function: initTransforms |
---|
740 | * Finalize the initialization of the Proj object |
---|
741 | * |
---|
742 | */ |
---|
743 | initTransforms: function() { |
---|
744 | Proj4js.extend(this, Proj4js.Proj[this.projName]); |
---|
745 | this.init(); |
---|
746 | this.readyToUse = true; |
---|
747 | if( this.queue ) { |
---|
748 | var item; |
---|
749 | while( (item = this.queue.shift()) ) { |
---|
750 | item.call( this, this ); |
---|
751 | } |
---|
752 | } |
---|
753 | }, |
---|
754 | |
---|
755 | /** |
---|
756 | * Function: parseWKT |
---|
757 | * Parses a WKT string to get initialization parameters |
---|
758 | * |
---|
759 | */ |
---|
760 | wktRE: /^(\w+)\[(.*)\]$/, |
---|
761 | parseWKT: function(wkt) { |
---|
762 | var wktMatch = wkt.match(this.wktRE); |
---|
763 | if (!wktMatch) return; |
---|
764 | var wktObject = wktMatch[1]; |
---|
765 | var wktContent = wktMatch[2]; |
---|
766 | var wktTemp = wktContent.split(","); |
---|
767 | var wktName; |
---|
768 | if (wktObject.toUpperCase() == "TOWGS84") { |
---|
769 | wktName = wktObject; //no name supplied for the TOWGS84 array |
---|
770 | } else { |
---|
771 | wktName = wktTemp.shift(); |
---|
772 | } |
---|
773 | wktName = wktName.replace(/^\"/,""); |
---|
774 | wktName = wktName.replace(/\"$/,""); |
---|
775 | |
---|
776 | /* |
---|
777 | wktContent = wktTemp.join(","); |
---|
778 | var wktArray = wktContent.split("],"); |
---|
779 | for (var i=0; i<wktArray.length-1; ++i) { |
---|
780 | wktArray[i] += "]"; |
---|
781 | } |
---|
782 | */ |
---|
783 | |
---|
784 | var wktArray = new Array(); |
---|
785 | var bkCount = 0; |
---|
786 | var obj = ""; |
---|
787 | for (var i=0; i<wktTemp.length; ++i) { |
---|
788 | var token = wktTemp[i]; |
---|
789 | for (var j=0; j<token.length; ++j) { |
---|
790 | if (token.charAt(j) == "[") ++bkCount; |
---|
791 | if (token.charAt(j) == "]") --bkCount; |
---|
792 | } |
---|
793 | obj += token; |
---|
794 | if (bkCount === 0) { |
---|
795 | wktArray.push(obj); |
---|
796 | obj = ""; |
---|
797 | } else { |
---|
798 | obj += ","; |
---|
799 | } |
---|
800 | } |
---|
801 | |
---|
802 | //do something based on the type of the wktObject being parsed |
---|
803 | //add in variations in the spelling as required |
---|
804 | switch (wktObject) { |
---|
805 | case 'LOCAL_CS': |
---|
806 | this.projName = 'identity' |
---|
807 | this.localCS = true; |
---|
808 | this.srsCode = wktName; |
---|
809 | break; |
---|
810 | case 'GEOGCS': |
---|
811 | this.projName = 'longlat' |
---|
812 | this.geocsCode = wktName; |
---|
813 | if (!this.srsCode) this.srsCode = wktName; |
---|
814 | break; |
---|
815 | case 'PROJCS': |
---|
816 | this.srsCode = wktName; |
---|
817 | break; |
---|
818 | case 'GEOCCS': |
---|
819 | break; |
---|
820 | case 'PROJECTION': |
---|
821 | this.projName = Proj4js.wktProjections[wktName] |
---|
822 | break; |
---|
823 | case 'DATUM': |
---|
824 | this.datumName = wktName; |
---|
825 | break; |
---|
826 | case 'LOCAL_DATUM': |
---|
827 | this.datumCode = 'none'; |
---|
828 | break; |
---|
829 | case 'SPHEROID': |
---|
830 | this.ellps = wktName; |
---|
831 | this.a = parseFloat(wktArray.shift()); |
---|
832 | this.rf = parseFloat(wktArray.shift()); |
---|
833 | break; |
---|
834 | case 'PRIMEM': |
---|
835 | this.from_greenwich = parseFloat(wktArray.shift()); //to radians? |
---|
836 | break; |
---|
837 | case 'UNIT': |
---|
838 | this.units = wktName; |
---|
839 | this.unitsPerMeter = parseFloat(wktArray.shift()); |
---|
840 | break; |
---|
841 | case 'PARAMETER': |
---|
842 | var name = wktName.toLowerCase(); |
---|
843 | var value = parseFloat(wktArray.shift()); |
---|
844 | //there may be many variations on the wktName values, add in case |
---|
845 | //statements as required |
---|
846 | switch (name) { |
---|
847 | case 'false_easting': |
---|
848 | this.x0 = value; |
---|
849 | break; |
---|
850 | case 'false_northing': |
---|
851 | this.y0 = value; |
---|
852 | break; |
---|
853 | case 'scale_factor': |
---|
854 | this.k0 = value; |
---|
855 | break; |
---|
856 | case 'central_meridian': |
---|
857 | this.long0 = value*Proj4js.common.D2R; |
---|
858 | break; |
---|
859 | case 'latitude_of_origin': |
---|
860 | this.lat0 = value*Proj4js.common.D2R; |
---|
861 | break; |
---|
862 | case 'more_here': |
---|
863 | break; |
---|
864 | default: |
---|
865 | break; |
---|
866 | } |
---|
867 | break; |
---|
868 | case 'TOWGS84': |
---|
869 | this.datum_params = wktArray; |
---|
870 | break; |
---|
871 | //DGR 2010-11-12: AXIS |
---|
872 | case 'AXIS': |
---|
873 | var name= wktName.toLowerCase(); |
---|
874 | var value= wktArray.shift(); |
---|
875 | switch (value) { |
---|
876 | case 'EAST' : value= 'e'; break; |
---|
877 | case 'WEST' : value= 'w'; break; |
---|
878 | case 'NORTH': value= 'n'; break; |
---|
879 | case 'SOUTH': value= 's'; break; |
---|
880 | case 'UP' : value= 'u'; break; |
---|
881 | case 'DOWN' : value= 'd'; break; |
---|
882 | case 'OTHER': |
---|
883 | default : value= ' '; break;//FIXME |
---|
884 | } |
---|
885 | if (!this.axis) { this.axis= "enu"; } |
---|
886 | switch(name) { |
---|
887 | case 'X': this.axis= value + this.axis.substr(1,2); break; |
---|
888 | case 'Y': this.axis= this.axis.substr(0,1) + value + this.axis.substr(2,1); break; |
---|
889 | case 'Z': this.axis= this.axis.substr(0,2) + value ; break; |
---|
890 | default : break; |
---|
891 | } |
---|
892 | case 'MORE_HERE': |
---|
893 | break; |
---|
894 | default: |
---|
895 | break; |
---|
896 | } |
---|
897 | for (var i=0; i<wktArray.length; ++i) { |
---|
898 | this.parseWKT(wktArray[i]); |
---|
899 | } |
---|
900 | }, |
---|
901 | |
---|
902 | /** |
---|
903 | * Function: parseDefs |
---|
904 | * Parses the PROJ.4 initialization string and sets the associated properties. |
---|
905 | * |
---|
906 | */ |
---|
907 | parseDefs: function() { |
---|
908 | this.defData = Proj4js.defs[this.srsCode]; |
---|
909 | var paramName, paramVal; |
---|
910 | if (!this.defData) { |
---|
911 | return; |
---|
912 | } |
---|
913 | var paramArray=this.defData.split("+"); |
---|
914 | |
---|
915 | for (var prop=0; prop<paramArray.length; prop++) { |
---|
916 | var property = paramArray[prop].split("="); |
---|
917 | paramName = property[0].toLowerCase(); |
---|
918 | paramVal = property[1]; |
---|
919 | |
---|
920 | switch (paramName.replace(/\s/gi,"")) { // trim out spaces |
---|
921 | case "": break; // throw away nameless parameter |
---|
922 | case "title": this.title = paramVal; break; |
---|
923 | case "proj": this.projName = paramVal.replace(/\s/gi,""); break; |
---|
924 | case "units": this.units = paramVal.replace(/\s/gi,""); break; |
---|
925 | case "datum": this.datumCode = paramVal.replace(/\s/gi,""); break; |
---|
926 | case "nadgrids": this.nagrids = paramVal.replace(/\s/gi,""); break; |
---|
927 | case "ellps": this.ellps = paramVal.replace(/\s/gi,""); break; |
---|
928 | case "a": this.a = parseFloat(paramVal); break; // semi-major radius |
---|
929 | case "b": this.b = parseFloat(paramVal); break; // semi-minor radius |
---|
930 | // DGR 2007-11-20 |
---|
931 | case "rf": this.rf = parseFloat(paramVal); break; // inverse flattening rf= a/(a-b) |
---|
932 | case "lat_0": this.lat0 = paramVal*Proj4js.common.D2R; break; // phi0, central latitude |
---|
933 | case "lat_1": this.lat1 = paramVal*Proj4js.common.D2R; break; //standard parallel 1 |
---|
934 | case "lat_2": this.lat2 = paramVal*Proj4js.common.D2R; break; //standard parallel 2 |
---|
935 | case "lat_ts": this.lat_ts = paramVal*Proj4js.common.D2R; break; // used in merc and eqc |
---|
936 | case "lon_0": this.long0 = paramVal*Proj4js.common.D2R; break; // lam0, central longitude |
---|
937 | case "alpha": this.alpha = parseFloat(paramVal)*Proj4js.common.D2R; break; //for somerc projection |
---|
938 | case "lonc": this.longc = paramVal*Proj4js.common.D2R; break; //for somerc projection |
---|
939 | case "x_0": this.x0 = parseFloat(paramVal); break; // false easting |
---|
940 | case "y_0": this.y0 = parseFloat(paramVal); break; // false northing |
---|
941 | case "k_0": this.k0 = parseFloat(paramVal); break; // projection scale factor |
---|
942 | case "k": this.k0 = parseFloat(paramVal); break; // both forms returned |
---|
943 | case "r_a": this.R_A = true; break; // sphere--area of ellipsoid |
---|
944 | case "zone": this.zone = parseInt(paramVal); break; // UTM Zone |
---|
945 | case "south": this.utmSouth = true; break; // UTM north/south |
---|
946 | case "towgs84":this.datum_params = paramVal.split(","); break; |
---|
947 | case "to_meter": this.to_meter = parseFloat(paramVal); break; // cartesian scaling |
---|
948 | case "from_greenwich": this.from_greenwich = paramVal*Proj4js.common.D2R; break; |
---|
949 | // DGR 2008-07-09 : if pm is not a well-known prime meridian take |
---|
950 | // the value instead of 0.0, then convert to radians |
---|
951 | case "pm": paramVal = paramVal.replace(/\s/gi,""); |
---|
952 | this.from_greenwich = Proj4js.PrimeMeridian[paramVal] ? |
---|
953 | Proj4js.PrimeMeridian[paramVal] : parseFloat(paramVal); |
---|
954 | this.from_greenwich *= Proj4js.common.D2R; |
---|
955 | break; |
---|
956 | // DGR 2010-11-12: axis |
---|
957 | case "axis": paramVal = paramVal.replace(/\s/gi,""); |
---|
958 | var legalAxis= "ewnsud"; |
---|
959 | if (paramVal.length==3 && |
---|
960 | legalAxis.indexOf(paramVal.substr(0,1))!=-1 && |
---|
961 | legalAxis.indexOf(paramVal.substr(1,1))!=-1 && |
---|
962 | legalAxis.indexOf(paramVal.substr(2,1))!=-1) { |
---|
963 | this.axis= paramVal; |
---|
964 | } //FIXME: be silent ? |
---|
965 | break |
---|
966 | case "no_defs": break; |
---|
967 | default: //alert("Unrecognized parameter: " + paramName); |
---|
968 | } // switch() |
---|
969 | } // for paramArray |
---|
970 | this.deriveConstants(); |
---|
971 | }, |
---|
972 | |
---|
973 | /** |
---|
974 | * Function: deriveConstants |
---|
975 | * Sets several derived constant values and initialization of datum and ellipse |
---|
976 | * parameters. |
---|
977 | * |
---|
978 | */ |
---|
979 | deriveConstants: function() { |
---|
980 | if (this.nagrids == '@null') this.datumCode = 'none'; |
---|
981 | if (this.datumCode && this.datumCode != 'none') { |
---|
982 | var datumDef = Proj4js.Datum[this.datumCode]; |
---|
983 | if (datumDef) { |
---|
984 | this.datum_params = datumDef.towgs84 ? datumDef.towgs84.split(',') : null; |
---|
985 | this.ellps = datumDef.ellipse; |
---|
986 | this.datumName = datumDef.datumName ? datumDef.datumName : this.datumCode; |
---|
987 | } |
---|
988 | } |
---|
989 | if (!this.a) { // do we have an ellipsoid? |
---|
990 | var ellipse = Proj4js.Ellipsoid[this.ellps] ? Proj4js.Ellipsoid[this.ellps] : Proj4js.Ellipsoid['WGS84']; |
---|
991 | Proj4js.extend(this, ellipse); |
---|
992 | } |
---|
993 | if (this.rf && !this.b) this.b = (1.0 - 1.0/this.rf) * this.a; |
---|
994 | if (Math.abs(this.a - this.b)<Proj4js.common.EPSLN) { |
---|
995 | this.sphere = true; |
---|
996 | this.b= this.a; |
---|
997 | } |
---|
998 | this.a2 = this.a * this.a; // used in geocentric |
---|
999 | this.b2 = this.b * this.b; // used in geocentric |
---|
1000 | this.es = (this.a2-this.b2)/this.a2; // e ^ 2 |
---|
1001 | this.e = Math.sqrt(this.es); // eccentricity |
---|
1002 | if (this.R_A) { |
---|
1003 | this.a *= 1. - this.es * (Proj4js.common.SIXTH + this.es * (Proj4js.common.RA4 + this.es * Proj4js.common.RA6)); |
---|
1004 | this.a2 = this.a * this.a; |
---|
1005 | this.b2 = this.b * this.b; |
---|
1006 | this.es = 0.; |
---|
1007 | } |
---|
1008 | this.ep2=(this.a2-this.b2)/this.b2; // used in geocentric |
---|
1009 | if (!this.k0) this.k0 = 1.0; //default value |
---|
1010 | //DGR 2010-11-12: axis |
---|
1011 | if (!this.axis) { this.axis= "enu"; } |
---|
1012 | |
---|
1013 | this.datum = new Proj4js.datum(this); |
---|
1014 | } |
---|
1015 | }); |
---|
1016 | |
---|
1017 | Proj4js.Proj.longlat = { |
---|
1018 | init: function() { |
---|
1019 | //no-op for longlat |
---|
1020 | }, |
---|
1021 | forward: function(pt) { |
---|
1022 | //identity transform |
---|
1023 | return pt; |
---|
1024 | }, |
---|
1025 | inverse: function(pt) { |
---|
1026 | //identity transform |
---|
1027 | return pt; |
---|
1028 | } |
---|
1029 | }; |
---|
1030 | Proj4js.Proj.identity = Proj4js.Proj.longlat; |
---|
1031 | |
---|
1032 | /** |
---|
1033 | Proj4js.defs is a collection of coordinate system definition objects in the |
---|
1034 | PROJ.4 command line format. |
---|
1035 | Generally a def is added by means of a separate .js file for example: |
---|
1036 | |
---|
1037 | <SCRIPT type="text/javascript" src="defs/EPSG26912.js"></SCRIPT> |
---|
1038 | |
---|
1039 | def is a CS definition in PROJ.4 WKT format, for example: |
---|
1040 | +proj="tmerc" //longlat, etc. |
---|
1041 | +a=majorRadius |
---|
1042 | +b=minorRadius |
---|
1043 | +lat0=somenumber |
---|
1044 | +long=somenumber |
---|
1045 | */ |
---|
1046 | Proj4js.defs = { |
---|
1047 | // These are so widely used, we'll go ahead and throw them in |
---|
1048 | // without requiring a separate .js file |
---|
1049 | 'WGS84': "+title=long/lat:WGS84 +proj=longlat +ellps=WGS84 +datum=WGS84 +units=degrees", |
---|
1050 | 'EPSG:4326': "+title=long/lat:WGS84 +proj=longlat +a=6378137.0 +b=6356752.31424518 +ellps=WGS84 +datum=WGS84 +units=degrees", |
---|
1051 | 'EPSG:4269': "+title=long/lat:NAD83 +proj=longlat +a=6378137.0 +b=6356752.31414036 +ellps=GRS80 +datum=NAD83 +units=degrees", |
---|
1052 | 'EPSG:3785': "+title= Google Mercator +proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +no_defs" |
---|
1053 | }; |
---|
1054 | Proj4js.defs['GOOGLE'] = Proj4js.defs['EPSG:3785']; |
---|
1055 | Proj4js.defs['EPSG:900913'] = Proj4js.defs['EPSG:3785']; |
---|
1056 | Proj4js.defs['EPSG:102113'] = Proj4js.defs['EPSG:3785']; |
---|
1057 | |
---|
1058 | Proj4js.common = { |
---|
1059 | PI : 3.141592653589793238, //Math.PI, |
---|
1060 | HALF_PI : 1.570796326794896619, //Math.PI*0.5, |
---|
1061 | TWO_PI : 6.283185307179586477, //Math.PI*2, |
---|
1062 | FORTPI : 0.78539816339744833, |
---|
1063 | R2D : 57.29577951308232088, |
---|
1064 | D2R : 0.01745329251994329577, |
---|
1065 | SEC_TO_RAD : 4.84813681109535993589914102357e-6, /* SEC_TO_RAD = Pi/180/3600 */ |
---|
1066 | EPSLN : 1.0e-10, |
---|
1067 | MAX_ITER : 20, |
---|
1068 | // following constants from geocent.c |
---|
1069 | COS_67P5 : 0.38268343236508977, /* cosine of 67.5 degrees */ |
---|
1070 | AD_C : 1.0026000, /* Toms region 1 constant */ |
---|
1071 | |
---|
1072 | /* datum_type values */ |
---|
1073 | PJD_UNKNOWN : 0, |
---|
1074 | PJD_3PARAM : 1, |
---|
1075 | PJD_7PARAM : 2, |
---|
1076 | PJD_GRIDSHIFT: 3, |
---|
1077 | PJD_WGS84 : 4, // WGS84 or equivalent |
---|
1078 | PJD_NODATUM : 5, // WGS84 or equivalent |
---|
1079 | SRS_WGS84_SEMIMAJOR : 6378137.0, // only used in grid shift transforms |
---|
1080 | |
---|
1081 | // ellipoid pj_set_ell.c |
---|
1082 | SIXTH : .1666666666666666667, /* 1/6 */ |
---|
1083 | RA4 : .04722222222222222222, /* 17/360 */ |
---|
1084 | RA6 : .02215608465608465608, /* 67/3024 */ |
---|
1085 | RV4 : .06944444444444444444, /* 5/72 */ |
---|
1086 | RV6 : .04243827160493827160, /* 55/1296 */ |
---|
1087 | |
---|
1088 | // Function to compute the constant small m which is the radius of |
---|
1089 | // a parallel of latitude, phi, divided by the semimajor axis. |
---|
1090 | // ----------------------------------------------------------------- |
---|
1091 | msfnz : function(eccent, sinphi, cosphi) { |
---|
1092 | var con = eccent * sinphi; |
---|
1093 | return cosphi/(Math.sqrt(1.0 - con * con)); |
---|
1094 | }, |
---|
1095 | |
---|
1096 | // Function to compute the constant small t for use in the forward |
---|
1097 | // computations in the Lambert Conformal Conic and the Polar |
---|
1098 | // Stereographic projections. |
---|
1099 | // ----------------------------------------------------------------- |
---|
1100 | tsfnz : function(eccent, phi, sinphi) { |
---|
1101 | var con = eccent * sinphi; |
---|
1102 | var com = .5 * eccent; |
---|
1103 | con = Math.pow(((1.0 - con) / (1.0 + con)), com); |
---|
1104 | return (Math.tan(.5 * (this.HALF_PI - phi))/con); |
---|
1105 | }, |
---|
1106 | |
---|
1107 | // Function to compute the latitude angle, phi2, for the inverse of the |
---|
1108 | // Lambert Conformal Conic and Polar Stereographic projections. |
---|
1109 | // ---------------------------------------------------------------- |
---|
1110 | phi2z : function(eccent, ts) { |
---|
1111 | var eccnth = .5 * eccent; |
---|
1112 | var con, dphi; |
---|
1113 | var phi = this.HALF_PI - 2 * Math.atan(ts); |
---|
1114 | for (var i = 0; i <= 15; i++) { |
---|
1115 | con = eccent * Math.sin(phi); |
---|
1116 | dphi = this.HALF_PI - 2 * Math.atan(ts *(Math.pow(((1.0 - con)/(1.0 + con)),eccnth))) - phi; |
---|
1117 | phi += dphi; |
---|
1118 | if (Math.abs(dphi) <= .0000000001) return phi; |
---|
1119 | } |
---|
1120 | alert("phi2z has NoConvergence"); |
---|
1121 | return (-9999); |
---|
1122 | }, |
---|
1123 | |
---|
1124 | /* Function to compute constant small q which is the radius of a |
---|
1125 | parallel of latitude, phi, divided by the semimajor axis. |
---|
1126 | ------------------------------------------------------------*/ |
---|
1127 | qsfnz : function(eccent,sinphi) { |
---|
1128 | var con; |
---|
1129 | if (eccent > 1.0e-7) { |
---|
1130 | con = eccent * sinphi; |
---|
1131 | return (( 1.0- eccent * eccent) * (sinphi /(1.0 - con * con) - (.5/eccent)*Math.log((1.0 - con)/(1.0 + con)))); |
---|
1132 | } else { |
---|
1133 | return(2.0 * sinphi); |
---|
1134 | } |
---|
1135 | }, |
---|
1136 | |
---|
1137 | /* Function to eliminate roundoff errors in asin |
---|
1138 | ----------------------------------------------*/ |
---|
1139 | asinz : function(x) { |
---|
1140 | if (Math.abs(x)>1.0) { |
---|
1141 | x=(x>1.0)?1.0:-1.0; |
---|
1142 | } |
---|
1143 | return Math.asin(x); |
---|
1144 | }, |
---|
1145 | |
---|
1146 | // following functions from gctpc cproj.c for transverse mercator projections |
---|
1147 | e0fn : function(x) {return(1.0-0.25*x*(1.0+x/16.0*(3.0+1.25*x)));}, |
---|
1148 | e1fn : function(x) {return(0.375*x*(1.0+0.25*x*(1.0+0.46875*x)));}, |
---|
1149 | e2fn : function(x) {return(0.05859375*x*x*(1.0+0.75*x));}, |
---|
1150 | e3fn : function(x) {return(x*x*x*(35.0/3072.0));}, |
---|
1151 | mlfn : function(e0,e1,e2,e3,phi) {return(e0*phi-e1*Math.sin(2.0*phi)+e2*Math.sin(4.0*phi)-e3*Math.sin(6.0*phi));}, |
---|
1152 | |
---|
1153 | srat : function(esinp, exp) { |
---|
1154 | return(Math.pow((1.0-esinp)/(1.0+esinp), exp)); |
---|
1155 | }, |
---|
1156 | |
---|
1157 | // Function to return the sign of an argument |
---|
1158 | sign : function(x) { if (x < 0.0) return(-1); else return(1);}, |
---|
1159 | |
---|
1160 | // Function to adjust longitude to -180 to 180; input in radians |
---|
1161 | adjust_lon : function(x) { |
---|
1162 | x = (Math.abs(x) < this.PI) ? x: (x - (this.sign(x)*this.TWO_PI) ); |
---|
1163 | return x; |
---|
1164 | }, |
---|
1165 | |
---|
1166 | // IGNF - DGR : algorithms used by IGN France |
---|
1167 | |
---|
1168 | // Function to adjust latitude to -90 to 90; input in radians |
---|
1169 | adjust_lat : function(x) { |
---|
1170 | x= (Math.abs(x) < this.HALF_PI) ? x: (x - (this.sign(x)*this.PI) ); |
---|
1171 | return x; |
---|
1172 | }, |
---|
1173 | |
---|
1174 | // Latitude Isometrique - close to tsfnz ... |
---|
1175 | latiso : function(eccent, phi, sinphi) { |
---|
1176 | if (Math.abs(phi) > this.HALF_PI) return +Number.NaN; |
---|
1177 | if (phi==this.HALF_PI) return Number.POSITIVE_INFINITY; |
---|
1178 | if (phi==-1.0*this.HALF_PI) return -1.0*Number.POSITIVE_INFINITY; |
---|
1179 | |
---|
1180 | var con= eccent*sinphi; |
---|
1181 | return Math.log(Math.tan((this.HALF_PI+phi)/2.0))+eccent*Math.log((1.0-con)/(1.0+con))/2.0; |
---|
1182 | }, |
---|
1183 | |
---|
1184 | fL : function(x,L) { |
---|
1185 | return 2.0*Math.atan(x*Math.exp(L)) - this.HALF_PI; |
---|
1186 | }, |
---|
1187 | |
---|
1188 | // Inverse Latitude Isometrique - close to ph2z |
---|
1189 | invlatiso : function(eccent, ts) { |
---|
1190 | var phi= this.fL(1.0,ts); |
---|
1191 | var Iphi= 0.0; |
---|
1192 | var con= 0.0; |
---|
1193 | do { |
---|
1194 | Iphi= phi; |
---|
1195 | con= eccent*Math.sin(Iphi); |
---|
1196 | phi= this.fL(Math.exp(eccent*Math.log((1.0+con)/(1.0-con))/2.0),ts) |
---|
1197 | } while (Math.abs(phi-Iphi)>1.0e-12); |
---|
1198 | return phi; |
---|
1199 | }, |
---|
1200 | |
---|
1201 | // Needed for Gauss Schreiber |
---|
1202 | // Original: Denis Makarov (info@binarythings.com) |
---|
1203 | // Web Site: http://www.binarythings.com |
---|
1204 | sinh : function(x) |
---|
1205 | { |
---|
1206 | var r= Math.exp(x); |
---|
1207 | r= (r-1.0/r)/2.0; |
---|
1208 | return r; |
---|
1209 | }, |
---|
1210 | |
---|
1211 | cosh : function(x) |
---|
1212 | { |
---|
1213 | var r= Math.exp(x); |
---|
1214 | r= (r+1.0/r)/2.0; |
---|
1215 | return r; |
---|
1216 | }, |
---|
1217 | |
---|
1218 | tanh : function(x) |
---|
1219 | { |
---|
1220 | var r= Math.exp(x); |
---|
1221 | r= (r-1.0/r)/(r+1.0/r); |
---|
1222 | return r; |
---|
1223 | }, |
---|
1224 | |
---|
1225 | asinh : function(x) |
---|
1226 | { |
---|
1227 | var s= (x>= 0? 1.0:-1.0); |
---|
1228 | return s*(Math.log( Math.abs(x) + Math.sqrt(x*x+1.0) )); |
---|
1229 | }, |
---|
1230 | |
---|
1231 | acosh : function(x) |
---|
1232 | { |
---|
1233 | return 2.0*Math.log(Math.sqrt((x+1.0)/2.0) + Math.sqrt((x-1.0)/2.0)); |
---|
1234 | }, |
---|
1235 | |
---|
1236 | atanh : function(x) |
---|
1237 | { |
---|
1238 | return Math.log((x-1.0)/(x+1.0))/2.0; |
---|
1239 | }, |
---|
1240 | |
---|
1241 | // Grande Normale |
---|
1242 | gN : function(a,e,sinphi) |
---|
1243 | { |
---|
1244 | var temp= e*sinphi; |
---|
1245 | return a/Math.sqrt(1.0 - temp*temp); |
---|
1246 | } |
---|
1247 | |
---|
1248 | }; |
---|
1249 | |
---|
1250 | /** datum object |
---|
1251 | */ |
---|
1252 | Proj4js.datum = Proj4js.Class({ |
---|
1253 | |
---|
1254 | initialize : function(proj) { |
---|
1255 | this.datum_type = Proj4js.common.PJD_WGS84; //default setting |
---|
1256 | if (proj.datumCode && proj.datumCode == 'none') { |
---|
1257 | this.datum_type = Proj4js.common.PJD_NODATUM; |
---|
1258 | } |
---|
1259 | if (proj && proj.datum_params) { |
---|
1260 | for (var i=0; i<proj.datum_params.length; i++) { |
---|
1261 | proj.datum_params[i]=parseFloat(proj.datum_params[i]); |
---|
1262 | } |
---|
1263 | if (proj.datum_params[0] != 0 || proj.datum_params[1] != 0 || proj.datum_params[2] != 0 ) { |
---|
1264 | this.datum_type = Proj4js.common.PJD_3PARAM; |
---|
1265 | } |
---|
1266 | if (proj.datum_params.length > 3) { |
---|
1267 | if (proj.datum_params[3] != 0 || proj.datum_params[4] != 0 || |
---|
1268 | proj.datum_params[5] != 0 || proj.datum_params[6] != 0 ) { |
---|
1269 | this.datum_type = Proj4js.common.PJD_7PARAM; |
---|
1270 | proj.datum_params[3] *= Proj4js.common.SEC_TO_RAD; |
---|
1271 | proj.datum_params[4] *= Proj4js.common.SEC_TO_RAD; |
---|
1272 | proj.datum_params[5] *= Proj4js.common.SEC_TO_RAD; |
---|
1273 | proj.datum_params[6] = (proj.datum_params[6]/1000000.0) + 1.0; |
---|
1274 | } |
---|
1275 | } |
---|
1276 | } |
---|
1277 | if (proj) { |
---|
1278 | this.a = proj.a; //datum object also uses these values |
---|
1279 | this.b = proj.b; |
---|
1280 | this.es = proj.es; |
---|
1281 | this.ep2 = proj.ep2; |
---|
1282 | this.datum_params = proj.datum_params; |
---|
1283 | } |
---|
1284 | }, |
---|
1285 | |
---|
1286 | /****************************************************************/ |
---|
1287 | // cs_compare_datums() |
---|
1288 | // Returns 1 (TRUE) if the two datums match, otherwise 0 (FALSE). |
---|
1289 | compare_datums : function( dest ) { |
---|
1290 | if( this.datum_type != dest.datum_type ) { |
---|
1291 | return false; // false, datums are not equal |
---|
1292 | } else if( this.a != dest.a || Math.abs(this.es-dest.es) > 0.000000000050 ) { |
---|
1293 | // the tolerence for es is to ensure that GRS80 and WGS84 |
---|
1294 | // are considered identical |
---|
1295 | return false; |
---|
1296 | } else if( this.datum_type == Proj4js.common.PJD_3PARAM ) { |
---|
1297 | return (this.datum_params[0] == dest.datum_params[0] |
---|
1298 | && this.datum_params[1] == dest.datum_params[1] |
---|
1299 | && this.datum_params[2] == dest.datum_params[2]); |
---|
1300 | } else if( this.datum_type == Proj4js.common.PJD_7PARAM ) { |
---|
1301 | return (this.datum_params[0] == dest.datum_params[0] |
---|
1302 | && this.datum_params[1] == dest.datum_params[1] |
---|
1303 | && this.datum_params[2] == dest.datum_params[2] |
---|
1304 | && this.datum_params[3] == dest.datum_params[3] |
---|
1305 | && this.datum_params[4] == dest.datum_params[4] |
---|
1306 | && this.datum_params[5] == dest.datum_params[5] |
---|
1307 | && this.datum_params[6] == dest.datum_params[6]); |
---|
1308 | } else if( this.datum_type == Proj4js.common.PJD_GRIDSHIFT ) { |
---|
1309 | return strcmp( pj_param(this.params,"snadgrids").s, |
---|
1310 | pj_param(dest.params,"snadgrids").s ) == 0; |
---|
1311 | } else { |
---|
1312 | return true; // datums are equal |
---|
1313 | } |
---|
1314 | }, // cs_compare_datums() |
---|
1315 | |
---|
1316 | /* |
---|
1317 | * The function Convert_Geodetic_To_Geocentric converts geodetic coordinates |
---|
1318 | * (latitude, longitude, and height) to geocentric coordinates (X, Y, Z), |
---|
1319 | * according to the current ellipsoid parameters. |
---|
1320 | * |
---|
1321 | * Latitude : Geodetic latitude in radians (input) |
---|
1322 | * Longitude : Geodetic longitude in radians (input) |
---|
1323 | * Height : Geodetic height, in meters (input) |
---|
1324 | * X : Calculated Geocentric X coordinate, in meters (output) |
---|
1325 | * Y : Calculated Geocentric Y coordinate, in meters (output) |
---|
1326 | * Z : Calculated Geocentric Z coordinate, in meters (output) |
---|
1327 | * |
---|
1328 | */ |
---|
1329 | geodetic_to_geocentric : function(p) { |
---|
1330 | var Longitude = p.x; |
---|
1331 | var Latitude = p.y; |
---|
1332 | var Height = p.z ? p.z : 0; //Z value not always supplied |
---|
1333 | var X; // output |
---|
1334 | var Y; |
---|
1335 | var Z; |
---|
1336 | |
---|
1337 | var Error_Code=0; // GEOCENT_NO_ERROR; |
---|
1338 | var Rn; /* Earth radius at location */ |
---|
1339 | var Sin_Lat; /* Math.sin(Latitude) */ |
---|
1340 | var Sin2_Lat; /* Square of Math.sin(Latitude) */ |
---|
1341 | var Cos_Lat; /* Math.cos(Latitude) */ |
---|
1342 | |
---|
1343 | /* |
---|
1344 | ** Don't blow up if Latitude is just a little out of the value |
---|
1345 | ** range as it may just be a rounding issue. Also removed longitude |
---|
1346 | ** test, it should be wrapped by Math.cos() and Math.sin(). NFW for PROJ.4, Sep/2001. |
---|
1347 | */ |
---|
1348 | if( Latitude < -Proj4js.common.HALF_PI && Latitude > -1.001 * Proj4js.common.HALF_PI ) { |
---|
1349 | Latitude = -Proj4js.common.HALF_PI; |
---|
1350 | } else if( Latitude > Proj4js.common.HALF_PI && Latitude < 1.001 * Proj4js.common.HALF_PI ) { |
---|
1351 | Latitude = Proj4js.common.HALF_PI; |
---|
1352 | } else if ((Latitude < -Proj4js.common.HALF_PI) || (Latitude > Proj4js.common.HALF_PI)) { |
---|
1353 | /* Latitude out of range */ |
---|
1354 | Proj4js.reportError('geocent:lat out of range:'+Latitude); |
---|
1355 | return null; |
---|
1356 | } |
---|
1357 | |
---|
1358 | if (Longitude > Proj4js.common.PI) Longitude -= (2*Proj4js.common.PI); |
---|
1359 | Sin_Lat = Math.sin(Latitude); |
---|
1360 | Cos_Lat = Math.cos(Latitude); |
---|
1361 | Sin2_Lat = Sin_Lat * Sin_Lat; |
---|
1362 | Rn = this.a / (Math.sqrt(1.0e0 - this.es * Sin2_Lat)); |
---|
1363 | X = (Rn + Height) * Cos_Lat * Math.cos(Longitude); |
---|
1364 | Y = (Rn + Height) * Cos_Lat * Math.sin(Longitude); |
---|
1365 | Z = ((Rn * (1 - this.es)) + Height) * Sin_Lat; |
---|
1366 | |
---|
1367 | p.x = X; |
---|
1368 | p.y = Y; |
---|
1369 | p.z = Z; |
---|
1370 | return Error_Code; |
---|
1371 | }, // cs_geodetic_to_geocentric() |
---|
1372 | |
---|
1373 | |
---|
1374 | geocentric_to_geodetic : function (p) { |
---|
1375 | /* local defintions and variables */ |
---|
1376 | /* end-criterium of loop, accuracy of sin(Latitude) */ |
---|
1377 | var genau = 1.E-12; |
---|
1378 | var genau2 = (genau*genau); |
---|
1379 | var maxiter = 30; |
---|
1380 | |
---|
1381 | var P; /* distance between semi-minor axis and location */ |
---|
1382 | var RR; /* distance between center and location */ |
---|
1383 | var CT; /* sin of geocentric latitude */ |
---|
1384 | var ST; /* cos of geocentric latitude */ |
---|
1385 | var RX; |
---|
1386 | var RK; |
---|
1387 | var RN; /* Earth radius at location */ |
---|
1388 | var CPHI0; /* cos of start or old geodetic latitude in iterations */ |
---|
1389 | var SPHI0; /* sin of start or old geodetic latitude in iterations */ |
---|
1390 | var CPHI; /* cos of searched geodetic latitude */ |
---|
1391 | var SPHI; /* sin of searched geodetic latitude */ |
---|
1392 | var SDPHI; /* end-criterium: addition-theorem of sin(Latitude(iter)-Latitude(iter-1)) */ |
---|
1393 | var At_Pole; /* indicates location is in polar region */ |
---|
1394 | var iter; /* # of continous iteration, max. 30 is always enough (s.a.) */ |
---|
1395 | |
---|
1396 | var X = p.x; |
---|
1397 | var Y = p.y; |
---|
1398 | var Z = p.z ? p.z : 0.0; //Z value not always supplied |
---|
1399 | var Longitude; |
---|
1400 | var Latitude; |
---|
1401 | var Height; |
---|
1402 | |
---|
1403 | At_Pole = false; |
---|
1404 | P = Math.sqrt(X*X+Y*Y); |
---|
1405 | RR = Math.sqrt(X*X+Y*Y+Z*Z); |
---|
1406 | |
---|
1407 | /* special cases for latitude and longitude */ |
---|
1408 | if (P/this.a < genau) { |
---|
1409 | |
---|
1410 | /* special case, if P=0. (X=0., Y=0.) */ |
---|
1411 | At_Pole = true; |
---|
1412 | Longitude = 0.0; |
---|
1413 | |
---|
1414 | /* if (X,Y,Z)=(0.,0.,0.) then Height becomes semi-minor axis |
---|
1415 | * of ellipsoid (=center of mass), Latitude becomes PI/2 */ |
---|
1416 | if (RR/this.a < genau) { |
---|
1417 | Latitude = Proj4js.common.HALF_PI; |
---|
1418 | Height = -this.b; |
---|
1419 | return; |
---|
1420 | } |
---|
1421 | } else { |
---|
1422 | /* ellipsoidal (geodetic) longitude |
---|
1423 | * interval: -PI < Longitude <= +PI */ |
---|
1424 | Longitude=Math.atan2(Y,X); |
---|
1425 | } |
---|
1426 | |
---|
1427 | /* -------------------------------------------------------------- |
---|
1428 | * Following iterative algorithm was developped by |
---|
1429 | * "Institut fï¿œr Erdmessung", University of Hannover, July 1988. |
---|
1430 | * Internet: www.ife.uni-hannover.de |
---|
1431 | * Iterative computation of CPHI,SPHI and Height. |
---|
1432 | * Iteration of CPHI and SPHI to 10**-12 radian resp. |
---|
1433 | * 2*10**-7 arcsec. |
---|
1434 | * -------------------------------------------------------------- |
---|
1435 | */ |
---|
1436 | CT = Z/RR; |
---|
1437 | ST = P/RR; |
---|
1438 | RX = 1.0/Math.sqrt(1.0-this.es*(2.0-this.es)*ST*ST); |
---|
1439 | CPHI0 = ST*(1.0-this.es)*RX; |
---|
1440 | SPHI0 = CT*RX; |
---|
1441 | iter = 0; |
---|
1442 | |
---|
1443 | /* loop to find sin(Latitude) resp. Latitude |
---|
1444 | * until |sin(Latitude(iter)-Latitude(iter-1))| < genau */ |
---|
1445 | do |
---|
1446 | { |
---|
1447 | iter++; |
---|
1448 | RN = this.a/Math.sqrt(1.0-this.es*SPHI0*SPHI0); |
---|
1449 | |
---|
1450 | /* ellipsoidal (geodetic) height */ |
---|
1451 | Height = P*CPHI0+Z*SPHI0-RN*(1.0-this.es*SPHI0*SPHI0); |
---|
1452 | |
---|
1453 | RK = this.es*RN/(RN+Height); |
---|
1454 | RX = 1.0/Math.sqrt(1.0-RK*(2.0-RK)*ST*ST); |
---|
1455 | CPHI = ST*(1.0-RK)*RX; |
---|
1456 | SPHI = CT*RX; |
---|
1457 | SDPHI = SPHI*CPHI0-CPHI*SPHI0; |
---|
1458 | CPHI0 = CPHI; |
---|
1459 | SPHI0 = SPHI; |
---|
1460 | } |
---|
1461 | while (SDPHI*SDPHI > genau2 && iter < maxiter); |
---|
1462 | |
---|
1463 | /* ellipsoidal (geodetic) latitude */ |
---|
1464 | Latitude=Math.atan(SPHI/Math.abs(CPHI)); |
---|
1465 | |
---|
1466 | p.x = Longitude; |
---|
1467 | p.y = Latitude; |
---|
1468 | p.z = Height; |
---|
1469 | return p; |
---|
1470 | }, // cs_geocentric_to_geodetic() |
---|
1471 | |
---|
1472 | /** Convert_Geocentric_To_Geodetic |
---|
1473 | * The method used here is derived from 'An Improved Algorithm for |
---|
1474 | * Geocentric to Geodetic Coordinate Conversion', by Ralph Toms, Feb 1996 |
---|
1475 | */ |
---|
1476 | geocentric_to_geodetic_noniter : function (p) { |
---|
1477 | var X = p.x; |
---|
1478 | var Y = p.y; |
---|
1479 | var Z = p.z ? p.z : 0; //Z value not always supplied |
---|
1480 | var Longitude; |
---|
1481 | var Latitude; |
---|
1482 | var Height; |
---|
1483 | |
---|
1484 | var W; /* distance from Z axis */ |
---|
1485 | var W2; /* square of distance from Z axis */ |
---|
1486 | var T0; /* initial estimate of vertical component */ |
---|
1487 | var T1; /* corrected estimate of vertical component */ |
---|
1488 | var S0; /* initial estimate of horizontal component */ |
---|
1489 | var S1; /* corrected estimate of horizontal component */ |
---|
1490 | var Sin_B0; /* Math.sin(B0), B0 is estimate of Bowring aux variable */ |
---|
1491 | var Sin3_B0; /* cube of Math.sin(B0) */ |
---|
1492 | var Cos_B0; /* Math.cos(B0) */ |
---|
1493 | var Sin_p1; /* Math.sin(phi1), phi1 is estimated latitude */ |
---|
1494 | var Cos_p1; /* Math.cos(phi1) */ |
---|
1495 | var Rn; /* Earth radius at location */ |
---|
1496 | var Sum; /* numerator of Math.cos(phi1) */ |
---|
1497 | var At_Pole; /* indicates location is in polar region */ |
---|
1498 | |
---|
1499 | X = parseFloat(X); // cast from string to float |
---|
1500 | Y = parseFloat(Y); |
---|
1501 | Z = parseFloat(Z); |
---|
1502 | |
---|
1503 | At_Pole = false; |
---|
1504 | if (X != 0.0) |
---|
1505 | { |
---|
1506 | Longitude = Math.atan2(Y,X); |
---|
1507 | } |
---|
1508 | else |
---|
1509 | { |
---|
1510 | if (Y > 0) |
---|
1511 | { |
---|
1512 | Longitude = Proj4js.common.HALF_PI; |
---|
1513 | } |
---|
1514 | else if (Y < 0) |
---|
1515 | { |
---|
1516 | Longitude = -Proj4js.common.HALF_PI; |
---|
1517 | } |
---|
1518 | else |
---|
1519 | { |
---|
1520 | At_Pole = true; |
---|
1521 | Longitude = 0.0; |
---|
1522 | if (Z > 0.0) |
---|
1523 | { /* north pole */ |
---|
1524 | Latitude = Proj4js.common.HALF_PI; |
---|
1525 | } |
---|
1526 | else if (Z < 0.0) |
---|
1527 | { /* south pole */ |
---|
1528 | Latitude = -Proj4js.common.HALF_PI; |
---|
1529 | } |
---|
1530 | else |
---|
1531 | { /* center of earth */ |
---|
1532 | Latitude = Proj4js.common.HALF_PI; |
---|
1533 | Height = -this.b; |
---|
1534 | return; |
---|
1535 | } |
---|
1536 | } |
---|
1537 | } |
---|
1538 | W2 = X*X + Y*Y; |
---|
1539 | W = Math.sqrt(W2); |
---|
1540 | T0 = Z * Proj4js.common.AD_C; |
---|
1541 | S0 = Math.sqrt(T0 * T0 + W2); |
---|
1542 | Sin_B0 = T0 / S0; |
---|
1543 | Cos_B0 = W / S0; |
---|
1544 | Sin3_B0 = Sin_B0 * Sin_B0 * Sin_B0; |
---|
1545 | T1 = Z + this.b * this.ep2 * Sin3_B0; |
---|
1546 | Sum = W - this.a * this.es * Cos_B0 * Cos_B0 * Cos_B0; |
---|
1547 | S1 = Math.sqrt(T1*T1 + Sum * Sum); |
---|
1548 | Sin_p1 = T1 / S1; |
---|
1549 | Cos_p1 = Sum / S1; |
---|
1550 | Rn = this.a / Math.sqrt(1.0 - this.es * Sin_p1 * Sin_p1); |
---|
1551 | if (Cos_p1 >= Proj4js.common.COS_67P5) |
---|
1552 | { |
---|
1553 | Height = W / Cos_p1 - Rn; |
---|
1554 | } |
---|
1555 | else if (Cos_p1 <= -Proj4js.common.COS_67P5) |
---|
1556 | { |
---|
1557 | Height = W / -Cos_p1 - Rn; |
---|
1558 | } |
---|
1559 | else |
---|
1560 | { |
---|
1561 | Height = Z / Sin_p1 + Rn * (this.es - 1.0); |
---|
1562 | } |
---|
1563 | if (At_Pole == false) |
---|
1564 | { |
---|
1565 | Latitude = Math.atan(Sin_p1 / Cos_p1); |
---|
1566 | } |
---|
1567 | |
---|
1568 | p.x = Longitude; |
---|
1569 | p.y = Latitude; |
---|
1570 | p.z = Height; |
---|
1571 | return p; |
---|
1572 | }, // geocentric_to_geodetic_noniter() |
---|
1573 | |
---|
1574 | /****************************************************************/ |
---|
1575 | // pj_geocentic_to_wgs84( p ) |
---|
1576 | // p = point to transform in geocentric coordinates (x,y,z) |
---|
1577 | geocentric_to_wgs84 : function ( p ) { |
---|
1578 | |
---|
1579 | if( this.datum_type == Proj4js.common.PJD_3PARAM ) |
---|
1580 | { |
---|
1581 | // if( x[io] == HUGE_VAL ) |
---|
1582 | // continue; |
---|
1583 | p.x += this.datum_params[0]; |
---|
1584 | p.y += this.datum_params[1]; |
---|
1585 | p.z += this.datum_params[2]; |
---|
1586 | |
---|
1587 | } |
---|
1588 | else if (this.datum_type == Proj4js.common.PJD_7PARAM) |
---|
1589 | { |
---|
1590 | var Dx_BF =this.datum_params[0]; |
---|
1591 | var Dy_BF =this.datum_params[1]; |
---|
1592 | var Dz_BF =this.datum_params[2]; |
---|
1593 | var Rx_BF =this.datum_params[3]; |
---|
1594 | var Ry_BF =this.datum_params[4]; |
---|
1595 | var Rz_BF =this.datum_params[5]; |
---|
1596 | var M_BF =this.datum_params[6]; |
---|
1597 | // if( x[io] == HUGE_VAL ) |
---|
1598 | // continue; |
---|
1599 | var x_out = M_BF*( p.x - Rz_BF*p.y + Ry_BF*p.z) + Dx_BF; |
---|
1600 | var y_out = M_BF*( Rz_BF*p.x + p.y - Rx_BF*p.z) + Dy_BF; |
---|
1601 | var z_out = M_BF*(-Ry_BF*p.x + Rx_BF*p.y + p.z) + Dz_BF; |
---|
1602 | p.x = x_out; |
---|
1603 | p.y = y_out; |
---|
1604 | p.z = z_out; |
---|
1605 | } |
---|
1606 | }, // cs_geocentric_to_wgs84 |
---|
1607 | |
---|
1608 | /****************************************************************/ |
---|
1609 | // pj_geocentic_from_wgs84() |
---|
1610 | // coordinate system definition, |
---|
1611 | // point to transform in geocentric coordinates (x,y,z) |
---|
1612 | geocentric_from_wgs84 : function( p ) { |
---|
1613 | |
---|
1614 | if( this.datum_type == Proj4js.common.PJD_3PARAM ) |
---|
1615 | { |
---|
1616 | //if( x[io] == HUGE_VAL ) |
---|
1617 | // continue; |
---|
1618 | p.x -= this.datum_params[0]; |
---|
1619 | p.y -= this.datum_params[1]; |
---|
1620 | p.z -= this.datum_params[2]; |
---|
1621 | |
---|
1622 | } |
---|
1623 | else if (this.datum_type == Proj4js.common.PJD_7PARAM) |
---|
1624 | { |
---|
1625 | var Dx_BF =this.datum_params[0]; |
---|
1626 | var Dy_BF =this.datum_params[1]; |
---|
1627 | var Dz_BF =this.datum_params[2]; |
---|
1628 | var Rx_BF =this.datum_params[3]; |
---|
1629 | var Ry_BF =this.datum_params[4]; |
---|
1630 | var Rz_BF =this.datum_params[5]; |
---|
1631 | var M_BF =this.datum_params[6]; |
---|
1632 | var x_tmp = (p.x - Dx_BF) / M_BF; |
---|
1633 | var y_tmp = (p.y - Dy_BF) / M_BF; |
---|
1634 | var z_tmp = (p.z - Dz_BF) / M_BF; |
---|
1635 | //if( x[io] == HUGE_VAL ) |
---|
1636 | // continue; |
---|
1637 | |
---|
1638 | p.x = x_tmp + Rz_BF*y_tmp - Ry_BF*z_tmp; |
---|
1639 | p.y = -Rz_BF*x_tmp + y_tmp + Rx_BF*z_tmp; |
---|
1640 | p.z = Ry_BF*x_tmp - Rx_BF*y_tmp + z_tmp; |
---|
1641 | } //cs_geocentric_from_wgs84() |
---|
1642 | } |
---|
1643 | }); |
---|
1644 | |
---|
1645 | /** point object, nothing fancy, just allows values to be |
---|
1646 | passed back and forth by reference rather than by value. |
---|
1647 | Other point classes may be used as long as they have |
---|
1648 | x and y properties, which will get modified in the transform method. |
---|
1649 | */ |
---|
1650 | Proj4js.Point = Proj4js.Class({ |
---|
1651 | |
---|
1652 | /** |
---|
1653 | * Constructor: Proj4js.Point |
---|
1654 | * |
---|
1655 | * Parameters: |
---|
1656 | * - x {float} or {Array} either the first coordinates component or |
---|
1657 | * the full coordinates |
---|
1658 | * - y {float} the second component |
---|
1659 | * - z {float} the third component, optional. |
---|
1660 | */ |
---|
1661 | initialize : function(x,y,z) { |
---|
1662 | if (typeof x == 'object') { |
---|
1663 | this.x = x[0]; |
---|
1664 | this.y = x[1]; |
---|
1665 | this.z = x[2] || 0.0; |
---|
1666 | } else if (typeof x == 'string' && typeof y == 'undefined') { |
---|
1667 | var coords = x.split(','); |
---|
1668 | this.x = parseFloat(coords[0]); |
---|
1669 | this.y = parseFloat(coords[1]); |
---|
1670 | this.z = parseFloat(coords[2]) || 0.0; |
---|
1671 | } else { |
---|
1672 | this.x = x; |
---|
1673 | this.y = y; |
---|
1674 | this.z = z || 0.0; |
---|
1675 | } |
---|
1676 | }, |
---|
1677 | |
---|
1678 | /** |
---|
1679 | * APIMethod: clone |
---|
1680 | * Build a copy of a Proj4js.Point object. |
---|
1681 | * |
---|
1682 | * Return: |
---|
1683 | * {Proj4js}.Point the cloned point. |
---|
1684 | */ |
---|
1685 | clone : function() { |
---|
1686 | return new Proj4js.Point(this.x, this.y, this.z); |
---|
1687 | }, |
---|
1688 | |
---|
1689 | /** |
---|
1690 | * APIMethod: toString |
---|
1691 | * Return a readable string version of the point |
---|
1692 | * |
---|
1693 | * Return: |
---|
1694 | * {String} String representation of Proj4js.Point object. |
---|
1695 | * (ex. <i>"x=5,y=42"</i>) |
---|
1696 | */ |
---|
1697 | toString : function() { |
---|
1698 | return ("x=" + this.x + ",y=" + this.y); |
---|
1699 | }, |
---|
1700 | |
---|
1701 | /** |
---|
1702 | * APIMethod: toShortString |
---|
1703 | * Return a short string version of the point. |
---|
1704 | * |
---|
1705 | * Return: |
---|
1706 | * {String} Shortened String representation of Proj4js.Point object. |
---|
1707 | * (ex. <i>"5, 42"</i>) |
---|
1708 | */ |
---|
1709 | toShortString : function() { |
---|
1710 | return (this.x + ", " + this.y); |
---|
1711 | } |
---|
1712 | }); |
---|
1713 | |
---|
1714 | Proj4js.PrimeMeridian = { |
---|
1715 | "greenwich": 0.0, //"0dE", |
---|
1716 | "lisbon": -9.131906111111, //"9d07'54.862\"W", |
---|
1717 | "paris": 2.337229166667, //"2d20'14.025\"E", |
---|
1718 | "bogota": -74.080916666667, //"74d04'51.3\"W", |
---|
1719 | "madrid": -3.687938888889, //"3d41'16.58\"W", |
---|
1720 | "rome": 12.452333333333, //"12d27'8.4\"E", |
---|
1721 | "bern": 7.439583333333, //"7d26'22.5\"E", |
---|
1722 | "jakarta": 106.807719444444, //"106d48'27.79\"E", |
---|
1723 | "ferro": -17.666666666667, //"17d40'W", |
---|
1724 | "brussels": 4.367975, //"4d22'4.71\"E", |
---|
1725 | "stockholm": 18.058277777778, //"18d3'29.8\"E", |
---|
1726 | "athens": 23.7163375, //"23d42'58.815\"E", |
---|
1727 | "oslo": 10.722916666667 //"10d43'22.5\"E" |
---|
1728 | }; |
---|
1729 | |
---|
1730 | Proj4js.Ellipsoid = { |
---|
1731 | "MERIT": {a:6378137.0, rf:298.257, ellipseName:"MERIT 1983"}, |
---|
1732 | "SGS85": {a:6378136.0, rf:298.257, ellipseName:"Soviet Geodetic System 85"}, |
---|
1733 | "GRS80": {a:6378137.0, rf:298.257222101, ellipseName:"GRS 1980(IUGG, 1980)"}, |
---|
1734 | "IAU76": {a:6378140.0, rf:298.257, ellipseName:"IAU 1976"}, |
---|
1735 | "airy": {a:6377563.396, b:6356256.910, ellipseName:"Airy 1830"}, |
---|
1736 | "APL4.": {a:6378137, rf:298.25, ellipseName:"Appl. Physics. 1965"}, |
---|
1737 | "NWL9D": {a:6378145.0, rf:298.25, ellipseName:"Naval Weapons Lab., 1965"}, |
---|
1738 | "mod_airy": {a:6377340.189, b:6356034.446, ellipseName:"Modified Airy"}, |
---|
1739 | "andrae": {a:6377104.43, rf:300.0, ellipseName:"Andrae 1876 (Den., Iclnd.)"}, |
---|
1740 | "aust_SA": {a:6378160.0, rf:298.25, ellipseName:"Australian Natl & S. Amer. 1969"}, |
---|
1741 | "GRS67": {a:6378160.0, rf:298.2471674270, ellipseName:"GRS 67(IUGG 1967)"}, |
---|
1742 | "bessel": {a:6377397.155, rf:299.1528128, ellipseName:"Bessel 1841"}, |
---|
1743 | "bess_nam": {a:6377483.865, rf:299.1528128, ellipseName:"Bessel 1841 (Namibia)"}, |
---|
1744 | "clrk66": {a:6378206.4, b:6356583.8, ellipseName:"Clarke 1866"}, |
---|
1745 | "clrk80": {a:6378249.145, rf:293.4663, ellipseName:"Clarke 1880 mod."}, |
---|
1746 | "CPM": {a:6375738.7, rf:334.29, ellipseName:"Comm. des Poids et Mesures 1799"}, |
---|
1747 | "delmbr": {a:6376428.0, rf:311.5, ellipseName:"Delambre 1810 (Belgium)"}, |
---|
1748 | "engelis": {a:6378136.05, rf:298.2566, ellipseName:"Engelis 1985"}, |
---|
1749 | "evrst30": {a:6377276.345, rf:300.8017, ellipseName:"Everest 1830"}, |
---|
1750 | "evrst48": {a:6377304.063, rf:300.8017, ellipseName:"Everest 1948"}, |
---|
1751 | "evrst56": {a:6377301.243, rf:300.8017, ellipseName:"Everest 1956"}, |
---|
1752 | "evrst69": {a:6377295.664, rf:300.8017, ellipseName:"Everest 1969"}, |
---|
1753 | "evrstSS": {a:6377298.556, rf:300.8017, ellipseName:"Everest (Sabah & Sarawak)"}, |
---|
1754 | "fschr60": {a:6378166.0, rf:298.3, ellipseName:"Fischer (Mercury Datum) 1960"}, |
---|
1755 | "fschr60m": {a:6378155.0, rf:298.3, ellipseName:"Fischer 1960"}, |
---|
1756 | "fschr68": {a:6378150.0, rf:298.3, ellipseName:"Fischer 1968"}, |
---|
1757 | "helmert": {a:6378200.0, rf:298.3, ellipseName:"Helmert 1906"}, |
---|
1758 | "hough": {a:6378270.0, rf:297.0, ellipseName:"Hough"}, |
---|
1759 | "intl": {a:6378388.0, rf:297.0, ellipseName:"International 1909 (Hayford)"}, |
---|
1760 | "kaula": {a:6378163.0, rf:298.24, ellipseName:"Kaula 1961"}, |
---|
1761 | "lerch": {a:6378139.0, rf:298.257, ellipseName:"Lerch 1979"}, |
---|
1762 | "mprts": {a:6397300.0, rf:191.0, ellipseName:"Maupertius 1738"}, |
---|
1763 | "new_intl": {a:6378157.5, b:6356772.2, ellipseName:"New International 1967"}, |
---|
1764 | "plessis": {a:6376523.0, rf:6355863.0, ellipseName:"Plessis 1817 (France)"}, |
---|
1765 | "krass": {a:6378245.0, rf:298.3, ellipseName:"Krassovsky, 1942"}, |
---|
1766 | "SEasia": {a:6378155.0, b:6356773.3205, ellipseName:"Southeast Asia"}, |
---|
1767 | "walbeck": {a:6376896.0, b:6355834.8467, ellipseName:"Walbeck"}, |
---|
1768 | "WGS60": {a:6378165.0, rf:298.3, ellipseName:"WGS 60"}, |
---|
1769 | "WGS66": {a:6378145.0, rf:298.25, ellipseName:"WGS 66"}, |
---|
1770 | "WGS72": {a:6378135.0, rf:298.26, ellipseName:"WGS 72"}, |
---|
1771 | "WGS84": {a:6378137.0, rf:298.257223563, ellipseName:"WGS 84"}, |
---|
1772 | "sphere": {a:6370997.0, b:6370997.0, ellipseName:"Normal Sphere (r=6370997)"} |
---|
1773 | }; |
---|
1774 | |
---|
1775 | Proj4js.Datum = { |
---|
1776 | "WGS84": {towgs84: "0,0,0", ellipse: "WGS84", datumName: "WGS84"}, |
---|
1777 | "GGRS87": {towgs84: "-199.87,74.79,246.62", ellipse: "GRS80", datumName: "Greek_Geodetic_Reference_System_1987"}, |
---|
1778 | "NAD83": {towgs84: "0,0,0", ellipse: "GRS80", datumName: "North_American_Datum_1983"}, |
---|
1779 | "NAD27": {nadgrids: "@conus,@alaska,@ntv2_0.gsb,@ntv1_can.dat", ellipse: "clrk66", datumName: "North_American_Datum_1927"}, |
---|
1780 | "potsdam": {towgs84: "606.0,23.0,413.0", ellipse: "bessel", datumName: "Potsdam Rauenberg 1950 DHDN"}, |
---|
1781 | "carthage": {towgs84: "-263.0,6.0,431.0", ellipse: "clark80", datumName: "Carthage 1934 Tunisia"}, |
---|
1782 | "hermannskogel": {towgs84: "653.0,-212.0,449.0", ellipse: "bessel", datumName: "Hermannskogel"}, |
---|
1783 | "ire65": {towgs84: "482.530,-130.596,564.557,-1.042,-0.214,-0.631,8.15", ellipse: "mod_airy", datumName: "Ireland 1965"}, |
---|
1784 | "nzgd49": {towgs84: "59.47,-5.04,187.44,0.47,-0.1,1.024,-4.5993", ellipse: "intl", datumName: "New Zealand Geodetic Datum 1949"}, |
---|
1785 | "OSGB36": {towgs84: "446.448,-125.157,542.060,0.1502,0.2470,0.8421,-20.4894", ellipse: "airy", datumName: "Airy 1830"} |
---|
1786 | }; |
---|
1787 | |
---|
1788 | Proj4js.WGS84 = new Proj4js.Proj('WGS84'); |
---|
1789 | Proj4js.Datum['OSB36'] = Proj4js.Datum['OSGB36']; //as returned from spatialreference.org |
---|
1790 | |
---|
1791 | //lookup table to go from the projection name in WKT to the Proj4js projection name |
---|
1792 | //build this out as required |
---|
1793 | Proj4js.wktProjections = { |
---|
1794 | "Lambert Tangential Conformal Conic Projection": "lcc", |
---|
1795 | "Mercator": "merc", |
---|
1796 | "Popular Visualisation Pseudo Mercator": "merc", |
---|
1797 | "Transverse_Mercator": "tmerc", |
---|
1798 | "Transverse Mercator": "tmerc", |
---|
1799 | "Lambert Azimuthal Equal Area": "laea", |
---|
1800 | "Universal Transverse Mercator System": "utm" |
---|
1801 | }; |
---|
1802 | |
---|
1803 | |
---|
1804 | /* ====================================================================== |
---|
1805 | projCode/aea.js |
---|
1806 | ====================================================================== */ |
---|
1807 | |
---|
1808 | /******************************************************************************* |
---|
1809 | NAME ALBERS CONICAL EQUAL AREA |
---|
1810 | |
---|
1811 | PURPOSE: Transforms input longitude and latitude to Easting and Northing |
---|
1812 | for the Albers Conical Equal Area projection. The longitude |
---|
1813 | and latitude must be in radians. The Easting and Northing |
---|
1814 | values will be returned in meters. |
---|
1815 | |
---|
1816 | PROGRAMMER DATE |
---|
1817 | ---------- ---- |
---|
1818 | T. Mittan, Feb, 1992 |
---|
1819 | |
---|
1820 | ALGORITHM REFERENCES |
---|
1821 | |
---|
1822 | 1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological |
---|
1823 | Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United |
---|
1824 | State Government Printing Office, Washington D.C., 1987. |
---|
1825 | |
---|
1826 | 2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections", |
---|
1827 | U.S. Geological Survey Professional Paper 1453 , United State Government |
---|
1828 | Printing Office, Washington D.C., 1989. |
---|
1829 | *******************************************************************************/ |
---|
1830 | |
---|
1831 | |
---|
1832 | Proj4js.Proj.aea = { |
---|
1833 | init : function() { |
---|
1834 | |
---|
1835 | if (Math.abs(this.lat1 + this.lat2) < Proj4js.common.EPSLN) { |
---|
1836 | Proj4js.reportError("aeaInitEqualLatitudes"); |
---|
1837 | return; |
---|
1838 | } |
---|
1839 | this.temp = this.b / this.a; |
---|
1840 | this.es = 1.0 - Math.pow(this.temp,2); |
---|
1841 | this.e3 = Math.sqrt(this.es); |
---|
1842 | |
---|
1843 | this.sin_po=Math.sin(this.lat1); |
---|
1844 | this.cos_po=Math.cos(this.lat1); |
---|
1845 | this.t1=this.sin_po; |
---|
1846 | this.con = this.sin_po; |
---|
1847 | this.ms1 = Proj4js.common.msfnz(this.e3,this.sin_po,this.cos_po); |
---|
1848 | this.qs1 = Proj4js.common.qsfnz(this.e3,this.sin_po,this.cos_po); |
---|
1849 | |
---|
1850 | this.sin_po=Math.sin(this.lat2); |
---|
1851 | this.cos_po=Math.cos(this.lat2); |
---|
1852 | this.t2=this.sin_po; |
---|
1853 | this.ms2 = Proj4js.common.msfnz(this.e3,this.sin_po,this.cos_po); |
---|
1854 | this.qs2 = Proj4js.common.qsfnz(this.e3,this.sin_po,this.cos_po); |
---|
1855 | |
---|
1856 | this.sin_po=Math.sin(this.lat0); |
---|
1857 | this.cos_po=Math.cos(this.lat0); |
---|
1858 | this.t3=this.sin_po; |
---|
1859 | this.qs0 = Proj4js.common.qsfnz(this.e3,this.sin_po,this.cos_po); |
---|
1860 | |
---|
1861 | if (Math.abs(this.lat1 - this.lat2) > Proj4js.common.EPSLN) { |
---|
1862 | this.ns0 = (this.ms1 * this.ms1 - this.ms2 *this.ms2)/ (this.qs2 - this.qs1); |
---|
1863 | } else { |
---|
1864 | this.ns0 = this.con; |
---|
1865 | } |
---|
1866 | this.c = this.ms1 * this.ms1 + this.ns0 * this.qs1; |
---|
1867 | this.rh = this.a * Math.sqrt(this.c - this.ns0 * this.qs0)/this.ns0; |
---|
1868 | }, |
---|
1869 | |
---|
1870 | /* Albers Conical Equal Area forward equations--mapping lat,long to x,y |
---|
1871 | -------------------------------------------------------------------*/ |
---|
1872 | forward: function(p){ |
---|
1873 | |
---|
1874 | var lon=p.x; |
---|
1875 | var lat=p.y; |
---|
1876 | |
---|
1877 | this.sin_phi=Math.sin(lat); |
---|
1878 | this.cos_phi=Math.cos(lat); |
---|
1879 | |
---|
1880 | var qs = Proj4js.common.qsfnz(this.e3,this.sin_phi,this.cos_phi); |
---|
1881 | var rh1 =this.a * Math.sqrt(this.c - this.ns0 * qs)/this.ns0; |
---|
1882 | var theta = this.ns0 * Proj4js.common.adjust_lon(lon - this.long0); |
---|
1883 | var x = rh1 * Math.sin(theta) + this.x0; |
---|
1884 | var y = this.rh - rh1 * Math.cos(theta) + this.y0; |
---|
1885 | |
---|
1886 | p.x = x; |
---|
1887 | p.y = y; |
---|
1888 | return p; |
---|
1889 | }, |
---|
1890 | |
---|
1891 | |
---|
1892 | inverse: function(p) { |
---|
1893 | var rh1,qs,con,theta,lon,lat; |
---|
1894 | |
---|
1895 | p.x -= this.x0; |
---|
1896 | p.y = this.rh - p.y + this.y0; |
---|
1897 | if (this.ns0 >= 0) { |
---|
1898 | rh1 = Math.sqrt(p.x *p.x + p.y * p.y); |
---|
1899 | con = 1.0; |
---|
1900 | } else { |
---|
1901 | rh1 = -Math.sqrt(p.x * p.x + p.y *p.y); |
---|
1902 | con = -1.0; |
---|
1903 | } |
---|
1904 | theta = 0.0; |
---|
1905 | if (rh1 != 0.0) { |
---|
1906 | theta = Math.atan2(con * p.x, con * p.y); |
---|
1907 | } |
---|
1908 | con = rh1 * this.ns0 / this.a; |
---|
1909 | qs = (this.c - con * con) / this.ns0; |
---|
1910 | if (this.e3 >= 1e-10) { |
---|
1911 | con = 1 - .5 * (1.0 -this.es) * Math.log((1.0 - this.e3) / (1.0 + this.e3))/this.e3; |
---|
1912 | if (Math.abs(Math.abs(con) - Math.abs(qs)) > .0000000001 ) { |
---|
1913 | lat = this.phi1z(this.e3,qs); |
---|
1914 | } else { |
---|
1915 | if (qs >= 0) { |
---|
1916 | lat = .5 * PI; |
---|
1917 | } else { |
---|
1918 | lat = -.5 * PI; |
---|
1919 | } |
---|
1920 | } |
---|
1921 | } else { |
---|
1922 | lat = this.phi1z(e3,qs); |
---|
1923 | } |
---|
1924 | |
---|
1925 | lon = Proj4js.common.adjust_lon(theta/this.ns0 + this.long0); |
---|
1926 | p.x = lon; |
---|
1927 | p.y = lat; |
---|
1928 | return p; |
---|
1929 | }, |
---|
1930 | |
---|
1931 | /* Function to compute phi1, the latitude for the inverse of the |
---|
1932 | Albers Conical Equal-Area projection. |
---|
1933 | -------------------------------------------*/ |
---|
1934 | phi1z: function (eccent,qs) { |
---|
1935 | var con, com, dphi; |
---|
1936 | var phi = Proj4js.common.asinz(.5 * qs); |
---|
1937 | if (eccent < Proj4js.common.EPSLN) return phi; |
---|
1938 | |
---|
1939 | var eccnts = eccent * eccent; |
---|
1940 | for (var i = 1; i <= 25; i++) { |
---|
1941 | sinphi = Math.sin(phi); |
---|
1942 | cosphi = Math.cos(phi); |
---|
1943 | con = eccent * sinphi; |
---|
1944 | com = 1.0 - con * con; |
---|
1945 | dphi = .5 * com * com / cosphi * (qs / (1.0 - eccnts) - sinphi / com + .5 / eccent * Math.log((1.0 - con) / (1.0 + con))); |
---|
1946 | phi = phi + dphi; |
---|
1947 | if (Math.abs(dphi) <= 1e-7) return phi; |
---|
1948 | } |
---|
1949 | Proj4js.reportError("aea:phi1z:Convergence error"); |
---|
1950 | return null; |
---|
1951 | } |
---|
1952 | |
---|
1953 | }; |
---|
1954 | |
---|
1955 | |
---|
1956 | |
---|
1957 | /* ====================================================================== |
---|
1958 | projCode/sterea.js |
---|
1959 | ====================================================================== */ |
---|
1960 | |
---|
1961 | |
---|
1962 | Proj4js.Proj.sterea = { |
---|
1963 | dependsOn : 'gauss', |
---|
1964 | |
---|
1965 | init : function() { |
---|
1966 | Proj4js.Proj['gauss'].init.apply(this); |
---|
1967 | if (!this.rc) { |
---|
1968 | Proj4js.reportError("sterea:init:E_ERROR_0"); |
---|
1969 | return; |
---|
1970 | } |
---|
1971 | this.sinc0 = Math.sin(this.phic0); |
---|
1972 | this.cosc0 = Math.cos(this.phic0); |
---|
1973 | this.R2 = 2.0 * this.rc; |
---|
1974 | if (!this.title) this.title = "Oblique Stereographic Alternative"; |
---|
1975 | }, |
---|
1976 | |
---|
1977 | forward : function(p) { |
---|
1978 | p.x = Proj4js.common.adjust_lon(p.x-this.long0); /* adjust del longitude */ |
---|
1979 | Proj4js.Proj['gauss'].forward.apply(this, [p]); |
---|
1980 | sinc = Math.sin(p.y); |
---|
1981 | cosc = Math.cos(p.y); |
---|
1982 | cosl = Math.cos(p.x); |
---|
1983 | k = this.k0 * this.R2 / (1.0 + this.sinc0 * sinc + this.cosc0 * cosc * cosl); |
---|
1984 | p.x = k * cosc * Math.sin(p.x); |
---|
1985 | p.y = k * (this.cosc0 * sinc - this.sinc0 * cosc * cosl); |
---|
1986 | p.x = this.a * p.x + this.x0; |
---|
1987 | p.y = this.a * p.y + this.y0; |
---|
1988 | return p; |
---|
1989 | }, |
---|
1990 | |
---|
1991 | inverse : function(p) { |
---|
1992 | var lon,lat; |
---|
1993 | p.x = (p.x - this.x0) / this.a; /* descale and de-offset */ |
---|
1994 | p.y = (p.y - this.y0) / this.a; |
---|
1995 | |
---|
1996 | p.x /= this.k0; |
---|
1997 | p.y /= this.k0; |
---|
1998 | if ( (rho = Math.sqrt(p.x*p.x + p.y*p.y)) ) { |
---|
1999 | c = 2.0 * Math.atan2(rho, this.R2); |
---|
2000 | sinc = Math.sin(c); |
---|
2001 | cosc = Math.cos(c); |
---|
2002 | lat = Math.asin(cosc * this.sinc0 + p.y * sinc * this.cosc0 / rho); |
---|
2003 | lon = Math.atan2(p.x * sinc, rho * this.cosc0 * cosc - p.y * this.sinc0 * sinc); |
---|
2004 | } else { |
---|
2005 | lat = this.phic0; |
---|
2006 | lon = 0.; |
---|
2007 | } |
---|
2008 | |
---|
2009 | p.x = lon; |
---|
2010 | p.y = lat; |
---|
2011 | Proj4js.Proj['gauss'].inverse.apply(this,[p]); |
---|
2012 | p.x = Proj4js.common.adjust_lon(p.x + this.long0); /* adjust longitude to CM */ |
---|
2013 | return p; |
---|
2014 | } |
---|
2015 | }; |
---|
2016 | |
---|
2017 | /* ====================================================================== |
---|
2018 | projCode/poly.js |
---|
2019 | ====================================================================== */ |
---|
2020 | |
---|
2021 | /* Function to compute, phi4, the latitude for the inverse of the |
---|
2022 | Polyconic projection. |
---|
2023 | ------------------------------------------------------------*/ |
---|
2024 | function phi4z (eccent,e0,e1,e2,e3,a,b,c,phi) { |
---|
2025 | var sinphi, sin2ph, tanph, ml, mlp, con1, con2, con3, dphi, i; |
---|
2026 | |
---|
2027 | phi = a; |
---|
2028 | for (i = 1; i <= 15; i++) { |
---|
2029 | sinphi = Math.sin(phi); |
---|
2030 | tanphi = Math.tan(phi); |
---|
2031 | c = tanphi * Math.sqrt (1.0 - eccent * sinphi * sinphi); |
---|
2032 | sin2ph = Math.sin (2.0 * phi); |
---|
2033 | /* |
---|
2034 | ml = e0 * *phi - e1 * sin2ph + e2 * sin (4.0 * *phi); |
---|
2035 | mlp = e0 - 2.0 * e1 * cos (2.0 * *phi) + 4.0 * e2 * cos (4.0 * *phi); |
---|
2036 | */ |
---|
2037 | ml = e0 * phi - e1 * sin2ph + e2 * Math.sin (4.0 * phi) - e3 * Math.sin (6.0 * phi); |
---|
2038 | mlp = e0 - 2.0 * e1 * Math.cos (2.0 * phi) + 4.0 * e2 * Math.cos (4.0 * phi) - 6.0 * e3 * Math.cos (6.0 * phi); |
---|
2039 | con1 = 2.0 * ml + c * (ml * ml + b) - 2.0 * a * (c * ml + 1.0); |
---|
2040 | con2 = eccent * sin2ph * (ml * ml + b - 2.0 * a * ml) / (2.0 *c); |
---|
2041 | con3 = 2.0 * (a - ml) * (c * mlp - 2.0 / sin2ph) - 2.0 * mlp; |
---|
2042 | dphi = con1 / (con2 + con3); |
---|
2043 | phi += dphi; |
---|
2044 | if (Math.abs(dphi) <= .0000000001 ) return(phi); |
---|
2045 | } |
---|
2046 | Proj4js.reportError("phi4z: No convergence"); |
---|
2047 | return null; |
---|
2048 | } |
---|
2049 | |
---|
2050 | |
---|
2051 | /* Function to compute the constant e4 from the input of the eccentricity |
---|
2052 | of the spheroid, x. This constant is used in the Polar Stereographic |
---|
2053 | projection. |
---|
2054 | --------------------------------------------------------------------*/ |
---|
2055 | function e4fn(x) { |
---|
2056 | var con, com; |
---|
2057 | con = 1.0 + x; |
---|
2058 | com = 1.0 - x; |
---|
2059 | return (Math.sqrt((Math.pow(con,con))*(Math.pow(com,com)))); |
---|
2060 | } |
---|
2061 | |
---|
2062 | |
---|
2063 | |
---|
2064 | |
---|
2065 | |
---|
2066 | /******************************************************************************* |
---|
2067 | NAME POLYCONIC |
---|
2068 | |
---|
2069 | PURPOSE: Transforms input longitude and latitude to Easting and |
---|
2070 | Northing for the Polyconic projection. The |
---|
2071 | longitude and latitude must be in radians. The Easting |
---|
2072 | and Northing values will be returned in meters. |
---|
2073 | |
---|
2074 | PROGRAMMER DATE |
---|
2075 | ---------- ---- |
---|
2076 | T. Mittan Mar, 1993 |
---|
2077 | |
---|
2078 | ALGORITHM REFERENCES |
---|
2079 | |
---|
2080 | 1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological |
---|
2081 | Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United |
---|
2082 | State Government Printing Office, Washington D.C., 1987. |
---|
2083 | |
---|
2084 | 2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections", |
---|
2085 | U.S. Geological Survey Professional Paper 1453 , United State Government |
---|
2086 | Printing Office, Washington D.C., 1989. |
---|
2087 | *******************************************************************************/ |
---|
2088 | |
---|
2089 | Proj4js.Proj.poly = { |
---|
2090 | |
---|
2091 | /* Initialize the POLYCONIC projection |
---|
2092 | ----------------------------------*/ |
---|
2093 | init: function() { |
---|
2094 | var temp; /* temporary variable */ |
---|
2095 | if (this.lat0=0) this.lat0=90;//this.lat0 ca |
---|
2096 | |
---|
2097 | /* Place parameters in static storage for common use |
---|
2098 | -------------------------------------------------*/ |
---|
2099 | this.temp = this.b / this.a; |
---|
2100 | this.es = 1.0 - Math.pow(this.temp,2);// devait etre dans tmerc.js mais n y est pas donc je commente sinon retour de valeurs nulles |
---|
2101 | this.e = Math.sqrt(this.es); |
---|
2102 | this.e0 = Proj4js.common.e0fn(this.es); |
---|
2103 | this.e1 = Proj4js.common.e1fn(this.es); |
---|
2104 | this.e2 = Proj4js.common.e2fn(this.es); |
---|
2105 | this.e3 = Proj4js.common.e3fn(this.es); |
---|
2106 | this.ml0 = Proj4js.common.mlfn(this.e0, this.e1,this.e2, this.e3, this.lat0);//si que des zeros le calcul ne se fait pas |
---|
2107 | //if (!this.ml0) {this.ml0=0;} |
---|
2108 | }, |
---|
2109 | |
---|
2110 | |
---|
2111 | /* Polyconic forward equations--mapping lat,long to x,y |
---|
2112 | ---------------------------------------------------*/ |
---|
2113 | forward: function(p) { |
---|
2114 | var sinphi, cosphi; /* sin and cos value */ |
---|
2115 | var al; /* temporary values */ |
---|
2116 | var c; /* temporary values */ |
---|
2117 | var con, ml; /* cone constant, small m */ |
---|
2118 | var ms; /* small m */ |
---|
2119 | var x,y; |
---|
2120 | |
---|
2121 | var lon=p.x; |
---|
2122 | var lat=p.y; |
---|
2123 | |
---|
2124 | con = Proj4js.common.adjust_lon(lon - this.long0); |
---|
2125 | if (Math.abs(lat) <= .0000001) { |
---|
2126 | x = this.x0 + this.a * con; |
---|
2127 | y = this.y0 - this.a * this.ml0; |
---|
2128 | } else { |
---|
2129 | sinphi = Math.sin(lat); |
---|
2130 | cosphi = Math.cos(lat); |
---|
2131 | |
---|
2132 | ml = Proj4js.common.mlfn(this.e0, this.e1, this.e2, this.e3, lat); |
---|
2133 | ms = Proj4js.common.msfnz(this.e,sinphi,cosphi); |
---|
2134 | con = sinphi; |
---|
2135 | x = this.x0 + this.a * ms * Math.sin(con)/sinphi; |
---|
2136 | y = this.y0 + this.a * (ml - this.ml0 + ms * (1.0 - Math.cos(con))/sinphi); |
---|
2137 | } |
---|
2138 | |
---|
2139 | p.x=x; |
---|
2140 | p.y=y; |
---|
2141 | return p; |
---|
2142 | }, |
---|
2143 | |
---|
2144 | |
---|
2145 | /* Inverse equations |
---|
2146 | -----------------*/ |
---|
2147 | inverse: function(p) { |
---|
2148 | var sin_phi, cos_phi; /* sin and cos value */ |
---|
2149 | var al; /* temporary values */ |
---|
2150 | var b; /* temporary values */ |
---|
2151 | var c; /* temporary values */ |
---|
2152 | var con, ml; /* cone constant, small m */ |
---|
2153 | var iflg; /* error flag */ |
---|
2154 | var lon,lat; |
---|
2155 | p.x -= this.x0; |
---|
2156 | p.y -= this.y0; |
---|
2157 | al = this.ml0 + p.y/this.a; |
---|
2158 | iflg = 0; |
---|
2159 | |
---|
2160 | if (Math.abs(al) <= .0000001) { |
---|
2161 | lon = p.x/this.a + this.long0; |
---|
2162 | lat = 0.0; |
---|
2163 | } else { |
---|
2164 | b = al * al + (p.x/this.a) * (p.x/this.a); |
---|
2165 | iflg = phi4z(this.es,this.e0,this.e1,this.e2,this.e3,this.al,b,c,lat); |
---|
2166 | if (iflg != 1) return(iflg); |
---|
2167 | lon = Proj4js.common.adjust_lon((Proj4js.common.asinz(p.x * c / this.a) / Math.sin(lat)) + this.long0); |
---|
2168 | } |
---|
2169 | |
---|
2170 | p.x=lon; |
---|
2171 | p.y=lat; |
---|
2172 | return p; |
---|
2173 | } |
---|
2174 | }; |
---|
2175 | |
---|
2176 | |
---|
2177 | |
---|
2178 | /* ====================================================================== |
---|
2179 | projCode/equi.js |
---|
2180 | ====================================================================== */ |
---|
2181 | |
---|
2182 | /******************************************************************************* |
---|
2183 | NAME EQUIRECTANGULAR |
---|
2184 | |
---|
2185 | PURPOSE: Transforms input longitude and latitude to Easting and |
---|
2186 | Northing for the Equirectangular projection. The |
---|
2187 | longitude and latitude must be in radians. The Easting |
---|
2188 | and Northing values will be returned in meters. |
---|
2189 | |
---|
2190 | PROGRAMMER DATE |
---|
2191 | ---------- ---- |
---|
2192 | T. Mittan Mar, 1993 |
---|
2193 | |
---|
2194 | ALGORITHM REFERENCES |
---|
2195 | |
---|
2196 | 1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological |
---|
2197 | Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United |
---|
2198 | State Government Printing Office, Washington D.C., 1987. |
---|
2199 | |
---|
2200 | 2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections", |
---|
2201 | U.S. Geological Survey Professional Paper 1453 , United State Government |
---|
2202 | Printing Office, Washington D.C., 1989. |
---|
2203 | *******************************************************************************/ |
---|
2204 | Proj4js.Proj.equi = { |
---|
2205 | |
---|
2206 | init: function() { |
---|
2207 | if(!this.x0) this.x0=0; |
---|
2208 | if(!this.y0) this.y0=0; |
---|
2209 | if(!this.lat0) this.lat0=0; |
---|
2210 | if(!this.long0) this.long0=0; |
---|
2211 | ///this.t2; |
---|
2212 | }, |
---|
2213 | |
---|
2214 | |
---|
2215 | |
---|
2216 | /* Equirectangular forward equations--mapping lat,long to x,y |
---|
2217 | ---------------------------------------------------------*/ |
---|
2218 | forward: function(p) { |
---|
2219 | |
---|
2220 | var lon=p.x; |
---|
2221 | var lat=p.y; |
---|
2222 | |
---|
2223 | var dlon = Proj4js.common.adjust_lon(lon - this.long0); |
---|
2224 | var x = this.x0 +this. a * dlon *Math.cos(this.lat0); |
---|
2225 | var y = this.y0 + this.a * lat; |
---|
2226 | |
---|
2227 | this.t1=x; |
---|
2228 | this.t2=Math.cos(this.lat0); |
---|
2229 | p.x=x; |
---|
2230 | p.y=y; |
---|
2231 | return p; |
---|
2232 | }, //equiFwd() |
---|
2233 | |
---|
2234 | |
---|
2235 | |
---|
2236 | /* Equirectangular inverse equations--mapping x,y to lat/long |
---|
2237 | ---------------------------------------------------------*/ |
---|
2238 | inverse: function(p) { |
---|
2239 | |
---|
2240 | p.x -= this.x0; |
---|
2241 | p.y -= this.y0; |
---|
2242 | var lat = p.y /this. a; |
---|
2243 | |
---|
2244 | if ( Math.abs(lat) > Proj4js.common.HALF_PI) { |
---|
2245 | Proj4js.reportError("equi:Inv:DataError"); |
---|
2246 | } |
---|
2247 | var lon = Proj4js.common.adjust_lon(this.long0 + p.x / (this.a * Math.cos(this.lat0))); |
---|
2248 | p.x=lon; |
---|
2249 | p.y=lat; |
---|
2250 | }//equiInv() |
---|
2251 | }; |
---|
2252 | |
---|
2253 | |
---|
2254 | /* ====================================================================== |
---|
2255 | projCode/merc.js |
---|
2256 | ====================================================================== */ |
---|
2257 | |
---|
2258 | /******************************************************************************* |
---|
2259 | NAME MERCATOR |
---|
2260 | |
---|
2261 | PURPOSE: Transforms input longitude and latitude to Easting and |
---|
2262 | Northing for the Mercator projection. The |
---|
2263 | longitude and latitude must be in radians. The Easting |
---|
2264 | and Northing values will be returned in meters. |
---|
2265 | |
---|
2266 | PROGRAMMER DATE |
---|
2267 | ---------- ---- |
---|
2268 | D. Steinwand, EROS Nov, 1991 |
---|
2269 | T. Mittan Mar, 1993 |
---|
2270 | |
---|
2271 | ALGORITHM REFERENCES |
---|
2272 | |
---|
2273 | 1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological |
---|
2274 | Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United |
---|
2275 | State Government Printing Office, Washington D.C., 1987. |
---|
2276 | |
---|
2277 | 2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections", |
---|
2278 | U.S. Geological Survey Professional Paper 1453 , United State Government |
---|
2279 | Printing Office, Washington D.C., 1989. |
---|
2280 | *******************************************************************************/ |
---|
2281 | |
---|
2282 | //static double r_major = a; /* major axis */ |
---|
2283 | //static double r_minor = b; /* minor axis */ |
---|
2284 | //static double lon_center = long0; /* Center longitude (projection center) */ |
---|
2285 | //static double lat_origin = lat0; /* center latitude */ |
---|
2286 | //static double e,es; /* eccentricity constants */ |
---|
2287 | //static double m1; /* small value m */ |
---|
2288 | //static double false_northing = y0; /* y offset in meters */ |
---|
2289 | //static double false_easting = x0; /* x offset in meters */ |
---|
2290 | //scale_fact = k0 |
---|
2291 | |
---|
2292 | Proj4js.Proj.merc = { |
---|
2293 | init : function() { |
---|
2294 | //?this.temp = this.r_minor / this.r_major; |
---|
2295 | //this.temp = this.b / this.a; |
---|
2296 | //this.es = 1.0 - Math.sqrt(this.temp); |
---|
2297 | //this.e = Math.sqrt( this.es ); |
---|
2298 | //?this.m1 = Math.cos(this.lat_origin) / (Math.sqrt( 1.0 - this.es * Math.sin(this.lat_origin) * Math.sin(this.lat_origin))); |
---|
2299 | //this.m1 = Math.cos(0.0) / (Math.sqrt( 1.0 - this.es * Math.sin(0.0) * Math.sin(0.0))); |
---|
2300 | if (this.lat_ts) { |
---|
2301 | if (this.sphere) { |
---|
2302 | this.k0 = Math.cos(this.lat_ts); |
---|
2303 | } else { |
---|
2304 | this.k0 = Proj4js.common.msfnz(this.es, Math.sin(this.lat_ts), Math.cos(this.lat_ts)); |
---|
2305 | } |
---|
2306 | } |
---|
2307 | }, |
---|
2308 | |
---|
2309 | /* Mercator forward equations--mapping lat,long to x,y |
---|
2310 | --------------------------------------------------*/ |
---|
2311 | |
---|
2312 | forward : function(p) { |
---|
2313 | //alert("ll2m coords : "+coords); |
---|
2314 | var lon = p.x; |
---|
2315 | var lat = p.y; |
---|
2316 | // convert to radians |
---|
2317 | if ( lat*Proj4js.common.R2D > 90.0 && |
---|
2318 | lat*Proj4js.common.R2D < -90.0 && |
---|
2319 | lon*Proj4js.common.R2D > 180.0 && |
---|
2320 | lon*Proj4js.common.R2D < -180.0) { |
---|
2321 | Proj4js.reportError("merc:forward: llInputOutOfRange: "+ lon +" : " + lat); |
---|
2322 | return null; |
---|
2323 | } |
---|
2324 | |
---|
2325 | var x,y; |
---|
2326 | if(Math.abs( Math.abs(lat) - Proj4js.common.HALF_PI) <= Proj4js.common.EPSLN) { |
---|
2327 | Proj4js.reportError("merc:forward: ll2mAtPoles"); |
---|
2328 | return null; |
---|
2329 | } else { |
---|
2330 | if (this.sphere) { |
---|
2331 | x = this.x0 + this.a * this.k0 * Proj4js.common.adjust_lon(lon - this.long0); |
---|
2332 | y = this.y0 + this.a * this.k0 * Math.log(Math.tan(Proj4js.common.FORTPI + 0.5*lat)); |
---|
2333 | } else { |
---|
2334 | var sinphi = Math.sin(lat); |
---|
2335 | var ts = Proj4js.common.tsfnz(this.e,lat,sinphi); |
---|
2336 | x = this.x0 + this.a * this.k0 * Proj4js.common.adjust_lon(lon - this.long0); |
---|
2337 | y = this.y0 - this.a * this.k0 * Math.log(ts); |
---|
2338 | } |
---|
2339 | p.x = x; |
---|
2340 | p.y = y; |
---|
2341 | return p; |
---|
2342 | } |
---|
2343 | }, |
---|
2344 | |
---|
2345 | |
---|
2346 | /* Mercator inverse equations--mapping x,y to lat/long |
---|
2347 | --------------------------------------------------*/ |
---|
2348 | inverse : function(p) { |
---|
2349 | |
---|
2350 | var x = p.x - this.x0; |
---|
2351 | var y = p.y - this.y0; |
---|
2352 | var lon,lat; |
---|
2353 | |
---|
2354 | if (this.sphere) { |
---|
2355 | lat = Proj4js.common.HALF_PI - 2.0 * Math.atan(Math.exp(-y / this.a * this.k0)); |
---|
2356 | } else { |
---|
2357 | var ts = Math.exp(-y / (this.a * this.k0)); |
---|
2358 | lat = Proj4js.common.phi2z(this.e,ts); |
---|
2359 | if(lat == -9999) { |
---|
2360 | Proj4js.reportError("merc:inverse: lat = -9999"); |
---|
2361 | return null; |
---|
2362 | } |
---|
2363 | } |
---|
2364 | lon = Proj4js.common.adjust_lon(this.long0+ x / (this.a * this.k0)); |
---|
2365 | |
---|
2366 | p.x = lon; |
---|
2367 | p.y = lat; |
---|
2368 | return p; |
---|
2369 | } |
---|
2370 | }; |
---|
2371 | |
---|
2372 | |
---|
2373 | /* ====================================================================== |
---|
2374 | projCode/utm.js |
---|
2375 | ====================================================================== */ |
---|
2376 | |
---|
2377 | /******************************************************************************* |
---|
2378 | NAME TRANSVERSE MERCATOR |
---|
2379 | |
---|
2380 | PURPOSE: Transforms input longitude and latitude to Easting and |
---|
2381 | Northing for the Transverse Mercator projection. The |
---|
2382 | longitude and latitude must be in radians. The Easting |
---|
2383 | and Northing values will be returned in meters. |
---|
2384 | |
---|
2385 | ALGORITHM REFERENCES |
---|
2386 | |
---|
2387 | 1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological |
---|
2388 | Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United |
---|
2389 | State Government Printing Office, Washington D.C., 1987. |
---|
2390 | |
---|
2391 | 2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections", |
---|
2392 | U.S. Geological Survey Professional Paper 1453 , United State Government |
---|
2393 | Printing Office, Washington D.C., 1989. |
---|
2394 | *******************************************************************************/ |
---|
2395 | |
---|
2396 | |
---|
2397 | /** |
---|
2398 | Initialize Transverse Mercator projection |
---|
2399 | */ |
---|
2400 | |
---|
2401 | Proj4js.Proj.utm = { |
---|
2402 | dependsOn : 'tmerc', |
---|
2403 | |
---|
2404 | init : function() { |
---|
2405 | if (!this.zone) { |
---|
2406 | Proj4js.reportError("utm:init: zone must be specified for UTM"); |
---|
2407 | return; |
---|
2408 | } |
---|
2409 | this.lat0 = 0.0; |
---|
2410 | this.long0 = ((6 * Math.abs(this.zone)) - 183) * Proj4js.common.D2R; |
---|
2411 | this.x0 = 500000.0; |
---|
2412 | this.y0 = this.utmSouth ? 10000000.0 : 0.0; |
---|
2413 | this.k0 = 0.9996; |
---|
2414 | |
---|
2415 | Proj4js.Proj['tmerc'].init.apply(this); |
---|
2416 | this.forward = Proj4js.Proj['tmerc'].forward; |
---|
2417 | this.inverse = Proj4js.Proj['tmerc'].inverse; |
---|
2418 | } |
---|
2419 | }; |
---|
2420 | /* ====================================================================== |
---|
2421 | projCode/eqdc.js |
---|
2422 | ====================================================================== */ |
---|
2423 | |
---|
2424 | /******************************************************************************* |
---|
2425 | NAME EQUIDISTANT CONIC |
---|
2426 | |
---|
2427 | PURPOSE: Transforms input longitude and latitude to Easting and Northing |
---|
2428 | for the Equidistant Conic projection. The longitude and |
---|
2429 | latitude must be in radians. The Easting and Northing values |
---|
2430 | will be returned in meters. |
---|
2431 | |
---|
2432 | PROGRAMMER DATE |
---|
2433 | ---------- ---- |
---|
2434 | T. Mittan Mar, 1993 |
---|
2435 | |
---|
2436 | ALGORITHM REFERENCES |
---|
2437 | |
---|
2438 | 1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological |
---|
2439 | Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United |
---|
2440 | State Government Printing Office, Washington D.C., 1987. |
---|
2441 | |
---|
2442 | 2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections", |
---|
2443 | U.S. Geological Survey Professional Paper 1453 , United State Government |
---|
2444 | Printing Office, Washington D.C., 1989. |
---|
2445 | *******************************************************************************/ |
---|
2446 | |
---|
2447 | /* Variables common to all subroutines in this code file |
---|
2448 | -----------------------------------------------------*/ |
---|
2449 | |
---|
2450 | Proj4js.Proj.eqdc = { |
---|
2451 | |
---|
2452 | /* Initialize the Equidistant Conic projection |
---|
2453 | ------------------------------------------*/ |
---|
2454 | init: function() { |
---|
2455 | |
---|
2456 | /* Place parameters in static storage for common use |
---|
2457 | -------------------------------------------------*/ |
---|
2458 | |
---|
2459 | if(!this.mode) this.mode=0;//chosen default mode |
---|
2460 | this.temp = this.b / this.a; |
---|
2461 | this.es = 1.0 - Math.pow(this.temp,2); |
---|
2462 | this.e = Math.sqrt(this.es); |
---|
2463 | this.e0 = Proj4js.common.e0fn(this.es); |
---|
2464 | this.e1 = Proj4js.common.e1fn(this.es); |
---|
2465 | this.e2 = Proj4js.common.e2fn(this.es); |
---|
2466 | this.e3 = Proj4js.common.e3fn(this.es); |
---|
2467 | |
---|
2468 | this.sinphi=Math.sin(this.lat1); |
---|
2469 | this.cosphi=Math.cos(this.lat1); |
---|
2470 | |
---|
2471 | this.ms1 = Proj4js.common.msfnz(this.e,this.sinphi,this.cosphi); |
---|
2472 | this.ml1 = Proj4js.common.mlfn(this.e0, this.e1, this.e2,this.e3, this.lat1); |
---|
2473 | |
---|
2474 | /* format B |
---|
2475 | ---------*/ |
---|
2476 | if (this.mode != 0) { |
---|
2477 | if (Math.abs(this.lat1 + this.lat2) < Proj4js.common.EPSLN) { |
---|
2478 | Proj4js.reportError("eqdc:Init:EqualLatitudes"); |
---|
2479 | //return(81); |
---|
2480 | } |
---|
2481 | this.sinphi=Math.sin(this.lat2); |
---|
2482 | this.cosphi=Math.cos(this.lat2); |
---|
2483 | |
---|
2484 | this.ms2 = Proj4js.common.msfnz(this.e,this.sinphi,this.cosphi); |
---|
2485 | this.ml2 = Proj4js.common.mlfn(this.e0, this.e1, this.e2, this.e3, this.lat2); |
---|
2486 | if (Math.abs(this.lat1 - this.lat2) >= Proj4js.common.EPSLN) { |
---|
2487 | this.ns = (this.ms1 - this.ms2) / (this.ml2 - this.ml1); |
---|
2488 | } else { |
---|
2489 | this.ns = this.sinphi; |
---|
2490 | } |
---|
2491 | } else { |
---|
2492 | this.ns = this.sinphi; |
---|
2493 | } |
---|
2494 | this.g = this.ml1 + this.ms1/this.ns; |
---|
2495 | this.ml0 = Proj4js.common.mlfn(this.e0, this.e1,this. e2, this.e3, this.lat0); |
---|
2496 | this.rh = this.a * (this.g - this.ml0); |
---|
2497 | }, |
---|
2498 | |
---|
2499 | |
---|
2500 | /* Equidistant Conic forward equations--mapping lat,long to x,y |
---|
2501 | -----------------------------------------------------------*/ |
---|
2502 | forward: function(p) { |
---|
2503 | var lon=p.x; |
---|
2504 | var lat=p.y; |
---|
2505 | |
---|
2506 | /* Forward equations |
---|
2507 | -----------------*/ |
---|
2508 | var ml = Proj4js.common.mlfn(this.e0, this.e1, this.e2, this.e3, lat); |
---|
2509 | var rh1 = this.a * (this.g - ml); |
---|
2510 | var theta = this.ns * Proj4js.common.adjust_lon(lon - this.long0); |
---|
2511 | |
---|
2512 | var x = this.x0 + rh1 * Math.sin(theta); |
---|
2513 | var y = this.y0 + this.rh - rh1 * Math.cos(theta); |
---|
2514 | p.x=x; |
---|
2515 | p.y=y; |
---|
2516 | return p; |
---|
2517 | }, |
---|
2518 | |
---|
2519 | /* Inverse equations |
---|
2520 | -----------------*/ |
---|
2521 | inverse: function(p) { |
---|
2522 | p.x -= this.x0; |
---|
2523 | p.y = this.rh - p.y + this.y0; |
---|
2524 | var con, rh1; |
---|
2525 | if (this.ns >= 0) { |
---|
2526 | var rh1 = Math.sqrt(p.x *p.x + p.y * p.y); |
---|
2527 | var con = 1.0; |
---|
2528 | } else { |
---|
2529 | rh1 = -Math.sqrt(p.x *p. x +p. y * p.y); |
---|
2530 | con = -1.0; |
---|
2531 | } |
---|
2532 | var theta = 0.0; |
---|
2533 | if (rh1 != 0.0) theta = Math.atan2(con *p.x, con *p.y); |
---|
2534 | var ml = this.g - rh1 /this.a; |
---|
2535 | var lat = this.phi3z(ml,this.e0,this.e1,this.e2,this.e3); |
---|
2536 | var lon = Proj4js.common.adjust_lon(this.long0 + theta / this.ns); |
---|
2537 | |
---|
2538 | p.x=lon; |
---|
2539 | p.y=lat; |
---|
2540 | return p; |
---|
2541 | }, |
---|
2542 | |
---|
2543 | /* Function to compute latitude, phi3, for the inverse of the Equidistant |
---|
2544 | Conic projection. |
---|
2545 | -----------------------------------------------------------------*/ |
---|
2546 | phi3z: function(ml,e0,e1,e2,e3) { |
---|
2547 | var phi; |
---|
2548 | var dphi; |
---|
2549 | |
---|
2550 | phi = ml; |
---|
2551 | for (var i = 0; i < 15; i++) { |
---|
2552 | dphi = (ml + e1 * Math.sin(2.0 * phi) - e2 * Math.sin(4.0 * phi) + e3 * Math.sin(6.0 * phi))/ e0 - phi; |
---|
2553 | phi += dphi; |
---|
2554 | if (Math.abs(dphi) <= .0000000001) { |
---|
2555 | return phi; |
---|
2556 | } |
---|
2557 | } |
---|
2558 | Proj4js.reportError("PHI3Z-CONV:Latitude failed to converge after 15 iterations"); |
---|
2559 | return null; |
---|
2560 | } |
---|
2561 | |
---|
2562 | |
---|
2563 | }; |
---|
2564 | /* ====================================================================== |
---|
2565 | projCode/tmerc.js |
---|
2566 | ====================================================================== */ |
---|
2567 | |
---|
2568 | /******************************************************************************* |
---|
2569 | NAME TRANSVERSE MERCATOR |
---|
2570 | |
---|
2571 | PURPOSE: Transforms input longitude and latitude to Easting and |
---|
2572 | Northing for the Transverse Mercator projection. The |
---|
2573 | longitude and latitude must be in radians. The Easting |
---|
2574 | and Northing values will be returned in meters. |
---|
2575 | |
---|
2576 | ALGORITHM REFERENCES |
---|
2577 | |
---|
2578 | 1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological |
---|
2579 | Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United |
---|
2580 | State Government Printing Office, Washington D.C., 1987. |
---|
2581 | |
---|
2582 | 2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections", |
---|
2583 | U.S. Geological Survey Professional Paper 1453 , United State Government |
---|
2584 | Printing Office, Washington D.C., 1989. |
---|
2585 | *******************************************************************************/ |
---|
2586 | |
---|
2587 | |
---|
2588 | /** |
---|
2589 | Initialize Transverse Mercator projection |
---|
2590 | */ |
---|
2591 | |
---|
2592 | Proj4js.Proj.tmerc = { |
---|
2593 | init : function() { |
---|
2594 | this.e0 = Proj4js.common.e0fn(this.es); |
---|
2595 | this.e1 = Proj4js.common.e1fn(this.es); |
---|
2596 | this.e2 = Proj4js.common.e2fn(this.es); |
---|
2597 | this.e3 = Proj4js.common.e3fn(this.es); |
---|
2598 | this.ml0 = this.a * Proj4js.common.mlfn(this.e0, this.e1, this.e2, this.e3, this.lat0); |
---|
2599 | }, |
---|
2600 | |
---|
2601 | /** |
---|
2602 | Transverse Mercator Forward - long/lat to x/y |
---|
2603 | long/lat in radians |
---|
2604 | */ |
---|
2605 | forward : function(p) { |
---|
2606 | var lon = p.x; |
---|
2607 | var lat = p.y; |
---|
2608 | |
---|
2609 | var delta_lon = Proj4js.common.adjust_lon(lon - this.long0); // Delta longitude |
---|
2610 | var con; // cone constant |
---|
2611 | var x, y; |
---|
2612 | var sin_phi=Math.sin(lat); |
---|
2613 | var cos_phi=Math.cos(lat); |
---|
2614 | |
---|
2615 | if (this.sphere) { /* spherical form */ |
---|
2616 | var b = cos_phi * Math.sin(delta_lon); |
---|
2617 | if ((Math.abs(Math.abs(b) - 1.0)) < .0000000001) { |
---|
2618 | Proj4js.reportError("tmerc:forward: Point projects into infinity"); |
---|
2619 | return(93); |
---|
2620 | } else { |
---|
2621 | x = .5 * this.a * this.k0 * Math.log((1.0 + b)/(1.0 - b)); |
---|
2622 | con = Math.acos(cos_phi * Math.cos(delta_lon)/Math.sqrt(1.0 - b*b)); |
---|
2623 | if (lat < 0) con = - con; |
---|
2624 | y = this.a * this.k0 * (con - this.lat0); |
---|
2625 | } |
---|
2626 | } else { |
---|
2627 | var al = cos_phi * delta_lon; |
---|
2628 | var als = Math.pow(al,2); |
---|
2629 | var c = this.ep2 * Math.pow(cos_phi,2); |
---|
2630 | var tq = Math.tan(lat); |
---|
2631 | var t = Math.pow(tq,2); |
---|
2632 | con = 1.0 - this.es * Math.pow(sin_phi,2); |
---|
2633 | var n = this.a / Math.sqrt(con); |
---|
2634 | var ml = this.a * Proj4js.common.mlfn(this.e0, this.e1, this.e2, this.e3, lat); |
---|
2635 | |
---|
2636 | x = this.k0 * n * al * (1.0 + als / 6.0 * (1.0 - t + c + als / 20.0 * (5.0 - 18.0 * t + Math.pow(t,2) + 72.0 * c - 58.0 * this.ep2))) + this.x0; |
---|
2637 | y = this.k0 * (ml - this.ml0 + n * tq * (als * (0.5 + als / 24.0 * (5.0 - t + 9.0 * c + 4.0 * Math.pow(c,2) + als / 30.0 * (61.0 - 58.0 * t + Math.pow(t,2) + 600.0 * c - 330.0 * this.ep2))))) + this.y0; |
---|
2638 | |
---|
2639 | } |
---|
2640 | p.x = x; p.y = y; |
---|
2641 | return p; |
---|
2642 | }, // tmercFwd() |
---|
2643 | |
---|
2644 | /** |
---|
2645 | Transverse Mercator Inverse - x/y to long/lat |
---|
2646 | */ |
---|
2647 | inverse : function(p) { |
---|
2648 | var con, phi; /* temporary angles */ |
---|
2649 | var delta_phi; /* difference between longitudes */ |
---|
2650 | var i; |
---|
2651 | var max_iter = 6; /* maximun number of iterations */ |
---|
2652 | var lat, lon; |
---|
2653 | |
---|
2654 | if (this.sphere) { /* spherical form */ |
---|
2655 | var f = Math.exp(p.x/(this.a * this.k0)); |
---|
2656 | var g = .5 * (f - 1/f); |
---|
2657 | var temp = this.lat0 + p.y/(this.a * this.k0); |
---|
2658 | var h = Math.cos(temp); |
---|
2659 | con = Math.sqrt((1.0 - h * h)/(1.0 + g * g)); |
---|
2660 | lat = Proj4js.common.asinz(con); |
---|
2661 | if (temp < 0) |
---|
2662 | lat = -lat; |
---|
2663 | if ((g == 0) && (h == 0)) { |
---|
2664 | lon = this.long0; |
---|
2665 | } else { |
---|
2666 | lon = Proj4js.common.adjust_lon(Math.atan2(g,h) + this.long0); |
---|
2667 | } |
---|
2668 | } else { // ellipsoidal form |
---|
2669 | var x = p.x - this.x0; |
---|
2670 | var y = p.y - this.y0; |
---|
2671 | |
---|
2672 | con = (this.ml0 + y / this.k0) / this.a; |
---|
2673 | phi = con; |
---|
2674 | for (i=0;true;i++) { |
---|
2675 | delta_phi=((con + this.e1 * Math.sin(2.0*phi) - this.e2 * Math.sin(4.0*phi) + this.e3 * Math.sin(6.0*phi)) / this.e0) - phi; |
---|
2676 | phi += delta_phi; |
---|
2677 | if (Math.abs(delta_phi) <= Proj4js.common.EPSLN) break; |
---|
2678 | if (i >= max_iter) { |
---|
2679 | Proj4js.reportError("tmerc:inverse: Latitude failed to converge"); |
---|
2680 | return(95); |
---|
2681 | } |
---|
2682 | } // for() |
---|
2683 | if (Math.abs(phi) < Proj4js.common.HALF_PI) { |
---|
2684 | // sincos(phi, &sin_phi, &cos_phi); |
---|
2685 | var sin_phi=Math.sin(phi); |
---|
2686 | var cos_phi=Math.cos(phi); |
---|
2687 | var tan_phi = Math.tan(phi); |
---|
2688 | var c = this.ep2 * Math.pow(cos_phi,2); |
---|
2689 | var cs = Math.pow(c,2); |
---|
2690 | var t = Math.pow(tan_phi,2); |
---|
2691 | var ts = Math.pow(t,2); |
---|
2692 | con = 1.0 - this.es * Math.pow(sin_phi,2); |
---|
2693 | var n = this.a / Math.sqrt(con); |
---|
2694 | var r = n * (1.0 - this.es) / con; |
---|
2695 | var d = x / (n * this.k0); |
---|
2696 | var ds = Math.pow(d,2); |
---|
2697 | lat = phi - (n * tan_phi * ds / r) * (0.5 - ds / 24.0 * (5.0 + 3.0 * t + 10.0 * c - 4.0 * cs - 9.0 * this.ep2 - ds / 30.0 * (61.0 + 90.0 * t + 298.0 * c + 45.0 * ts - 252.0 * this.ep2 - 3.0 * cs))); |
---|
2698 | lon = Proj4js.common.adjust_lon(this.long0 + (d * (1.0 - ds / 6.0 * (1.0 + 2.0 * t + c - ds / 20.0 * (5.0 - 2.0 * c + 28.0 * t - 3.0 * cs + 8.0 * this.ep2 + 24.0 * ts))) / cos_phi)); |
---|
2699 | } else { |
---|
2700 | lat = Proj4js.common.HALF_PI * Proj4js.common.sign(y); |
---|
2701 | lon = this.long0; |
---|
2702 | } |
---|
2703 | } |
---|
2704 | p.x = lon; |
---|
2705 | p.y = lat; |
---|
2706 | return p; |
---|
2707 | } // tmercInv() |
---|
2708 | }; |
---|
2709 | /* ====================================================================== |
---|
2710 | defs/GOOGLE.js |
---|
2711 | ====================================================================== */ |
---|
2712 | |
---|
2713 | Proj4js.defs["GOOGLE"]="+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +no_defs"; |
---|
2714 | Proj4js.defs["EPSG:900913"]=Proj4js.defs["GOOGLE"]; |
---|
2715 | /* ====================================================================== |
---|
2716 | projCode/gstmerc.js |
---|
2717 | ====================================================================== */ |
---|
2718 | |
---|
2719 | Proj4js.Proj.gstmerc = { |
---|
2720 | init : function() { |
---|
2721 | |
---|
2722 | // array of: a, b, lon0, lat0, k0, x0, y0 |
---|
2723 | var temp= this.b / this.a; |
---|
2724 | this.e= Math.sqrt(1.0 - temp*temp); |
---|
2725 | this.lc= this.long0; |
---|
2726 | this.rs= Math.sqrt(1.0+this.e*this.e*Math.pow(Math.cos(this.lat0),4.0)/(1.0-this.e*this.e)); |
---|
2727 | var sinz= Math.sin(this.lat0); |
---|
2728 | var pc= Math.asin(sinz/this.rs); |
---|
2729 | var sinzpc= Math.sin(pc); |
---|
2730 | this.cp= Proj4js.common.latiso(0.0,pc,sinzpc)-this.rs*Proj4js.common.latiso(this.e,this.lat0,sinz); |
---|
2731 | this.n2= this.k0*this.a*Math.sqrt(1.0-this.e*this.e)/(1.0-this.e*this.e*sinz*sinz); |
---|
2732 | this.xs= this.x0; |
---|
2733 | this.ys= this.y0-this.n2*pc; |
---|
2734 | |
---|
2735 | if (!this.title) this.title = "Gauss Schreiber transverse mercator"; |
---|
2736 | }, |
---|
2737 | |
---|
2738 | |
---|
2739 | // forward equations--mapping lat,long to x,y |
---|
2740 | // ----------------------------------------------------------------- |
---|
2741 | forward : function(p) { |
---|
2742 | |
---|
2743 | var lon= p.x; |
---|
2744 | var lat= p.y; |
---|
2745 | |
---|
2746 | var L= this.rs*(lon-this.lc); |
---|
2747 | var Ls= this.cp+(this.rs*Proj4js.common.latiso(this.e,lat,Math.sin(lat))); |
---|
2748 | var lat1= Math.asin(Math.sin(L)/Proj4js.common.cosh(Ls)); |
---|
2749 | var Ls1= Proj4js.common.latiso(0.0,lat1,Math.sin(lat1)); |
---|
2750 | p.x= this.xs+(this.n2*Ls1); |
---|
2751 | p.y= this.ys+(this.n2*Math.atan(Proj4js.common.sinh(Ls)/Math.cos(L))); |
---|
2752 | return p; |
---|
2753 | }, |
---|
2754 | |
---|
2755 | // inverse equations--mapping x,y to lat/long |
---|
2756 | // ----------------------------------------------------------------- |
---|
2757 | inverse : function(p) { |
---|
2758 | |
---|
2759 | var x= p.x; |
---|
2760 | var y= p.y; |
---|
2761 | |
---|
2762 | var L= Math.atan(Proj4js.common.sinh((x-this.xs)/this.n2)/Math.cos((y-this.ys)/this.n2)); |
---|
2763 | var lat1= Math.asin(Math.sin((y-this.ys)/this.n2)/Proj4js.common.cosh((x-this.xs)/this.n2)); |
---|
2764 | var LC= Proj4js.common.latiso(0.0,lat1,Math.sin(lat1)); |
---|
2765 | p.x= this.lc+L/this.rs; |
---|
2766 | p.y= Proj4js.common.invlatiso(this.e,(LC-this.cp)/this.rs); |
---|
2767 | return p; |
---|
2768 | } |
---|
2769 | |
---|
2770 | }; |
---|
2771 | /* ====================================================================== |
---|
2772 | projCode/ortho.js |
---|
2773 | ====================================================================== */ |
---|
2774 | |
---|
2775 | /******************************************************************************* |
---|
2776 | NAME ORTHOGRAPHIC |
---|
2777 | |
---|
2778 | PURPOSE: Transforms input longitude and latitude to Easting and |
---|
2779 | Northing for the Orthographic projection. The |
---|
2780 | longitude and latitude must be in radians. The Easting |
---|
2781 | and Northing values will be returned in meters. |
---|
2782 | |
---|
2783 | PROGRAMMER DATE |
---|
2784 | ---------- ---- |
---|
2785 | T. Mittan Mar, 1993 |
---|
2786 | |
---|
2787 | ALGORITHM REFERENCES |
---|
2788 | |
---|
2789 | 1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological |
---|
2790 | Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United |
---|
2791 | State Government Printing Office, Washington D.C., 1987. |
---|
2792 | |
---|
2793 | 2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections", |
---|
2794 | U.S. Geological Survey Professional Paper 1453 , United State Government |
---|
2795 | Printing Office, Washington D.C., 1989. |
---|
2796 | *******************************************************************************/ |
---|
2797 | |
---|
2798 | Proj4js.Proj.ortho = { |
---|
2799 | |
---|
2800 | /* Initialize the Orthographic projection |
---|
2801 | -------------------------------------*/ |
---|
2802 | init: function(def) { |
---|
2803 | //double temp; /* temporary variable */ |
---|
2804 | |
---|
2805 | /* Place parameters in static storage for common use |
---|
2806 | -------------------------------------------------*/; |
---|
2807 | this.sin_p14=Math.sin(this.lat0); |
---|
2808 | this.cos_p14=Math.cos(this.lat0); |
---|
2809 | }, |
---|
2810 | |
---|
2811 | |
---|
2812 | /* Orthographic forward equations--mapping lat,long to x,y |
---|
2813 | ---------------------------------------------------*/ |
---|
2814 | forward: function(p) { |
---|
2815 | var sinphi, cosphi; /* sin and cos value */ |
---|
2816 | var dlon; /* delta longitude value */ |
---|
2817 | var coslon; /* cos of longitude */ |
---|
2818 | var ksp; /* scale factor */ |
---|
2819 | var g; |
---|
2820 | var lon=p.x; |
---|
2821 | var lat=p.y; |
---|
2822 | /* Forward equations |
---|
2823 | -----------------*/ |
---|
2824 | dlon = Proj4js.common.adjust_lon(lon - this.long0); |
---|
2825 | |
---|
2826 | sinphi=Math.sin(lat); |
---|
2827 | cosphi=Math.cos(lat); |
---|
2828 | |
---|
2829 | coslon = Math.cos(dlon); |
---|
2830 | g = this.sin_p14 * sinphi + this.cos_p14 * cosphi * coslon; |
---|
2831 | ksp = 1.0; |
---|
2832 | if ((g > 0) || (Math.abs(g) <= Proj4js.common.EPSLN)) { |
---|
2833 | var x = this.a * ksp * cosphi * Math.sin(dlon); |
---|
2834 | var y = this.y0 + this.a * ksp * (this.cos_p14 * sinphi - this.sin_p14 * cosphi * coslon); |
---|
2835 | } else { |
---|
2836 | Proj4js.reportError("orthoFwdPointError"); |
---|
2837 | } |
---|
2838 | p.x=x; |
---|
2839 | p.y=y; |
---|
2840 | return p; |
---|
2841 | }, |
---|
2842 | |
---|
2843 | |
---|
2844 | inverse: function(p) { |
---|
2845 | var rh; /* height above ellipsoid */ |
---|
2846 | var z; /* angle */ |
---|
2847 | var sinz,cosz; /* sin of z and cos of z */ |
---|
2848 | var temp; |
---|
2849 | var con; |
---|
2850 | var lon , lat; |
---|
2851 | /* Inverse equations |
---|
2852 | -----------------*/ |
---|
2853 | p.x -= this.x0; |
---|
2854 | p.y -= this.y0; |
---|
2855 | rh = Math.sqrt(p.x * p.x + p.y * p.y); |
---|
2856 | if (rh > this.a + .0000001) { |
---|
2857 | Proj4js.reportError("orthoInvDataError"); |
---|
2858 | } |
---|
2859 | z = Proj4js.common.asinz(rh / this.a); |
---|
2860 | |
---|
2861 | sinz=Math.sin(z); |
---|
2862 | cosz=Math.cos(z); |
---|
2863 | |
---|
2864 | lon = this.long0; |
---|
2865 | if (Math.abs(rh) <= Proj4js.common.EPSLN) { |
---|
2866 | lat = this.lat0; |
---|
2867 | } |
---|
2868 | lat = Proj4js.common.asinz(cosz * this.sin_p14 + (p.y * sinz * this.cos_p14)/rh); |
---|
2869 | con = Math.abs(this.lat0) - Proj4js.common.HALF_PI; |
---|
2870 | if (Math.abs(con) <= Proj4js.common.EPSLN) { |
---|
2871 | if (this.lat0 >= 0) { |
---|
2872 | lon = Proj4js.common.adjust_lon(this.long0 + Math.atan2(p.x, -p.y)); |
---|
2873 | } else { |
---|
2874 | lon = Proj4js.common.adjust_lon(this.long0 -Math.atan2(-p.x, p.y)); |
---|
2875 | } |
---|
2876 | } |
---|
2877 | con = cosz - this.sin_p14 * Math.sin(lat); |
---|
2878 | p.x=lon; |
---|
2879 | p.y=lat; |
---|
2880 | return p; |
---|
2881 | } |
---|
2882 | }; |
---|
2883 | |
---|
2884 | |
---|
2885 | /* ====================================================================== |
---|
2886 | projCode/somerc.js |
---|
2887 | ====================================================================== */ |
---|
2888 | |
---|
2889 | /******************************************************************************* |
---|
2890 | NAME SWISS OBLIQUE MERCATOR |
---|
2891 | |
---|
2892 | PURPOSE: Swiss projection. |
---|
2893 | WARNING: X and Y are inverted (weird) in the swiss coordinate system. Not |
---|
2894 | here, since we want X to be horizontal and Y vertical. |
---|
2895 | |
---|
2896 | ALGORITHM REFERENCES |
---|
2897 | 1. "Formules et constantes pour le Calcul pour la |
---|
2898 | projection cylindrique conforme à axe oblique et pour la transformation entre |
---|
2899 | des systÚmes de référence". |
---|
2900 | http://www.swisstopo.admin.ch/internet/swisstopo/fr/home/topics/survey/sys/refsys/switzerland.parsysrelated1.31216.downloadList.77004.DownloadFile.tmp/swissprojectionfr.pdf |
---|
2901 | |
---|
2902 | *******************************************************************************/ |
---|
2903 | |
---|
2904 | Proj4js.Proj.somerc = { |
---|
2905 | |
---|
2906 | init: function() { |
---|
2907 | var phy0 = this.lat0; |
---|
2908 | this.lambda0 = this.long0; |
---|
2909 | var sinPhy0 = Math.sin(phy0); |
---|
2910 | var semiMajorAxis = this.a; |
---|
2911 | var invF = this.rf; |
---|
2912 | var flattening = 1 / invF; |
---|
2913 | var e2 = 2 * flattening - Math.pow(flattening, 2); |
---|
2914 | var e = this.e = Math.sqrt(e2); |
---|
2915 | this.R = this.k0 * semiMajorAxis * Math.sqrt(1 - e2) / (1 - e2 * Math.pow(sinPhy0, 2.0)); |
---|
2916 | this.alpha = Math.sqrt(1 + e2 / (1 - e2) * Math.pow(Math.cos(phy0), 4.0)); |
---|
2917 | this.b0 = Math.asin(sinPhy0 / this.alpha); |
---|
2918 | this.K = Math.log(Math.tan(Math.PI / 4.0 + this.b0 / 2.0)) |
---|
2919 | - this.alpha |
---|
2920 | * Math.log(Math.tan(Math.PI / 4.0 + phy0 / 2.0)) |
---|
2921 | + this.alpha |
---|
2922 | * e / 2 |
---|
2923 | * Math.log((1 + e * sinPhy0) |
---|
2924 | / (1 - e * sinPhy0)); |
---|
2925 | }, |
---|
2926 | |
---|
2927 | |
---|
2928 | forward: function(p) { |
---|
2929 | var Sa1 = Math.log(Math.tan(Math.PI / 4.0 - p.y / 2.0)); |
---|
2930 | var Sa2 = this.e / 2.0 |
---|
2931 | * Math.log((1 + this.e * Math.sin(p.y)) |
---|
2932 | / (1 - this.e * Math.sin(p.y))); |
---|
2933 | var S = -this.alpha * (Sa1 + Sa2) + this.K; |
---|
2934 | |
---|
2935 | // spheric latitude |
---|
2936 | var b = 2.0 * (Math.atan(Math.exp(S)) - Math.PI / 4.0); |
---|
2937 | |
---|
2938 | // spheric longitude |
---|
2939 | var I = this.alpha * (p.x - this.lambda0); |
---|
2940 | |
---|
2941 | // psoeudo equatorial rotation |
---|
2942 | var rotI = Math.atan(Math.sin(I) |
---|
2943 | / (Math.sin(this.b0) * Math.tan(b) + |
---|
2944 | Math.cos(this.b0) * Math.cos(I))); |
---|
2945 | |
---|
2946 | var rotB = Math.asin(Math.cos(this.b0) * Math.sin(b) - |
---|
2947 | Math.sin(this.b0) * Math.cos(b) * Math.cos(I)); |
---|
2948 | |
---|
2949 | p.y = this.R / 2.0 |
---|
2950 | * Math.log((1 + Math.sin(rotB)) / (1 - Math.sin(rotB))) |
---|
2951 | + this.y0; |
---|
2952 | p.x = this.R * rotI + this.x0; |
---|
2953 | return p; |
---|
2954 | }, |
---|
2955 | |
---|
2956 | inverse: function(p) { |
---|
2957 | var Y = p.x - this.x0; |
---|
2958 | var X = p.y - this.y0; |
---|
2959 | |
---|
2960 | var rotI = Y / this.R; |
---|
2961 | var rotB = 2 * (Math.atan(Math.exp(X / this.R)) - Math.PI / 4.0); |
---|
2962 | |
---|
2963 | var b = Math.asin(Math.cos(this.b0) * Math.sin(rotB) |
---|
2964 | + Math.sin(this.b0) * Math.cos(rotB) * Math.cos(rotI)); |
---|
2965 | var I = Math.atan(Math.sin(rotI) |
---|
2966 | / (Math.cos(this.b0) * Math.cos(rotI) - Math.sin(this.b0) |
---|
2967 | * Math.tan(rotB))); |
---|
2968 | |
---|
2969 | var lambda = this.lambda0 + I / this.alpha; |
---|
2970 | |
---|
2971 | var S = 0.0; |
---|
2972 | var phy = b; |
---|
2973 | var prevPhy = -1000.0; |
---|
2974 | var iteration = 0; |
---|
2975 | while (Math.abs(phy - prevPhy) > 0.0000001) |
---|
2976 | { |
---|
2977 | if (++iteration > 20) |
---|
2978 | { |
---|
2979 | Proj4js.reportError("omercFwdInfinity"); |
---|
2980 | return; |
---|
2981 | } |
---|
2982 | //S = Math.log(Math.tan(Math.PI / 4.0 + phy / 2.0)); |
---|
2983 | S = 1.0 |
---|
2984 | / this.alpha |
---|
2985 | * (Math.log(Math.tan(Math.PI / 4.0 + b / 2.0)) - this.K) |
---|
2986 | + this.e |
---|
2987 | * Math.log(Math.tan(Math.PI / 4.0 |
---|
2988 | + Math.asin(this.e * Math.sin(phy)) |
---|
2989 | / 2.0)); |
---|
2990 | prevPhy = phy; |
---|
2991 | phy = 2.0 * Math.atan(Math.exp(S)) - Math.PI / 2.0; |
---|
2992 | } |
---|
2993 | |
---|
2994 | p.x = lambda; |
---|
2995 | p.y = phy; |
---|
2996 | return p; |
---|
2997 | } |
---|
2998 | }; |
---|
2999 | /* ====================================================================== |
---|
3000 | projCode/stere.js |
---|
3001 | ====================================================================== */ |
---|
3002 | |
---|
3003 | |
---|
3004 | // Initialize the Stereographic projection |
---|
3005 | |
---|
3006 | Proj4js.Proj.stere = { |
---|
3007 | ssfn_: function(phit, sinphi, eccen) { |
---|
3008 | sinphi *= eccen; |
---|
3009 | return (Math.tan (.5 * (Proj4js.common.HALF_PI + phit)) * Math.pow((1. - sinphi) / (1. + sinphi), .5 * eccen)); |
---|
3010 | }, |
---|
3011 | TOL: 1.e-8, |
---|
3012 | NITER: 8, |
---|
3013 | CONV: 1.e-10, |
---|
3014 | S_POLE: 0, |
---|
3015 | N_POLE: 1, |
---|
3016 | OBLIQ: 2, |
---|
3017 | EQUIT: 3, |
---|
3018 | |
---|
3019 | init : function() { |
---|
3020 | this.phits = this.lat_ts ? this.lat_ts : Proj4js.common.HALF_PI; |
---|
3021 | var t = Math.abs(this.lat0); |
---|
3022 | if ((Math.abs(t) - Proj4js.common.HALF_PI) < Proj4js.common.EPSLN) { |
---|
3023 | this.mode = this.lat0 < 0. ? this.S_POLE : this.N_POLE; |
---|
3024 | } else { |
---|
3025 | this.mode = t > Proj4js.common.EPSLN ? this.OBLIQ : this.EQUIT; |
---|
3026 | } |
---|
3027 | this.phits = Math.abs(this.phits); |
---|
3028 | if (this.es) { |
---|
3029 | var X; |
---|
3030 | |
---|
3031 | switch (this.mode) { |
---|
3032 | case this.N_POLE: |
---|
3033 | case this.S_POLE: |
---|
3034 | if (Math.abs(this.phits - Proj4js.common.HALF_PI) < Proj4js.common.EPSLN) { |
---|
3035 | this.akm1 = 2. * this.k0 / Math.sqrt(Math.pow(1+this.e,1+this.e)*Math.pow(1-this.e,1-this.e)); |
---|
3036 | } else { |
---|
3037 | t = Math.sin(this.phits); |
---|
3038 | this.akm1 = Math.cos(this.phits) / Proj4js.common.tsfnz(this.e, this.phits, t); |
---|
3039 | t *= this.e; |
---|
3040 | this.akm1 /= Math.sqrt(1. - t * t); |
---|
3041 | } |
---|
3042 | break; |
---|
3043 | case this.EQUIT: |
---|
3044 | this.akm1 = 2. * this.k0; |
---|
3045 | break; |
---|
3046 | case this.OBLIQ: |
---|
3047 | t = Math.sin(this.lat0); |
---|
3048 | X = 2. * Math.atan(this.ssfn_(this.lat0, t, this.e)) - Proj4js.common.HALF_PI; |
---|
3049 | t *= this.e; |
---|
3050 | this.akm1 = 2. * this.k0 * Math.cos(this.lat0) / Math.sqrt(1. - t * t); |
---|
3051 | this.sinX1 = Math.sin(X); |
---|
3052 | this.cosX1 = Math.cos(X); |
---|
3053 | break; |
---|
3054 | } |
---|
3055 | } else { |
---|
3056 | switch (this.mode) { |
---|
3057 | case this.OBLIQ: |
---|
3058 | this.sinph0 = Math.sin(this.lat0); |
---|
3059 | this.cosph0 = Math.cos(this.lat0); |
---|
3060 | case this.EQUIT: |
---|
3061 | this.akm1 = 2. * this.k0; |
---|
3062 | break; |
---|
3063 | case this.S_POLE: |
---|
3064 | case this.N_POLE: |
---|
3065 | this.akm1 = Math.abs(this.phits - Proj4js.common.HALF_PI) >= Proj4js.common.EPSLN ? |
---|
3066 | Math.cos(this.phits) / Math.tan(Proj4js.common.FORTPI - .5 * this.phits) : |
---|
3067 | 2. * this.k0 ; |
---|
3068 | break; |
---|
3069 | } |
---|
3070 | } |
---|
3071 | }, |
---|
3072 | |
---|
3073 | // Stereographic forward equations--mapping lat,long to x,y |
---|
3074 | forward: function(p) { |
---|
3075 | var lon = p.x; |
---|
3076 | lon = Proj4js.common.adjust_lon(lon - this.long0); |
---|
3077 | var lat = p.y; |
---|
3078 | var x, y; |
---|
3079 | |
---|
3080 | if (this.sphere) { |
---|
3081 | var sinphi, cosphi, coslam, sinlam; |
---|
3082 | |
---|
3083 | sinphi = Math.sin(lat); |
---|
3084 | cosphi = Math.cos(lat); |
---|
3085 | coslam = Math.cos(lon); |
---|
3086 | sinlam = Math.sin(lon); |
---|
3087 | switch (this.mode) { |
---|
3088 | case this.EQUIT: |
---|
3089 | y = 1. + cosphi * coslam; |
---|
3090 | if (y <= Proj4js.common.EPSLN) { |
---|
3091 | F_ERROR; |
---|
3092 | } |
---|
3093 | y = this.akm1 / y; |
---|
3094 | x = y * cosphi * sinlam; |
---|
3095 | y *= sinphi; |
---|
3096 | break; |
---|
3097 | case this.OBLIQ: |
---|
3098 | y = 1. + this.sinph0 * sinphi + this.cosph0 * cosphi * coslam; |
---|
3099 | if (y <= Proj4js.common.EPSLN) { |
---|
3100 | F_ERROR; |
---|
3101 | } |
---|
3102 | y = this.akm1 / y; |
---|
3103 | x = y * cosphi * sinlam; |
---|
3104 | y *= this.cosph0 * sinphi - this.sinph0 * cosphi * coslam; |
---|
3105 | break; |
---|
3106 | case this.N_POLE: |
---|
3107 | coslam = -coslam; |
---|
3108 | lat = -lat; |
---|
3109 | //Note no break here so it conitnues through S_POLE |
---|
3110 | case this.S_POLE: |
---|
3111 | if (Math.abs(lat - Proj4js.common.HALF_PI) < this.TOL) { |
---|
3112 | F_ERROR; |
---|
3113 | } |
---|
3114 | y = this.akm1 * Math.tan(Proj4js.common.FORTPI + .5 * lat); |
---|
3115 | x = sinlam * y; |
---|
3116 | y *= coslam; |
---|
3117 | break; |
---|
3118 | } |
---|
3119 | } else { |
---|
3120 | coslam = Math.cos(lon); |
---|
3121 | sinlam = Math.sin(lon); |
---|
3122 | sinphi = Math.sin(lat); |
---|
3123 | if (this.mode == this.OBLIQ || this.mode == this.EQUIT) { |
---|
3124 | X = 2. * Math.atan(this.ssfn_(lat, sinphi, this.e)); |
---|
3125 | sinX = Math.sin(X - Proj4js.common.HALF_PI); |
---|
3126 | cosX = Math.cos(X); |
---|
3127 | } |
---|
3128 | switch (this.mode) { |
---|
3129 | case this.OBLIQ: |
---|
3130 | A = this.akm1 / (this.cosX1 * (1. + this.sinX1 * sinX + this.cosX1 * cosX * coslam)); |
---|
3131 | y = A * (this.cosX1 * sinX - this.sinX1 * cosX * coslam); |
---|
3132 | x = A * cosX; |
---|
3133 | break; |
---|
3134 | case this.EQUIT: |
---|
3135 | A = 2. * this.akm1 / (1. + cosX * coslam); |
---|
3136 | y = A * sinX; |
---|
3137 | x = A * cosX; |
---|
3138 | break; |
---|
3139 | case this.S_POLE: |
---|
3140 | lat = -lat; |
---|
3141 | coslam = - coslam; |
---|
3142 | sinphi = -sinphi; |
---|
3143 | case this.N_POLE: |
---|
3144 | x = this.akm1 * Proj4js.common.tsfnz(this.e, lat, sinphi); |
---|
3145 | y = - x * coslam; |
---|
3146 | break; |
---|
3147 | } |
---|
3148 | x = x * sinlam; |
---|
3149 | } |
---|
3150 | p.x = x*this.a + this.x0; |
---|
3151 | p.y = y*this.a + this.y0; |
---|
3152 | return p; |
---|
3153 | }, |
---|
3154 | |
---|
3155 | |
---|
3156 | //* Stereographic inverse equations--mapping x,y to lat/long |
---|
3157 | inverse: function(p) { |
---|
3158 | var x = (p.x - this.x0)/this.a; /* descale and de-offset */ |
---|
3159 | var y = (p.y - this.y0)/this.a; |
---|
3160 | var lon, lat; |
---|
3161 | |
---|
3162 | var cosphi, sinphi, tp=0.0, phi_l=0.0, rho, halfe=0.0, pi2=0.0; |
---|
3163 | var i; |
---|
3164 | |
---|
3165 | if (this.sphere) { |
---|
3166 | var c, rh, sinc, cosc; |
---|
3167 | |
---|
3168 | rh = Math.sqrt(x*x + y*y); |
---|
3169 | c = 2. * Math.atan(rh / this.akm1); |
---|
3170 | sinc = Math.sin(c); |
---|
3171 | cosc = Math.cos(c); |
---|
3172 | lon = 0.; |
---|
3173 | switch (this.mode) { |
---|
3174 | case this.EQUIT: |
---|
3175 | if (Math.abs(rh) <= Proj4js.common.EPSLN) { |
---|
3176 | lat = 0.; |
---|
3177 | } else { |
---|
3178 | lat = Math.asin(y * sinc / rh); |
---|
3179 | } |
---|
3180 | if (cosc != 0. || x != 0.) lon = Math.atan2(x * sinc, cosc * rh); |
---|
3181 | break; |
---|
3182 | case this.OBLIQ: |
---|
3183 | if (Math.abs(rh) <= Proj4js.common.EPSLN) { |
---|
3184 | lat = this.phi0; |
---|
3185 | } else { |
---|
3186 | lat = Math.asin(cosc * sinph0 + y * sinc * cosph0 / rh); |
---|
3187 | } |
---|
3188 | c = cosc - sinph0 * Math.sin(lat); |
---|
3189 | if (c != 0. || x != 0.) { |
---|
3190 | lon = Math.atan2(x * sinc * cosph0, c * rh); |
---|
3191 | } |
---|
3192 | break; |
---|
3193 | case this.N_POLE: |
---|
3194 | y = -y; |
---|
3195 | case this.S_POLE: |
---|
3196 | if (Math.abs(rh) <= Proj4js.common.EPSLN) { |
---|
3197 | lat = this.phi0; |
---|
3198 | } else { |
---|
3199 | lat = Math.asin(this.mode == this.S_POLE ? -cosc : cosc); |
---|
3200 | } |
---|
3201 | lon = (x == 0. && y == 0.) ? 0. : Math.atan2(x, y); |
---|
3202 | break; |
---|
3203 | } |
---|
3204 | p.x = Proj4js.common.adjust_lon(lon + this.long0); |
---|
3205 | p.y = lat; |
---|
3206 | } else { |
---|
3207 | rho = Math.sqrt(x*x + y*y); |
---|
3208 | switch (this.mode) { |
---|
3209 | case this.OBLIQ: |
---|
3210 | case this.EQUIT: |
---|
3211 | tp = 2. * Math.atan2(rho * this.cosX1 , this.akm1); |
---|
3212 | cosphi = Math.cos(tp); |
---|
3213 | sinphi = Math.sin(tp); |
---|
3214 | if( rho == 0.0 ) { |
---|
3215 | phi_l = Math.asin(cosphi * this.sinX1); |
---|
3216 | } else { |
---|
3217 | phi_l = Math.asin(cosphi * this.sinX1 + (y * sinphi * this.cosX1 / rho)); |
---|
3218 | } |
---|
3219 | |
---|
3220 | tp = Math.tan(.5 * (Proj4js.common.HALF_PI + phi_l)); |
---|
3221 | x *= sinphi; |
---|
3222 | y = rho * this.cosX1 * cosphi - y * this.sinX1* sinphi; |
---|
3223 | pi2 = Proj4js.common.HALF_PI; |
---|
3224 | halfe = .5 * this.e; |
---|
3225 | break; |
---|
3226 | case this.N_POLE: |
---|
3227 | y = -y; |
---|
3228 | case this.S_POLE: |
---|
3229 | tp = - rho / this.akm1; |
---|
3230 | phi_l = Proj4js.common.HALF_PI - 2. * Math.atan(tp); |
---|
3231 | pi2 = -Proj4js.common.HALF_PI; |
---|
3232 | halfe = -.5 * this.e; |
---|
3233 | break; |
---|
3234 | } |
---|
3235 | for (i = this.NITER; i--; phi_l = lat) { //check this |
---|
3236 | sinphi = this.e * Math.sin(phi_l); |
---|
3237 | lat = 2. * Math.atan(tp * Math.pow((1.+sinphi)/(1.-sinphi), halfe)) - pi2; |
---|
3238 | if (Math.abs(phi_l - lat) < this.CONV) { |
---|
3239 | if (this.mode == this.S_POLE) lat = -lat; |
---|
3240 | lon = (x == 0. && y == 0.) ? 0. : Math.atan2(x, y); |
---|
3241 | p.x = Proj4js.common.adjust_lon(lon + this.long0); |
---|
3242 | p.y = lat; |
---|
3243 | return p; |
---|
3244 | } |
---|
3245 | } |
---|
3246 | } |
---|
3247 | } |
---|
3248 | }; |
---|
3249 | /* ====================================================================== |
---|
3250 | projCode/nzmg.js |
---|
3251 | ====================================================================== */ |
---|
3252 | |
---|
3253 | /******************************************************************************* |
---|
3254 | NAME NEW ZEALAND MAP GRID |
---|
3255 | |
---|
3256 | PURPOSE: Transforms input longitude and latitude to Easting and |
---|
3257 | Northing for the New Zealand Map Grid projection. The |
---|
3258 | longitude and latitude must be in radians. The Easting |
---|
3259 | and Northing values will be returned in meters. |
---|
3260 | |
---|
3261 | |
---|
3262 | ALGORITHM REFERENCES |
---|
3263 | |
---|
3264 | 1. Department of Land and Survey Technical Circular 1973/32 |
---|
3265 | http://www.linz.govt.nz/docs/miscellaneous/nz-map-definition.pdf |
---|
3266 | |
---|
3267 | 2. OSG Technical Report 4.1 |
---|
3268 | http://www.linz.govt.nz/docs/miscellaneous/nzmg.pdf |
---|
3269 | |
---|
3270 | |
---|
3271 | IMPLEMENTATION NOTES |
---|
3272 | |
---|
3273 | The two references use different symbols for the calculated values. This |
---|
3274 | implementation uses the variable names similar to the symbols in reference [1]. |
---|
3275 | |
---|
3276 | The alogrithm uses different units for delta latitude and delta longitude. |
---|
3277 | The delta latitude is assumed to be in units of seconds of arc x 10^-5. |
---|
3278 | The delta longitude is the usual radians. Look out for these conversions. |
---|
3279 | |
---|
3280 | The algorithm is described using complex arithmetic. There were three |
---|
3281 | options: |
---|
3282 | * find and use a Javascript library for complex arithmetic |
---|
3283 | * write my own complex library |
---|
3284 | * expand the complex arithmetic by hand to simple arithmetic |
---|
3285 | |
---|
3286 | This implementation has expanded the complex multiplication operations |
---|
3287 | into parallel simple arithmetic operations for the real and imaginary parts. |
---|
3288 | The imaginary part is way over to the right of the display; this probably |
---|
3289 | violates every coding standard in the world, but, to me, it makes it much |
---|
3290 | more obvious what is going on. |
---|
3291 | |
---|
3292 | The following complex operations are used: |
---|
3293 | - addition |
---|
3294 | - multiplication |
---|
3295 | - division |
---|
3296 | - complex number raised to integer power |
---|
3297 | - summation |
---|
3298 | |
---|
3299 | A summary of complex arithmetic operations: |
---|
3300 | (from http://en.wikipedia.org/wiki/Complex_arithmetic) |
---|
3301 | addition: (a + bi) + (c + di) = (a + c) + (b + d)i |
---|
3302 | subtraction: (a + bi) - (c + di) = (a - c) + (b - d)i |
---|
3303 | multiplication: (a + bi) x (c + di) = (ac - bd) + (bc + ad)i |
---|
3304 | division: (a + bi) / (c + di) = [(ac + bd)/(cc + dd)] + [(bc - ad)/(cc + dd)]i |
---|
3305 | |
---|
3306 | The algorithm needs to calculate summations of simple and complex numbers. This is |
---|
3307 | implemented using a for-loop, pre-loading the summed value to zero. |
---|
3308 | |
---|
3309 | The algorithm needs to calculate theta^2, theta^3, etc while doing a summation. |
---|
3310 | There are three possible implementations: |
---|
3311 | - use Math.pow in the summation loop - except for complex numbers |
---|
3312 | - precalculate the values before running the loop |
---|
3313 | - calculate theta^n = theta^(n-1) * theta during the loop |
---|
3314 | This implementation uses the third option for both real and complex arithmetic. |
---|
3315 | |
---|
3316 | For example |
---|
3317 | psi_n = 1; |
---|
3318 | sum = 0; |
---|
3319 | for (n = 1; n <=6; n++) { |
---|
3320 | psi_n1 = psi_n * psi; // calculate psi^(n+1) |
---|
3321 | psi_n = psi_n1; |
---|
3322 | sum = sum + A[n] * psi_n; |
---|
3323 | } |
---|
3324 | |
---|
3325 | |
---|
3326 | TEST VECTORS |
---|
3327 | |
---|
3328 | NZMG E, N: 2487100.638 6751049.719 metres |
---|
3329 | NZGD49 long, lat: 172.739194 -34.444066 degrees |
---|
3330 | |
---|
3331 | NZMG E, N: 2486533.395 6077263.661 metres |
---|
3332 | NZGD49 long, lat: 172.723106 -40.512409 degrees |
---|
3333 | |
---|
3334 | NZMG E, N: 2216746.425 5388508.765 metres |
---|
3335 | NZGD49 long, lat: 169.172062 -46.651295 degrees |
---|
3336 | |
---|
3337 | Note that these test vectors convert from NZMG metres to lat/long referenced |
---|
3338 | to NZGD49, not the more usual WGS84. The difference is about 70m N/S and about |
---|
3339 | 10m E/W. |
---|
3340 | |
---|
3341 | These test vectors are provided in reference [1]. Many more test |
---|
3342 | vectors are available in |
---|
3343 | http://www.linz.govt.nz/docs/topography/topographicdata/placenamesdatabase/nznamesmar08.zip |
---|
3344 | which is a catalog of names on the 260-series maps. |
---|
3345 | |
---|
3346 | |
---|
3347 | EPSG CODES |
---|
3348 | |
---|
3349 | NZMG EPSG:27200 |
---|
3350 | NZGD49 EPSG:4272 |
---|
3351 | |
---|
3352 | http://spatialreference.org/ defines these as |
---|
3353 | Proj4js.defs["EPSG:4272"] = "+proj=longlat +ellps=intl +datum=nzgd49 +no_defs "; |
---|
3354 | Proj4js.defs["EPSG:27200"] = "+proj=nzmg +lat_0=-41 +lon_0=173 +x_0=2510000 +y_0=6023150 +ellps=intl +datum=nzgd49 +units=m +no_defs "; |
---|
3355 | |
---|
3356 | |
---|
3357 | LICENSE |
---|
3358 | Copyright: Stephen Irons 2008 |
---|
3359 | Released under terms of the LGPL as per: http://www.gnu.org/copyleft/lesser.html |
---|
3360 | |
---|
3361 | *******************************************************************************/ |
---|
3362 | |
---|
3363 | |
---|
3364 | /** |
---|
3365 | Initialize New Zealand Map Grip projection |
---|
3366 | */ |
---|
3367 | |
---|
3368 | Proj4js.Proj.nzmg = { |
---|
3369 | |
---|
3370 | /** |
---|
3371 | * iterations: Number of iterations to refine inverse transform. |
---|
3372 | * 0 -> km accuracy |
---|
3373 | * 1 -> m accuracy -- suitable for most mapping applications |
---|
3374 | * 2 -> mm accuracy |
---|
3375 | */ |
---|
3376 | iterations: 1, |
---|
3377 | |
---|
3378 | init : function() { |
---|
3379 | this.A = new Array(); |
---|
3380 | this.A[1] = +0.6399175073; |
---|
3381 | this.A[2] = -0.1358797613; |
---|
3382 | this.A[3] = +0.063294409; |
---|
3383 | this.A[4] = -0.02526853; |
---|
3384 | this.A[5] = +0.0117879; |
---|
3385 | this.A[6] = -0.0055161; |
---|
3386 | this.A[7] = +0.0026906; |
---|
3387 | this.A[8] = -0.001333; |
---|
3388 | this.A[9] = +0.00067; |
---|
3389 | this.A[10] = -0.00034; |
---|
3390 | |
---|
3391 | this.B_re = new Array(); this.B_im = new Array(); |
---|
3392 | this.B_re[1] = +0.7557853228; this.B_im[1] = 0.0; |
---|
3393 | this.B_re[2] = +0.249204646; this.B_im[2] = +0.003371507; |
---|
3394 | this.B_re[3] = -0.001541739; this.B_im[3] = +0.041058560; |
---|
3395 | this.B_re[4] = -0.10162907; this.B_im[4] = +0.01727609; |
---|
3396 | this.B_re[5] = -0.26623489; this.B_im[5] = -0.36249218; |
---|
3397 | this.B_re[6] = -0.6870983; this.B_im[6] = -1.1651967; |
---|
3398 | |
---|
3399 | this.C_re = new Array(); this.C_im = new Array(); |
---|
3400 | this.C_re[1] = +1.3231270439; this.C_im[1] = 0.0; |
---|
3401 | this.C_re[2] = -0.577245789; this.C_im[2] = -0.007809598; |
---|
3402 | this.C_re[3] = +0.508307513; this.C_im[3] = -0.112208952; |
---|
3403 | this.C_re[4] = -0.15094762; this.C_im[4] = +0.18200602; |
---|
3404 | this.C_re[5] = +1.01418179; this.C_im[5] = +1.64497696; |
---|
3405 | this.C_re[6] = +1.9660549; this.C_im[6] = +2.5127645; |
---|
3406 | |
---|
3407 | this.D = new Array(); |
---|
3408 | this.D[1] = +1.5627014243; |
---|
3409 | this.D[2] = +0.5185406398; |
---|
3410 | this.D[3] = -0.03333098; |
---|
3411 | this.D[4] = -0.1052906; |
---|
3412 | this.D[5] = -0.0368594; |
---|
3413 | this.D[6] = +0.007317; |
---|
3414 | this.D[7] = +0.01220; |
---|
3415 | this.D[8] = +0.00394; |
---|
3416 | this.D[9] = -0.0013; |
---|
3417 | }, |
---|
3418 | |
---|
3419 | /** |
---|
3420 | New Zealand Map Grid Forward - long/lat to x/y |
---|
3421 | long/lat in radians |
---|
3422 | */ |
---|
3423 | forward : function(p) { |
---|
3424 | var lon = p.x; |
---|
3425 | var lat = p.y; |
---|
3426 | |
---|
3427 | var delta_lat = lat - this.lat0; |
---|
3428 | var delta_lon = lon - this.long0; |
---|
3429 | |
---|
3430 | // 1. Calculate d_phi and d_psi ... // and d_lambda |
---|
3431 | // For this algorithm, delta_latitude is in seconds of arc x 10-5, so we need to scale to those units. Longitude is radians. |
---|
3432 | var d_phi = delta_lat / Proj4js.common.SEC_TO_RAD * 1E-5; var d_lambda = delta_lon; |
---|
3433 | var d_phi_n = 1; // d_phi^0 |
---|
3434 | |
---|
3435 | var d_psi = 0; |
---|
3436 | for (n = 1; n <= 10; n++) { |
---|
3437 | d_phi_n = d_phi_n * d_phi; |
---|
3438 | d_psi = d_psi + this.A[n] * d_phi_n; |
---|
3439 | } |
---|
3440 | |
---|
3441 | // 2. Calculate theta |
---|
3442 | var th_re = d_psi; var th_im = d_lambda; |
---|
3443 | |
---|
3444 | // 3. Calculate z |
---|
3445 | var th_n_re = 1; var th_n_im = 0; // theta^0 |
---|
3446 | var th_n_re1; var th_n_im1; |
---|
3447 | |
---|
3448 | var z_re = 0; var z_im = 0; |
---|
3449 | for (n = 1; n <= 6; n++) { |
---|
3450 | th_n_re1 = th_n_re*th_re - th_n_im*th_im; th_n_im1 = th_n_im*th_re + th_n_re*th_im; |
---|
3451 | th_n_re = th_n_re1; th_n_im = th_n_im1; |
---|
3452 | z_re = z_re + this.B_re[n]*th_n_re - this.B_im[n]*th_n_im; z_im = z_im + this.B_im[n]*th_n_re + this.B_re[n]*th_n_im; |
---|
3453 | } |
---|
3454 | |
---|
3455 | // 4. Calculate easting and northing |
---|
3456 | x = (z_im * this.a) + this.x0; |
---|
3457 | y = (z_re * this.a) + this.y0; |
---|
3458 | |
---|
3459 | p.x = x; p.y = y; |
---|
3460 | |
---|
3461 | return p; |
---|
3462 | }, |
---|
3463 | |
---|
3464 | |
---|
3465 | /** |
---|
3466 | New Zealand Map Grid Inverse - x/y to long/lat |
---|
3467 | */ |
---|
3468 | inverse : function(p) { |
---|
3469 | |
---|
3470 | var x = p.x; |
---|
3471 | var y = p.y; |
---|
3472 | |
---|
3473 | var delta_x = x - this.x0; |
---|
3474 | var delta_y = y - this.y0; |
---|
3475 | |
---|
3476 | // 1. Calculate z |
---|
3477 | var z_re = delta_y / this.a; var z_im = delta_x / this.a; |
---|
3478 | |
---|
3479 | // 2a. Calculate theta - first approximation gives km accuracy |
---|
3480 | var z_n_re = 1; var z_n_im = 0; // z^0 |
---|
3481 | var z_n_re1; var z_n_im1; |
---|
3482 | |
---|
3483 | var th_re = 0; var th_im = 0; |
---|
3484 | for (n = 1; n <= 6; n++) { |
---|
3485 | z_n_re1 = z_n_re*z_re - z_n_im*z_im; z_n_im1 = z_n_im*z_re + z_n_re*z_im; |
---|
3486 | z_n_re = z_n_re1; z_n_im = z_n_im1; |
---|
3487 | th_re = th_re + this.C_re[n]*z_n_re - this.C_im[n]*z_n_im; th_im = th_im + this.C_im[n]*z_n_re + this.C_re[n]*z_n_im; |
---|
3488 | } |
---|
3489 | |
---|
3490 | // 2b. Iterate to refine the accuracy of the calculation |
---|
3491 | // 0 iterations gives km accuracy |
---|
3492 | // 1 iteration gives m accuracy -- good enough for most mapping applications |
---|
3493 | // 2 iterations bives mm accuracy |
---|
3494 | for (i = 0; i < this.iterations; i++) { |
---|
3495 | var th_n_re = th_re; var th_n_im = th_im; |
---|
3496 | var th_n_re1; var th_n_im1; |
---|
3497 | |
---|
3498 | var num_re = z_re; var num_im = z_im; |
---|
3499 | for (n = 2; n <= 6; n++) { |
---|
3500 | th_n_re1 = th_n_re*th_re - th_n_im*th_im; th_n_im1 = th_n_im*th_re + th_n_re*th_im; |
---|
3501 | th_n_re = th_n_re1; th_n_im = th_n_im1; |
---|
3502 | num_re = num_re + (n-1)*(this.B_re[n]*th_n_re - this.B_im[n]*th_n_im); num_im = num_im + (n-1)*(this.B_im[n]*th_n_re + this.B_re[n]*th_n_im); |
---|
3503 | } |
---|
3504 | |
---|
3505 | th_n_re = 1; th_n_im = 0; |
---|
3506 | var den_re = this.B_re[1]; var den_im = this.B_im[1]; |
---|
3507 | for (n = 2; n <= 6; n++) { |
---|
3508 | th_n_re1 = th_n_re*th_re - th_n_im*th_im; th_n_im1 = th_n_im*th_re + th_n_re*th_im; |
---|
3509 | th_n_re = th_n_re1; th_n_im = th_n_im1; |
---|
3510 | den_re = den_re + n * (this.B_re[n]*th_n_re - this.B_im[n]*th_n_im); den_im = den_im + n * (this.B_im[n]*th_n_re + this.B_re[n]*th_n_im); |
---|
3511 | } |
---|
3512 | |
---|
3513 | // Complex division |
---|
3514 | var den2 = den_re*den_re + den_im*den_im; |
---|
3515 | th_re = (num_re*den_re + num_im*den_im) / den2; th_im = (num_im*den_re - num_re*den_im) / den2; |
---|
3516 | } |
---|
3517 | |
---|
3518 | // 3. Calculate d_phi ... // and d_lambda |
---|
3519 | var d_psi = th_re; var d_lambda = th_im; |
---|
3520 | var d_psi_n = 1; // d_psi^0 |
---|
3521 | |
---|
3522 | var d_phi = 0; |
---|
3523 | for (n = 1; n <= 9; n++) { |
---|
3524 | d_psi_n = d_psi_n * d_psi; |
---|
3525 | d_phi = d_phi + this.D[n] * d_psi_n; |
---|
3526 | } |
---|
3527 | |
---|
3528 | // 4. Calculate latitude and longitude |
---|
3529 | // d_phi is calcuated in second of arc * 10^-5, so we need to scale back to radians. d_lambda is in radians. |
---|
3530 | var lat = this.lat0 + (d_phi * Proj4js.common.SEC_TO_RAD * 1E5); |
---|
3531 | var lon = this.long0 + d_lambda; |
---|
3532 | |
---|
3533 | p.x = lon; |
---|
3534 | p.y = lat; |
---|
3535 | |
---|
3536 | return p; |
---|
3537 | } |
---|
3538 | }; |
---|
3539 | /* ====================================================================== |
---|
3540 | projCode/mill.js |
---|
3541 | ====================================================================== */ |
---|
3542 | |
---|
3543 | /******************************************************************************* |
---|
3544 | NAME MILLER CYLINDRICAL |
---|
3545 | |
---|
3546 | PURPOSE: Transforms input longitude and latitude to Easting and |
---|
3547 | Northing for the Miller Cylindrical projection. The |
---|
3548 | longitude and latitude must be in radians. The Easting |
---|
3549 | and Northing values will be returned in meters. |
---|
3550 | |
---|
3551 | PROGRAMMER DATE |
---|
3552 | ---------- ---- |
---|
3553 | T. Mittan March, 1993 |
---|
3554 | |
---|
3555 | This function was adapted from the Lambert Azimuthal Equal Area projection |
---|
3556 | code (FORTRAN) in the General Cartographic Transformation Package software |
---|
3557 | which is available from the U.S. Geological Survey National Mapping Division. |
---|
3558 | |
---|
3559 | ALGORITHM REFERENCES |
---|
3560 | |
---|
3561 | 1. "New Equal-Area Map Projections for Noncircular Regions", John P. Snyder, |
---|
3562 | The American Cartographer, Vol 15, No. 4, October 1988, pp. 341-355. |
---|
3563 | |
---|
3564 | 2. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological |
---|
3565 | Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United |
---|
3566 | State Government Printing Office, Washington D.C., 1987. |
---|
3567 | |
---|
3568 | 3. "Software Documentation for GCTP General Cartographic Transformation |
---|
3569 | Package", U.S. Geological Survey National Mapping Division, May 1982. |
---|
3570 | *******************************************************************************/ |
---|
3571 | |
---|
3572 | Proj4js.Proj.mill = { |
---|
3573 | |
---|
3574 | /* Initialize the Miller Cylindrical projection |
---|
3575 | -------------------------------------------*/ |
---|
3576 | init: function() { |
---|
3577 | //no-op |
---|
3578 | }, |
---|
3579 | |
---|
3580 | |
---|
3581 | /* Miller Cylindrical forward equations--mapping lat,long to x,y |
---|
3582 | ------------------------------------------------------------*/ |
---|
3583 | forward: function(p) { |
---|
3584 | var lon=p.x; |
---|
3585 | var lat=p.y; |
---|
3586 | /* Forward equations |
---|
3587 | -----------------*/ |
---|
3588 | var dlon = Proj4js.common.adjust_lon(lon -this.long0); |
---|
3589 | var x = this.x0 + this.a * dlon; |
---|
3590 | var y = this.y0 + this.a * Math.log(Math.tan((Proj4js.common.PI / 4.0) + (lat / 2.5))) * 1.25; |
---|
3591 | |
---|
3592 | p.x=x; |
---|
3593 | p.y=y; |
---|
3594 | return p; |
---|
3595 | },//millFwd() |
---|
3596 | |
---|
3597 | /* Miller Cylindrical inverse equations--mapping x,y to lat/long |
---|
3598 | ------------------------------------------------------------*/ |
---|
3599 | inverse: function(p) { |
---|
3600 | p.x -= this.x0; |
---|
3601 | p.y -= this.y0; |
---|
3602 | |
---|
3603 | var lon = Proj4js.common.adjust_lon(this.long0 + p.x /this.a); |
---|
3604 | var lat = 2.5 * (Math.atan(Math.exp(0.8*p.y/this.a)) - Proj4js.common.PI / 4.0); |
---|
3605 | |
---|
3606 | p.x=lon; |
---|
3607 | p.y=lat; |
---|
3608 | return p; |
---|
3609 | }//millInv() |
---|
3610 | }; |
---|
3611 | /* ====================================================================== |
---|
3612 | projCode/gnom.js |
---|
3613 | ====================================================================== */ |
---|
3614 | |
---|
3615 | /***************************************************************************** |
---|
3616 | NAME GNOMONIC |
---|
3617 | |
---|
3618 | PURPOSE: Transforms input longitude and latitude to Easting and |
---|
3619 | Northing for the Gnomonic Projection. |
---|
3620 | Implementation based on the existing sterea and ortho |
---|
3621 | implementations. |
---|
3622 | |
---|
3623 | PROGRAMMER DATE |
---|
3624 | ---------- ---- |
---|
3625 | Richard Marsden November 2009 |
---|
3626 | |
---|
3627 | ALGORITHM REFERENCES |
---|
3628 | |
---|
3629 | 1. Snyder, John P., "Flattening the Earth - Two Thousand Years of Map |
---|
3630 | Projections", University of Chicago Press 1993 |
---|
3631 | |
---|
3632 | 2. Wolfram Mathworld "Gnomonic Projection" |
---|
3633 | http://mathworld.wolfram.com/GnomonicProjection.html |
---|
3634 | Accessed: 12th November 2009 |
---|
3635 | ******************************************************************************/ |
---|
3636 | |
---|
3637 | Proj4js.Proj.gnom = { |
---|
3638 | |
---|
3639 | /* Initialize the Gnomonic projection |
---|
3640 | -------------------------------------*/ |
---|
3641 | init: function(def) { |
---|
3642 | |
---|
3643 | /* Place parameters in static storage for common use |
---|
3644 | -------------------------------------------------*/ |
---|
3645 | this.sin_p14=Math.sin(this.lat0); |
---|
3646 | this.cos_p14=Math.cos(this.lat0); |
---|
3647 | // Approximation for projecting points to the horizon (infinity) |
---|
3648 | this.infinity_dist = 1000 * this.a; |
---|
3649 | this.rc = 1; |
---|
3650 | }, |
---|
3651 | |
---|
3652 | |
---|
3653 | /* Gnomonic forward equations--mapping lat,long to x,y |
---|
3654 | ---------------------------------------------------*/ |
---|
3655 | forward: function(p) { |
---|
3656 | var sinphi, cosphi; /* sin and cos value */ |
---|
3657 | var dlon; /* delta longitude value */ |
---|
3658 | var coslon; /* cos of longitude */ |
---|
3659 | var ksp; /* scale factor */ |
---|
3660 | var g; |
---|
3661 | var lon=p.x; |
---|
3662 | var lat=p.y; |
---|
3663 | /* Forward equations |
---|
3664 | -----------------*/ |
---|
3665 | dlon = Proj4js.common.adjust_lon(lon - this.long0); |
---|
3666 | |
---|
3667 | sinphi=Math.sin(lat); |
---|
3668 | cosphi=Math.cos(lat); |
---|
3669 | |
---|
3670 | coslon = Math.cos(dlon); |
---|
3671 | g = this.sin_p14 * sinphi + this.cos_p14 * cosphi * coslon; |
---|
3672 | ksp = 1.0; |
---|
3673 | if ((g > 0) || (Math.abs(g) <= Proj4js.common.EPSLN)) { |
---|
3674 | x = this.x0 + this.a * ksp * cosphi * Math.sin(dlon) / g; |
---|
3675 | y = this.y0 + this.a * ksp * (this.cos_p14 * sinphi - this.sin_p14 * cosphi * coslon) / g; |
---|
3676 | } else { |
---|
3677 | Proj4js.reportError("orthoFwdPointError"); |
---|
3678 | |
---|
3679 | // Point is in the opposing hemisphere and is unprojectable |
---|
3680 | // We still need to return a reasonable point, so we project |
---|
3681 | // to infinity, on a bearing |
---|
3682 | // equivalent to the northern hemisphere equivalent |
---|
3683 | // This is a reasonable approximation for short shapes and lines that |
---|
3684 | // straddle the horizon. |
---|
3685 | |
---|
3686 | x = this.x0 + this.infinity_dist * cosphi * Math.sin(dlon); |
---|
3687 | y = this.y0 + this.infinity_dist * (this.cos_p14 * sinphi - this.sin_p14 * cosphi * coslon); |
---|
3688 | |
---|
3689 | } |
---|
3690 | p.x=x; |
---|
3691 | p.y=y; |
---|
3692 | return p; |
---|
3693 | }, |
---|
3694 | |
---|
3695 | |
---|
3696 | inverse: function(p) { |
---|
3697 | var rh; /* Rho */ |
---|
3698 | var z; /* angle */ |
---|
3699 | var sinc, cosc; |
---|
3700 | var c; |
---|
3701 | var lon , lat; |
---|
3702 | |
---|
3703 | /* Inverse equations |
---|
3704 | -----------------*/ |
---|
3705 | p.x = (p.x - this.x0) / this.a; |
---|
3706 | p.y = (p.y - this.y0) / this.a; |
---|
3707 | |
---|
3708 | p.x /= this.k0; |
---|
3709 | p.y /= this.k0; |
---|
3710 | |
---|
3711 | if ( (rh = Math.sqrt(p.x * p.x + p.y * p.y)) ) { |
---|
3712 | c = Math.atan2(rh, this.rc); |
---|
3713 | sinc = Math.sin(c); |
---|
3714 | cosc = Math.cos(c); |
---|
3715 | |
---|
3716 | lat = Proj4js.common.asinz(cosc*this.sin_p14 + (p.y*sinc*this.cos_p14) / rh); |
---|
3717 | lon = Math.atan2(p.x*sinc, rh*this.cos_p14*cosc - p.y*this.sin_p14*sinc); |
---|
3718 | lon = Proj4js.common.adjust_lon(this.long0+lon); |
---|
3719 | } else { |
---|
3720 | lat = this.phic0; |
---|
3721 | lon = 0.0; |
---|
3722 | } |
---|
3723 | |
---|
3724 | p.x=lon; |
---|
3725 | p.y=lat; |
---|
3726 | return p; |
---|
3727 | } |
---|
3728 | }; |
---|
3729 | |
---|
3730 | |
---|
3731 | /* ====================================================================== |
---|
3732 | projCode/sinu.js |
---|
3733 | ====================================================================== */ |
---|
3734 | |
---|
3735 | /******************************************************************************* |
---|
3736 | NAME SINUSOIDAL |
---|
3737 | |
---|
3738 | PURPOSE: Transforms input longitude and latitude to Easting and |
---|
3739 | Northing for the Sinusoidal projection. The |
---|
3740 | longitude and latitude must be in radians. The Easting |
---|
3741 | and Northing values will be returned in meters. |
---|
3742 | |
---|
3743 | PROGRAMMER DATE |
---|
3744 | ---------- ---- |
---|
3745 | D. Steinwand, EROS May, 1991 |
---|
3746 | |
---|
3747 | This function was adapted from the Sinusoidal projection code (FORTRAN) in the |
---|
3748 | General Cartographic Transformation Package software which is available from |
---|
3749 | the U.S. Geological Survey National Mapping Division. |
---|
3750 | |
---|
3751 | ALGORITHM REFERENCES |
---|
3752 | |
---|
3753 | 1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological |
---|
3754 | Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United |
---|
3755 | State Government Printing Office, Washington D.C., 1987. |
---|
3756 | |
---|
3757 | 2. "Software Documentation for GCTP General Cartographic Transformation |
---|
3758 | Package", U.S. Geological Survey National Mapping Division, May 1982. |
---|
3759 | *******************************************************************************/ |
---|
3760 | |
---|
3761 | Proj4js.Proj.sinu = { |
---|
3762 | |
---|
3763 | /* Initialize the Sinusoidal projection |
---|
3764 | ------------------------------------*/ |
---|
3765 | init: function() { |
---|
3766 | /* Place parameters in static storage for common use |
---|
3767 | -------------------------------------------------*/ |
---|
3768 | this.R = 6370997.0; //Radius of earth |
---|
3769 | }, |
---|
3770 | |
---|
3771 | /* Sinusoidal forward equations--mapping lat,long to x,y |
---|
3772 | -----------------------------------------------------*/ |
---|
3773 | forward: function(p) { |
---|
3774 | var x,y,delta_lon; |
---|
3775 | var lon=p.x; |
---|
3776 | var lat=p.y; |
---|
3777 | /* Forward equations |
---|
3778 | -----------------*/ |
---|
3779 | delta_lon = Proj4js.common.adjust_lon(lon - this.long0); |
---|
3780 | x = this.R * delta_lon * Math.cos(lat) + this.x0; |
---|
3781 | y = this.R * lat + this.y0; |
---|
3782 | |
---|
3783 | p.x=x; |
---|
3784 | p.y=y; |
---|
3785 | return p; |
---|
3786 | }, |
---|
3787 | |
---|
3788 | inverse: function(p) { |
---|
3789 | var lat,temp,lon; |
---|
3790 | |
---|
3791 | /* Inverse equations |
---|
3792 | -----------------*/ |
---|
3793 | p.x -= this.x0; |
---|
3794 | p.y -= this.y0; |
---|
3795 | lat = p.y / this.R; |
---|
3796 | if (Math.abs(lat) > Proj4js.common.HALF_PI) { |
---|
3797 | Proj4js.reportError("sinu:Inv:DataError"); |
---|
3798 | } |
---|
3799 | temp = Math.abs(lat) - Proj4js.common.HALF_PI; |
---|
3800 | if (Math.abs(temp) > Proj4js.common.EPSLN) { |
---|
3801 | temp = this.long0+ p.x / (this.R *Math.cos(lat)); |
---|
3802 | lon = Proj4js.common.adjust_lon(temp); |
---|
3803 | } else { |
---|
3804 | lon = this.long0; |
---|
3805 | } |
---|
3806 | |
---|
3807 | p.x=lon; |
---|
3808 | p.y=lat; |
---|
3809 | return p; |
---|
3810 | } |
---|
3811 | }; |
---|
3812 | |
---|
3813 | |
---|
3814 | /* ====================================================================== |
---|
3815 | projCode/vandg.js |
---|
3816 | ====================================================================== */ |
---|
3817 | |
---|
3818 | /******************************************************************************* |
---|
3819 | NAME VAN DER GRINTEN |
---|
3820 | |
---|
3821 | PURPOSE: Transforms input Easting and Northing to longitude and |
---|
3822 | latitude for the Van der Grinten projection. The |
---|
3823 | Easting and Northing must be in meters. The longitude |
---|
3824 | and latitude values will be returned in radians. |
---|
3825 | |
---|
3826 | PROGRAMMER DATE |
---|
3827 | ---------- ---- |
---|
3828 | T. Mittan March, 1993 |
---|
3829 | |
---|
3830 | This function was adapted from the Van Der Grinten projection code |
---|
3831 | (FORTRAN) in the General Cartographic Transformation Package software |
---|
3832 | which is available from the U.S. Geological Survey National Mapping Division. |
---|
3833 | |
---|
3834 | ALGORITHM REFERENCES |
---|
3835 | |
---|
3836 | 1. "New Equal-Area Map Projections for Noncircular Regions", John P. Snyder, |
---|
3837 | The American Cartographer, Vol 15, No. 4, October 1988, pp. 341-355. |
---|
3838 | |
---|
3839 | 2. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological |
---|
3840 | Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United |
---|
3841 | State Government Printing Office, Washington D.C., 1987. |
---|
3842 | |
---|
3843 | 3. "Software Documentation for GCTP General Cartographic Transformation |
---|
3844 | Package", U.S. Geological Survey National Mapping Division, May 1982. |
---|
3845 | *******************************************************************************/ |
---|
3846 | |
---|
3847 | Proj4js.Proj.vandg = { |
---|
3848 | |
---|
3849 | /* Initialize the Van Der Grinten projection |
---|
3850 | ----------------------------------------*/ |
---|
3851 | init: function() { |
---|
3852 | this.R = 6370997.0; //Radius of earth |
---|
3853 | }, |
---|
3854 | |
---|
3855 | forward: function(p) { |
---|
3856 | |
---|
3857 | var lon=p.x; |
---|
3858 | var lat=p.y; |
---|
3859 | |
---|
3860 | /* Forward equations |
---|
3861 | -----------------*/ |
---|
3862 | var dlon = Proj4js.common.adjust_lon(lon - this.long0); |
---|
3863 | var x,y; |
---|
3864 | |
---|
3865 | if (Math.abs(lat) <= Proj4js.common.EPSLN) { |
---|
3866 | x = this.x0 + this.R * dlon; |
---|
3867 | y = this.y0; |
---|
3868 | } |
---|
3869 | var theta = Proj4js.common.asinz(2.0 * Math.abs(lat / Proj4js.common.PI)); |
---|
3870 | if ((Math.abs(dlon) <= Proj4js.common.EPSLN) || (Math.abs(Math.abs(lat) - Proj4js.common.HALF_PI) <= Proj4js.common.EPSLN)) { |
---|
3871 | x = this.x0; |
---|
3872 | if (lat >= 0) { |
---|
3873 | y = this.y0 + Proj4js.common.PI * this.R * Math.tan(.5 * theta); |
---|
3874 | } else { |
---|
3875 | y = this.y0 + Proj4js.common.PI * this.R * - Math.tan(.5 * theta); |
---|
3876 | } |
---|
3877 | // return(OK); |
---|
3878 | } |
---|
3879 | var al = .5 * Math.abs((Proj4js.common.PI / dlon) - (dlon / Proj4js.common.PI)); |
---|
3880 | var asq = al * al; |
---|
3881 | var sinth = Math.sin(theta); |
---|
3882 | var costh = Math.cos(theta); |
---|
3883 | |
---|
3884 | var g = costh / (sinth + costh - 1.0); |
---|
3885 | var gsq = g * g; |
---|
3886 | var m = g * (2.0 / sinth - 1.0); |
---|
3887 | var msq = m * m; |
---|
3888 | var con = Proj4js.common.PI * this.R * (al * (g - msq) + Math.sqrt(asq * (g - msq) * (g - msq) - (msq + asq) * (gsq - msq))) / (msq + asq); |
---|
3889 | if (dlon < 0) { |
---|
3890 | con = -con; |
---|
3891 | } |
---|
3892 | x = this.x0 + con; |
---|
3893 | con = Math.abs(con / (Proj4js.common.PI * this.R)); |
---|
3894 | if (lat >= 0) { |
---|
3895 | y = this.y0 + Proj4js.common.PI * this.R * Math.sqrt(1.0 - con * con - 2.0 * al * con); |
---|
3896 | } else { |
---|
3897 | y = this.y0 - Proj4js.common.PI * this.R * Math.sqrt(1.0 - con * con - 2.0 * al * con); |
---|
3898 | } |
---|
3899 | p.x = x; |
---|
3900 | p.y = y; |
---|
3901 | return p; |
---|
3902 | }, |
---|
3903 | |
---|
3904 | /* Van Der Grinten inverse equations--mapping x,y to lat/long |
---|
3905 | ---------------------------------------------------------*/ |
---|
3906 | inverse: function(p) { |
---|
3907 | var dlon; |
---|
3908 | var xx,yy,xys,c1,c2,c3; |
---|
3909 | var al,asq; |
---|
3910 | var a1; |
---|
3911 | var m1; |
---|
3912 | var con; |
---|
3913 | var th1; |
---|
3914 | var d; |
---|
3915 | |
---|
3916 | /* inverse equations |
---|
3917 | -----------------*/ |
---|
3918 | p.x -= this.x0; |
---|
3919 | p.y -= this.y0; |
---|
3920 | con = Proj4js.common.PI * this.R; |
---|
3921 | xx = p.x / con; |
---|
3922 | yy =p.y / con; |
---|
3923 | xys = xx * xx + yy * yy; |
---|
3924 | c1 = -Math.abs(yy) * (1.0 + xys); |
---|
3925 | c2 = c1 - 2.0 * yy * yy + xx * xx; |
---|
3926 | c3 = -2.0 * c1 + 1.0 + 2.0 * yy * yy + xys * xys; |
---|
3927 | d = yy * yy / c3 + (2.0 * c2 * c2 * c2 / c3 / c3 / c3 - 9.0 * c1 * c2 / c3 /c3) / 27.0; |
---|
3928 | a1 = (c1 - c2 * c2 / 3.0 / c3) / c3; |
---|
3929 | m1 = 2.0 * Math.sqrt( -a1 / 3.0); |
---|
3930 | con = ((3.0 * d) / a1) / m1; |
---|
3931 | if (Math.abs(con) > 1.0) { |
---|
3932 | if (con >= 0.0) { |
---|
3933 | con = 1.0; |
---|
3934 | } else { |
---|
3935 | con = -1.0; |
---|
3936 | } |
---|
3937 | } |
---|
3938 | th1 = Math.acos(con) / 3.0; |
---|
3939 | if (p.y >= 0) { |
---|
3940 | lat = (-m1 *Math.cos(th1 + Proj4js.common.PI / 3.0) - c2 / 3.0 / c3) * Proj4js.common.PI; |
---|
3941 | } else { |
---|
3942 | lat = -(-m1 * Math.cos(th1 + Proj4js.common.PI / 3.0) - c2 / 3.0 / c3) * Proj4js.common.PI; |
---|
3943 | } |
---|
3944 | |
---|
3945 | if (Math.abs(xx) < Proj4js.common.EPSLN) { |
---|
3946 | lon = this.long0; |
---|
3947 | } |
---|
3948 | lon = Proj4js.common.adjust_lon(this.long0 + Proj4js.common.PI * (xys - 1.0 + Math.sqrt(1.0 + 2.0 * (xx * xx - yy * yy) + xys * xys)) / 2.0 / xx); |
---|
3949 | |
---|
3950 | p.x=lon; |
---|
3951 | p.y=lat; |
---|
3952 | return p; |
---|
3953 | } |
---|
3954 | }; |
---|
3955 | /* ====================================================================== |
---|
3956 | projCode/cea.js |
---|
3957 | ====================================================================== */ |
---|
3958 | |
---|
3959 | /******************************************************************************* |
---|
3960 | NAME LAMBERT CYLINDRICAL EQUAL AREA |
---|
3961 | |
---|
3962 | PURPOSE: Transforms input longitude and latitude to Easting and |
---|
3963 | Northing for the Lambert Cylindrical Equal Area projection. |
---|
3964 | This class of projection includes the Behrmann and |
---|
3965 | Gall-Peters Projections. The |
---|
3966 | longitude and latitude must be in radians. The Easting |
---|
3967 | and Northing values will be returned in meters. |
---|
3968 | |
---|
3969 | PROGRAMMER DATE |
---|
3970 | ---------- ---- |
---|
3971 | R. Marsden August 2009 |
---|
3972 | Winwaed Software Tech LLC, http://www.winwaed.com |
---|
3973 | |
---|
3974 | This function was adapted from the Miller Cylindrical Projection in the Proj4JS |
---|
3975 | library. |
---|
3976 | |
---|
3977 | Note: This implementation assumes a Spherical Earth. The (commented) code |
---|
3978 | has been included for the ellipsoidal forward transform, but derivation of |
---|
3979 | the ellispoidal inverse transform is beyond me. Note that most of the |
---|
3980 | Proj4JS implementations do NOT currently support ellipsoidal figures. |
---|
3981 | Therefore this is not seen as a problem - especially this lack of support |
---|
3982 | is explicitly stated here. |
---|
3983 | |
---|
3984 | ALGORITHM REFERENCES |
---|
3985 | |
---|
3986 | 1. "Cartographic Projection Procedures for the UNIX Environment - |
---|
3987 | A User's Manual" by Gerald I. Evenden, USGS Open File Report 90-284 |
---|
3988 | and Release 4 Interim Reports (2003) |
---|
3989 | |
---|
3990 | 2. Snyder, John P., "Flattening the Earth - Two Thousand Years of Map |
---|
3991 | Projections", Univ. Chicago Press, 1993 |
---|
3992 | *******************************************************************************/ |
---|
3993 | |
---|
3994 | Proj4js.Proj.cea = { |
---|
3995 | |
---|
3996 | /* Initialize the Cylindrical Equal Area projection |
---|
3997 | -------------------------------------------*/ |
---|
3998 | init: function() { |
---|
3999 | //no-op |
---|
4000 | }, |
---|
4001 | |
---|
4002 | |
---|
4003 | /* Cylindrical Equal Area forward equations--mapping lat,long to x,y |
---|
4004 | ------------------------------------------------------------*/ |
---|
4005 | forward: function(p) { |
---|
4006 | var lon=p.x; |
---|
4007 | var lat=p.y; |
---|
4008 | /* Forward equations |
---|
4009 | -----------------*/ |
---|
4010 | dlon = Proj4js.common.adjust_lon(lon -this.long0); |
---|
4011 | var x = this.x0 + this.a * dlon * Math.cos(this.lat_ts); |
---|
4012 | var y = this.y0 + this.a * Math.sin(lat) / Math.cos(this.lat_ts); |
---|
4013 | /* Elliptical Forward Transform |
---|
4014 | Not implemented due to a lack of a matchign inverse function |
---|
4015 | { |
---|
4016 | var Sin_Lat = Math.sin(lat); |
---|
4017 | var Rn = this.a * (Math.sqrt(1.0e0 - this.es * Sin_Lat * Sin_Lat )); |
---|
4018 | x = this.x0 + this.a * dlon * Math.cos(this.lat_ts); |
---|
4019 | y = this.y0 + Rn * Math.sin(lat) / Math.cos(this.lat_ts); |
---|
4020 | } |
---|
4021 | */ |
---|
4022 | |
---|
4023 | |
---|
4024 | p.x=x; |
---|
4025 | p.y=y; |
---|
4026 | return p; |
---|
4027 | },//ceaFwd() |
---|
4028 | |
---|
4029 | /* Cylindrical Equal Area inverse equations--mapping x,y to lat/long |
---|
4030 | ------------------------------------------------------------*/ |
---|
4031 | inverse: function(p) { |
---|
4032 | p.x -= this.x0; |
---|
4033 | p.y -= this.y0; |
---|
4034 | |
---|
4035 | var lon = Proj4js.common.adjust_lon( this.long0 + (p.x / this.a) / Math.cos(this.lat_ts) ); |
---|
4036 | |
---|
4037 | var lat = Math.asin( (p.y/this.a) * Math.cos(this.lat_ts) ); |
---|
4038 | |
---|
4039 | p.x=lon; |
---|
4040 | p.y=lat; |
---|
4041 | return p; |
---|
4042 | }//ceaInv() |
---|
4043 | }; |
---|
4044 | /* ====================================================================== |
---|
4045 | projCode/eqc.js |
---|
4046 | ====================================================================== */ |
---|
4047 | |
---|
4048 | /* similar to equi.js FIXME proj4 uses eqc */ |
---|
4049 | Proj4js.Proj.eqc = { |
---|
4050 | init : function() { |
---|
4051 | |
---|
4052 | if(!this.x0) this.x0=0; |
---|
4053 | if(!this.y0) this.y0=0; |
---|
4054 | if(!this.lat0) this.lat0=0; |
---|
4055 | if(!this.long0) this.long0=0; |
---|
4056 | if(!this.lat_ts) this.lat_ts=0; |
---|
4057 | if (!this.title) this.title = "Equidistant Cylindrical (Plate Carre)"; |
---|
4058 | |
---|
4059 | this.rc= Math.cos(this.lat_ts); |
---|
4060 | }, |
---|
4061 | |
---|
4062 | |
---|
4063 | // forward equations--mapping lat,long to x,y |
---|
4064 | // ----------------------------------------------------------------- |
---|
4065 | forward : function(p) { |
---|
4066 | |
---|
4067 | var lon= p.x; |
---|
4068 | var lat= p.y; |
---|
4069 | |
---|
4070 | var dlon = Proj4js.common.adjust_lon(lon - this.long0); |
---|
4071 | var dlat = Proj4js.common.adjust_lat(lat - this.lat0 ); |
---|
4072 | p.x= this.x0 + (this.a*dlon*this.rc); |
---|
4073 | p.y= this.y0 + (this.a*dlat ); |
---|
4074 | return p; |
---|
4075 | }, |
---|
4076 | |
---|
4077 | // inverse equations--mapping x,y to lat/long |
---|
4078 | // ----------------------------------------------------------------- |
---|
4079 | inverse : function(p) { |
---|
4080 | |
---|
4081 | var x= p.x; |
---|
4082 | var y= p.y; |
---|
4083 | |
---|
4084 | p.x= Proj4js.common.adjust_lon(this.long0 + ((x - this.x0)/(this.a*this.rc))); |
---|
4085 | p.y= Proj4js.common.adjust_lat(this.lat0 + ((y - this.y0)/(this.a ))); |
---|
4086 | return p; |
---|
4087 | } |
---|
4088 | |
---|
4089 | }; |
---|
4090 | /* ====================================================================== |
---|
4091 | projCode/cass.js |
---|
4092 | ====================================================================== */ |
---|
4093 | |
---|
4094 | /******************************************************************************* |
---|
4095 | NAME CASSINI |
---|
4096 | |
---|
4097 | PURPOSE: Transforms input longitude and latitude to Easting and |
---|
4098 | Northing for the Cassini projection. The |
---|
4099 | longitude and latitude must be in radians. The Easting |
---|
4100 | and Northing values will be returned in meters. |
---|
4101 | Ported from PROJ.4. |
---|
4102 | |
---|
4103 | |
---|
4104 | ALGORITHM REFERENCES |
---|
4105 | |
---|
4106 | 1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological |
---|
4107 | Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United |
---|
4108 | State Government Printing Office, Washington D.C., 1987. |
---|
4109 | |
---|
4110 | 2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections", |
---|
4111 | U.S. Geological Survey Professional Paper 1453 , United State Government |
---|
4112 | *******************************************************************************/ |
---|
4113 | |
---|
4114 | |
---|
4115 | //Proj4js.defs["EPSG:28191"] = "+proj=cass +lat_0=31.73409694444445 +lon_0=35.21208055555556 +x_0=170251.555 +y_0=126867.909 +a=6378300.789 +b=6356566.435 +towgs84=-275.722,94.7824,340.894,-8.001,-4.42,-11.821,1 +units=m +no_defs"; |
---|
4116 | |
---|
4117 | // Initialize the Cassini projection |
---|
4118 | // ----------------------------------------------------------------- |
---|
4119 | |
---|
4120 | Proj4js.Proj.cass = { |
---|
4121 | init : function() { |
---|
4122 | if (!this.sphere) { |
---|
4123 | this.en = this.pj_enfn(this.es) |
---|
4124 | this.m0 = this.pj_mlfn(this.lat0, Math.sin(this.lat0), Math.cos(this.lat0), this.en); |
---|
4125 | } |
---|
4126 | }, |
---|
4127 | |
---|
4128 | C1: .16666666666666666666, |
---|
4129 | C2: .00833333333333333333, |
---|
4130 | C3: .04166666666666666666, |
---|
4131 | C4: .33333333333333333333, |
---|
4132 | C5: .06666666666666666666, |
---|
4133 | |
---|
4134 | |
---|
4135 | /* Cassini forward equations--mapping lat,long to x,y |
---|
4136 | -----------------------------------------------------------------------*/ |
---|
4137 | forward: function(p) { |
---|
4138 | |
---|
4139 | /* Forward equations |
---|
4140 | -----------------*/ |
---|
4141 | var x,y; |
---|
4142 | var lam=p.x; |
---|
4143 | var phi=p.y; |
---|
4144 | lam = Proj4js.common.adjust_lon(lam - this.long0); |
---|
4145 | |
---|
4146 | if (this.sphere) { |
---|
4147 | x = Math.asin(Math.cos(phi) * Math.sin(lam)); |
---|
4148 | y = Math.atan2(Math.tan(phi) , Math.cos(lam)) - this.phi0; |
---|
4149 | } else { |
---|
4150 | //ellipsoid |
---|
4151 | this.n = Math.sin(phi); |
---|
4152 | this.c = Math.cos(phi); |
---|
4153 | y = this.pj_mlfn(phi, this.n, this.c, this.en); |
---|
4154 | this.n = 1./Math.sqrt(1. - this.es * this.n * this.n); |
---|
4155 | this.tn = Math.tan(phi); |
---|
4156 | this.t = this.tn * this.tn; |
---|
4157 | this.a1 = lam * this.c; |
---|
4158 | this.c *= this.es * this.c / (1 - this.es); |
---|
4159 | this.a2 = this.a1 * this.a1; |
---|
4160 | x = this.n * this.a1 * (1. - this.a2 * this.t * (this.C1 - (8. - this.t + 8. * this.c) * this.a2 * this.C2)); |
---|
4161 | y -= this.m0 - this.n * this.tn * this.a2 * (.5 + (5. - this.t + 6. * this.c) * this.a2 * this.C3); |
---|
4162 | } |
---|
4163 | |
---|
4164 | p.x = this.a*x + this.x0; |
---|
4165 | p.y = this.a*y + this.y0; |
---|
4166 | return p; |
---|
4167 | },//cassFwd() |
---|
4168 | |
---|
4169 | /* Inverse equations |
---|
4170 | -----------------*/ |
---|
4171 | inverse: function(p) { |
---|
4172 | p.x -= this.x0; |
---|
4173 | p.y -= this.y0; |
---|
4174 | var x = p.x/this.a; |
---|
4175 | var y = p.y/this.a; |
---|
4176 | |
---|
4177 | if (this.sphere) { |
---|
4178 | this.dd = y + this.lat0; |
---|
4179 | phi = Math.asin(Math.sin(this.dd) * Math.cos(x)); |
---|
4180 | lam = Math.atan2(Math.tan(x), Math.cos(this.dd)); |
---|
4181 | } else { |
---|
4182 | /* ellipsoid */ |
---|
4183 | ph1 = this.pj_inv_mlfn(this.m0 + y, this.es, this.en); |
---|
4184 | this.tn = Math.tan(ph1); |
---|
4185 | this.t = this.tn * this.tn; |
---|
4186 | this.n = Math.sin(ph1); |
---|
4187 | this.r = 1. / (1. - this.es * this.n * this.n); |
---|
4188 | this.n = Math.sqrt(this.r); |
---|
4189 | this.r *= (1. - this.es) * this.n; |
---|
4190 | this.dd = x / this.n; |
---|
4191 | this.d2 = this.dd * this.dd; |
---|
4192 | phi = ph1 - (this.n * this.tn / this.r) * this.d2 * (.5 - (1. + 3. * this.t) * this.d2 * this.C3); |
---|
4193 | lam = this.dd * (1. + this.t * this.d2 * (-this.C4 + (1. + 3. * this.t) * this.d2 * this.C5)) / Math.cos(ph1); |
---|
4194 | } |
---|
4195 | p.x = Proj4js.common.adjust_lon(this.long0+lam); |
---|
4196 | p.y = phi; |
---|
4197 | return p; |
---|
4198 | },//lamazInv() |
---|
4199 | |
---|
4200 | |
---|
4201 | //code from the PROJ.4 pj_mlfn.c file; this may be useful for other projections |
---|
4202 | pj_enfn: function(es) { |
---|
4203 | en = new Array(); |
---|
4204 | en[0] = this.C00 - es * (this.C02 + es * (this.C04 + es * (this.C06 + es * this.C08))); |
---|
4205 | en[1] = es * (this.C22 - es * (this.C04 + es * (this.C06 + es * this.C08))); |
---|
4206 | var t = es * es; |
---|
4207 | en[2] = t * (this.C44 - es * (this.C46 + es * this.C48)); |
---|
4208 | t *= es; |
---|
4209 | en[3] = t * (this.C66 - es * this.C68); |
---|
4210 | en[4] = t * es * this.C88; |
---|
4211 | return en; |
---|
4212 | }, |
---|
4213 | |
---|
4214 | pj_mlfn: function(phi, sphi, cphi, en) { |
---|
4215 | cphi *= sphi; |
---|
4216 | sphi *= sphi; |
---|
4217 | return(en[0] * phi - cphi * (en[1] + sphi*(en[2]+ sphi*(en[3] + sphi*en[4])))); |
---|
4218 | }, |
---|
4219 | |
---|
4220 | pj_inv_mlfn: function(arg, es, en) { |
---|
4221 | k = 1./(1.-es); |
---|
4222 | phi = arg; |
---|
4223 | for (i = Proj4js.common.MAX_ITER; i ; --i) { /* rarely goes over 2 iterations */ |
---|
4224 | s = Math.sin(phi); |
---|
4225 | t = 1. - es * s * s; |
---|
4226 | //t = this.pj_mlfn(phi, s, Math.cos(phi), en) - arg; |
---|
4227 | //phi -= t * (t * Math.sqrt(t)) * k; |
---|
4228 | t = (this.pj_mlfn(phi, s, Math.cos(phi), en) - arg) * (t * Math.sqrt(t)) * k; |
---|
4229 | phi -= t; |
---|
4230 | if (Math.abs(t) < Proj4js.common.EPSLN) |
---|
4231 | return phi; |
---|
4232 | } |
---|
4233 | Proj4js.reportError("cass:pj_inv_mlfn: Convergence error"); |
---|
4234 | return phi; |
---|
4235 | }, |
---|
4236 | |
---|
4237 | /* meridinal distance for ellipsoid and inverse |
---|
4238 | ** 8th degree - accurate to < 1e-5 meters when used in conjuction |
---|
4239 | ** with typical major axis values. |
---|
4240 | ** Inverse determines phi to EPS (1e-11) radians, about 1e-6 seconds. |
---|
4241 | */ |
---|
4242 | C00: 1.0, |
---|
4243 | C02: .25, |
---|
4244 | C04: .046875, |
---|
4245 | C06: .01953125, |
---|
4246 | C08: .01068115234375, |
---|
4247 | C22: .75, |
---|
4248 | C44: .46875, |
---|
4249 | C46: .01302083333333333333, |
---|
4250 | C48: .00712076822916666666, |
---|
4251 | C66: .36458333333333333333, |
---|
4252 | C68: .00569661458333333333, |
---|
4253 | C88: .3076171875 |
---|
4254 | |
---|
4255 | } |
---|
4256 | /* ====================================================================== |
---|
4257 | projCode/gauss.js |
---|
4258 | ====================================================================== */ |
---|
4259 | |
---|
4260 | |
---|
4261 | Proj4js.Proj.gauss = { |
---|
4262 | |
---|
4263 | init : function() { |
---|
4264 | sphi = Math.sin(this.lat0); |
---|
4265 | cphi = Math.cos(this.lat0); |
---|
4266 | cphi *= cphi; |
---|
4267 | this.rc = Math.sqrt(1.0 - this.es) / (1.0 - this.es * sphi * sphi); |
---|
4268 | this.C = Math.sqrt(1.0 + this.es * cphi * cphi / (1.0 - this.es)); |
---|
4269 | this.phic0 = Math.asin(sphi / this.C); |
---|
4270 | this.ratexp = 0.5 * this.C * this.e; |
---|
4271 | this.K = Math.tan(0.5 * this.phic0 + Proj4js.common.FORTPI) / (Math.pow(Math.tan(0.5*this.lat0 + Proj4js.common.FORTPI), this.C) * Proj4js.common.srat(this.e*sphi, this.ratexp)); |
---|
4272 | }, |
---|
4273 | |
---|
4274 | forward : function(p) { |
---|
4275 | var lon = p.x; |
---|
4276 | var lat = p.y; |
---|
4277 | |
---|
4278 | p.y = 2.0 * Math.atan( this.K * Math.pow(Math.tan(0.5 * lat + Proj4js.common.FORTPI), this.C) * Proj4js.common.srat(this.e * Math.sin(lat), this.ratexp) ) - Proj4js.common.HALF_PI; |
---|
4279 | p.x = this.C * lon; |
---|
4280 | return p; |
---|
4281 | }, |
---|
4282 | |
---|
4283 | inverse : function(p) { |
---|
4284 | var DEL_TOL = 1e-14; |
---|
4285 | var lon = p.x / this.C; |
---|
4286 | var lat = p.y; |
---|
4287 | num = Math.pow(Math.tan(0.5 * lat + Proj4js.common.FORTPI)/this.K, 1./this.C); |
---|
4288 | for (var i = Proj4js.common.MAX_ITER; i>0; --i) { |
---|
4289 | lat = 2.0 * Math.atan(num * Proj4js.common.srat(this.e * Math.sin(p.y), -0.5 * this.e)) - Proj4js.common.HALF_PI; |
---|
4290 | if (Math.abs(lat - p.y) < DEL_TOL) break; |
---|
4291 | p.y = lat; |
---|
4292 | } |
---|
4293 | /* convergence failed */ |
---|
4294 | if (!i) { |
---|
4295 | Proj4js.reportError("gauss:inverse:convergence failed"); |
---|
4296 | return null; |
---|
4297 | } |
---|
4298 | p.x = lon; |
---|
4299 | p.y = lat; |
---|
4300 | return p; |
---|
4301 | } |
---|
4302 | }; |
---|
4303 | |
---|
4304 | /* ====================================================================== |
---|
4305 | projCode/omerc.js |
---|
4306 | ====================================================================== */ |
---|
4307 | |
---|
4308 | /******************************************************************************* |
---|
4309 | NAME OBLIQUE MERCATOR (HOTINE) |
---|
4310 | |
---|
4311 | PURPOSE: Transforms input longitude and latitude to Easting and |
---|
4312 | Northing for the Oblique Mercator projection. The |
---|
4313 | longitude and latitude must be in radians. The Easting |
---|
4314 | and Northing values will be returned in meters. |
---|
4315 | |
---|
4316 | PROGRAMMER DATE |
---|
4317 | ---------- ---- |
---|
4318 | T. Mittan Mar, 1993 |
---|
4319 | |
---|
4320 | ALGORITHM REFERENCES |
---|
4321 | |
---|
4322 | 1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological |
---|
4323 | Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United |
---|
4324 | State Government Printing Office, Washington D.C., 1987. |
---|
4325 | |
---|
4326 | 2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections", |
---|
4327 | U.S. Geological Survey Professional Paper 1453 , United State Government |
---|
4328 | Printing Office, Washington D.C., 1989. |
---|
4329 | *******************************************************************************/ |
---|
4330 | |
---|
4331 | Proj4js.Proj.omerc = { |
---|
4332 | |
---|
4333 | /* Initialize the Oblique Mercator projection |
---|
4334 | ------------------------------------------*/ |
---|
4335 | init: function() { |
---|
4336 | if (!this.mode) this.mode=0; |
---|
4337 | if (!this.lon1) {this.lon1=0;this.mode=1;} |
---|
4338 | if (!this.lon2) this.lon2=0; |
---|
4339 | if (!this.lat2) this.lat2=0; |
---|
4340 | |
---|
4341 | /* Place parameters in static storage for common use |
---|
4342 | -------------------------------------------------*/ |
---|
4343 | var temp = this.b/ this.a; |
---|
4344 | var es = 1.0 - Math.pow(temp,2); |
---|
4345 | var e = Math.sqrt(es); |
---|
4346 | |
---|
4347 | this.sin_p20=Math.sin(this.lat0); |
---|
4348 | this.cos_p20=Math.cos(this.lat0); |
---|
4349 | |
---|
4350 | this.con = 1.0 - this.es * this.sin_p20 * this.sin_p20; |
---|
4351 | this.com = Math.sqrt(1.0 - es); |
---|
4352 | this.bl = Math.sqrt(1.0 + this.es * Math.pow(this.cos_p20,4.0)/(1.0 - es)); |
---|
4353 | this.al = this.a * this.bl * this.k0 * this.com / this.con; |
---|
4354 | if (Math.abs(this.lat0) < Proj4js.common.EPSLN) { |
---|
4355 | this.ts = 1.0; |
---|
4356 | this.d = 1.0; |
---|
4357 | this.el = 1.0; |
---|
4358 | } else { |
---|
4359 | this.ts = Proj4js.common.tsfnz(this.e,this.lat0,this.sin_p20); |
---|
4360 | this.con = Math.sqrt(this.con); |
---|
4361 | this.d = this.bl * this.com / (this.cos_p20 * this.con); |
---|
4362 | if ((this.d * this.d - 1.0) > 0.0) { |
---|
4363 | if (this.lat0 >= 0.0) { |
---|
4364 | this.f = this.d + Math.sqrt(this.d * this.d - 1.0); |
---|
4365 | } else { |
---|
4366 | this.f = this.d - Math.sqrt(this.d * this.d - 1.0); |
---|
4367 | } |
---|
4368 | } else { |
---|
4369 | this.f = this.d; |
---|
4370 | } |
---|
4371 | this.el = this.f * Math.pow(this.ts,this.bl); |
---|
4372 | } |
---|
4373 | |
---|
4374 | //this.longc=52.60353916666667; |
---|
4375 | |
---|
4376 | if (this.mode != 0) { |
---|
4377 | this.g = .5 * (this.f - 1.0/this.f); |
---|
4378 | this.gama = Proj4js.common.asinz(Math.sin(this.alpha) / this.d); |
---|
4379 | this.longc= this.longc - Proj4js.common.asinz(this.g * Math.tan(this.gama))/this.bl; |
---|
4380 | |
---|
4381 | /* Report parameters common to format B |
---|
4382 | -------------------------------------*/ |
---|
4383 | //genrpt(azimuth * R2D,"Azimuth of Central Line: "); |
---|
4384 | //cenlon(lon_origin); |
---|
4385 | // cenlat(lat_origin); |
---|
4386 | |
---|
4387 | this.con = Math.abs(this.lat0); |
---|
4388 | if ((this.con > Proj4js.common.EPSLN) && (Math.abs(this.con - Proj4js.common.HALF_PI) > Proj4js.common.EPSLN)) { |
---|
4389 | this.singam=Math.sin(this.gama); |
---|
4390 | this.cosgam=Math.cos(this.gama); |
---|
4391 | |
---|
4392 | this.sinaz=Math.sin(this.alpha); |
---|
4393 | this.cosaz=Math.cos(this.alpha); |
---|
4394 | |
---|
4395 | if (this.lat0>= 0) { |
---|
4396 | this.u = (this.al / this.bl) * Math.atan(Math.sqrt(this.d*this.d - 1.0)/this.cosaz); |
---|
4397 | } else { |
---|
4398 | this.u = -(this.al / this.bl) *Math.atan(Math.sqrt(this.d*this.d - 1.0)/this.cosaz); |
---|
4399 | } |
---|
4400 | } else { |
---|
4401 | Proj4js.reportError("omerc:Init:DataError"); |
---|
4402 | } |
---|
4403 | } else { |
---|
4404 | this.sinphi =Math. sin(this.at1); |
---|
4405 | this.ts1 = Proj4js.common.tsfnz(this.e,this.lat1,this.sinphi); |
---|
4406 | this.sinphi = Math.sin(this.lat2); |
---|
4407 | this.ts2 = Proj4js.common.tsfnz(this.e,this.lat2,this.sinphi); |
---|
4408 | this.h = Math.pow(this.ts1,this.bl); |
---|
4409 | this.l = Math.pow(this.ts2,this.bl); |
---|
4410 | this.f = this.el/this.h; |
---|
4411 | this.g = .5 * (this.f - 1.0/this.f); |
---|
4412 | this.j = (this.el * this.el - this.l * this.h)/(this.el * this.el + this.l * this.h); |
---|
4413 | this.p = (this.l - this.h) / (this.l + this.h); |
---|
4414 | this.dlon = this.lon1 - this.lon2; |
---|
4415 | if (this.dlon < -Proj4js.common.PI) this.lon2 = this.lon2 - 2.0 * Proj4js.common.PI; |
---|
4416 | if (this.dlon > Proj4js.common.PI) this.lon2 = this.lon2 + 2.0 * Proj4js.common.PI; |
---|
4417 | this.dlon = this.lon1 - this.lon2; |
---|
4418 | this.longc = .5 * (this.lon1 + this.lon2) -Math.atan(this.j * Math.tan(.5 * this.bl * this.dlon)/this.p)/this.bl; |
---|
4419 | this.dlon = Proj4js.common.adjust_lon(this.lon1 - this.longc); |
---|
4420 | this.gama = Math.atan(Math.sin(this.bl * this.dlon)/this.g); |
---|
4421 | this.alpha = Proj4js.common.asinz(this.d * Math.sin(this.gama)); |
---|
4422 | |
---|
4423 | /* Report parameters common to format A |
---|
4424 | -------------------------------------*/ |
---|
4425 | |
---|
4426 | if (Math.abs(this.lat1 - this.lat2) <= Proj4js.common.EPSLN) { |
---|
4427 | Proj4js.reportError("omercInitDataError"); |
---|
4428 | //return(202); |
---|
4429 | } else { |
---|
4430 | this.con = Math.abs(this.lat1); |
---|
4431 | } |
---|
4432 | if ((this.con <= Proj4js.common.EPSLN) || (Math.abs(this.con - HALF_PI) <= Proj4js.common.EPSLN)) { |
---|
4433 | Proj4js.reportError("omercInitDataError"); |
---|
4434 | //return(202); |
---|
4435 | } else { |
---|
4436 | if (Math.abs(Math.abs(this.lat0) - Proj4js.common.HALF_PI) <= Proj4js.common.EPSLN) { |
---|
4437 | Proj4js.reportError("omercInitDataError"); |
---|
4438 | //return(202); |
---|
4439 | } |
---|
4440 | } |
---|
4441 | |
---|
4442 | this.singam=Math.sin(this.gam); |
---|
4443 | this.cosgam=Math.cos(this.gam); |
---|
4444 | |
---|
4445 | this.sinaz=Math.sin(this.alpha); |
---|
4446 | this.cosaz=Math.cos(this.alpha); |
---|
4447 | |
---|
4448 | |
---|
4449 | if (this.lat0 >= 0) { |
---|
4450 | this.u = (this.al/this.bl) * Math.atan(Math.sqrt(this.d * this.d - 1.0)/this.cosaz); |
---|
4451 | } else { |
---|
4452 | this.u = -(this.al/this.bl) * Math.atan(Math.sqrt(this.d * this.d - 1.0)/this.cosaz); |
---|
4453 | } |
---|
4454 | } |
---|
4455 | }, |
---|
4456 | |
---|
4457 | |
---|
4458 | /* Oblique Mercator forward equations--mapping lat,long to x,y |
---|
4459 | ----------------------------------------------------------*/ |
---|
4460 | forward: function(p) { |
---|
4461 | var theta; /* angle */ |
---|
4462 | var sin_phi, cos_phi;/* sin and cos value */ |
---|
4463 | var b; /* temporary values */ |
---|
4464 | var c, t, tq; /* temporary values */ |
---|
4465 | var con, n, ml; /* cone constant, small m */ |
---|
4466 | var q,us,vl; |
---|
4467 | var ul,vs; |
---|
4468 | var s; |
---|
4469 | var dlon; |
---|
4470 | var ts1; |
---|
4471 | |
---|
4472 | var lon=p.x; |
---|
4473 | var lat=p.y; |
---|
4474 | /* Forward equations |
---|
4475 | -----------------*/ |
---|
4476 | sin_phi = Math.sin(lat); |
---|
4477 | dlon = Proj4js.common.adjust_lon(lon - this.longc); |
---|
4478 | vl = Math.sin(this.bl * dlon); |
---|
4479 | if (Math.abs(Math.abs(lat) - Proj4js.common.HALF_PI) > Proj4js.common.EPSLN) { |
---|
4480 | ts1 = Proj4js.common.tsfnz(this.e,lat,sin_phi); |
---|
4481 | q = this.el / (Math.pow(ts1,this.bl)); |
---|
4482 | s = .5 * (q - 1.0 / q); |
---|
4483 | t = .5 * (q + 1.0/ q); |
---|
4484 | ul = (s * this.singam - vl * this.cosgam) / t; |
---|
4485 | con = Math.cos(this.bl * dlon); |
---|
4486 | if (Math.abs(con) < .0000001) { |
---|
4487 | us = this.al * this.bl * dlon; |
---|
4488 | } else { |
---|
4489 | us = this.al * Math.atan((s * this.cosgam + vl * this.singam) / con)/this.bl; |
---|
4490 | if (con < 0) us = us + Proj4js.common.PI * this.al / this.bl; |
---|
4491 | } |
---|
4492 | } else { |
---|
4493 | if (lat >= 0) { |
---|
4494 | ul = this.singam; |
---|
4495 | } else { |
---|
4496 | ul = -this.singam; |
---|
4497 | } |
---|
4498 | us = this.al * lat / this.bl; |
---|
4499 | } |
---|
4500 | if (Math.abs(Math.abs(ul) - 1.0) <= Proj4js.common.EPSLN) { |
---|
4501 | //alert("Point projects into infinity","omer-for"); |
---|
4502 | Proj4js.reportError("omercFwdInfinity"); |
---|
4503 | //return(205); |
---|
4504 | } |
---|
4505 | vs = .5 * this.al * Math.log((1.0 - ul)/(1.0 + ul)) / this.bl; |
---|
4506 | us = us - this.u; |
---|
4507 | var x = this.x0 + vs * this.cosaz + us * this.sinaz; |
---|
4508 | var y = this.y0 + us * this.cosaz - vs * this.sinaz; |
---|
4509 | |
---|
4510 | p.x=x; |
---|
4511 | p.y=y; |
---|
4512 | return p; |
---|
4513 | }, |
---|
4514 | |
---|
4515 | inverse: function(p) { |
---|
4516 | var delta_lon; /* Delta longitude (Given longitude - center */ |
---|
4517 | var theta; /* angle */ |
---|
4518 | var delta_theta; /* adjusted longitude */ |
---|
4519 | var sin_phi, cos_phi;/* sin and cos value */ |
---|
4520 | var b; /* temporary values */ |
---|
4521 | var c, t, tq; /* temporary values */ |
---|
4522 | var con, n, ml; /* cone constant, small m */ |
---|
4523 | var vs,us,q,s,ts1; |
---|
4524 | var vl,ul,bs; |
---|
4525 | var dlon; |
---|
4526 | var flag; |
---|
4527 | |
---|
4528 | /* Inverse equations |
---|
4529 | -----------------*/ |
---|
4530 | p.x -= this.x0; |
---|
4531 | p.y -= this.y0; |
---|
4532 | flag = 0; |
---|
4533 | vs = p.x * this.cosaz - p.y * this.sinaz; |
---|
4534 | us = p.y * this.cosaz + p.x * this.sinaz; |
---|
4535 | us = us + this.u; |
---|
4536 | q = Math.exp(-this.bl * vs / this.al); |
---|
4537 | s = .5 * (q - 1.0/q); |
---|
4538 | t = .5 * (q + 1.0/q); |
---|
4539 | vl = Math.sin(this.bl * us / this.al); |
---|
4540 | ul = (vl * this.cosgam + s * this.singam)/t; |
---|
4541 | if (Math.abs(Math.abs(ul) - 1.0) <= Proj4js.common.EPSLN) |
---|
4542 | { |
---|
4543 | lon = this.longc; |
---|
4544 | if (ul >= 0.0) { |
---|
4545 | lat = Proj4js.common.HALF_PI; |
---|
4546 | } else { |
---|
4547 | lat = -Proj4js.common.HALF_PI; |
---|
4548 | } |
---|
4549 | } else { |
---|
4550 | con = 1.0 / this.bl; |
---|
4551 | ts1 =Math.pow((this.el / Math.sqrt((1.0 + ul) / (1.0 - ul))),con); |
---|
4552 | lat = Proj4js.common.phi2z(this.e,ts1); |
---|
4553 | //if (flag != 0) |
---|
4554 | //return(flag); |
---|
4555 | //~ con = Math.cos(this.bl * us /al); |
---|
4556 | theta = this.longc - Math.atan2((s * this.cosgam - vl * this.singam) , con)/this.bl; |
---|
4557 | lon = Proj4js.common.adjust_lon(theta); |
---|
4558 | } |
---|
4559 | p.x=lon; |
---|
4560 | p.y=lat; |
---|
4561 | return p; |
---|
4562 | } |
---|
4563 | }; |
---|
4564 | /* ====================================================================== |
---|
4565 | projCode/lcc.js |
---|
4566 | ====================================================================== */ |
---|
4567 | |
---|
4568 | /******************************************************************************* |
---|
4569 | NAME LAMBERT CONFORMAL CONIC |
---|
4570 | |
---|
4571 | PURPOSE: Transforms input longitude and latitude to Easting and |
---|
4572 | Northing for the Lambert Conformal Conic projection. The |
---|
4573 | longitude and latitude must be in radians. The Easting |
---|
4574 | and Northing values will be returned in meters. |
---|
4575 | |
---|
4576 | |
---|
4577 | ALGORITHM REFERENCES |
---|
4578 | |
---|
4579 | 1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological |
---|
4580 | Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United |
---|
4581 | State Government Printing Office, Washington D.C., 1987. |
---|
4582 | |
---|
4583 | 2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections", |
---|
4584 | U.S. Geological Survey Professional Paper 1453 , United State Government |
---|
4585 | *******************************************************************************/ |
---|
4586 | |
---|
4587 | |
---|
4588 | //<2104> +proj=lcc +lat_1=10.16666666666667 +lat_0=10.16666666666667 +lon_0=-71.60561777777777 +k_0=1 +x0=-17044 +x0=-23139.97 +ellps=intl +units=m +no_defs no_defs |
---|
4589 | |
---|
4590 | // Initialize the Lambert Conformal conic projection |
---|
4591 | // ----------------------------------------------------------------- |
---|
4592 | |
---|
4593 | //Proj4js.Proj.lcc = Class.create(); |
---|
4594 | Proj4js.Proj.lcc = { |
---|
4595 | init : function() { |
---|
4596 | |
---|
4597 | // array of: r_maj,r_min,lat1,lat2,c_lon,c_lat,false_east,false_north |
---|
4598 | //double c_lat; /* center latitude */ |
---|
4599 | //double c_lon; /* center longitude */ |
---|
4600 | //double lat1; /* first standard parallel */ |
---|
4601 | //double lat2; /* second standard parallel */ |
---|
4602 | //double r_maj; /* major axis */ |
---|
4603 | //double r_min; /* minor axis */ |
---|
4604 | //double false_east; /* x offset in meters */ |
---|
4605 | //double false_north; /* y offset in meters */ |
---|
4606 | |
---|
4607 | if (!this.lat2){this.lat2=this.lat0;}//if lat2 is not defined |
---|
4608 | if (!this.k0) this.k0 = 1.0; |
---|
4609 | |
---|
4610 | // Standard Parallels cannot be equal and on opposite sides of the equator |
---|
4611 | if (Math.abs(this.lat1+this.lat2) < Proj4js.common.EPSLN) { |
---|
4612 | Proj4js.reportError("lcc:init: Equal Latitudes"); |
---|
4613 | return; |
---|
4614 | } |
---|
4615 | |
---|
4616 | var temp = this.b / this.a; |
---|
4617 | this.e = Math.sqrt(1.0 - temp*temp); |
---|
4618 | |
---|
4619 | var sin1 = Math.sin(this.lat1); |
---|
4620 | var cos1 = Math.cos(this.lat1); |
---|
4621 | var ms1 = Proj4js.common.msfnz(this.e, sin1, cos1); |
---|
4622 | var ts1 = Proj4js.common.tsfnz(this.e, this.lat1, sin1); |
---|
4623 | |
---|
4624 | var sin2 = Math.sin(this.lat2); |
---|
4625 | var cos2 = Math.cos(this.lat2); |
---|
4626 | var ms2 = Proj4js.common.msfnz(this.e, sin2, cos2); |
---|
4627 | var ts2 = Proj4js.common.tsfnz(this.e, this.lat2, sin2); |
---|
4628 | |
---|
4629 | var ts0 = Proj4js.common.tsfnz(this.e, this.lat0, Math.sin(this.lat0)); |
---|
4630 | |
---|
4631 | if (Math.abs(this.lat1 - this.lat2) > Proj4js.common.EPSLN) { |
---|
4632 | this.ns = Math.log(ms1/ms2)/Math.log(ts1/ts2); |
---|
4633 | } else { |
---|
4634 | this.ns = sin1; |
---|
4635 | } |
---|
4636 | this.f0 = ms1 / (this.ns * Math.pow(ts1, this.ns)); |
---|
4637 | this.rh = this.a * this.f0 * Math.pow(ts0, this.ns); |
---|
4638 | if (!this.title) this.title = "Lambert Conformal Conic"; |
---|
4639 | }, |
---|
4640 | |
---|
4641 | |
---|
4642 | // Lambert Conformal conic forward equations--mapping lat,long to x,y |
---|
4643 | // ----------------------------------------------------------------- |
---|
4644 | forward : function(p) { |
---|
4645 | |
---|
4646 | var lon = p.x; |
---|
4647 | var lat = p.y; |
---|
4648 | |
---|
4649 | // convert to radians |
---|
4650 | if ( lat <= 90.0 && lat >= -90.0 && lon <= 180.0 && lon >= -180.0) { |
---|
4651 | //lon = lon * Proj4js.common.D2R; |
---|
4652 | //lat = lat * Proj4js.common.D2R; |
---|
4653 | } else { |
---|
4654 | Proj4js.reportError("lcc:forward: llInputOutOfRange: "+ lon +" : " + lat); |
---|
4655 | return null; |
---|
4656 | } |
---|
4657 | |
---|
4658 | var con = Math.abs( Math.abs(lat) - Proj4js.common.HALF_PI); |
---|
4659 | var ts, rh1; |
---|
4660 | if (con > Proj4js.common.EPSLN) { |
---|
4661 | ts = Proj4js.common.tsfnz(this.e, lat, Math.sin(lat) ); |
---|
4662 | rh1 = this.a * this.f0 * Math.pow(ts, this.ns); |
---|
4663 | } else { |
---|
4664 | con = lat * this.ns; |
---|
4665 | if (con <= 0) { |
---|
4666 | Proj4js.reportError("lcc:forward: No Projection"); |
---|
4667 | return null; |
---|
4668 | } |
---|
4669 | rh1 = 0; |
---|
4670 | } |
---|
4671 | var theta = this.ns * Proj4js.common.adjust_lon(lon - this.long0); |
---|
4672 | p.x = this.k0 * (rh1 * Math.sin(theta)) + this.x0; |
---|
4673 | p.y = this.k0 * (this.rh - rh1 * Math.cos(theta)) + this.y0; |
---|
4674 | |
---|
4675 | return p; |
---|
4676 | }, |
---|
4677 | |
---|
4678 | // Lambert Conformal Conic inverse equations--mapping x,y to lat/long |
---|
4679 | // ----------------------------------------------------------------- |
---|
4680 | inverse : function(p) { |
---|
4681 | |
---|
4682 | var rh1, con, ts; |
---|
4683 | var lat, lon; |
---|
4684 | var x = (p.x - this.x0)/this.k0; |
---|
4685 | var y = (this.rh - (p.y - this.y0)/this.k0); |
---|
4686 | if (this.ns > 0) { |
---|
4687 | rh1 = Math.sqrt (x * x + y * y); |
---|
4688 | con = 1.0; |
---|
4689 | } else { |
---|
4690 | rh1 = -Math.sqrt (x * x + y * y); |
---|
4691 | con = -1.0; |
---|
4692 | } |
---|
4693 | var theta = 0.0; |
---|
4694 | if (rh1 != 0) { |
---|
4695 | theta = Math.atan2((con * x),(con * y)); |
---|
4696 | } |
---|
4697 | if ((rh1 != 0) || (this.ns > 0.0)) { |
---|
4698 | con = 1.0/this.ns; |
---|
4699 | ts = Math.pow((rh1/(this.a * this.f0)), con); |
---|
4700 | lat = Proj4js.common.phi2z(this.e, ts); |
---|
4701 | if (lat == -9999) return null; |
---|
4702 | } else { |
---|
4703 | lat = -Proj4js.common.HALF_PI; |
---|
4704 | } |
---|
4705 | lon = Proj4js.common.adjust_lon(theta/this.ns + this.long0); |
---|
4706 | |
---|
4707 | p.x = lon; |
---|
4708 | p.y = lat; |
---|
4709 | return p; |
---|
4710 | } |
---|
4711 | }; |
---|
4712 | |
---|
4713 | |
---|
4714 | |
---|
4715 | |
---|
4716 | /* ====================================================================== |
---|
4717 | projCode/laea.js |
---|
4718 | ====================================================================== */ |
---|
4719 | |
---|
4720 | /******************************************************************************* |
---|
4721 | NAME LAMBERT AZIMUTHAL EQUAL-AREA |
---|
4722 | |
---|
4723 | PURPOSE: Transforms input longitude and latitude to Easting and |
---|
4724 | Northing for the Lambert Azimuthal Equal-Area projection. The |
---|
4725 | longitude and latitude must be in radians. The Easting |
---|
4726 | and Northing values will be returned in meters. |
---|
4727 | |
---|
4728 | PROGRAMMER DATE |
---|
4729 | ---------- ---- |
---|
4730 | D. Steinwand, EROS March, 1991 |
---|
4731 | |
---|
4732 | This function was adapted from the Lambert Azimuthal Equal Area projection |
---|
4733 | code (FORTRAN) in the General Cartographic Transformation Package software |
---|
4734 | which is available from the U.S. Geological Survey National Mapping Division. |
---|
4735 | |
---|
4736 | ALGORITHM REFERENCES |
---|
4737 | |
---|
4738 | 1. "New Equal-Area Map Projections for Noncircular Regions", John P. Snyder, |
---|
4739 | The American Cartographer, Vol 15, No. 4, October 1988, pp. 341-355. |
---|
4740 | |
---|
4741 | 2. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological |
---|
4742 | Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United |
---|
4743 | State Government Printing Office, Washington D.C., 1987. |
---|
4744 | |
---|
4745 | 3. "Software Documentation for GCTP General Cartographic Transformation |
---|
4746 | Package", U.S. Geological Survey National Mapping Division, May 1982. |
---|
4747 | *******************************************************************************/ |
---|
4748 | |
---|
4749 | Proj4js.Proj.laea = { |
---|
4750 | S_POLE: 1, |
---|
4751 | N_POLE: 2, |
---|
4752 | EQUIT: 3, |
---|
4753 | OBLIQ: 4, |
---|
4754 | |
---|
4755 | |
---|
4756 | /* Initialize the Lambert Azimuthal Equal Area projection |
---|
4757 | ------------------------------------------------------*/ |
---|
4758 | init: function() { |
---|
4759 | var t = Math.abs(this.lat0); |
---|
4760 | if (Math.abs(t - Proj4js.common.HALF_PI) < Proj4js.common.EPSLN) { |
---|
4761 | this.mode = this.lat0 < 0. ? this.S_POLE : this.N_POLE; |
---|
4762 | } else if (Math.abs(t) < Proj4js.common.EPSLN) { |
---|
4763 | this.mode = this.EQUIT; |
---|
4764 | } else { |
---|
4765 | this.mode = this.OBLIQ; |
---|
4766 | } |
---|
4767 | if (this.es > 0) { |
---|
4768 | var sinphi; |
---|
4769 | |
---|
4770 | this.qp = Proj4js.common.qsfnz(this.e, 1.0); |
---|
4771 | this.mmf = .5 / (1. - this.es); |
---|
4772 | this.apa = this.authset(this.es); |
---|
4773 | switch (this.mode) { |
---|
4774 | case this.N_POLE: |
---|
4775 | case this.S_POLE: |
---|
4776 | this.dd = 1.; |
---|
4777 | break; |
---|
4778 | case this.EQUIT: |
---|
4779 | this.rq = Math.sqrt(.5 * this.qp); |
---|
4780 | this.dd = 1. / this.rq; |
---|
4781 | this.xmf = 1.; |
---|
4782 | this.ymf = .5 * this.qp; |
---|
4783 | break; |
---|
4784 | case this.OBLIQ: |
---|
4785 | this.rq = Math.sqrt(.5 * this.qp); |
---|
4786 | sinphi = Math.sin(this.lat0); |
---|
4787 | this.sinb1 = Proj4js.common.qsfnz(this.e, sinphi) / this.qp; |
---|
4788 | this.cosb1 = Math.sqrt(1. - this.sinb1 * this.sinb1); |
---|
4789 | this.dd = Math.cos(this.lat0) / (Math.sqrt(1. - this.es * sinphi * sinphi) * this.rq * this.cosb1); |
---|
4790 | this.ymf = (this.xmf = this.rq) / this.dd; |
---|
4791 | this.xmf *= this.dd; |
---|
4792 | break; |
---|
4793 | } |
---|
4794 | } else { |
---|
4795 | if (this.mode == this.OBLIQ) { |
---|
4796 | this.sinph0 = Math.sin(this.lat0); |
---|
4797 | this.cosph0 = Math.cos(this.lat0); |
---|
4798 | } |
---|
4799 | } |
---|
4800 | }, |
---|
4801 | |
---|
4802 | /* Lambert Azimuthal Equal Area forward equations--mapping lat,long to x,y |
---|
4803 | -----------------------------------------------------------------------*/ |
---|
4804 | forward: function(p) { |
---|
4805 | |
---|
4806 | /* Forward equations |
---|
4807 | -----------------*/ |
---|
4808 | var x,y; |
---|
4809 | var lam=p.x; |
---|
4810 | var phi=p.y; |
---|
4811 | lam = Proj4js.common.adjust_lon(lam - this.long0); |
---|
4812 | |
---|
4813 | if (this.sphere) { |
---|
4814 | var coslam, cosphi, sinphi; |
---|
4815 | |
---|
4816 | sinphi = Math.sin(phi); |
---|
4817 | cosphi = Math.cos(phi); |
---|
4818 | coslam = Math.cos(lam); |
---|
4819 | switch (this.mode) { |
---|
4820 | case this.OBLIQ: |
---|
4821 | case this.EQUIT: |
---|
4822 | y = (this.mode == this.EQUIT) ? 1. + cosphi * coslam : 1. + this.sinph0 * sinphi + this.cosph0 * cosphi * coslam; |
---|
4823 | if (y <= Proj4js.common.EPSLN) { |
---|
4824 | Proj4js.reportError("laea:fwd:y less than eps"); |
---|
4825 | return null; |
---|
4826 | } |
---|
4827 | y = Math.sqrt(2. / y); |
---|
4828 | x = y * cosphi * Math.sin(lam); |
---|
4829 | y *= (this.mode == this.EQUIT) ? sinphi : this.cosph0 * sinphi - this.sinph0 * cosphi * coslam; |
---|
4830 | break; |
---|
4831 | case this.N_POLE: |
---|
4832 | coslam = -coslam; |
---|
4833 | case this.S_POLE: |
---|
4834 | if (Math.abs(phi + this.phi0) < Proj4js.common.EPSLN) { |
---|
4835 | Proj4js.reportError("laea:fwd:phi < eps"); |
---|
4836 | return null; |
---|
4837 | } |
---|
4838 | y = Proj4js.common.FORTPI - phi * .5; |
---|
4839 | y = 2. * ((this.mode == this.S_POLE) ? Math.cos(y) : Math.sin(y)); |
---|
4840 | x = y * Math.sin(lam); |
---|
4841 | y *= coslam; |
---|
4842 | break; |
---|
4843 | } |
---|
4844 | } else { |
---|
4845 | var coslam, sinlam, sinphi, q, sinb=0.0, cosb=0.0, b=0.0; |
---|
4846 | |
---|
4847 | coslam = Math.cos(lam); |
---|
4848 | sinlam = Math.sin(lam); |
---|
4849 | sinphi = Math.sin(phi); |
---|
4850 | q = Proj4js.common.qsfnz(this.e, sinphi); |
---|
4851 | if (this.mode == this.OBLIQ || this.mode == this.EQUIT) { |
---|
4852 | sinb = q / this.qp; |
---|
4853 | cosb = Math.sqrt(1. - sinb * sinb); |
---|
4854 | } |
---|
4855 | switch (this.mode) { |
---|
4856 | case this.OBLIQ: |
---|
4857 | b = 1. + this.sinb1 * sinb + this.cosb1 * cosb * coslam; |
---|
4858 | break; |
---|
4859 | case this.EQUIT: |
---|
4860 | b = 1. + cosb * coslam; |
---|
4861 | break; |
---|
4862 | case this.N_POLE: |
---|
4863 | b = Proj4js.common.HALF_PI + phi; |
---|
4864 | q = this.qp - q; |
---|
4865 | break; |
---|
4866 | case this.S_POLE: |
---|
4867 | b = phi - Proj4js.common.HALF_PI; |
---|
4868 | q = this.qp + q; |
---|
4869 | break; |
---|
4870 | } |
---|
4871 | if (Math.abs(b) < Proj4js.common.EPSLN) { |
---|
4872 | Proj4js.reportError("laea:fwd:b < eps"); |
---|
4873 | return null; |
---|
4874 | } |
---|
4875 | switch (this.mode) { |
---|
4876 | case this.OBLIQ: |
---|
4877 | case this.EQUIT: |
---|
4878 | b = Math.sqrt(2. / b); |
---|
4879 | if (this.mode == this.OBLIQ) { |
---|
4880 | y = this.ymf * b * (this.cosb1 * sinb - this.sinb1 * cosb * coslam); |
---|
4881 | } else { |
---|
4882 | y = (b = Math.sqrt(2. / (1. + cosb * coslam))) * sinb * this.ymf; |
---|
4883 | } |
---|
4884 | x = this.xmf * b * cosb * sinlam; |
---|
4885 | break; |
---|
4886 | case this.N_POLE: |
---|
4887 | case this.S_POLE: |
---|
4888 | if (q >= 0.) { |
---|
4889 | x = (b = Math.sqrt(q)) * sinlam; |
---|
4890 | y = coslam * ((this.mode == this.S_POLE) ? b : -b); |
---|
4891 | } else { |
---|
4892 | x = y = 0.; |
---|
4893 | } |
---|
4894 | break; |
---|
4895 | } |
---|
4896 | } |
---|
4897 | |
---|
4898 | //v 1.0 |
---|
4899 | /* |
---|
4900 | var sin_lat=Math.sin(lat); |
---|
4901 | var cos_lat=Math.cos(lat); |
---|
4902 | |
---|
4903 | var sin_delta_lon=Math.sin(delta_lon); |
---|
4904 | var cos_delta_lon=Math.cos(delta_lon); |
---|
4905 | |
---|
4906 | var g =this.sin_lat_o * sin_lat +this.cos_lat_o * cos_lat * cos_delta_lon; |
---|
4907 | if (g == -1.0) { |
---|
4908 | Proj4js.reportError("laea:fwd:Point projects to a circle of radius "+ 2.0 * R); |
---|
4909 | return null; |
---|
4910 | } |
---|
4911 | var ksp = this.a * Math.sqrt(2.0 / (1.0 + g)); |
---|
4912 | var x = ksp * cos_lat * sin_delta_lon + this.x0; |
---|
4913 | var y = ksp * (this.cos_lat_o * sin_lat - this.sin_lat_o * cos_lat * cos_delta_lon) + this.y0; |
---|
4914 | */ |
---|
4915 | p.x = this.a*x + this.x0; |
---|
4916 | p.y = this.a*y + this.y0; |
---|
4917 | return p; |
---|
4918 | },//lamazFwd() |
---|
4919 | |
---|
4920 | /* Inverse equations |
---|
4921 | -----------------*/ |
---|
4922 | inverse: function(p) { |
---|
4923 | p.x -= this.x0; |
---|
4924 | p.y -= this.y0; |
---|
4925 | var x = p.x/this.a; |
---|
4926 | var y = p.y/this.a; |
---|
4927 | |
---|
4928 | if (this.sphere) { |
---|
4929 | var cosz=0.0, rh, sinz=0.0; |
---|
4930 | |
---|
4931 | rh = Math.sqrt(x*x + y*y); |
---|
4932 | var phi = rh * .5; |
---|
4933 | if (phi > 1.) { |
---|
4934 | Proj4js.reportError("laea:Inv:DataError"); |
---|
4935 | return null; |
---|
4936 | } |
---|
4937 | phi = 2. * Math.asin(phi); |
---|
4938 | if (this.mode == this.OBLIQ || this.mode == this.EQUIT) { |
---|
4939 | sinz = Math.sin(phi); |
---|
4940 | cosz = Math.cos(phi); |
---|
4941 | } |
---|
4942 | switch (this.mode) { |
---|
4943 | case this.EQUIT: |
---|
4944 | phi = (Math.abs(rh) <= Proj4js.common.EPSLN) ? 0. : Math.asin(y * sinz / rh); |
---|
4945 | x *= sinz; |
---|
4946 | y = cosz * rh; |
---|
4947 | break; |
---|
4948 | case this.OBLIQ: |
---|
4949 | phi = (Math.abs(rh) <= Proj4js.common.EPSLN) ? this.phi0 : Math.asin(cosz * sinph0 + y * sinz * cosph0 / rh); |
---|
4950 | x *= sinz * cosph0; |
---|
4951 | y = (cosz - Math.sin(phi) * sinph0) * rh; |
---|
4952 | break; |
---|
4953 | case this.N_POLE: |
---|
4954 | y = -y; |
---|
4955 | phi = Proj4js.common.HALF_PI - phi; |
---|
4956 | break; |
---|
4957 | case this.S_POLE: |
---|
4958 | phi -= Proj4js.common.HALF_PI; |
---|
4959 | break; |
---|
4960 | } |
---|
4961 | lam = (y == 0. && (this.mode == this.EQUIT || this.mode == this.OBLIQ)) ? 0. : Math.atan2(x, y); |
---|
4962 | } else { |
---|
4963 | var cCe, sCe, q, rho, ab=0.0; |
---|
4964 | |
---|
4965 | switch (this.mode) { |
---|
4966 | case this.EQUIT: |
---|
4967 | case this.OBLIQ: |
---|
4968 | x /= this.dd; |
---|
4969 | y *= this.dd; |
---|
4970 | rho = Math.sqrt(x*x + y*y); |
---|
4971 | if (rho < Proj4js.common.EPSLN) { |
---|
4972 | p.x = 0.; |
---|
4973 | p.y = this.phi0; |
---|
4974 | return p; |
---|
4975 | } |
---|
4976 | sCe = 2. * Math.asin(.5 * rho / this.rq); |
---|
4977 | cCe = Math.cos(sCe); |
---|
4978 | x *= (sCe = Math.sin(sCe)); |
---|
4979 | if (this.mode == this.OBLIQ) { |
---|
4980 | ab = cCe * this.sinb1 + y * sCe * this.cosb1 / rho |
---|
4981 | q = this.qp * ab; |
---|
4982 | y = rho * this.cosb1 * cCe - y * this.sinb1 * sCe; |
---|
4983 | } else { |
---|
4984 | ab = y * sCe / rho; |
---|
4985 | q = this.qp * ab; |
---|
4986 | y = rho * cCe; |
---|
4987 | } |
---|
4988 | break; |
---|
4989 | case this.N_POLE: |
---|
4990 | y = -y; |
---|
4991 | case this.S_POLE: |
---|
4992 | q = (x * x + y * y); |
---|
4993 | if (!q ) { |
---|
4994 | p.x = 0.; |
---|
4995 | p.y = this.phi0; |
---|
4996 | return p; |
---|
4997 | } |
---|
4998 | /* |
---|
4999 | q = this.qp - q; |
---|
5000 | */ |
---|
5001 | ab = 1. - q / this.qp; |
---|
5002 | if (this.mode == this.S_POLE) { |
---|
5003 | ab = - ab; |
---|
5004 | } |
---|
5005 | break; |
---|
5006 | } |
---|
5007 | lam = Math.atan2(x, y); |
---|
5008 | phi = this.authlat(Math.asin(ab), this.apa); |
---|
5009 | } |
---|
5010 | |
---|
5011 | /* |
---|
5012 | var Rh = Math.Math.sqrt(p.x *p.x +p.y * p.y); |
---|
5013 | var temp = Rh / (2.0 * this.a); |
---|
5014 | |
---|
5015 | if (temp > 1) { |
---|
5016 | Proj4js.reportError("laea:Inv:DataError"); |
---|
5017 | return null; |
---|
5018 | } |
---|
5019 | |
---|
5020 | var z = 2.0 * Proj4js.common.asinz(temp); |
---|
5021 | var sin_z=Math.sin(z); |
---|
5022 | var cos_z=Math.cos(z); |
---|
5023 | |
---|
5024 | var lon =this.long0; |
---|
5025 | if (Math.abs(Rh) > Proj4js.common.EPSLN) { |
---|
5026 | var lat = Proj4js.common.asinz(this.sin_lat_o * cos_z +this. cos_lat_o * sin_z *p.y / Rh); |
---|
5027 | var temp =Math.abs(this.lat0) - Proj4js.common.HALF_PI; |
---|
5028 | if (Math.abs(temp) > Proj4js.common.EPSLN) { |
---|
5029 | temp = cos_z -this.sin_lat_o * Math.sin(lat); |
---|
5030 | if(temp!=0.0) lon=Proj4js.common.adjust_lon(this.long0+Math.atan2(p.x*sin_z*this.cos_lat_o,temp*Rh)); |
---|
5031 | } else if (this.lat0 < 0.0) { |
---|
5032 | lon = Proj4js.common.adjust_lon(this.long0 - Math.atan2(-p.x,p.y)); |
---|
5033 | } else { |
---|
5034 | lon = Proj4js.common.adjust_lon(this.long0 + Math.atan2(p.x, -p.y)); |
---|
5035 | } |
---|
5036 | } else { |
---|
5037 | lat = this.lat0; |
---|
5038 | } |
---|
5039 | */ |
---|
5040 | //return(OK); |
---|
5041 | p.x = Proj4js.common.adjust_lon(this.long0+lam); |
---|
5042 | p.y = phi; |
---|
5043 | return p; |
---|
5044 | },//lamazInv() |
---|
5045 | |
---|
5046 | /* determine latitude from authalic latitude */ |
---|
5047 | P00: .33333333333333333333, |
---|
5048 | P01: .17222222222222222222, |
---|
5049 | P02: .10257936507936507936, |
---|
5050 | P10: .06388888888888888888, |
---|
5051 | P11: .06640211640211640211, |
---|
5052 | P20: .01641501294219154443, |
---|
5053 | |
---|
5054 | authset: function(es) { |
---|
5055 | var t; |
---|
5056 | var APA = new Array(); |
---|
5057 | APA[0] = es * this.P00; |
---|
5058 | t = es * es; |
---|
5059 | APA[0] += t * this.P01; |
---|
5060 | APA[1] = t * this.P10; |
---|
5061 | t *= es; |
---|
5062 | APA[0] += t * this.P02; |
---|
5063 | APA[1] += t * this.P11; |
---|
5064 | APA[2] = t * this.P20; |
---|
5065 | return APA; |
---|
5066 | }, |
---|
5067 | |
---|
5068 | authlat: function(beta, APA) { |
---|
5069 | var t = beta+beta; |
---|
5070 | return(beta + APA[0] * Math.sin(t) + APA[1] * Math.sin(t+t) + APA[2] * Math.sin(t+t+t)); |
---|
5071 | } |
---|
5072 | |
---|
5073 | }; |
---|
5074 | |
---|
5075 | |
---|
5076 | |
---|
5077 | /* ====================================================================== |
---|
5078 | projCode/aeqd.js |
---|
5079 | ====================================================================== */ |
---|
5080 | |
---|
5081 | Proj4js.Proj.aeqd = { |
---|
5082 | |
---|
5083 | init : function() { |
---|
5084 | this.sin_p12=Math.sin(this.lat0); |
---|
5085 | this.cos_p12=Math.cos(this.lat0); |
---|
5086 | }, |
---|
5087 | |
---|
5088 | forward: function(p) { |
---|
5089 | var lon=p.x; |
---|
5090 | var lat=p.y; |
---|
5091 | var ksp; |
---|
5092 | |
---|
5093 | var sinphi=Math.sin(p.y); |
---|
5094 | var cosphi=Math.cos(p.y); |
---|
5095 | var dlon = Proj4js.common.adjust_lon(lon - this.long0); |
---|
5096 | var coslon = Math.cos(dlon); |
---|
5097 | var g = this.sin_p12 * sinphi + this.cos_p12 * cosphi * coslon; |
---|
5098 | if (Math.abs(Math.abs(g) - 1.0) < Proj4js.common.EPSLN) { |
---|
5099 | ksp = 1.0; |
---|
5100 | if (g < 0.0) { |
---|
5101 | Proj4js.reportError("aeqd:Fwd:PointError"); |
---|
5102 | return; |
---|
5103 | } |
---|
5104 | } else { |
---|
5105 | var z = Math.acos(g); |
---|
5106 | ksp = z/Math.sin(z); |
---|
5107 | } |
---|
5108 | p.x = this.x0 + this.a * ksp * cosphi * Math.sin(dlon); |
---|
5109 | p.y = this.y0 + this.a * ksp * (this.cos_p12 * sinphi - this.sin_p12 * cosphi * coslon); |
---|
5110 | return p; |
---|
5111 | }, |
---|
5112 | |
---|
5113 | inverse: function(p){ |
---|
5114 | p.x -= this.x0; |
---|
5115 | p.y -= this.y0; |
---|
5116 | |
---|
5117 | var rh = Math.sqrt(p.x * p.x + p.y *p.y); |
---|
5118 | if (rh > (2.0 * Proj4js.common.HALF_PI * this.a)) { |
---|
5119 | Proj4js.reportError("aeqdInvDataError"); |
---|
5120 | return; |
---|
5121 | } |
---|
5122 | var z = rh / this.a; |
---|
5123 | |
---|
5124 | var sinz=Math.sin(z); |
---|
5125 | var cosz=Math.cos(z); |
---|
5126 | |
---|
5127 | var lon = this.long0; |
---|
5128 | var lat; |
---|
5129 | if (Math.abs(rh) <= Proj4js.common.EPSLN) { |
---|
5130 | lat = this.lat0; |
---|
5131 | } else { |
---|
5132 | lat = Proj4js.common.asinz(cosz * this.sin_p12 + (p.y * sinz * this.cos_p12) / rh); |
---|
5133 | var con = Math.abs(this.lat0) - Proj4js.common.HALF_PI; |
---|
5134 | if (Math.abs(con) <= Proj4js.common.EPSLN) { |
---|
5135 | if (lat0 >= 0.0) { |
---|
5136 | lon = Proj4js.common.adjust_lon(this.long0 + Math.atan2(p.x , -p.y)); |
---|
5137 | } else { |
---|
5138 | lon = Proj4js.common.adjust_lon(this.long0 - Math.atan2(-p.x , p.y)); |
---|
5139 | } |
---|
5140 | } else { |
---|
5141 | con = cosz - this.sin_p12 * Math.sin(lat); |
---|
5142 | if ((Math.abs(con) < Proj4js.common.EPSLN) && (Math.abs(p.x) < Proj4js.common.EPSLN)) { |
---|
5143 | //no-op, just keep the lon value as is |
---|
5144 | } else { |
---|
5145 | var temp = Math.atan2((p.x * sinz * this.cos_p12), (con * rh)); |
---|
5146 | lon = Proj4js.common.adjust_lon(this.long0 + Math.atan2((p.x * sinz * this.cos_p12), (con * rh))); |
---|
5147 | } |
---|
5148 | } |
---|
5149 | } |
---|
5150 | |
---|
5151 | p.x = lon; |
---|
5152 | p.y = lat; |
---|
5153 | return p; |
---|
5154 | } |
---|
5155 | }; |
---|
5156 | /* ====================================================================== |
---|
5157 | projCode/moll.js |
---|
5158 | ====================================================================== */ |
---|
5159 | |
---|
5160 | /******************************************************************************* |
---|
5161 | NAME MOLLWEIDE |
---|
5162 | |
---|
5163 | PURPOSE: Transforms input longitude and latitude to Easting and |
---|
5164 | Northing for the MOllweide projection. The |
---|
5165 | longitude and latitude must be in radians. The Easting |
---|
5166 | and Northing values will be returned in meters. |
---|
5167 | |
---|
5168 | PROGRAMMER DATE |
---|
5169 | ---------- ---- |
---|
5170 | D. Steinwand, EROS May, 1991; Updated Sept, 1992; Updated Feb, 1993 |
---|
5171 | S. Nelson, EDC Jun, 2993; Made corrections in precision and |
---|
5172 | number of iterations. |
---|
5173 | |
---|
5174 | ALGORITHM REFERENCES |
---|
5175 | |
---|
5176 | 1. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections", |
---|
5177 | U.S. Geological Survey Professional Paper 1453 , United State Government |
---|
5178 | Printing Office, Washington D.C., 1989. |
---|
5179 | |
---|
5180 | 2. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological |
---|
5181 | Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United |
---|
5182 | State Government Printing Office, Washington D.C., 1987. |
---|
5183 | *******************************************************************************/ |
---|
5184 | |
---|
5185 | Proj4js.Proj.moll = { |
---|
5186 | |
---|
5187 | /* Initialize the Mollweide projection |
---|
5188 | ------------------------------------*/ |
---|
5189 | init: function(){ |
---|
5190 | //no-op |
---|
5191 | }, |
---|
5192 | |
---|
5193 | /* Mollweide forward equations--mapping lat,long to x,y |
---|
5194 | ----------------------------------------------------*/ |
---|
5195 | forward: function(p) { |
---|
5196 | |
---|
5197 | /* Forward equations |
---|
5198 | -----------------*/ |
---|
5199 | var lon=p.x; |
---|
5200 | var lat=p.y; |
---|
5201 | |
---|
5202 | var delta_lon = Proj4js.common.adjust_lon(lon - this.long0); |
---|
5203 | var theta = lat; |
---|
5204 | var con = Proj4js.common.PI * Math.sin(lat); |
---|
5205 | |
---|
5206 | /* Iterate using the Newton-Raphson method to find theta |
---|
5207 | -----------------------------------------------------*/ |
---|
5208 | for (var i=0;true;i++) { |
---|
5209 | var delta_theta = -(theta + Math.sin(theta) - con)/ (1.0 + Math.cos(theta)); |
---|
5210 | theta += delta_theta; |
---|
5211 | if (Math.abs(delta_theta) < Proj4js.common.EPSLN) break; |
---|
5212 | if (i >= 50) { |
---|
5213 | Proj4js.reportError("moll:Fwd:IterationError"); |
---|
5214 | //return(241); |
---|
5215 | } |
---|
5216 | } |
---|
5217 | theta /= 2.0; |
---|
5218 | |
---|
5219 | /* If the latitude is 90 deg, force the x coordinate to be "0 + false easting" |
---|
5220 | this is done here because of precision problems with "cos(theta)" |
---|
5221 | --------------------------------------------------------------------------*/ |
---|
5222 | if (Proj4js.common.PI/2 - Math.abs(lat) < Proj4js.common.EPSLN) delta_lon =0; |
---|
5223 | var x = 0.900316316158 * this.a * delta_lon * Math.cos(theta) + this.x0; |
---|
5224 | var y = 1.4142135623731 * this.a * Math.sin(theta) + this.y0; |
---|
5225 | |
---|
5226 | p.x=x; |
---|
5227 | p.y=y; |
---|
5228 | return p; |
---|
5229 | }, |
---|
5230 | |
---|
5231 | inverse: function(p){ |
---|
5232 | var theta; |
---|
5233 | var arg; |
---|
5234 | |
---|
5235 | /* Inverse equations |
---|
5236 | -----------------*/ |
---|
5237 | p.x-= this.x0; |
---|
5238 | //~ p.y -= this.y0; |
---|
5239 | var arg = p.y / (1.4142135623731 * this.a); |
---|
5240 | |
---|
5241 | /* Because of division by zero problems, 'arg' can not be 1.0. Therefore |
---|
5242 | a number very close to one is used instead. |
---|
5243 | -------------------------------------------------------------------*/ |
---|
5244 | if(Math.abs(arg) > 0.999999999999) arg=0.999999999999; |
---|
5245 | var theta =Math.asin(arg); |
---|
5246 | var lon = Proj4js.common.adjust_lon(this.long0 + (p.x / (0.900316316158 * this.a * Math.cos(theta)))); |
---|
5247 | if(lon < (-Proj4js.common.PI)) lon= -Proj4js.common.PI; |
---|
5248 | if(lon > Proj4js.common.PI) lon= Proj4js.common.PI; |
---|
5249 | arg = (2.0 * theta + Math.sin(2.0 * theta)) / Proj4js.common.PI; |
---|
5250 | if(Math.abs(arg) > 1.0)arg=1.0; |
---|
5251 | var lat = Math.asin(arg); |
---|
5252 | //return(OK); |
---|
5253 | |
---|
5254 | p.x=lon; |
---|
5255 | p.y=lat; |
---|
5256 | return p; |
---|
5257 | } |
---|
5258 | }; |
---|
5259 | |
---|