|
1
|
/*
|
|
2
|
proj4js.js -- Javascript reprojection library.
|
|
3
|
|
|
4
|
Author: Mike Adair madairATdmsolutions.ca
|
|
5
|
Richard Greenwood rich@greenwoodmap.com
|
|
6
|
License: LGPL as per: http://www.gnu.org/copyleft/lesser.html
|
|
7
|
Note: This program is an almost direct port of the C library
|
|
8
|
Proj4.
|
|
9
|
*/
|
|
10
|
/* ======================================================================
|
|
11
|
proj4js.js
|
|
12
|
====================================================================== */
|
|
13
|
|
|
14
|
/*
|
|
15
|
Author: Mike Adair madairATdmsolutions.ca
|
|
16
|
Richard Greenwood rich@greenwoodmap.com
|
|
17
|
License: LGPL as per: http://www.gnu.org/copyleft/lesser.html
|
|
18
|
Note: This program is an almost direct port of the C library Proj4.
|
|
19
|
$Id: Proj.js 2956 2007-07-09 12:17:52Z steven $
|
|
20
|
*/
|
|
21
|
|
|
22
|
/**
|
|
23
|
* Provides methods for coordinate transformations between map projections and
|
|
24
|
* longitude/latitude, including datum transformations.
|
|
25
|
*
|
|
26
|
* Initialization of Proj objects is with a projection code, usually EPSG codes.
|
|
27
|
* The code passed in will be stripped of colons (':') and converted to uppercase
|
|
28
|
* for internal use.
|
|
29
|
* If you know what map projections your application will be dealing with, the
|
|
30
|
* definition for the projections can be included with the script tag when the
|
|
31
|
* application is being coded. Otherwise, practically any projection definition
|
|
32
|
* can be loaded dynamically at run-time with an AJAX request to a lookup service
|
|
33
|
* such as spatialreference.org.
|
|
34
|
* The actual code supporting the forward and inverse tansformations for each
|
|
35
|
* projection class is loaded dynamically at run-time. These may also be
|
|
36
|
* specified when the application is coded if the projections to be used are known
|
|
37
|
* beforehand.
|
|
38
|
* A projection object has properties for units and title strings.
|
|
39
|
* All coordinates are handled as points which is a 2 element array where x is
|
|
40
|
* the first element and y is the second.
|
|
41
|
* For the transform() method pass in mapXY and a destination projection object
|
|
42
|
* and it returns a map XY coordinate in the other projection
|
|
43
|
*/
|
|
44
|
|
|
45
|
/**
|
|
46
|
* Global namespace object for Proj4js library to use
|
|
47
|
*/
|
|
48
|
Proj4js = {
|
|
49
|
|
|
50
|
/**
|
|
51
|
* Property: defaultDatum
|
|
52
|
* The datum to use when no others a specified
|
|
53
|
*/
|
|
54
|
defaultDatum: 'WGS84', //default datum
|
|
55
|
|
|
56
|
/**
|
|
57
|
* Property: proxyScript
|
|
58
|
* A proxy script to execute AJAX requests in other domains.
|
|
59
|
*/
|
|
60
|
proxyScript: null, //TBD: customize this for spatialreference.org output
|
|
61
|
|
|
62
|
/**
|
|
63
|
* Property: defsLookupService
|
|
64
|
* AJAX service to retreive projection definition parameters from
|
|
65
|
*/
|
|
66
|
defsLookupService: 'http://spatialreference.org/ref',
|
|
67
|
|
|
68
|
/**
|
|
69
|
* Property: libPath
|
|
70
|
* internal: http server path to library code.
|
|
71
|
* TBD figure this out automatically
|
|
72
|
*/
|
|
73
|
libPath: '../lib/',
|
|
74
|
|
|
75
|
/**
|
|
76
|
* Method: transform(source, dest, point)
|
|
77
|
* Transform a point coordinate from one map projection to another.
|
|
78
|
*
|
|
79
|
* Parameters:
|
|
80
|
* source - {Proj4js.Proj} source map projection for the transformation
|
|
81
|
* dest - {Proj4js.Proj} destination map projection for the transformation
|
|
82
|
* point - {Object} point to transform, may be geodetic (long, lat) or
|
|
83
|
* projected Cartesian (x,y), but should always have x,y properties.
|
|
84
|
*/
|
|
85
|
transform : function(source, dest, point) {
|
|
86
|
if (!source.readyToUse || !dest.readyToUse) {
|
|
87
|
this.reportError("Proj4js initialization for "+source.srsCode+" not yet complete");
|
|
88
|
return;
|
|
89
|
}
|
|
90
|
|
|
91
|
if (point.transformed) {
|
|
92
|
this.log("point already transformed");
|
|
93
|
return;
|
|
94
|
}
|
|
95
|
|
|
96
|
// Workaround for Spherical Mercator
|
|
97
|
if ((source.srsProjNumber =="900913" && dest.datumCode != "WGS84") ||
|
|
98
|
(dest.srsProjNumber == "900913" && source.datumCode != "WGS84")) {
|
|
99
|
var wgs84 = Proj4js.WGS84;
|
|
100
|
this.transform(source, wgs84, point);
|
|
101
|
point.transformed = false;
|
|
102
|
source = wgs84;
|
|
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
|
point.transformed = true;
|
|
142
|
return point;
|
|
143
|
}, // transform()
|
|
144
|
|
|
145
|
/** datum_transform()
|
|
146
|
source coordinate system definition,
|
|
147
|
destination coordinate system definition,
|
|
148
|
point to transform in geodetic coordinates (long, lat, height)
|
|
149
|
*/
|
|
150
|
datum_transform : function( source, dest, point ) {
|
|
151
|
|
|
152
|
// Short cut if the datums are identical.
|
|
153
|
if( source.compare_datums( dest ) ) {
|
|
154
|
return point; // in this case, zero is sucess,
|
|
155
|
// whereas cs_compare_datums returns 1 to indicate TRUE
|
|
156
|
// confusing, should fix this
|
|
157
|
}
|
|
158
|
|
|
159
|
// Explicitly skip datum transform by setting 'datum=none' as parameter for either source or dest
|
|
160
|
if( source.datum_type == Proj4js.common.PJD_NODATUM
|
|
161
|
|| dest.datum_type == Proj4js.common.PJD_NODATUM) {
|
|
162
|
return point;
|
|
163
|
}
|
|
164
|
|
|
165
|
// If this datum requires grid shifts, then apply it to geodetic coordinates.
|
|
166
|
if( source.datum_type == Proj4js.common.PJD_GRIDSHIFT )
|
|
167
|
{
|
|
168
|
alert("ERROR: Grid shift transformations are not implemented yet.");
|
|
169
|
/*
|
|
170
|
pj_apply_gridshift( pj_param(source.params,"snadgrids").s, 0,
|
|
171
|
point_count, point_offset, x, y, z );
|
|
172
|
CHECK_RETURN;
|
|
173
|
|
|
174
|
src_a = SRS_WGS84_SEMIMAJOR;
|
|
175
|
src_es = 0.006694379990;
|
|
176
|
*/
|
|
177
|
}
|
|
178
|
|
|
179
|
if( dest.datum_type == Proj4js.common.PJD_GRIDSHIFT )
|
|
180
|
{
|
|
181
|
alert("ERROR: Grid shift transformations are not implemented yet.");
|
|
182
|
/*
|
|
183
|
dst_a = ;
|
|
184
|
dst_es = 0.006694379990;
|
|
185
|
*/
|
|
186
|
}
|
|
187
|
|
|
188
|
// Do we need to go through geocentric coordinates?
|
|
189
|
if( source.es != dest.es || source.a != dest.a
|
|
190
|
|| source.datum_type == Proj4js.common.PJD_3PARAM
|
|
191
|
|| source.datum_type == Proj4js.common.PJD_7PARAM
|
|
192
|
|| dest.datum_type == Proj4js.common.PJD_3PARAM
|
|
193
|
|| dest.datum_type == Proj4js.common.PJD_7PARAM)
|
|
194
|
{
|
|
195
|
|
|
196
|
// Convert to geocentric coordinates.
|
|
197
|
source.geodetic_to_geocentric( point );
|
|
198
|
// CHECK_RETURN;
|
|
199
|
|
|
200
|
// Convert between datums
|
|
201
|
if( source.datum_type == Proj4js.common.PJD_3PARAM || source.datum_type == Proj4js.common.PJD_7PARAM ) {
|
|
202
|
source.geocentric_to_wgs84(point);
|
|
203
|
// CHECK_RETURN;
|
|
204
|
}
|
|
205
|
|
|
206
|
if( dest.datum_type == Proj4js.common.PJD_3PARAM || dest.datum_type == Proj4js.common.PJD_7PARAM ) {
|
|
207
|
dest.geocentric_from_wgs84(point);
|
|
208
|
// CHECK_RETURN;
|
|
209
|
}
|
|
210
|
|
|
211
|
// Convert back to geodetic coordinates
|
|
212
|
dest.geocentric_to_geodetic( point );
|
|
213
|
// CHECK_RETURN;
|
|
214
|
}
|
|
215
|
|
|
216
|
// Apply grid shift to destination if required
|
|
217
|
if( dest.datum_type == Proj4js.common.PJD_GRIDSHIFT )
|
|
218
|
{
|
|
219
|
alert("ERROR: Grid shift transformations are not implemented yet.");
|
|
220
|
// pj_apply_gridshift( pj_param(dest.params,"snadgrids").s, 1, point);
|
|
221
|
// CHECK_RETURN;
|
|
222
|
}
|
|
223
|
return point;
|
|
224
|
}, // cs_datum_transform
|
|
225
|
|
|
226
|
/**
|
|
227
|
* Function: reportError
|
|
228
|
* An internal method to report errors back to user. Should be overridden
|
|
229
|
* by applications to deliver error messages.
|
|
230
|
*/
|
|
231
|
reportError: function(msg) {
|
|
232
|
},
|
|
233
|
|
|
234
|
/**
|
|
235
|
* Function: log
|
|
236
|
* An internal method to log events.
|
|
237
|
*/
|
|
238
|
log: function(msg) {
|
|
239
|
},
|
|
240
|
|
|
241
|
loadProjDefinition : function(proj) {
|
|
242
|
|
|
243
|
//check in memory
|
|
244
|
if (this.defs[proj.srsCode]) return this.defs[proj.srsCode];
|
|
245
|
|
|
246
|
//set AJAX options
|
|
247
|
var options = {
|
|
248
|
method: 'get',
|
|
249
|
asynchronous: false, //need to wait until defs are loaded before proceeding
|
|
250
|
onSuccess: this.defsLoadedFromDisk.bind(this,proj.srsCode)
|
|
251
|
}
|
|
252
|
|
|
253
|
//else check for def on the server
|
|
254
|
var url = this.libPath + 'defs/' + proj.srsAuth.toUpperCase() + proj.srsProjNumber + '.js';
|
|
255
|
new OpenLayers.Ajax.Request(url, options);
|
|
256
|
if ( this.defs[proj.srsCode] ) return this.defs[proj.srsCode];
|
|
257
|
|
|
258
|
//else load from web service via AJAX request
|
|
259
|
if (this.proxyScript) {
|
|
260
|
var url = this.proxyScript + this.defsLookupService +'/' + proj.srsAuth +'/'+ proj.srsProjNumber + '/proj4';
|
|
261
|
options.onSuccess = this.defsLoadedFromService.bind(this,proj.srsCode)
|
|
262
|
options.onFailure = this.defsFailed.bind(this,proj.srsCode);
|
|
263
|
new OpenLayers.Ajax.Request(url, options);
|
|
264
|
}
|
|
265
|
|
|
266
|
//may return null here if the defs are not found
|
|
267
|
return this.defs[proj.srsCode];
|
|
268
|
},
|
|
269
|
|
|
270
|
defsLoadedFromDisk: function(srsCode, transport) {
|
|
271
|
eval(transport.responseText);
|
|
272
|
},
|
|
273
|
|
|
274
|
defsLoadedFromService: function(srsCode, transport) {
|
|
275
|
this.defs[srsCode] = transport.responseText;
|
|
276
|
// save this also in the prototype, so we don't need to fetch it again
|
|
277
|
Proj4js.defs[srsCode] = transport.responseText;
|
|
278
|
},
|
|
279
|
|
|
280
|
defsFailed: function(srsCode) {
|
|
281
|
this.reportError('failed to load projection definition for: '+srsCode);
|
|
282
|
OpenLayers.Util.extend(this.defs[srsCode], this.defs['WGS84']); //set it to something so it can at least continue
|
|
283
|
},
|
|
284
|
|
|
285
|
loadProjCode : function(projName) {
|
|
286
|
if (this.Proj[projName]) return;
|
|
287
|
|
|
288
|
//set AJAX options
|
|
289
|
var options = {
|
|
290
|
method: 'get',
|
|
291
|
asynchronous: false, //need to wait until defs are loaded before proceeding
|
|
292
|
onSuccess: this.loadProjCodeSuccess.bind(this, projName),
|
|
293
|
onFailure: this.loadProjCodeFailure.bind(this, projName)
|
|
294
|
};
|
|
295
|
|
|
296
|
//load the projection class
|
|
297
|
var url = this.libPath + 'projCode/' + projName + '.js';
|
|
298
|
new OpenLayers.Ajax.Request(url, options);
|
|
299
|
},
|
|
300
|
|
|
301
|
loadProjCodeSuccess : function(projName, transport) {
|
|
302
|
eval(transport.responseText);
|
|
303
|
if (this.Proj[projName].dependsOn){
|
|
304
|
this.loadProjCode(this.Proj[projName].dependsOn);
|
|
305
|
}
|
|
306
|
},
|
|
307
|
|
|
308
|
loadProjCodeFailure : function(projName) {
|
|
309
|
Proj4js.reportError("failed to find projection file for: " + projName);
|
|
310
|
//TBD initialize with identity transforms so proj will still work
|
|
311
|
}
|
|
312
|
|
|
313
|
};
|
|
314
|
|
|
315
|
/**
|
|
316
|
* Class: Proj4js.Proj
|
|
317
|
* Projection objects provide coordinate transformation methods for point coordinates
|
|
318
|
* once they have been initialized with a projection code.
|
|
319
|
*/
|
|
320
|
Proj4js.Proj = OpenLayers.Class({
|
|
321
|
|
|
322
|
/**
|
|
323
|
* Property: readyToUse
|
|
324
|
* Flag to indicate if initialization is complete for this Proj object
|
|
325
|
*/
|
|
326
|
readyToUse : false,
|
|
327
|
|
|
328
|
/**
|
|
329
|
* Property: title
|
|
330
|
* The title to describe the projection
|
|
331
|
*/
|
|
332
|
title: null,
|
|
333
|
|
|
334
|
/**
|
|
335
|
* Property: projName
|
|
336
|
* The projection class for this projection, e.g. lcc (lambert conformal conic,
|
|
337
|
* or merc for mercator. These are exactly equicvalent to their Proj4
|
|
338
|
* counterparts.
|
|
339
|
*/
|
|
340
|
projName: null,
|
|
341
|
/**
|
|
342
|
* Property: units
|
|
343
|
* The units of the projection. Values include 'm' and 'degrees'
|
|
344
|
*/
|
|
345
|
units: null,
|
|
346
|
/**
|
|
347
|
* Property: datum
|
|
348
|
* The datum specified for the projection
|
|
349
|
*/
|
|
350
|
datum: null,
|
|
351
|
|
|
352
|
/**
|
|
353
|
* Constructor: initialize
|
|
354
|
* Constructor for Proj4js.Proj objects
|
|
355
|
*
|
|
356
|
* Parameters:
|
|
357
|
* srsCode - a code for map projection definition parameters. These are usually
|
|
358
|
* (but not always) EPSG codes.
|
|
359
|
*/
|
|
360
|
initialize : function(srsCode) {
|
|
361
|
this.srsCode = srsCode.toUpperCase();
|
|
362
|
if (this.srsCode.indexOf("EPSG") == 0) {
|
|
363
|
this.srsCode = this.srsCode;
|
|
364
|
this.srsAuth = 'epsg';
|
|
365
|
this.srsProjNumber = this.srsCode.substring(5);
|
|
366
|
} else {
|
|
367
|
this.srsAuth = '';
|
|
368
|
this.srsProjNumber = this.srsCode;
|
|
369
|
}
|
|
370
|
|
|
371
|
var defs = Proj4js.loadProjDefinition(this);
|
|
372
|
if (defs) {
|
|
373
|
this.parseDefs(defs);
|
|
374
|
Proj4js.loadProjCode(this.projName);
|
|
375
|
this.callInit();
|
|
376
|
}
|
|
377
|
|
|
378
|
},
|
|
379
|
|
|
380
|
callInit : function() {
|
|
381
|
Proj4js.log('projection script loaded for:' + this.projName);
|
|
382
|
OpenLayers.Util.extend(this, Proj4js.Proj[this.projName]);
|
|
383
|
this.init();
|
|
384
|
this.mapXYToLonLat = this.inverse;
|
|
385
|
this.lonLatToMapXY = this.forward;
|
|
386
|
this.readyToUse = true;
|
|
387
|
},
|
|
388
|
|
|
389
|
parseDefs : function(proj4opts) {
|
|
390
|
this.defData = proj4opts;
|
|
391
|
var paramName, paramVal;
|
|
392
|
var paramArray=this.defData.split("+");
|
|
393
|
|
|
394
|
for (var prop=0; prop<paramArray.length; prop++) {
|
|
395
|
var property = paramArray[prop].split("=");
|
|
396
|
paramName = property[0].toLowerCase();
|
|
397
|
paramVal = property[1];
|
|
398
|
|
|
399
|
switch (paramName.replace(/\s/gi,"")) { // trim out spaces
|
|
400
|
case "": break; // throw away nameless parameter
|
|
401
|
case "title": this.title = paramVal; break;
|
|
402
|
case "proj": this.projName = paramVal.replace(/\s/gi,""); break;
|
|
403
|
case "units": this.units = paramVal.replace(/\s/gi,""); break;
|
|
404
|
case "datum": this.datumCode = paramVal.replace(/\s/gi,""); break;
|
|
405
|
case "nadgrids": this.nagrids = paramVal.replace(/\s/gi,""); break;
|
|
406
|
case "ellps": this.ellps = paramVal.replace(/\s/gi,""); break;
|
|
407
|
case "a": this.a = parseFloat(paramVal); break; // semi-major radius
|
|
408
|
case "b": this.b = parseFloat(paramVal); break; // semi-minor radius
|
|
409
|
case "lat_0": this.lat0 = paramVal*Proj4js.common.D2R; break; // phi0, central latitude
|
|
410
|
case "lat_1": this.lat1 = paramVal*Proj4js.common.D2R; break; //standard parallel 1
|
|
411
|
case "lat_2": this.lat2 = paramVal*Proj4js.common.D2R; break; //standard parallel 2
|
|
412
|
case "lat_ts": this.lat_ts = paramVal*Proj4js.common.D2R; break; //used in merc
|
|
413
|
case "lon_0": this.long0 = paramVal*Proj4js.common.D2R; break; // lam0, central longitude
|
|
414
|
case "x_0": this.x0 = parseFloat(paramVal); break; // false easting
|
|
415
|
case "y_0": this.y0 = parseFloat(paramVal); break; // false northing
|
|
416
|
case "k_0": this.k0 = parseFloat(paramVal); break; // projection scale factor
|
|
417
|
case "k": this.k0 = parseFloat(paramVal); break; // both forms returned
|
|
418
|
case "R_A": this.R = parseFloat(paramVal); break; //Spheroid radius
|
|
419
|
case "zone": this.zone = parseInt(paramVal); break; // UTM Zone
|
|
420
|
case "south": this.utmSouth = true; break; // UTM north/south
|
|
421
|
case "towgs84":this.datum_params = paramVal.split(","); break;
|
|
422
|
case "to_meter": this.to_meter = parseFloat(paramVal); break; // cartesian scaling
|
|
423
|
case "from_greenwich": this.from_greenwich = paramVal*Proj4js.common.D2R; break;
|
|
424
|
case "pm": paramVal = paramVal.replace(/\s/gi,"");
|
|
425
|
this.from_greenwich = Proj4js.PrimeMeridian[paramVal] ?
|
|
426
|
Proj4js.PrimeMeridian[paramVal]*Proj4js.common.D2R : 0.0; break;
|
|
427
|
case "no_defs": break;
|
|
428
|
default: Proj4js.log("Unrecognized parameter: " + paramName);
|
|
429
|
} // switch()
|
|
430
|
} // for paramArray
|
|
431
|
this.deriveConstants();
|
|
432
|
},
|
|
433
|
|
|
434
|
deriveConstants : function() {
|
|
435
|
if (this.nagrids == '@null') this.datumCode = 'none';
|
|
436
|
if (this.datumCode && this.datumCode != 'none') {
|
|
437
|
var datumDef = Proj4js.Datum[this.datumCode];
|
|
438
|
if (datumDef) {
|
|
439
|
this.datum_params = datumDef.towgs84.split(',');
|
|
440
|
this.ellps = datumDef.ellipse;
|
|
441
|
this.datumName = datumDef.datumName;
|
|
442
|
}
|
|
443
|
}
|
|
444
|
if (!this.a) { // do we have an ellipsoid?
|
|
445
|
var ellipse = Proj4js.Ellipsoid[this.ellps] ? Proj4js.Ellipsoid[this.ellps] : Proj4js.Ellipsoid['WGS84'];
|
|
446
|
OpenLayers.Util.extend(this, ellipse);
|
|
447
|
}
|
|
448
|
if (this.rf && !this.b) this.b = (1.0 - 1.0/this.rf) * this.a;
|
|
449
|
if (Math.abs(this.a - this.b)<Proj4js.common.EPSLN) this.sphere = true;
|
|
450
|
this.a2 = this.a * this.a; // used in geocentric
|
|
451
|
this.b2 = this.b * this.b; // used in geocentric
|
|
452
|
this.es = (this.a2-this.b2)/this.a2; // e ^ 2
|
|
453
|
//this.es=1-(Math.pow(this.b,2)/Math.pow(this.a,2));
|
|
454
|
this.e = Math.sqrt(this.es); // eccentricity
|
|
455
|
this.ep2=(this.a2-this.b2)/this.b2; // used in geocentric
|
|
456
|
if (!this.k0) this.k0 = 1.0; //default value
|
|
457
|
|
|
458
|
this.datum = new Proj4js.datum(this);
|
|
459
|
}
|
|
460
|
});
|
|
461
|
|
|
462
|
Proj4js.Proj.longlat = {
|
|
463
|
init : function() {
|
|
464
|
//no-op for longlat
|
|
465
|
},
|
|
466
|
forward : function(pt) {
|
|
467
|
//identity transform
|
|
468
|
return pt;
|
|
469
|
},
|
|
470
|
inverse : function(pt) {
|
|
471
|
//identity transform
|
|
472
|
return pt;
|
|
473
|
}
|
|
474
|
};
|
|
475
|
|
|
476
|
/**
|
|
477
|
Proj4js.defs is a collection of coordinate system definition objects in the
|
|
478
|
Proj4 command line format.
|
|
479
|
Generally a def is added by means of a separate .js file for example:
|
|
480
|
|
|
481
|
<SCRIPT type="text/javascript" src="defs/EPSG26912.js"></SCRIPT>
|
|
482
|
|
|
483
|
def is a CS definition in PROJ.4 WKT format, for example:
|
|
484
|
+proj="tmerc" //longlat, etc.
|
|
485
|
+a=majorRadius
|
|
486
|
+b=minorRadius
|
|
487
|
+lat0=somenumber
|
|
488
|
+long=somenumber
|
|
489
|
*/
|
|
490
|
Proj4js.defs = {
|
|
491
|
// These are so widely used, we'll go ahead and throw them in
|
|
492
|
// without requiring a separate .js file
|
|
493
|
'WGS84': "+title=long/lat:WGS84 +proj=longlat +ellps=WGS84 +datum=WGS84",
|
|
494
|
'EPSG:4326': "+title=long/lat:WGS84 +proj=longlat +a=6378137.0 +b=6356752.31424518 +ellps=WGS84 +datum=WGS84",
|
|
495
|
'EPSG:4269': "+title=long/lat:NAD83 +proj=longlat +a=6378137.0 +b=6356752.31414036 +ellps=GRS80 +datum=NAD83"
|
|
496
|
};
|
|
497
|
//+a=6378137.0 +b=6356752.31424518 +ellps=WGS84 +datum=WGS84",
|
|
498
|
Proj4js.common = {
|
|
499
|
PI : Math.PI,
|
|
500
|
HALF_PI : Math.PI*0.5,
|
|
501
|
TWO_PI : Math.PI*2,
|
|
502
|
FORTPI : 0.78539816339744833,
|
|
503
|
R2D : 57.2957795131,
|
|
504
|
D2R : 0.0174532925199,
|
|
505
|
SEC_TO_RAD : 4.84813681109535993589914102357e-6, /* SEC_TO_RAD = Pi/180/3600 */
|
|
506
|
EPSLN : 1.0e-10,
|
|
507
|
MAX_ITER : 20,
|
|
508
|
// following constants from geocent.c
|
|
509
|
COS_67P5 : 0.38268343236508977, /* cosine of 67.5 degrees */
|
|
510
|
AD_C : 1.0026000, /* Toms region 1 constant */
|
|
511
|
|
|
512
|
/* datum_type values */
|
|
513
|
PJD_UNKNOWN : 0,
|
|
514
|
PJD_3PARAM : 1,
|
|
515
|
PJD_7PARAM : 2,
|
|
516
|
PJD_GRIDSHIFT: 3,
|
|
517
|
PJD_WGS84 : 4, // WGS84 or equivalent
|
|
518
|
PJD_NODATUM : 5, // WGS84 or equivalent
|
|
519
|
SRS_WGS84_SEMIMAJOR : 6378137.0, // only used in grid shift transforms
|
|
520
|
|
|
521
|
// Function to compute the constant small m which is the radius of
|
|
522
|
// a parallel of latitude, phi, divided by the semimajor axis.
|
|
523
|
// -----------------------------------------------------------------
|
|
524
|
msfnz : function(eccent, sinphi, cosphi) {
|
|
525
|
var con = eccent * sinphi;
|
|
526
|
return cosphi/(Math.sqrt(1.0 - con * con));
|
|
527
|
},
|
|
528
|
|
|
529
|
// Function to compute the constant small t for use in the forward
|
|
530
|
// computations in the Lambert Conformal Conic and the Polar
|
|
531
|
// Stereographic projections.
|
|
532
|
// -----------------------------------------------------------------
|
|
533
|
tsfnz : function(eccent, phi, sinphi) {
|
|
534
|
var con = eccent * sinphi;
|
|
535
|
var com = .5 * eccent;
|
|
536
|
con = Math.pow(((1.0 - con) / (1.0 + con)), com);
|
|
537
|
return (Math.tan(.5 * (this.HALF_PI - phi))/con);
|
|
538
|
},
|
|
539
|
|
|
540
|
// Function to compute the latitude angle, phi2, for the inverse of the
|
|
541
|
// Lambert Conformal Conic and Polar Stereographic projections.
|
|
542
|
// ----------------------------------------------------------------
|
|
543
|
phi2z : function(eccent, ts) {
|
|
544
|
var eccnth = .5 * eccent;
|
|
545
|
var con, dphi;
|
|
546
|
var phi = this.HALF_PI - 2 * Math.atan(ts);
|
|
547
|
for (i = 0; i <= 15; i++) {
|
|
548
|
con = eccent * Math.sin(phi);
|
|
549
|
dphi = this.HALF_PI - 2 * Math.atan(ts *(Math.pow(((1.0 - con)/(1.0 + con)),eccnth))) - phi;
|
|
550
|
phi += dphi;
|
|
551
|
if (Math.abs(dphi) <= .0000000001) return phi;
|
|
552
|
}
|
|
553
|
alert("phi2z has NoConvergence");
|
|
554
|
return (-9999);
|
|
555
|
},
|
|
556
|
|
|
557
|
/* Function to compute constant small q which is the radius of a
|
|
558
|
parallel of latitude, phi, divided by the semimajor axis.
|
|
559
|
------------------------------------------------------------*/
|
|
560
|
qsfnz : function(eccent,sinphi,cosphi) {
|
|
561
|
var con;
|
|
562
|
if (eccent > 1.0e-7) {
|
|
563
|
con = eccent * sinphi;
|
|
564
|
return (( 1.0- eccent * eccent) * (sinphi /(1.0 - con * con) - (.5/eccent)*Math.log((1.0 - con)/(1.0 + con))));
|
|
565
|
} else {
|
|
566
|
return(2.0 * sinphi);
|
|
567
|
}
|
|
568
|
},
|
|
569
|
|
|
570
|
/* Function to eliminate roundoff errors in asin
|
|
571
|
----------------------------------------------*/
|
|
572
|
asinz : function(x) {
|
|
573
|
if (Math.abs(x)>1.0) {
|
|
574
|
x=(x>1.0)?1.0:-1.0;
|
|
575
|
}
|
|
576
|
return Math.asin(x);
|
|
577
|
},
|
|
578
|
|
|
579
|
// following functions from gctpc cproj.c for transverse mercator projections
|
|
580
|
e0fn : function(x) {return(1.0-0.25*x*(1.0+x/16.0*(3.0+1.25*x)));},
|
|
581
|
e1fn : function(x) {return(0.375*x*(1.0+0.25*x*(1.0+0.46875*x)));},
|
|
582
|
e2fn : function(x) {return(0.05859375*x*x*(1.0+0.75*x));},
|
|
583
|
e3fn : function(x) {return(x*x*x*(35.0/3072.0));},
|
|
584
|
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));},
|
|
585
|
|
|
586
|
srat : function(esinp, exp) {
|
|
587
|
return(Math.pow((1.0-esinp)/(1.0+esinp), exp));
|
|
588
|
},
|
|
589
|
|
|
590
|
// Function to return the sign of an argument
|
|
591
|
sign : function(x) { if (x < 0.0) return(-1); else return(1);},
|
|
592
|
|
|
593
|
// Function to adjust longitude to -180 to 180; input in radians
|
|
594
|
adjust_lon : function(x) {
|
|
595
|
x = (Math.abs(x) < this.PI) ? x: (x - (this.sign(x)*this.TWO_PI) );
|
|
596
|
return x;
|
|
597
|
}
|
|
598
|
|
|
599
|
};
|
|
600
|
|
|
601
|
/** datum object
|
|
602
|
*/
|
|
603
|
Proj4js.datum = OpenLayers.Class({
|
|
604
|
|
|
605
|
initialize : function(proj) {
|
|
606
|
this.datum_type = Proj4js.common.PJD_WGS84; //default setting
|
|
607
|
if (proj.datumCode && proj.datumCode == 'none') {
|
|
608
|
this.datum_type = Proj4js.common.PJD_NODATUM;
|
|
609
|
}
|
|
610
|
if (proj && proj.datum_params) {
|
|
611
|
for (var i=0; i<proj.datum_params.length; i++) {
|
|
612
|
proj.datum_params[i]=parseFloat(proj.datum_params[i]);
|
|
613
|
}
|
|
614
|
if (proj.datum_params[0] != 0 || proj.datum_params[1] != 0 || proj.datum_params[2] != 0 ) {
|
|
615
|
this.datum_type = Proj4js.common.PJD_3PARAM;
|
|
616
|
}
|
|
617
|
if (proj.datum_params.length > 3) {
|
|
618
|
if (proj.datum_params[3] != 0 || proj.datum_params[4] != 0 ||
|
|
619
|
proj.datum_params[5] != 0 || proj.datum_params[6] != 0 ) {
|
|
620
|
this.datum_type = Proj4js.common.PJD_7PARAM;
|
|
621
|
proj.datum_params[3] *= Proj4js.common.SEC_TO_RAD;
|
|
622
|
proj.datum_params[4] *= Proj4js.common.SEC_TO_RAD;
|
|
623
|
proj.datum_params[5] *= Proj4js.common.SEC_TO_RAD;
|
|
624
|
proj.datum_params[6] = (proj.datum_params[6]/1000000.0) + 1.0;
|
|
625
|
}
|
|
626
|
}
|
|
627
|
}
|
|
628
|
if (proj) {
|
|
629
|
this.a = proj.a; //datum object also uses these values
|
|
630
|
this.b = proj.b;
|
|
631
|
this.es = proj.es;
|
|
632
|
this.ep2 = proj.ep2;
|
|
633
|
this.datum_params = proj.datum_params;
|
|
634
|
}
|
|
635
|
},
|
|
636
|
|
|
637
|
/****************************************************************/
|
|
638
|
// cs_compare_datums()
|
|
639
|
// Returns 1 (TRUE) if the two datums match, otherwise 0 (FALSE).
|
|
640
|
compare_datums : function( dest ) {
|
|
641
|
if( this.datum_type != dest.datum_type ) {
|
|
642
|
return false; // false, datums are not equal
|
|
643
|
} else if (this.a != dest.a || Math.abs(this.es-dest.es) > 0.000000000050) {
|
|
644
|
// the tolerence for es is to ensure that GRS80 and WGS84
|
|
645
|
// are considered identical
|
|
646
|
return false;
|
|
647
|
} else if( this.datum_type == Proj4js.common.PJD_3PARAM ) {
|
|
648
|
return (this.datum_params[0] == dest.datum_params[0]
|
|
649
|
&& this.datum_params[1] == dest.datum_params[1]
|
|
650
|
&& this.datum_params[2] == dest.datum_params[2]);
|
|
651
|
} else if( this.datum_type == Proj4js.common.PJD_7PARAM ) {
|
|
652
|
return (this.datum_params[0] == dest.datum_params[0]
|
|
653
|
&& this.datum_params[1] == dest.datum_params[1]
|
|
654
|
&& this.datum_params[2] == dest.datum_params[2]
|
|
655
|
&& this.datum_params[3] == dest.datum_params[3]
|
|
656
|
&& this.datum_params[4] == dest.datum_params[4]
|
|
657
|
&& this.datum_params[5] == dest.datum_params[5]
|
|
658
|
&& this.datum_params[6] == dest.datum_params[6]);
|
|
659
|
} else if( this.datum_type == Proj4js.common.PJD_GRIDSHIFT ) {
|
|
660
|
return strcmp( pj_param(this.params,"snadgrids").s,
|
|
661
|
pj_param(dest.params,"snadgrids").s ) == 0;
|
|
662
|
} else {
|
|
663
|
return true; // datums are equal
|
|
664
|
}
|
|
665
|
}, // cs_compare_datums()
|
|
666
|
|
|
667
|
/*
|
|
668
|
* The function Convert_Geodetic_To_Geocentric converts geodetic coordinates
|
|
669
|
* (latitude, longitude, and height) to geocentric coordinates (X, Y, Z),
|
|
670
|
* according to the current ellipsoid parameters.
|
|
671
|
*
|
|
672
|
* Latitude : Geodetic latitude in radians (input)
|
|
673
|
* Longitude : Geodetic longitude in radians (input)
|
|
674
|
* Height : Geodetic height, in meters (input)
|
|
675
|
* X : Calculated Geocentric X coordinate, in meters (output)
|
|
676
|
* Y : Calculated Geocentric Y coordinate, in meters (output)
|
|
677
|
* Z : Calculated Geocentric Z coordinate, in meters (output)
|
|
678
|
*
|
|
679
|
*/
|
|
680
|
geodetic_to_geocentric : function(p) {
|
|
681
|
var Longitude = p.x;
|
|
682
|
var Latitude = p.y;
|
|
683
|
var Height = p.z ? p.z : 0; //Z value not always supplied
|
|
684
|
var X; // output
|
|
685
|
var Y;
|
|
686
|
var Z;
|
|
687
|
|
|
688
|
var Error_Code=0; // GEOCENT_NO_ERROR;
|
|
689
|
var Rn; /* Earth radius at location */
|
|
690
|
var Sin_Lat; /* Math.sin(Latitude) */
|
|
691
|
var Sin2_Lat; /* Square of Math.sin(Latitude) */
|
|
692
|
var Cos_Lat; /* Math.cos(Latitude) */
|
|
693
|
|
|
694
|
/*
|
|
695
|
** Don't blow up if Latitude is just a little out of the value
|
|
696
|
** range as it may just be a rounding issue. Also removed longitude
|
|
697
|
** test, it should be wrapped by Math.cos() and Math.sin(). NFW for PROJ.4, Sep/2001.
|
|
698
|
*/
|
|
699
|
if (Latitude < -Proj4js.common.HALF_PI && Latitude > -1.001 * Proj4js.common.HALF_PI ) {
|
|
700
|
Latitude = -Proj4js.common.HALF_PI;
|
|
701
|
} else if (Latitude > Proj4js.common.HALF_PI && Latitude < 1.001 * Proj4js.common.HALF_PI ) {
|
|
702
|
Latitude = Proj4js.common.HALF_PI;
|
|
703
|
} else if ((Latitude < -Proj4js.common.HALF_PI) || (Latitude > Proj4js.common.HALF_PI)) {
|
|
704
|
/* Latitude out of range */
|
|
705
|
Proj4js.reportError('geocent:lat out of range:'+Latitude);
|
|
706
|
return null;
|
|
707
|
}
|
|
708
|
|
|
709
|
if (Longitude > Proj4js.common.PI) Longitude -= (2*Proj4js.common.PI);
|
|
710
|
Sin_Lat = Math.sin(Latitude);
|
|
711
|
Cos_Lat = Math.cos(Latitude);
|
|
712
|
Sin2_Lat = Sin_Lat * Sin_Lat;
|
|
713
|
Rn = this.a / (Math.sqrt(1.0e0 - this.es * Sin2_Lat));
|
|
714
|
X = (Rn + Height) * Cos_Lat * Math.cos(Longitude);
|
|
715
|
Y = (Rn + Height) * Cos_Lat * Math.sin(Longitude);
|
|
716
|
Z = ((Rn * (1 - this.es)) + Height) * Sin_Lat;
|
|
717
|
|
|
718
|
p.x = X;
|
|
719
|
p.y = Y;
|
|
720
|
p.z = Z;
|
|
721
|
return Error_Code;
|
|
722
|
}, // cs_geodetic_to_geocentric()
|
|
723
|
|
|
724
|
geocentric_to_geodetic : function (p) {
|
|
725
|
/* local defintions and variables */
|
|
726
|
/* end-criterium of loop, accuracy of sin(Latitude) */
|
|
727
|
var genau = 1.E-12;
|
|
728
|
var genau2 = (genau*genau);
|
|
729
|
var maxiter = 30;
|
|
730
|
|
|
731
|
var P; /* distance between semi-minor axis and location */
|
|
732
|
var RR; /* distance between center and location */
|
|
733
|
var CT; /* sin of geocentric latitude */
|
|
734
|
var ST; /* cos of geocentric latitude */
|
|
735
|
var RX;
|
|
736
|
var RK;
|
|
737
|
var RN; /* Earth radius at location */
|
|
738
|
var CPHI0; /* cos of start or old geodetic latitude in iterations */
|
|
739
|
var SPHI0; /* sin of start or old geodetic latitude in iterations */
|
|
740
|
var CPHI; /* cos of searched geodetic latitude */
|
|
741
|
var SPHI; /* sin of searched geodetic latitude */
|
|
742
|
var SDPHI; /* end-criterium: addition-theorem of sin(Latitude(iter)-Latitude(iter-1)) */
|
|
743
|
var At_Pole; /* indicates location is in polar region */
|
|
744
|
var iter; /* # of continous iteration, max. 30 is always enough (s.a.) */
|
|
745
|
|
|
746
|
var X =p.x;
|
|
747
|
var Y = p.y;
|
|
748
|
var Z = p.z ? p.z : 0.0; //Z value not always supplied
|
|
749
|
var Longitude;
|
|
750
|
var Latitude;
|
|
751
|
var Height;
|
|
752
|
|
|
753
|
At_Pole = false;
|
|
754
|
P = Math.sqrt(X*X+Y*Y);
|
|
755
|
RR = Math.sqrt(X*X+Y*Y+Z*Z);
|
|
756
|
|
|
757
|
/* special cases for latitude and longitude */
|
|
758
|
if (P/this.a < genau) {
|
|
759
|
|
|
760
|
/* special case, if P=0. (X=0., Y=0.) */
|
|
761
|
At_Pole = true;
|
|
762
|
Longitude = 0.0;
|
|
763
|
|
|
764
|
/* if (X,Y,Z)=(0.,0.,0.) then Height becomes semi-minor axis
|
|
765
|
* of ellipsoid (=center of mass), Latitude becomes PI/2 */
|
|
766
|
if (RR/this.a < genau) {
|
|
767
|
Latitude = Proj4js.common.HALF_PI;
|
|
768
|
Height = -this.b;
|
|
769
|
return;
|
|
770
|
}
|
|
771
|
} else {
|
|
772
|
/* ellipsoidal (geodetic) longitude
|
|
773
|
* interval: -PI < Longitude <= +PI */
|
|
774
|
Longitude=Math.atan2(Y,X);
|
|
775
|
}
|
|
776
|
|
|
777
|
/* --------------------------------------------------------------
|
|
778
|
* Following iterative algorithm was developped by
|
|
779
|
* "Institut f?r Erdmessung", University of Hannover, July 1988.
|
|
780
|
* Internet: www.ife.uni-hannover.de
|
|
781
|
* Iterative computation of CPHI,SPHI and Height.
|
|
782
|
* Iteration of CPHI and SPHI to 10**-12 radian resp.
|
|
783
|
* 2*10**-7 arcsec.
|
|
784
|
* --------------------------------------------------------------
|
|
785
|
*/
|
|
786
|
CT = Z/RR;
|
|
787
|
ST = P/RR;
|
|
788
|
RX = 1.0/Math.sqrt(1.0-this.es*(2.0-this.es)*ST*ST);
|
|
789
|
CPHI0 = ST*(1.0-this.es)*RX;
|
|
790
|
SPHI0 = CT*RX;
|
|
791
|
iter = 0;
|
|
792
|
|
|
793
|
/* loop to find sin(Latitude) resp. Latitude
|
|
794
|
* until |sin(Latitude(iter)-Latitude(iter-1))| < genau */
|
|
795
|
do
|
|
796
|
{
|
|
797
|
iter++;
|
|
798
|
RN = this.a/Math.sqrt(1.0-this.es*SPHI0*SPHI0);
|
|
799
|
|
|
800
|
/* ellipsoidal (geodetic) height */
|
|
801
|
Height = P*CPHI0+Z*SPHI0-RN*(1.0-this.es*SPHI0*SPHI0);
|
|
802
|
|
|
803
|
RK = this.es*RN/(RN+Height);
|
|
804
|
RX = 1.0/Math.sqrt(1.0-RK*(2.0-RK)*ST*ST);
|
|
805
|
CPHI = ST*(1.0-RK)*RX;
|
|
806
|
SPHI = CT*RX;
|
|
807
|
SDPHI = SPHI*CPHI0-CPHI*SPHI0;
|
|
808
|
CPHI0 = CPHI;
|
|
809
|
SPHI0 = SPHI;
|
|
810
|
}
|
|
811
|
while (SDPHI*SDPHI > genau2 && iter < maxiter);
|
|
812
|
|
|
813
|
/* ellipsoidal (geodetic) latitude */
|
|
814
|
Latitude=Math.atan(SPHI/Math.abs(CPHI));
|
|
815
|
|
|
816
|
p.x = Longitude;
|
|
817
|
p.y =Latitude;
|
|
818
|
p.z = Height;
|
|
819
|
return p;
|
|
820
|
},
|
|
821
|
|
|
822
|
/** Convert_Geocentric_To_Geodetic
|
|
823
|
* The method used here is derived from 'An Improved Algorithm for
|
|
824
|
* Geocentric to Geodetic Coordinate Conversion', by Ralph Toms, Feb 1996
|
|
825
|
*/
|
|
826
|
geocentric_to_geodetic_noniter : function (p) {
|
|
827
|
var X =p.x;
|
|
828
|
var Y = p.y;
|
|
829
|
var Z = p.z ? p.z : 0; //Z value not always supplied
|
|
830
|
var Longitude;
|
|
831
|
var Latitude;
|
|
832
|
var Height;
|
|
833
|
|
|
834
|
var W; /* distance from Z axis */
|
|
835
|
var W2; /* square of distance from Z axis */
|
|
836
|
var T0; /* initial estimate of vertical component */
|
|
837
|
var T1; /* corrected estimate of vertical component */
|
|
838
|
var S0; /* initial estimate of horizontal component */
|
|
839
|
var S1; /* corrected estimate of horizontal component */
|
|
840
|
var Sin_B0; /* Math.sin(B0), B0 is estimate of Bowring aux variable */
|
|
841
|
var Sin3_B0; /* cube of Math.sin(B0) */
|
|
842
|
var Cos_B0; /* Math.cos(B0) */
|
|
843
|
var Sin_p1; /* Math.sin(phi1), phi1 is estimated latitude */
|
|
844
|
var Cos_p1; /* Math.cos(phi1) */
|
|
845
|
var Rn; /* Earth radius at location */
|
|
846
|
var Sum; /* numerator of Math.cos(phi1) */
|
|
847
|
var At_Pole; /* indicates location is in polar region */
|
|
848
|
|
|
849
|
X = parseFloat(X); // cast from string to float
|
|
850
|
Y = parseFloat(Y);
|
|
851
|
Z = parseFloat(Z);
|
|
852
|
|
|
853
|
At_Pole = false;
|
|
854
|
if (X != 0.0)
|
|
855
|
{
|
|
856
|
Longitude = Math.atan2(Y,X);
|
|
857
|
}
|
|
858
|
else
|
|
859
|
{
|
|
860
|
if (Y > 0)
|
|
861
|
{
|
|
862
|
Longitude = Proj4js.common.HALF_PI;
|
|
863
|
}
|
|
864
|
else if (Y < 0)
|
|
865
|
{
|
|
866
|
Longitude = -Proj4js.common.HALF_PI;
|
|
867
|
}
|
|
868
|
else
|
|
869
|
{
|
|
870
|
At_Pole = true;
|
|
871
|
Longitude = 0.0;
|
|
872
|
if (Z > 0.0)
|
|
873
|
{ /* north pole */
|
|
874
|
Latitude = Proj4js.common.HALF_PI;
|
|
875
|
}
|
|
876
|
else if (Z < 0.0)
|
|
877
|
{ /* south pole */
|
|
878
|
Latitude = -Proj4js.common.HALF_PI;
|
|
879
|
}
|
|
880
|
else
|
|
881
|
{ /* center of earth */
|
|
882
|
Latitude = Proj4js.common.HALF_PI;
|
|
883
|
Height = -this.b;
|
|
884
|
return;
|
|
885
|
}
|
|
886
|
}
|
|
887
|
}
|
|
888
|
W2 = X*X + Y*Y;
|
|
889
|
W = Math.sqrt(W2);
|
|
890
|
T0 = Z * Proj4js.common.AD_C;
|
|
891
|
S0 = Math.sqrt(T0 * T0 + W2);
|
|
892
|
Sin_B0 = T0 / S0;
|
|
893
|
Cos_B0 = W / S0;
|
|
894
|
Sin3_B0 = Sin_B0 * Sin_B0 * Sin_B0;
|
|
895
|
T1 = Z + this.b * this.ep2 * Sin3_B0;
|
|
896
|
Sum = W - this.a * this.es * Cos_B0 * Cos_B0 * Cos_B0;
|
|
897
|
S1 = Math.sqrt(T1*T1 + Sum * Sum);
|
|
898
|
Sin_p1 = T1 / S1;
|
|
899
|
Cos_p1 = Sum / S1;
|
|
900
|
Rn = this.a / Math.sqrt(1.0 - this.es * Sin_p1 * Sin_p1);
|
|
901
|
if (Cos_p1 >= Proj4js.common.COS_67P5)
|
|
902
|
{
|
|
903
|
Height = W / Cos_p1 - Rn;
|
|
904
|
}
|
|
905
|
else if (Cos_p1 <= -Proj4js.common.COS_67P5)
|
|
906
|
{
|
|
907
|
Height = W / -Cos_p1 - Rn;
|
|
908
|
}
|
|
909
|
else
|
|
910
|
{
|
|
911
|
Height = Z / Sin_p1 + Rn * (this.es - 1.0);
|
|
912
|
}
|
|
913
|
if (At_Pole == false)
|
|
914
|
{
|
|
915
|
Latitude = Math.atan(Sin_p1 / Cos_p1);
|
|
916
|
}
|
|
917
|
|
|
918
|
p.x = Longitude;
|
|
919
|
p.y =Latitude;
|
|
920
|
p.z = Height;
|
|
921
|
return p;
|
|
922
|
}, // cs_geocentric_to_geodetic()
|
|
923
|
|
|
924
|
/****************************************************************/
|
|
925
|
// pj_geocentic_to_wgs84( p )
|
|
926
|
// p = point to transform in geocentric coordinates (x,y,z)
|
|
927
|
geocentric_to_wgs84 : function ( p ) {
|
|
928
|
|
|
929
|
if( this.datum_type == Proj4js.common.PJD_3PARAM )
|
|
930
|
{
|
|
931
|
// if( x[io] == HUGE_VAL )
|
|
932
|
// continue;
|
|
933
|
p.x += this.datum_params[0];
|
|
934
|
p.y += this.datum_params[1];
|
|
935
|
p.z += this.datum_params[2];
|
|
936
|
|
|
937
|
}
|
|
938
|
else // if( this.datum_type == Proj4js.common.PJD_7PARAM )
|
|
939
|
{
|
|
940
|
var Dx_BF =this.datum_params[0];
|
|
941
|
var Dy_BF =this.datum_params[1];
|
|
942
|
var Dz_BF =this.datum_params[2];
|
|
943
|
var Rx_BF =this.datum_params[3];
|
|
944
|
var Ry_BF =this.datum_params[4];
|
|
945
|
var Rz_BF =this.datum_params[5];
|
|
946
|
var M_BF =this.datum_params[6];
|
|
947
|
// if( x[io] == HUGE_VAL )
|
|
948
|
// continue;
|
|
949
|
var x_out = M_BF*( p.x - Rz_BF*p.y + Ry_BF*p.z) + Dx_BF;
|
|
950
|
var y_out = M_BF*( Rz_BF*p.x + p.y - Rx_BF*p.z) + Dy_BF;
|
|
951
|
var z_out = M_BF*(-Ry_BF*p.x + Rx_BF*p.y + p.z) + Dz_BF;
|
|
952
|
p.x = x_out;
|
|
953
|
p.y = y_out;
|
|
954
|
p.z = z_out;
|
|
955
|
}
|
|
956
|
}, // cs_geocentric_to_wgs84
|
|
957
|
|
|
958
|
/****************************************************************/
|
|
959
|
// pj_geocentic_from_wgs84()
|
|
960
|
// coordinate system definition,
|
|
961
|
// point to transform in geocentric coordinates (x,y,z)
|
|
962
|
geocentric_from_wgs84 : function( p ) {
|
|
963
|
|
|
964
|
if( this.datum_type == Proj4js.common.PJD_3PARAM )
|
|
965
|
{
|
|
966
|
//if( x[io] == HUGE_VAL )
|
|
967
|
// continue;
|
|
968
|
p.x -= this.datum_params[0];
|
|
969
|
p.y -= this.datum_params[1];
|
|
970
|
p.z -= this.datum_params[2];
|
|
971
|
|
|
972
|
}
|
|
973
|
else // if( this.datum_type == Proj4js.common.PJD_7PARAM )
|
|
974
|
{
|
|
975
|
var Dx_BF =this.datum_params[0];
|
|
976
|
var Dy_BF =this.datum_params[1];
|
|
977
|
var Dz_BF =this.datum_params[2];
|
|
978
|
var Rx_BF =this.datum_params[3];
|
|
979
|
var Ry_BF =this.datum_params[4];
|
|
980
|
var Rz_BF =this.datum_params[5];
|
|
981
|
var M_BF =this.datum_params[6];
|
|
982
|
var x_tmp = (p.x - Dx_BF) / M_BF;
|
|
983
|
var y_tmp = (p.y - Dy_BF) / M_BF;
|
|
984
|
var z_tmp = (p.z - Dz_BF) / M_BF;
|
|
985
|
//if( x[io] == HUGE_VAL )
|
|
986
|
// continue;
|
|
987
|
|
|
988
|
p.x = x_tmp + Rz_BF*y_tmp - Ry_BF*z_tmp;
|
|
989
|
p.y = -Rz_BF*x_tmp + y_tmp + Rx_BF*z_tmp;
|
|
990
|
p.z = Ry_BF*x_tmp - Rx_BF*y_tmp + z_tmp;
|
|
991
|
} //cs_geocentric_from_wgs84()
|
|
992
|
}
|
|
993
|
});
|
|
994
|
|
|
995
|
/** point object, nothing fancy, just allows values to be
|
|
996
|
passed back and forth by reference rather than by value.
|
|
997
|
Other point classes may be used as long as they have
|
|
998
|
x and y properties, which will get modified in the transform method.
|
|
999
|
*/
|
|
1000
|
Proj4js.Point = OpenLayers.Class({
|
|
1001
|
|
|
1002
|
initialize : function(x,y,z) {
|
|
1003
|
if (typeof x == 'object') {
|
|
1004
|
this.x = x[0];
|
|
1005
|
this.y = x[1];
|
|
1006
|
this.z = x[2] || 0.0;
|
|
1007
|
} else {
|
|
1008
|
this.x = x;
|
|
1009
|
this.y = y;
|
|
1010
|
this.z = z || 0.0;
|
|
1011
|
}
|
|
1012
|
},
|
|
1013
|
|
|
1014
|
clone : function() {
|
|
1015
|
return new Proj4js.Point(this.x, this.y, this.z);
|
|
1016
|
},
|
|
1017
|
|
|
1018
|
/**
|
|
1019
|
* Method: toString
|
|
1020
|
* Return a readable string version of the lonlat
|
|
1021
|
*
|
|
1022
|
* Return:
|
|
1023
|
* {String} String representation of Proj4js.Point object.
|
|
1024
|
* (ex. <i>"x=5,y=42"</i>)
|
|
1025
|
*/
|
|
1026
|
toString:function() {
|
|
1027
|
return ("x=" + this.x + ",y=" + this.y);
|
|
1028
|
},
|
|
1029
|
|
|
1030
|
/**
|
|
1031
|
* APIMethod: toShortString
|
|
1032
|
*
|
|
1033
|
* Return:
|
|
1034
|
* {String} Shortened String representation of Proj4js.Point object.
|
|
1035
|
* (ex. <i>"5, 42"</i>)
|
|
1036
|
*/
|
|
1037
|
toShortString:function() {
|
|
1038
|
return (this.x + ", " + this.y);
|
|
1039
|
}
|
|
1040
|
});
|
|
1041
|
|
|
1042
|
Proj4js.PrimeMeridian = {
|
|
1043
|
"greenwich": 0.0, //"0dE",
|
|
1044
|
"lisbon": -9.131906111111, //"9d07'54.862\"W",
|
|
1045
|
"paris": 2.337229166667, //"2d20'14.025\"E",
|
|
1046
|
"bogota": -74.080916666667, //"74d04'51.3\"W",
|
|
1047
|
"madrid": -3.687938888889, //"3d41'16.58\"W",
|
|
1048
|
"rome": 12.452333333333, //"12d27'8.4\"E",
|
|
1049
|
"bern": 7.439583333333, //"7d26'22.5\"E",
|
|
1050
|
"jakarta": 106.807719444444, //"106d48'27.79\"E",
|
|
1051
|
"ferro": -17.666666666667, //"17d40'W",
|
|
1052
|
"brussels": 4.367975, //"4d22'4.71\"E",
|
|
1053
|
"stockholm": 18.058277777778, //"18d3'29.8\"E",
|
|
1054
|
"athens": 23.7163375, //"23d42'58.815\"E",
|
|
1055
|
"oslo": 10.722916666667 //"10d43'22.5\"E"
|
|
1056
|
};
|
|
1057
|
|
|
1058
|
Proj4js.Ellipsoid = {
|
|
1059
|
"MERIT": {a:6378137.0, rf:298.257, ellipseName:"MERIT 1983"},
|
|
1060
|
"SGS85": {a:6378136.0, rf:298.257, ellipseName:"Soviet Geodetic System 85"},
|
|
1061
|
"GRS80": {a:6378137.0, rf:298.257222101, ellipseName:"GRS 1980(IUGG, 1980)"},
|
|
1062
|
"IAU76": {a:6378140.0, rf:298.257, ellipseName:"IAU 1976"},
|
|
1063
|
"airy": {a:6377563.396, b:6356256.910, ellipseName:"Airy 1830"},
|
|
1064
|
"APL4.": {a:6378137, rf:298.25, ellipseName:"Appl. Physics. 1965"},
|
|
1065
|
"NWL9D": {a:6378145.0, rf:298.25, ellipseName:"Naval Weapons Lab., 1965"},
|
|
1066
|
"mod_airy": {a:6377340.189, b:6356034.446, ellipseName:"Modified Airy"},
|
|
1067
|
"andrae": {a:6377104.43, rf:300.0, ellipseName:"Andrae 1876 (Den., Iclnd.)"},
|
|
1068
|
"aust_SA": {a:6378160.0, rf:298.25, ellipseName:"Australian Natl & S. Amer. 1969"},
|
|
1069
|
"GRS67": {a:6378160.0, rf:298.2471674270, ellipseName:"GRS 67(IUGG 1967)"},
|
|
1070
|
"bessel": {a:6377397.155, rf:299.1528128, ellipseName:"Bessel 1841"},
|
|
1071
|
"bess_nam": {a:6377483.865, rf:299.1528128, ellipseName:"Bessel 1841 (Namibia)"},
|
|
1072
|
"clrk66": {a:6378206.4, b:6356583.8, ellipseName:"Clarke 1866"},
|
|
1073
|
"clrk80": {a:6378249.145, rf:293.4663, ellipseName:"Clarke 1880 mod."},
|
|
1074
|
"CPM": {a:6375738.7, rf:334.29, ellipseName:"Comm. des Poids et Mesures 1799"},
|
|
1075
|
"delmbr": {a:6376428.0, rf:311.5, ellipseName:"Delambre 1810 (Belgium)"},
|
|
1076
|
"engelis": {a:6378136.05, rf:298.2566, ellipseName:"Engelis 1985"},
|
|
1077
|
"evrst30": {a:6377276.345, rf:300.8017, ellipseName:"Everest 1830"},
|
|
1078
|
"evrst48": {a:6377304.063, rf:300.8017, ellipseName:"Everest 1948"},
|
|
1079
|
"evrst56": {a:6377301.243, rf:300.8017, ellipseName:"Everest 1956"},
|
|
1080
|
"evrst69": {a:6377295.664, rf:300.8017, ellipseName:"Everest 1969"},
|
|
1081
|
"evrstSS": {a:6377298.556, rf:300.8017, ellipseName:"Everest (Sabah & Sarawak)"},
|
|
1082
|
"fschr60": {a:6378166.0, rf:298.3, ellipseName:"Fischer (Mercury Datum) 1960"},
|
|
1083
|
"fschr60m": {a:6378155.0, rf:298.3, ellipseName:"Fischer 1960"},
|
|
1084
|
"fschr68": {a:6378150.0, rf:298.3, ellipseName:"Fischer 1968"},
|
|
1085
|
"helmert": {a:6378200.0, rf:298.3, ellipseName:"Helmert 1906"},
|
|
1086
|
"hough": {a:6378270.0, rf:297.0, ellipseName:"Hough"},
|
|
1087
|
"intl": {a:6378388.0, rf:297.0, ellipseName:"International 1909 (Hayford)"},
|
|
1088
|
"kaula": {a:6378163.0, rf:298.24, ellipseName:"Kaula 1961"},
|
|
1089
|
"lerch": {a:6378139.0, rf:298.257, ellipseName:"Lerch 1979"},
|
|
1090
|
"mprts": {a:6397300.0, rf:191.0, ellipseName:"Maupertius 1738"},
|
|
1091
|
"new_intl": {a:6378157.5, b:6356772.2, ellipseName:"New International 1967"},
|
|
1092
|
"plessis": {a:6376523.0, rf:6355863.0, ellipseName:"Plessis 1817 (France)"},
|
|
1093
|
"krass": {a:6378245.0, rf:298.3, ellipseName:"Krassovsky, 1942"},
|
|
1094
|
"SEasia": {a:6378155.0, b:6356773.3205, ellipseName:"Southeast Asia"},
|
|
1095
|
"walbeck": {a:6376896.0, b:6355834.8467, ellipseName:"Walbeck"},
|
|
1096
|
"WGS60": {a:6378165.0, rf:298.3, ellipseName:"WGS 60"},
|
|
1097
|
"WGS66": {a:6378145.0, rf:298.25, ellipseName:"WGS 66"},
|
|
1098
|
"WGS72": {a:6378135.0, rf:298.26, ellipseName:"WGS 72"},
|
|
1099
|
"WGS84": {a:6378137.0, rf:298.257223563, ellipseName:"WGS 84"},
|
|
1100
|
"sphere": {a:6370997.0, b:6370997.0, ellipseName:"Normal Sphere (r=6370997)"}
|
|
1101
|
};
|
|
1102
|
|
|
1103
|
Proj4js.Datum = {
|
|
1104
|
"WGS84": {towgs84: "0,0,0", ellipse: "WGS84", datumName: ""},
|
|
1105
|
"GGRS87": {towgs84: "-199.87,74.79,246.62", ellipse: "GRS80", datumName: "Greek_Geodetic_Reference_System_1987"},
|
|
1106
|
"NAD83": {towgs84: "0,0,0", ellipse: "GRS80", datumName: "North_American_Datum_1983"},
|
|
1107
|
"NAD27": {nadgrids: "@conus,@alaska,@ntv2_0.gsb,@ntv1_can.dat", ellipse: "clrk66", datumName: "North_American_Datum_1927"},
|
|
1108
|
"potsdam": {towgs84: "606.0,23.0,413.0", ellipse: "bessel", datumName: "Potsdam Rauenberg 1950 DHDN"},
|
|
1109
|
"carthage": {towgs84: "-263.0,6.0,431.0", ellipse: "clark80", datumName: "Carthage 1934 Tunisia"},
|
|
1110
|
"hermannskogel": {towgs84: "653.0,-212.0,449.0", ellipse: "bessel", datumName: "Hermannskogel"},
|
|
1111
|
"ire65": {towgs84: "482.530,-130.596,564.557,-1.042,-0.214,-0.631,8.15", ellipse: "mod_airy", datumName: "Ireland 1965"},
|
|
1112
|
"nzgd49": {towgs84: "59.47,-5.04,187.44,0.47,-0.1,1.024,-4.5993", ellipse: "intl", datumName: "New Zealand Geodetic Datum 1949"},
|
|
1113
|
"OSGB36": {towgs84: "446.448,-125.157,542.060,0.1502,0.2470,0.8421,-20.4894", ellipse: "airy", datumName: "Airy 1830"}
|
|
1114
|
};
|
|
1115
|
|
|
1116
|
Proj4js.WGS84 = new Proj4js.Proj('WGS84');
|
|
1117
|
Proj4js.Datum['OSB36'] = Proj4js.Datum['OSGB36']; //as returned from spatialreference.org
|
|
1118
|
/* ======================================================================
|
|
1119
|
projCode/sterea.js
|
|
1120
|
====================================================================== */
|
|
1121
|
|
|
1122
|
|
|
1123
|
Proj4js.Proj.sterea = {
|
|
1124
|
dependsOn : 'gauss',
|
|
1125
|
|
|
1126
|
init : function() {
|
|
1127
|
Proj4js.Proj['gauss'].init.apply(this);
|
|
1128
|
if (!this.rc) {
|
|
1129
|
Proj4js.reportError("sterea:init:E_ERROR_0");
|
|
1130
|
return;
|
|
1131
|
}
|
|
1132
|
this.sinc0 = Math.sin(this.phic0);
|
|
1133
|
this.cosc0 = Math.cos(this.phic0);
|
|
1134
|
this.R2 = 2.0 * this.rc;
|
|
1135
|
if (!this.title) this.title = "Oblique Stereographic Alternative";
|
|
1136
|
},
|
|
1137
|
|
|
1138
|
forward : function(p) {
|
|
1139
|
p.x = Proj4js.common.adjust_lon(p.x-this.long0); /* adjust del longitude */
|
|
1140
|
Proj4js.Proj['gauss'].forward.apply(this, [p]);
|
|
1141
|
sinc = Math.sin(p.y);
|
|
1142
|
cosc = Math.cos(p.y);
|
|
1143
|
cosl = Math.cos(p.x);
|
|
1144
|
k = this.k0 * this.R2 / (1.0 + this.sinc0 * sinc + this.cosc0 * cosc * cosl);
|
|
1145
|
p.x = k * cosc * Math.sin(p.x);
|
|
1146
|
p.y = k * (this.cosc0 * sinc - this.sinc0 * cosc * cosl);
|
|
1147
|
p.x = this.a * p.x + this.x0;
|
|
1148
|
p.y = this.a * p.y + this.y0;
|
|
1149
|
return p;
|
|
1150
|
},
|
|
1151
|
|
|
1152
|
inverse : function(p) {
|
|
1153
|
var lon,lat;
|
|
1154
|
p.x = (p.x - this.x0) / this.a; /* descale and de-offset */
|
|
1155
|
p.y = (p.y - this.y0) / this.a;
|
|
1156
|
|
|
1157
|
p.x /= this.k0;
|
|
1158
|
p.y /= this.k0;
|
|
1159
|
if ( (rho = Math.sqrt(p.x*p.x + p.y*p.y)) ) {
|
|
1160
|
c = 2.0 * Math.atan2(rho, this.R2);
|
|
1161
|
sinc = Math.sin(c);
|
|
1162
|
cosc = Math.cos(c);
|
|
1163
|
lat = Math.asin(cosc * this.sinc0 + p.y * sinc * this.cosc0 / rho);
|
|
1164
|
lon = Math.atan2(p.x * sinc, rho * this.cosc0 * cosc - p.y * this.sinc0 * sinc);
|
|
1165
|
} else {
|
|
1166
|
lat = this.phic0;
|
|
1167
|
lon = 0.;
|
|
1168
|
}
|
|
1169
|
|
|
1170
|
p.x = lon;
|
|
1171
|
p.y = lat;
|
|
1172
|
Proj4js.Proj['gauss'].inverse.apply(this,[p]);
|
|
1173
|
p.x = Proj4js.common.adjust_lon(p.x + this.long0); /* adjust longitude to CM */
|
|
1174
|
return p;
|
|
1175
|
}
|
|
1176
|
};
|
|
1177
|
|
|
1178
|
/* ======================================================================
|
|
1179
|
projCode/aea.js
|
|
1180
|
====================================================================== */
|
|
1181
|
|
|
1182
|
/*******************************************************************************
|
|
1183
|
NAME ALBERS CONICAL EQUAL AREA
|
|
1184
|
|
|
1185
|
PURPOSE: Transforms input longitude and latitude to Easting and Northing
|
|
1186
|
for the Albers Conical Equal Area projection. The longitude
|
|
1187
|
and latitude must be in radians. The Easting and Northing
|
|
1188
|
values will be returned in meters.
|
|
1189
|
|
|
1190
|
PROGRAMMER DATE
|
|
1191
|
---------- ----
|
|
1192
|
T. Mittan, Feb, 1992
|
|
1193
|
|
|
1194
|
ALGORITHM REFERENCES
|
|
1195
|
|
|
1196
|
1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
|
|
1197
|
Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
|
|
1198
|
State Government Printing Office, Washington D.C., 1987.
|
|
1199
|
|
|
1200
|
2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections",
|
|
1201
|
U.S. Geological Survey Professional Paper 1453 , United State Government
|
|
1202
|
Printing Office, Washington D.C., 1989.
|
|
1203
|
*******************************************************************************/
|
|
1204
|
|
|
1205
|
|
|
1206
|
Proj4js.Proj.aea = {
|
|
1207
|
init : function() {
|
|
1208
|
|
|
1209
|
if (Math.abs(this.lat1 + this.lat2) < Proj4js.common.EPSLN) {
|
|
1210
|
Proj4js.reportError("aeaInitEqualLatitudes");
|
|
1211
|
return;
|
|
1212
|
}
|
|
1213
|
this.temp = this.b / this.a;
|
|
1214
|
this.es = 1.0 - Math.pow(this.temp,2);
|
|
1215
|
this.e3 = Math.sqrt(this.es);
|
|
1216
|
|
|
1217
|
this.sin_po=Math.sin(this.lat1);
|
|
1218
|
this.cos_po=Math.cos(this.lat1);
|
|
1219
|
this.t1=this.sin_po
|
|
1220
|
this.con = this.sin_po;
|
|
1221
|
this.ms1 = Proj4js.common.msfnz(this.e3,this.sin_po,this.cos_po);
|
|
1222
|
this.qs1 = Proj4js.common.qsfnz(this.e3,this.sin_po,this.cos_po);
|
|
1223
|
|
|
1224
|
this.sin_po=Math.sin(this.lat2);
|
|
1225
|
this.cos_po=Math.cos(this.lat2);
|
|
1226
|
this.t2=this.sin_po;
|
|
1227
|
this.ms2 = Proj4js.common.msfnz(this.e3,this.sin_po,this.cos_po);
|
|
1228
|
this.qs2 = Proj4js.common.qsfnz(this.e3,this.sin_po,this.cos_po);
|
|
1229
|
|
|
1230
|
this.sin_po=Math.sin(this.lat0);
|
|
1231
|
this.cos_po=Math.cos(this.lat0);
|
|
1232
|
this.t3=this.sin_po;
|
|
1233
|
this.qs0 = Proj4js.common.qsfnz(this.e3,this.sin_po,this.cos_po);
|
|
1234
|
|
|
1235
|
if (Math.abs(this.lat1 - this.lat2) > Proj4js.common.EPSLN) {
|
|
1236
|
this.ns0 = (this.ms1 * this.ms1 - this.ms2 *this.ms2)/ (this.qs2 - this.qs1);
|
|
1237
|
} else {
|
|
1238
|
this.ns0 = this.con;
|
|
1239
|
}
|
|
1240
|
this.c = this.ms1 * this.ms1 + this.ns0 * this.qs1;
|
|
1241
|
this.rh = this.a * Math.sqrt(this.c - this.ns0 * this.qs0)/this.ns0;
|
|
1242
|
},
|
|
1243
|
|
|
1244
|
/* Albers Conical Equal Area forward equations--mapping lat,long to x,y
|
|
1245
|
-------------------------------------------------------------------*/
|
|
1246
|
forward: function(p){
|
|
1247
|
|
|
1248
|
var lon=p.x;
|
|
1249
|
var lat=p.y;
|
|
1250
|
|
|
1251
|
this.sin_phi=Math.sin(lat);
|
|
1252
|
this.cos_phi=Math.cos(lat);
|
|
1253
|
|
|
1254
|
var qs = Proj4js.common.qsfnz(this.e3,this.sin_phi,this.cos_phi);
|
|
1255
|
var rh1 =this.a * Math.sqrt(this.c - this.ns0 * qs)/this.ns0;
|
|
1256
|
var theta = this.ns0 * Proj4js.common.adjust_lon(lon - this.long0);
|
|
1257
|
var x = rh1 * Math.sin(theta) + this.x0;
|
|
1258
|
var y = this.rh - rh1 * Math.cos(theta) + this.y0;
|
|
1259
|
|
|
1260
|
p.x = x;
|
|
1261
|
p.y = y;
|
|
1262
|
return p;
|
|
1263
|
},
|
|
1264
|
|
|
1265
|
|
|
1266
|
inverse: function(p) {
|
|
1267
|
var rh1,qs,con,theta,lon,lat;
|
|
1268
|
|
|
1269
|
p.x -= this.x0;
|
|
1270
|
p.y = this.rh - p.y + this.y0;
|
|
1271
|
if (this.ns0 >= 0) {
|
|
1272
|
rh1 = Math.sqrt(p.x *p.x + p.y * p.y);
|
|
1273
|
con = 1.0;
|
|
1274
|
} else {
|
|
1275
|
rh1 = -Math.sqrt(p.x * p.x + p.y *p.y);
|
|
1276
|
con = -1.0;
|
|
1277
|
}
|
|
1278
|
theta = 0.0;
|
|
1279
|
if (rh1 != 0.0) {
|
|
1280
|
theta = Math.atan2(con * p.x, con * p.y);
|
|
1281
|
}
|
|
1282
|
con = rh1 * this.ns0 / this.a;
|
|
1283
|
qs = (this.c - con * con) / this.ns0;
|
|
1284
|
if (this.e3 >= 1e-10) {
|
|
1285
|
con = 1 - .5 * (1.0 -this.es) * Math.log((1.0 - this.e3) / (1.0 + this.e3))/this.e3;
|
|
1286
|
if (Math.abs(Math.abs(con) - Math.abs(qs)) > .0000000001 ) {
|
|
1287
|
lat = this.phi1z(this.e3,qs);
|
|
1288
|
} else {
|
|
1289
|
if (qs >= 0) {
|
|
1290
|
lat = .5 * PI;
|
|
1291
|
} else {
|
|
1292
|
lat = -.5 * PI;
|
|
1293
|
}
|
|
1294
|
}
|
|
1295
|
} else {
|
|
1296
|
lat = this.phi1z(e3,qs);
|
|
1297
|
}
|
|
1298
|
|
|
1299
|
lon = Proj4js.common.adjust_lon(theta/this.ns0 + this.long0);
|
|
1300
|
p.x = lon;
|
|
1301
|
p.y = lat;
|
|
1302
|
return p;
|
|
1303
|
},
|
|
1304
|
|
|
1305
|
/* Function to compute phi1, the latitude for the inverse of the
|
|
1306
|
Albers Conical Equal-Area projection.
|
|
1307
|
-------------------------------------------*/
|
|
1308
|
phi1z: function (eccent,qs) {
|
|
1309
|
var con, com, dphi;
|
|
1310
|
var phi = Proj4js.common.asinz(.5 * qs);
|
|
1311
|
if (eccent < Proj4js.common.EPSLN) return phi;
|
|
1312
|
|
|
1313
|
var eccnts = eccent * eccent;
|
|
1314
|
for (var i = 1; i <= 25; i++) {
|
|
1315
|
sinphi = Math.sin(phi);
|
|
1316
|
cosphi = Math.cos(phi);
|
|
1317
|
con = eccent * sinphi;
|
|
1318
|
com = 1.0 - con * con;
|
|
1319
|
dphi = .5 * com * com / cosphi * (qs / (1.0 - eccnts) - sinphi / com + .5 / eccent * Math.log((1.0 - con) / (1.0 + con)));
|
|
1320
|
phi = phi + dphi;
|
|
1321
|
if (Math.abs(dphi) <= 1e-7) return phi;
|
|
1322
|
}
|
|
1323
|
Proj4js.reportError("aea:phi1z:Convergence error");
|
|
1324
|
return null;
|
|
1325
|
}
|
|
1326
|
|
|
1327
|
};
|
|
1328
|
|
|
1329
|
|
|
1330
|
|
|
1331
|
/* ======================================================================
|
|
1332
|
projCode/poly.js
|
|
1333
|
====================================================================== */
|
|
1334
|
|
|
1335
|
/* Function to compute, phi4, the latitude for the inverse of the
|
|
1336
|
Polyconic projection.
|
|
1337
|
------------------------------------------------------------*/
|
|
1338
|
function phi4z (eccent,e0,e1,e2,e3,a,b,c,phi) {
|
|
1339
|
var sinphi, sin2ph, tanph, ml, mlp, con1, con2, con3, dphi, i;
|
|
1340
|
|
|
1341
|
phi = a;
|
|
1342
|
for (i = 1; i <= 15; i++) {
|
|
1343
|
sinphi = Math.sin(phi);
|
|
1344
|
tanphi = Math.tan(phi);
|
|
1345
|
c = tanphi * Math.sqrt (1.0 - eccent * sinphi * sinphi);
|
|
1346
|
sin2ph = Math.sin (2.0 * phi);
|
|
1347
|
/*
|
|
1348
|
ml = e0 * *phi - e1 * sin2ph + e2 * sin (4.0 * *phi);
|
|
1349
|
mlp = e0 - 2.0 * e1 * cos (2.0 * *phi) + 4.0 * e2 * cos (4.0 * *phi);
|
|
1350
|
*/
|
|
1351
|
ml = e0 * phi - e1 * sin2ph + e2 * Math.sin (4.0 * phi) - e3 * Math.sin (6.0 * phi);
|
|
1352
|
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);
|
|
1353
|
con1 = 2.0 * ml + c * (ml * ml + b) - 2.0 * a * (c * ml + 1.0);
|
|
1354
|
con2 = eccent * sin2ph * (ml * ml + b - 2.0 * a * ml) / (2.0 *c);
|
|
1355
|
con3 = 2.0 * (a - ml) * (c * mlp - 2.0 / sin2ph) - 2.0 * mlp;
|
|
1356
|
dphi = con1 / (con2 + con3);
|
|
1357
|
phi += dphi;
|
|
1358
|
if (Math.abs(dphi) <= .0000000001 ) return(phi);
|
|
1359
|
}
|
|
1360
|
Proj4js.reportError("phi4z: No convergence");
|
|
1361
|
return null;
|
|
1362
|
}
|
|
1363
|
|
|
1364
|
|
|
1365
|
/* Function to compute the constant e4 from the input of the eccentricity
|
|
1366
|
of the spheroid, x. This constant is used in the Polar Stereographic
|
|
1367
|
projection.
|
|
1368
|
--------------------------------------------------------------------*/
|
|
1369
|
function e4fn(x) {
|
|
1370
|
var con, com;
|
|
1371
|
con = 1.0 + x;
|
|
1372
|
com = 1.0 - x;
|
|
1373
|
return (Math.sqrt((Math.pow(con,con))*(Math.pow(com,com))));
|
|
1374
|
}
|
|
1375
|
|
|
1376
|
|
|
1377
|
|
|
1378
|
|
|
1379
|
|
|
1380
|
/*******************************************************************************
|
|
1381
|
NAME POLYCONIC
|
|
1382
|
|
|
1383
|
PURPOSE: Transforms input longitude and latitude to Easting and
|
|
1384
|
Northing for the Polyconic projection. The
|
|
1385
|
longitude and latitude must be in radians. The Easting
|
|
1386
|
and Northing values will be returned in meters.
|
|
1387
|
|
|
1388
|
PROGRAMMER DATE
|
|
1389
|
---------- ----
|
|
1390
|
T. Mittan Mar, 1993
|
|
1391
|
|
|
1392
|
ALGORITHM REFERENCES
|
|
1393
|
|
|
1394
|
1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
|
|
1395
|
Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
|
|
1396
|
State Government Printing Office, Washington D.C., 1987.
|
|
1397
|
|
|
1398
|
2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections",
|
|
1399
|
U.S. Geological Survey Professional Paper 1453 , United State Government
|
|
1400
|
Printing Office, Washington D.C., 1989.
|
|
1401
|
*******************************************************************************/
|
|
1402
|
|
|
1403
|
Proj4js.Proj.poly = {
|
|
1404
|
|
|
1405
|
/* Initialize the POLYCONIC projection
|
|
1406
|
----------------------------------*/
|
|
1407
|
init: function() {
|
|
1408
|
var temp; /* temporary variable */
|
|
1409
|
if (this.lat0=0) this.lat0=90;//this.lat0 ca
|
|
1410
|
|
|
1411
|
/* Place parameters in static storage for common use
|
|
1412
|
-------------------------------------------------*/
|
|
1413
|
this.temp = this.b / this.a;
|
|
1414
|
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
|
|
1415
|
this.e = Math.sqrt(this.es);
|
|
1416
|
this.e0 = Proj4js.common.e0fn(this.es);
|
|
1417
|
this.e1 = Proj4js.common.e1fn(this.es);
|
|
1418
|
this.e2 = Proj4js.common.e2fn(this.es);
|
|
1419
|
this.e3 = Proj4js.common.e3fn(this.es);
|
|
1420
|
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
|
|
1421
|
//if (!this.ml0) {this.ml0=0;}
|
|
1422
|
},
|
|
1423
|
|
|
1424
|
|
|
1425
|
/* Polyconic forward equations--mapping lat,long to x,y
|
|
1426
|
---------------------------------------------------*/
|
|
1427
|
forward: function(p) {
|
|
1428
|
var sinphi, cosphi; /* sin and cos value */
|
|
1429
|
var al; /* temporary values */
|
|
1430
|
var c; /* temporary values */
|
|
1431
|
var con, ml; /* cone constant, small m */
|
|
1432
|
var ms; /* small m */
|
|
1433
|
var x,y;
|
|
1434
|
|
|
1435
|
var lon=p.x;
|
|
1436
|
var lat=p.y;
|
|
1437
|
|
|
1438
|
con = Proj4js.common.adjust_lon(lon - this.long0);
|
|
1439
|
if (Math.abs(lat) <= .0000001) {
|
|
1440
|
x = this.x0 + this.a * con;
|
|
1441
|
y = this.y0 - this.a * this.ml0;
|
|
1442
|
} else {
|
|
1443
|
sinphi = Math.sin(lat);
|
|
1444
|
cosphi = Math.cos(lat);
|
|
1445
|
|
|
1446
|
ml = Proj4js.common.mlfn(this.e0, this.e1, this.e2, this.e3, lat);
|
|
1447
|
ms = Proj4js.common.msfnz(this.e,sinphi,cosphi);
|
|
1448
|
con = sinphi;
|
|
1449
|
x = this.x0 + this.a * ms * Math.sin(con)/sinphi;
|
|
1450
|
y = this.y0 + this.a * (ml - this.ml0 + ms * (1.0 - Math.cos(con))/sinphi);
|
|
1451
|
}
|
|
1452
|
|
|
1453
|
p.x=x;
|
|
1454
|
p.y=y;
|
|
1455
|
return p;
|
|
1456
|
},
|
|
1457
|
|
|
1458
|
|
|
1459
|
/* Inverse equations
|
|
1460
|
-----------------*/
|
|
1461
|
inverse: function(p) {
|
|
1462
|
var sin_phi, cos_phi; /* sin and cos value */
|
|
1463
|
var al; /* temporary values */
|
|
1464
|
var b; /* temporary values */
|
|
1465
|
var c; /* temporary values */
|
|
1466
|
var con, ml; /* cone constant, small m */
|
|
1467
|
var iflg; /* error flag */
|
|
1468
|
var lon,lat;
|
|
1469
|
p.x -= this.x0;
|
|
1470
|
p.y -= this.y0;
|
|
1471
|
al = this.ml0 + p.y/this.a;
|
|
1472
|
iflg = 0;
|
|
1473
|
|
|
1474
|
if (Math.abs(al) <= .0000001) {
|
|
1475
|
lon = p.x/this.a + this.long0;
|
|
1476
|
lat = 0.0;
|
|
1477
|
} else {
|
|
1478
|
b = al * al + (p.x/this.a) * (p.x/this.a);
|
|
1479
|
iflg = phi4z(this.es,this.e0,this.e1,this.e2,this.e3,this.al,b,c,lat);
|
|
1480
|
if (iflg != 1) return(iflg);
|
|
1481
|
lon = Proj4js.common.adjust_lon((asinz(p.x * c / this.a) / Math.sin(lat)) + this.long0);
|
|
1482
|
}
|
|
1483
|
|
|
1484
|
p.x=lon;
|
|
1485
|
p.y=lat;
|
|
1486
|
return p;
|
|
1487
|
}
|
|
1488
|
};
|
|
1489
|
|
|
1490
|
|
|
1491
|
|
|
1492
|
/* ======================================================================
|
|
1493
|
projCode/equi.js
|
|
1494
|
====================================================================== */
|
|
1495
|
|
|
1496
|
/*******************************************************************************
|
|
1497
|
NAME EQUIRECTANGULAR
|
|
1498
|
|
|
1499
|
PURPOSE: Transforms input longitude and latitude to Easting and
|
|
1500
|
Northing for the Equirectangular projection. The
|
|
1501
|
longitude and latitude must be in radians. The Easting
|
|
1502
|
and Northing values will be returned in meters.
|
|
1503
|
|
|
1504
|
PROGRAMMER DATE
|
|
1505
|
---------- ----
|
|
1506
|
T. Mittan Mar, 1993
|
|
1507
|
|
|
1508
|
ALGORITHM REFERENCES
|
|
1509
|
|
|
1510
|
1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
|
|
1511
|
Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
|
|
1512
|
State Government Printing Office, Washington D.C., 1987.
|
|
1513
|
|
|
1514
|
2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections",
|
|
1515
|
U.S. Geological Survey Professional Paper 1453 , United State Government
|
|
1516
|
Printing Office, Washington D.C., 1989.
|
|
1517
|
*******************************************************************************/
|
|
1518
|
Proj4js.Proj.equi = {
|
|
1519
|
|
|
1520
|
init: function() {
|
|
1521
|
if(!this.x0) this.x0=0;
|
|
1522
|
if(!this.y0) this.y0=0;
|
|
1523
|
if(!this.lat0) this.lat0=0;
|
|
1524
|
if(!this.long0) this.long0=0;
|
|
1525
|
///this.t2;
|
|
1526
|
},
|
|
1527
|
|
|
1528
|
|
|
1529
|
|
|
1530
|
/* Equirectangular forward equations--mapping lat,long to x,y
|
|
1531
|
---------------------------------------------------------*/
|
|
1532
|
forward: function(p) {
|
|
1533
|
|
|
1534
|
var lon=p.x;
|
|
1535
|
var lat=p.y;
|
|
1536
|
|
|
1537
|
var dlon = Proj4js.common.adjust_lon(lon - this.long0);
|
|
1538
|
var x = this.x0 +this. a * dlon *Math.cos(this.lat0);
|
|
1539
|
var y = this.y0 + this.a * lat;
|
|
1540
|
|
|
1541
|
this.t1=x;
|
|
1542
|
this.t2=Math.cos(this.lat0);
|
|
1543
|
p.x=x;
|
|
1544
|
p.y=y;
|
|
1545
|
return p;
|
|
1546
|
}, //equiFwd()
|
|
1547
|
|
|
1548
|
|
|
1549
|
|
|
1550
|
/* Equirectangular inverse equations--mapping x,y to lat/long
|
|
1551
|
---------------------------------------------------------*/
|
|
1552
|
inverse: function(p) {
|
|
1553
|
|
|
1554
|
p.x -= this.x0;
|
|
1555
|
p.y -= this.y0;
|
|
1556
|
var lat = p.y /this. a;
|
|
1557
|
|
|
1558
|
if ( Math.abs(lat) > Proj4js.common.HALF_PI) {
|
|
1559
|
Proj4js.reportError("equi:Inv:DataError");
|
|
1560
|
}
|
|
1561
|
var lon = Proj4js.common.adjust_lon(this.long0 + p.x / (this.a * Math.cos(this.lat0)));
|
|
1562
|
p.x=lon;
|
|
1563
|
p.y=lat;
|
|
1564
|
}//equiInv()
|
|
1565
|
};
|
|
1566
|
|
|
1567
|
|
|
1568
|
/* ======================================================================
|
|
1569
|
projCode/merc.js
|
|
1570
|
====================================================================== */
|
|
1571
|
|
|
1572
|
/*******************************************************************************
|
|
1573
|
NAME MERCATOR
|
|
1574
|
|
|
1575
|
PURPOSE: Transforms input longitude and latitude to Easting and
|
|
1576
|
Northing for the Mercator projection. The
|
|
1577
|
longitude and latitude must be in radians. The Easting
|
|
1578
|
and Northing values will be returned in meters.
|
|
1579
|
|
|
1580
|
PROGRAMMER DATE
|
|
1581
|
---------- ----
|
|
1582
|
D. Steinwand, EROS Nov, 1991
|
|
1583
|
T. Mittan Mar, 1993
|
|
1584
|
|
|
1585
|
ALGORITHM REFERENCES
|
|
1586
|
|
|
1587
|
1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
|
|
1588
|
Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
|
|
1589
|
State Government Printing Office, Washington D.C., 1987.
|
|
1590
|
|
|
1591
|
2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections",
|
|
1592
|
U.S. Geological Survey Professional Paper 1453 , United State Government
|
|
1593
|
Printing Office, Washington D.C., 1989.
|
|
1594
|
*******************************************************************************/
|
|
1595
|
|
|
1596
|
//static double r_major = a; /* major axis */
|
|
1597
|
//static double r_minor = b; /* minor axis */
|
|
1598
|
//static double lon_center = long0; /* Center longitude (projection center) */
|
|
1599
|
//static double lat_origin = lat0; /* center latitude */
|
|
1600
|
//static double e,es; /* eccentricity constants */
|
|
1601
|
//static double m1; /* small value m */
|
|
1602
|
//static double false_northing = y0; /* y offset in meters */
|
|
1603
|
//static double false_easting = x0; /* x offset in meters */
|
|
1604
|
//scale_fact = k0
|
|
1605
|
|
|
1606
|
Proj4js.Proj.merc = {
|
|
1607
|
init : function() {
|
|
1608
|
//?this.temp = this.r_minor / this.r_major;
|
|
1609
|
//this.temp = this.b / this.a;
|
|
1610
|
//this.es = 1.0 - Math.sqrt(this.temp);
|
|
1611
|
//this.e = Math.sqrt( this.es );
|
|
1612
|
//?this.m1 = Math.cos(this.lat_origin) / (Math.sqrt( 1.0 - this.es * Math.sin(this.lat_origin) * Math.sin(this.lat_origin)));
|
|
1613
|
//this.m1 = Math.cos(0.0) / (Math.sqrt( 1.0 - this.es * Math.sin(0.0) * Math.sin(0.0)));
|
|
1614
|
if (this.lat_ts) {
|
|
1615
|
if (this.sphere) {
|
|
1616
|
this.k0 = Math.cos(this.lat_ts);
|
|
1617
|
} else {
|
|
1618
|
this.k0 = Proj4js.common.msfnz(this.es, Math.sin(this.lat_ts), Math.cos(this.lat_ts));
|
|
1619
|
}
|
|
1620
|
}
|
|
1621
|
},
|
|
1622
|
|
|
1623
|
/* Mercator forward equations--mapping lat,long to x,y
|
|
1624
|
--------------------------------------------------*/
|
|
1625
|
|
|
1626
|
forward : function(p) {
|
|
1627
|
//alert("ll2m coords : "+coords);
|
|
1628
|
var lon = p.x;
|
|
1629
|
var lat = p.y;
|
|
1630
|
// convert to radians
|
|
1631
|
if ( lat*Proj4js.common.R2D > 90.0 &&
|
|
1632
|
lat*Proj4js.common.R2D < -90.0 &&
|
|
1633
|
lon*Proj4js.common.R2D > 180.0 &&
|
|
1634
|
lon*Proj4js.common.R2D < -180.0) {
|
|
1635
|
Proj4js.reportError("merc:forward: llInputOutOfRange: "+ lon +" : " + lat);
|
|
1636
|
return null;
|
|
1637
|
}
|
|
1638
|
|
|
1639
|
var x,y;
|
|
1640
|
if(Math.abs( Math.abs(lat) - Proj4js.common.HALF_PI) <= Proj4js.common.EPSLN) {
|
|
1641
|
Proj4js.reportError("merc:forward: ll2mAtPoles");
|
|
1642
|
return null;
|
|
1643
|
} else {
|
|
1644
|
if (this.sphere) {
|
|
1645
|
x = this.x0 + this.a * this.k0 * Proj4js.common.adjust_lon(lon - this.long0);
|
|
1646
|
y = this.y0 + this.a * this.k0 * Math.log(Math.tan(Proj4js.common.FORTPI + 0.5*lat));
|
|
1647
|
} else {
|
|
1648
|
var sinphi = Math.sin(lat);
|
|
1649
|
var ts = Proj4js.common.tsfnz(this.e,lat,sinphi);
|
|
1650
|
x = this.x0 + this.a * this.k0 * Proj4js.common.adjust_lon(lon - this.long0);
|
|
1651
|
y = this.y0 - this.a * this.k0 * Math.log(ts);
|
|
1652
|
}
|
|
1653
|
p.x = x;
|
|
1654
|
p.y = y;
|
|
1655
|
return p;
|
|
1656
|
}
|
|
1657
|
},
|
|
1658
|
|
|
1659
|
|
|
1660
|
/* Mercator inverse equations--mapping x,y to lat/long
|
|
1661
|
--------------------------------------------------*/
|
|
1662
|
inverse : function(p) {
|
|
1663
|
|
|
1664
|
var x = p.x - this.x0;
|
|
1665
|
var y = p.y - this.y0;
|
|
1666
|
var lon,lat;
|
|
1667
|
|
|
1668
|
if (this.sphere) {
|
|
1669
|
lat = Proj4js.common.HALF_PI - 2.0 * Math.atan(Math.exp(-y / this.a * this.k0));
|
|
1670
|
} else {
|
|
1671
|
var ts = Math.exp(-y / (this.a * this.k0));
|
|
1672
|
lat = Proj4js.common.phi2z(this.e,ts);
|
|
1673
|
if(lat == -9999) {
|
|
1674
|
Proj4js.reportError("merc:inverse: lat = -9999");
|
|
1675
|
return null;
|
|
1676
|
}
|
|
1677
|
}
|
|
1678
|
lon = Proj4js.common.adjust_lon(this.long0+ x / (this.a * this.k0));
|
|
1679
|
|
|
1680
|
p.x = lon;
|
|
1681
|
p.y = lat;
|
|
1682
|
return p;
|
|
1683
|
}
|
|
1684
|
};
|
|
1685
|
|
|
1686
|
|
|
1687
|
/* ======================================================================
|
|
1688
|
projCode/utm.js
|
|
1689
|
====================================================================== */
|
|
1690
|
|
|
1691
|
/*******************************************************************************
|
|
1692
|
NAME TRANSVERSE MERCATOR
|
|
1693
|
|
|
1694
|
PURPOSE: Transforms input longitude and latitude to Easting and
|
|
1695
|
Northing for the Transverse Mercator projection. The
|
|
1696
|
longitude and latitude must be in radians. The Easting
|
|
1697
|
and Northing values will be returned in meters.
|
|
1698
|
|
|
1699
|
ALGORITHM REFERENCES
|
|
1700
|
|
|
1701
|
1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
|
|
1702
|
Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
|
|
1703
|
State Government Printing Office, Washington D.C., 1987.
|
|
1704
|
|
|
1705
|
2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections",
|
|
1706
|
U.S. Geological Survey Professional Paper 1453 , United State Government
|
|
1707
|
Printing Office, Washington D.C., 1989.
|
|
1708
|
*******************************************************************************/
|
|
1709
|
|
|
1710
|
|
|
1711
|
/**
|
|
1712
|
Initialize Transverse Mercator projection
|
|
1713
|
*/
|
|
1714
|
|
|
1715
|
Proj4js.Proj.utm = {
|
|
1716
|
dependsOn : 'tmerc',
|
|
1717
|
|
|
1718
|
init : function() {
|
|
1719
|
if (!this.zone) {
|
|
1720
|
Proj4js.reportError("utm:init: zone must be specified for UTM");
|
|
1721
|
return;
|
|
1722
|
}
|
|
1723
|
this.lat0 = 0.0;
|
|
1724
|
this.long0 = ((6 * Math.abs(this.zone)) - 183) * Proj4js.common.D2R;
|
|
1725
|
this.x0 = 500000.0;
|
|
1726
|
this.y0 = this.utmSouth ? 10000000.0 : 0.0;
|
|
1727
|
this.k0 = 0.9996;
|
|
1728
|
|
|
1729
|
Proj4js.Proj['tmerc'].init.apply(this);
|
|
1730
|
this.forward = Proj4js.Proj['tmerc'].forward;
|
|
1731
|
this.inverse = Proj4js.Proj['tmerc'].inverse;
|
|
1732
|
}
|
|
1733
|
};
|
|
1734
|
/* ======================================================================
|
|
1735
|
projCode/eqdc.js
|
|
1736
|
====================================================================== */
|
|
1737
|
|
|
1738
|
/*******************************************************************************
|
|
1739
|
NAME EQUIDISTANT CONIC
|
|
1740
|
|
|
1741
|
PURPOSE: Transforms input longitude and latitude to Easting and Northing
|
|
1742
|
for the Equidistant Conic projection. The longitude and
|
|
1743
|
latitude must be in radians. The Easting and Northing values
|
|
1744
|
will be returned in meters.
|
|
1745
|
|
|
1746
|
PROGRAMMER DATE
|
|
1747
|
---------- ----
|
|
1748
|
T. Mittan Mar, 1993
|
|
1749
|
|
|
1750
|
ALGORITHM REFERENCES
|
|
1751
|
|
|
1752
|
1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
|
|
1753
|
Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
|
|
1754
|
State Government Printing Office, Washington D.C., 1987.
|
|
1755
|
|
|
1756
|
2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections",
|
|
1757
|
U.S. Geological Survey Professional Paper 1453 , United State Government
|
|
1758
|
Printing Office, Washington D.C., 1989.
|
|
1759
|
*******************************************************************************/
|
|
1760
|
|
|
1761
|
/* Variables common to all subroutines in this code file
|
|
1762
|
-----------------------------------------------------*/
|
|
1763
|
|
|
1764
|
Proj4js.Proj.eqdc = {
|
|
1765
|
|
|
1766
|
/* Initialize the Equidistant Conic projection
|
|
1767
|
------------------------------------------*/
|
|
1768
|
init: function() {
|
|
1769
|
|
|
1770
|
/* Place parameters in static storage for common use
|
|
1771
|
-------------------------------------------------*/
|
|
1772
|
|
|
1773
|
if(!this.mode) this.mode=0;//chosen default mode
|
|
1774
|
this.temp = this.b / this.a;
|
|
1775
|
this.es = 1.0 - Math.pow(this.temp,2);
|
|
1776
|
this.e = Math.sqrt(this.es);
|
|
1777
|
this.e0 = Proj4js.common.e0fn(this.es);
|
|
1778
|
this.e1 = Proj4js.common.e1fn(this.es);
|
|
1779
|
this.e2 = Proj4js.common.e2fn(this.es);
|
|
1780
|
this.e3 = Proj4js.common.e3fn(this.es);
|
|
1781
|
|
|
1782
|
this.sinphi=Math.sin(this.lat1);
|
|
1783
|
this.cosphi=Math.cos(this.lat1);
|
|
1784
|
|
|
1785
|
this.ms1 = Proj4js.common.msfnz(this.e,this.sinphi,this.cosphi);
|
|
1786
|
this.ml1 = Proj4js.common.mlfn(this.e0, this.e1, this.e2,this.e3, this.lat1);
|
|
1787
|
|
|
1788
|
/* format B
|
|
1789
|
---------*/
|
|
1790
|
if (this.mode != 0) {
|
|
1791
|
if (Math.abs(this.lat1 + this.lat2) < Proj4js.common.EPSLN) {
|
|
1792
|
Proj4js.reportError("eqdc:Init:EqualLatitudes");
|
|
1793
|
//return(81);
|
|
1794
|
}
|
|
1795
|
this.sinphi=Math.sin(this.lat2);
|
|
1796
|
this.cosphi=Math.cos(this.lat2);
|
|
1797
|
|
|
1798
|
this.ms2 = Proj4js.common.msfnz(this.e,this.sinphi,this.cosphi);
|
|
1799
|
this.ml2 = Proj4js.common.mlfn(this.e0, this.e1, this.e2, this.e3, this.lat2);
|
|
1800
|
if (Math.abs(this.lat1 - this.lat2) >= Proj4js.common.EPSLN) {
|
|
1801
|
this.ns = (this.ms1 - this.ms2) / (this.ml2 - this.ml1);
|
|
1802
|
} else {
|
|
1803
|
this.ns = this.sinphi;
|
|
1804
|
}
|
|
1805
|
} else {
|
|
1806
|
this.ns = this.sinphi;
|
|
1807
|
}
|
|
1808
|
this.g = this.ml1 + this.ms1/this.ns;
|
|
1809
|
this.ml0 = Proj4js.common.mlfn(this.e0, this.e1,this. e2, this.e3, this.lat0);
|
|
1810
|
this.rh = this.a * (this.g - this.ml0);
|
|
1811
|
},
|
|
1812
|
|
|
1813
|
|
|
1814
|
/* Equidistant Conic forward equations--mapping lat,long to x,y
|
|
1815
|
-----------------------------------------------------------*/
|
|
1816
|
forward: function(p) {
|
|
1817
|
var lon=p.x;
|
|
1818
|
var lat=p.y;
|
|
1819
|
|
|
1820
|
/* Forward equations
|
|
1821
|
-----------------*/
|
|
1822
|
var ml = Proj4js.common.mlfn(this.e0, this.e1, this.e2, this.e3, lat);
|
|
1823
|
var rh1 = this.a * (this.g - ml);
|
|
1824
|
var theta = this.ns * Proj4js.common.adjust_lon(lon - this.long0);
|
|
1825
|
|
|
1826
|
var x = this.x0 + rh1 * Math.sin(theta);
|
|
1827
|
var y = this.y0 + this.rh - rh1 * Math.cos(theta);
|
|
1828
|
p.x=x;
|
|
1829
|
p.y=y;
|
|
1830
|
return p;
|
|
1831
|
},
|
|
1832
|
|
|
1833
|
/* Inverse equations
|
|
1834
|
-----------------*/
|
|
1835
|
inverse: function(p) {
|
|
1836
|
p.x -= this.x0;
|
|
1837
|
p.y = this.rh - p.y + this.y0;
|
|
1838
|
var con, rh1;
|
|
1839
|
if (this.ns >= 0) {
|
|
1840
|
var rh1 = Math.sqrt(p.x *p.x + p.y * p.y);
|
|
1841
|
var con = 1.0;
|
|
1842
|
} else {
|
|
1843
|
rh1 = -Math.sqrt(p.x *p. x +p. y * p.y);
|
|
1844
|
con = -1.0;
|
|
1845
|
}
|
|
1846
|
var theta = 0.0;
|
|
1847
|
if (rh1 != 0.0) theta = Math.atan2(con *p.x, con *p.y);
|
|
1848
|
var ml = this.g - rh1 /this.a;
|
|
1849
|
var lat = this.phi3z(this.ml,this.e0,this.e1,this.e2,this.e3);
|
|
1850
|
var lon = Proj4js.common.adjust_lon(this.long0 + theta / this.ns);
|
|
1851
|
|
|
1852
|
p.x=lon;
|
|
1853
|
p.y=lat;
|
|
1854
|
return p;
|
|
1855
|
},
|
|
1856
|
|
|
1857
|
/* Function to compute latitude, phi3, for the inverse of the Equidistant
|
|
1858
|
Conic projection.
|
|
1859
|
-----------------------------------------------------------------*/
|
|
1860
|
phi3z: function(ml,e0,e1,e2,e3) {
|
|
1861
|
var phi;
|
|
1862
|
var dphi;
|
|
1863
|
|
|
1864
|
phi = ml;
|
|
1865
|
for (var i = 0; i < 15; i++) {
|
|
1866
|
dphi = (ml + e1 * Math.sin(2.0 * phi) - e2 * Math.sin(4.0 * phi) + e3 * Math.sin(6.0 * phi))/ e0 - phi;
|
|
1867
|
phi += dphi;
|
|
1868
|
if (Math.abs(dphi) <= .0000000001) {
|
|
1869
|
return phi;
|
|
1870
|
}
|
|
1871
|
}
|
|
1872
|
Proj4js.reportError("PHI3Z-CONV:Latitude failed to converge after 15 iterations");
|
|
1873
|
return null;
|
|
1874
|
}
|
|
1875
|
|
|
1876
|
|
|
1877
|
};
|
|
1878
|
/* ======================================================================
|
|
1879
|
projCode/tmerc.js
|
|
1880
|
====================================================================== */
|
|
1881
|
|
|
1882
|
/*******************************************************************************
|
|
1883
|
NAME TRANSVERSE MERCATOR
|
|
1884
|
|
|
1885
|
PURPOSE: Transforms input longitude and latitude to Easting and
|
|
1886
|
Northing for the Transverse Mercator projection. The
|
|
1887
|
longitude and latitude must be in radians. The Easting
|
|
1888
|
and Northing values will be returned in meters.
|
|
1889
|
|
|
1890
|
ALGORITHM REFERENCES
|
|
1891
|
|
|
1892
|
1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
|
|
1893
|
Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
|
|
1894
|
State Government Printing Office, Washington D.C., 1987.
|
|
1895
|
|
|
1896
|
2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections",
|
|
1897
|
U.S. Geological Survey Professional Paper 1453 , United State Government
|
|
1898
|
Printing Office, Washington D.C., 1989.
|
|
1899
|
*******************************************************************************/
|
|
1900
|
|
|
1901
|
|
|
1902
|
/**
|
|
1903
|
Initialize Transverse Mercator projection
|
|
1904
|
*/
|
|
1905
|
|
|
1906
|
Proj4js.Proj.tmerc = {
|
|
1907
|
init : function() {
|
|
1908
|
this.e0 = Proj4js.common.e0fn(this.es);
|
|
1909
|
this.e1 = Proj4js.common.e1fn(this.es);
|
|
1910
|
this.e2 = Proj4js.common.e2fn(this.es);
|
|
1911
|
this.e3 = Proj4js.common.e3fn(this.es);
|
|
1912
|
this.ml0 = this.a * Proj4js.common.mlfn(this.e0, this.e1, this.e2, this.e3, this.lat0);
|
|
1913
|
},
|
|
1914
|
|
|
1915
|
/**
|
|
1916
|
Transverse Mercator Forward - long/lat to x/y
|
|
1917
|
long/lat in radians
|
|
1918
|
*/
|
|
1919
|
forward : function(p) {
|
|
1920
|
var lon = p.x;
|
|
1921
|
var lat = p.y;
|
|
1922
|
|
|
1923
|
var delta_lon = Proj4js.common.adjust_lon(lon - this.long0); // Delta longitude
|
|
1924
|
var con; // cone constant
|
|
1925
|
var x, y;
|
|
1926
|
var sin_phi=Math.sin(lat);
|
|
1927
|
var cos_phi=Math.cos(lat);
|
|
1928
|
|
|
1929
|
if (this.sphere) { /* spherical form */
|
|
1930
|
var b = cos_phi * Math.sin(delta_lon);
|
|
1931
|
if ((Math.abs(Math.abs(b) - 1.0)) < .0000000001) {
|
|
1932
|
Proj4js.reportError("tmerc:forward: Point projects into infinity");
|
|
1933
|
return(93);
|
|
1934
|
} else {
|
|
1935
|
x = .5 * this.a * this.k0 * Math.log((1.0 + b)/(1.0 - b));
|
|
1936
|
con = Math.acos(cos_phi * Math.cos(delta_lon)/Math.sqrt(1.0 - b*b));
|
|
1937
|
if (lat < 0) con = - con;
|
|
1938
|
y = this.a * this.k0 * (con - this.lat0);
|
|
1939
|
}
|
|
1940
|
} else {
|
|
1941
|
var al = cos_phi * delta_lon;
|
|
1942
|
var als = Math.pow(al,2);
|
|
1943
|
var c = this.ep2 * Math.pow(cos_phi,2);
|
|
1944
|
var tq = Math.tan(lat);
|
|
1945
|
var t = Math.pow(tq,2);
|
|
1946
|
con = 1.0 - this.es * Math.pow(sin_phi,2);
|
|
1947
|
var n = this.a / Math.sqrt(con);
|
|
1948
|
var ml = this.a * Proj4js.common.mlfn(this.e0, this.e1, this.e2, this.e3, lat);
|
|
1949
|
|
|
1950
|
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;
|
|
1951
|
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;
|
|
1952
|
|
|
1953
|
}
|
|
1954
|
p.x = x; p.y = y;
|
|
1955
|
return p;
|
|
1956
|
}, // tmercFwd()
|
|
1957
|
|
|
1958
|
/**
|
|
1959
|
Transverse Mercator Inverse - x/y to long/lat
|
|
1960
|
*/
|
|
1961
|
inverse : function(p) {
|
|
1962
|
var con, phi; /* temporary angles */
|
|
1963
|
var delta_phi; /* difference between longitudes */
|
|
1964
|
var i;
|
|
1965
|
var max_iter = 6; /* maximun number of iterations */
|
|
1966
|
var lat, lon;
|
|
1967
|
|
|
1968
|
if (this.sphere) { /* spherical form */
|
|
1969
|
var f = Math.exp(p.x/(this.a * this.k0));
|
|
1970
|
var g = .5 * (f - 1/f);
|
|
1971
|
var temp = this.lat0 + p.y/(this.a * this.k0);
|
|
1972
|
var h = Math.cos(temp);
|
|
1973
|
con = Math.sqrt((1.0 - h * h)/(1.0 + g * g));
|
|
1974
|
lat = Math.asinz(con);
|
|
1975
|
if (temp < 0)
|
|
1976
|
lat = -lat;
|
|
1977
|
if ((g == 0) && (h == 0)) {
|
|
1978
|
lon = this.long0;
|
|
1979
|
} else {
|
|
1980
|
lon = Proj4js.common.adjust_lon(Math.atan2(g,h) + this.long0);
|
|
1981
|
}
|
|
1982
|
} else { // ellipsoidal form
|
|
1983
|
var x = p.x - this.x0;
|
|
1984
|
var y = p.y - this.y0;
|
|
1985
|
|
|
1986
|
con = (this.ml0 + y / this.k0) / this.a;
|
|
1987
|
phi = con;
|
|
1988
|
for (i=0;;i++) {
|
|
1989
|
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;
|
|
1990
|
phi += delta_phi;
|
|
1991
|
if (Math.abs(delta_phi) <= Proj4js.common.EPSLN) break;
|
|
1992
|
if (i >= max_iter) {
|
|
1993
|
Proj4js.reportError("tmerc:inverse: Latitude failed to converge");
|
|
1994
|
return(95);
|
|
1995
|
}
|
|
1996
|
} // for()
|
|
1997
|
if (Math.abs(phi) < Proj4js.common.HALF_PI) {
|
|
1998
|
// sincos(phi, &sin_phi, &cos_phi);
|
|
1999
|
var sin_phi=Math.sin(phi);
|
|
2000
|
var cos_phi=Math.cos(phi);
|
|
2001
|
var tan_phi = Math.tan(phi);
|
|
2002
|
var c = this.ep2 * Math.pow(cos_phi,2);
|
|
2003
|
var cs = Math.pow(c,2);
|
|
2004
|
var t = Math.pow(tan_phi,2);
|
|
2005
|
var ts = Math.pow(t,2);
|
|
2006
|
con = 1.0 - this.es * Math.pow(sin_phi,2);
|
|
2007
|
var n = this.a / Math.sqrt(con);
|
|
2008
|
var r = n * (1.0 - this.es) / con;
|
|
2009
|
var d = x / (n * this.k0);
|
|
2010
|
var ds = Math.pow(d,2);
|
|
2011
|
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)));
|
|
2012
|
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));
|
|
2013
|
} else {
|
|
2014
|
lat = Proj4js.common.HALF_PI * Proj4js.common.sign(y);
|
|
2015
|
lon = this.long0;
|
|
2016
|
}
|
|
2017
|
}
|
|
2018
|
p.x = lon;
|
|
2019
|
p.y = lat;
|
|
2020
|
return p;
|
|
2021
|
} // tmercInv()
|
|
2022
|
};
|
|
2023
|
/* ======================================================================
|
|
2024
|
defs/GOOGLE.js
|
|
2025
|
====================================================================== */
|
|
2026
|
|
|
2027
|
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";
|
|
2028
|
Proj4js.defs["EPSG:900913"]=Proj4js.defs["GOOGLE"];
|
|
2029
|
/* ======================================================================
|
|
2030
|
projCode/ortho.js
|
|
2031
|
====================================================================== */
|
|
2032
|
|
|
2033
|
/*******************************************************************************
|
|
2034
|
NAME ORTHOGRAPHIC
|
|
2035
|
|
|
2036
|
PURPOSE: Transforms input longitude and latitude to Easting and
|
|
2037
|
Northing for the Orthographic projection. The
|
|
2038
|
longitude and latitude must be in radians. The Easting
|
|
2039
|
and Northing values will be returned in meters.
|
|
2040
|
|
|
2041
|
PROGRAMMER DATE
|
|
2042
|
---------- ----
|
|
2043
|
T. Mittan Mar, 1993
|
|
2044
|
|
|
2045
|
ALGORITHM REFERENCES
|
|
2046
|
|
|
2047
|
1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
|
|
2048
|
Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
|
|
2049
|
State Government Printing Office, Washington D.C., 1987.
|
|
2050
|
|
|
2051
|
2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections",
|
|
2052
|
U.S. Geological Survey Professional Paper 1453 , United State Government
|
|
2053
|
Printing Office, Washington D.C., 1989.
|
|
2054
|
*******************************************************************************/
|
|
2055
|
|
|
2056
|
Proj4js.Proj.ortho = {
|
|
2057
|
|
|
2058
|
/* Initialize the Orthographic projection
|
|
2059
|
-------------------------------------*/
|
|
2060
|
init: function(def) {
|
|
2061
|
//double temp; /* temporary variable */
|
|
2062
|
|
|
2063
|
/* Place parameters in static storage for common use
|
|
2064
|
-------------------------------------------------*/;
|
|
2065
|
this.sin_p14=Math.sin(this.lat0);
|
|
2066
|
this.cos_p14=Math.cos(this.lat0);
|
|
2067
|
},
|
|
2068
|
|
|
2069
|
|
|
2070
|
/* Orthographic forward equations--mapping lat,long to x,y
|
|
2071
|
---------------------------------------------------*/
|
|
2072
|
forward: function(p) {
|
|
2073
|
var sinphi, cosphi; /* sin and cos value */
|
|
2074
|
var dlon; /* delta longitude value */
|
|
2075
|
var coslon; /* cos of longitude */
|
|
2076
|
var ksp; /* scale factor */
|
|
2077
|
var g;
|
|
2078
|
var lon=p.x;
|
|
2079
|
var lat=p.y;
|
|
2080
|
/* Forward equations
|
|
2081
|
-----------------*/
|
|
2082
|
dlon = Proj4js.common.adjust_lon(lon - this.long0);
|
|
2083
|
|
|
2084
|
sinphi=Math.sin(lat);
|
|
2085
|
cosphi=Math.cos(lat);
|
|
2086
|
|
|
2087
|
coslon = Math.cos(dlon);
|
|
2088
|
g = this.sin_p14 * sinphi + this.cos_p14 * cosphi * coslon;
|
|
2089
|
ksp = 1.0;
|
|
2090
|
if ((g > 0) || (Math.abs(g) <= Proj4js.common.EPSLN)) {
|
|
2091
|
var x = this.a * ksp * cosphi * Math.sin(dlon);
|
|
2092
|
var y = this.y0 + this.a * ksp * (this.cos_p14 * sinphi - this.sin_p14 * cosphi * coslon);
|
|
2093
|
} else {
|
|
2094
|
Proj4js.reportError("orthoFwdPointError");
|
|
2095
|
}
|
|
2096
|
p.x=x;
|
|
2097
|
p.y=y;
|
|
2098
|
return p;
|
|
2099
|
},
|
|
2100
|
|
|
2101
|
|
|
2102
|
inverse: function(p) {
|
|
2103
|
var rh; /* height above ellipsoid */
|
|
2104
|
var z; /* angle */
|
|
2105
|
var sinz,cosz; /* sin of z and cos of z */
|
|
2106
|
var temp;
|
|
2107
|
var con;
|
|
2108
|
var lon , lat;
|
|
2109
|
/* Inverse equations
|
|
2110
|
-----------------*/
|
|
2111
|
p.x -= this.x0;
|
|
2112
|
p.y -= this.y0;
|
|
2113
|
rh = Math.sqrt(p.x * p.x + p.y * p.y);
|
|
2114
|
if (rh > this.a + .0000001) {
|
|
2115
|
Proj4js.reportError("orthoInvDataError");
|
|
2116
|
}
|
|
2117
|
z = Proj4js.common.asinz(rh / this.a);
|
|
2118
|
|
|
2119
|
sinz=Math.sin(z);
|
|
2120
|
cosi=Math.cos(z);
|
|
2121
|
|
|
2122
|
lon = this.long0;
|
|
2123
|
if (Math.abs(rh) <= Proj4js.common.EPSLN) {
|
|
2124
|
lat = this.lat0;
|
|
2125
|
}
|
|
2126
|
lat = Proj4js.common.asinz(cosz * this.sin_p14 + (y * sinz * this.cos_p14)/rh);
|
|
2127
|
con = Math.abs(lat0) - Proj4js.common.HALF_PI;
|
|
2128
|
if (Math.abs(con) <= Proj4js.common.EPSLN) {
|
|
2129
|
if (this.lat0 >= 0) {
|
|
2130
|
lon = Proj4js.common.adjust_lon(this.long0 + Math.atan2(p.x, -p.y));
|
|
2131
|
} else {
|
|
2132
|
lon = Proj4js.common.adjust_lon(this.long0 -Math.atan2(-p.x, p.y));
|
|
2133
|
}
|
|
2134
|
}
|
|
2135
|
con = cosz - this.sin_p14 * Math.sin(lat);
|
|
2136
|
if ((Math.abs(con) >= Proj4js.common.EPSLN) || (Math.abs(x) >= Proj4js.common.EPSLN)) {
|
|
2137
|
lon = Proj4js.common.adjust_lon(this.long0 + Math.atan2((p.x * sinz * this.cos_p14), (con * rh)));
|
|
2138
|
}
|
|
2139
|
p.x=lon;
|
|
2140
|
p.y=lat;
|
|
2141
|
return p;
|
|
2142
|
}
|
|
2143
|
};
|
|
2144
|
|
|
2145
|
|
|
2146
|
/* ======================================================================
|
|
2147
|
projCode/stere.js
|
|
2148
|
====================================================================== */
|
|
2149
|
|
|
2150
|
|
|
2151
|
// Initialize the Stereographic projection
|
|
2152
|
|
|
2153
|
Proj4js.Proj.stere = {
|
|
2154
|
ssfn_: function(phit, sinphi, eccen) {
|
|
2155
|
sinphi *= eccen;
|
|
2156
|
return (Math.tan (.5 * (Proj4js.common.HALF_PI + phit)) * Math.pow((1. - sinphi) / (1. + sinphi), .5 * eccen));
|
|
2157
|
},
|
|
2158
|
TOL: 1.e-8,
|
|
2159
|
NITER: 8,
|
|
2160
|
CONV: 1.e-10,
|
|
2161
|
S_POLE: 0,
|
|
2162
|
N_POLE: 1,
|
|
2163
|
OBLIQ: 2,
|
|
2164
|
EQUIT: 3,
|
|
2165
|
|
|
2166
|
init : function() {
|
|
2167
|
this.phits = this.lat_ts ? this.lat_ts : Proj4js.common.HALF_PI;
|
|
2168
|
var t = Math.abs(this.lat0);
|
|
2169
|
if ((Math.abs(t) - Proj4js.common.HALF_PI) < Proj4js.common.EPSLN) {
|
|
2170
|
this.mode = this.lat0 < 0. ? this.S_POLE : this.N_POLE;
|
|
2171
|
} else {
|
|
2172
|
this.mode = t > Proj4js.common.EPSLN ? this.OBLIQ : this.EQUIT;
|
|
2173
|
}
|
|
2174
|
this.phits = Math.abs(this.phits);
|
|
2175
|
if (this.es) {
|
|
2176
|
var X;
|
|
2177
|
|
|
2178
|
switch (this.mode) {
|
|
2179
|
case this.N_POLE:
|
|
2180
|
case this.S_POLE:
|
|
2181
|
if (Math.abs(this.phits - Proj4js.common.HALF_PI) < Proj4js.common.EPSLN) {
|
|
2182
|
this.akm1 = 2. * this.k0 / Math.sqrt(Math.pow(1+this.e,1+this.e)*Math.pow(1-this.e,1-this.e));
|
|
2183
|
} else {
|
|
2184
|
t = Math.sin(this.phits);
|
|
2185
|
this.akm1 = Math.cos(this.phits) / Proj4js.common.tsfnz(this.e, this.phits, t);
|
|
2186
|
t *= this.e;
|
|
2187
|
this.akm1 /= Math.sqrt(1. - t * t);
|
|
2188
|
}
|
|
2189
|
break;
|
|
2190
|
case this.EQUIT:
|
|
2191
|
this.akm1 = 2. * this.k0;
|
|
2192
|
break;
|
|
2193
|
case this.OBLIQ:
|
|
2194
|
t = Math.sin(this.lat0);
|
|
2195
|
X = 2. * Math.atan(this.ssfn_(this.lat0, t, this.e)) - Proj4js.common.HALF_PI;
|
|
2196
|
t *= this.e;
|
|
2197
|
this.akm1 = 2. * this.k0 * Math.cos(this.lat0) / Math.sqrt(1. - t * t);
|
|
2198
|
this.sinX1 = Math.sin(X);
|
|
2199
|
this.cosX1 = Math.cos(X);
|
|
2200
|
break;
|
|
2201
|
}
|
|
2202
|
} else {
|
|
2203
|
switch (this.mode) {
|
|
2204
|
case this.OBLIQ:
|
|
2205
|
this.sinph0 = Math.sin(this.lat0);
|
|
2206
|
this.cosph0 = Math.cos(this.lat0);
|
|
2207
|
case this.EQUIT:
|
|
2208
|
this.akm1 = 2. * this.k0;
|
|
2209
|
break;
|
|
2210
|
case this.S_POLE:
|
|
2211
|
case this.N_POLE:
|
|
2212
|
this.akm1 = Math.abs(this.phits - Proj4js.common.HALF_PI) >= Proj4js.common.EPSLN ?
|
|
2213
|
Math.cos(this.phits) / Math.tan(Proj4js.common.FORTPI - .5 * this.phits) :
|
|
2214
|
2. * this.k0 ;
|
|
2215
|
break;
|
|
2216
|
}
|
|
2217
|
}
|
|
2218
|
},
|
|
2219
|
|
|
2220
|
// Stereographic forward equations--mapping lat,long to x,y
|
|
2221
|
forward: function(p) {
|
|
2222
|
var lon = p.x;
|
|
2223
|
var lat = p.y;
|
|
2224
|
var x, y
|
|
2225
|
|
|
2226
|
if (this.sphere) {
|
|
2227
|
var sinphi, cosphi, coslam, sinlam;
|
|
2228
|
|
|
2229
|
sinphi = Math.sin(lat);
|
|
2230
|
cosphi = Math.cos(lat);
|
|
2231
|
coslam = Math.cos(lon);
|
|
2232
|
sinlam = Math.sin(lon);
|
|
2233
|
switch (this.mode) {
|
|
2234
|
case this.EQUIT:
|
|
2235
|
y = 1. + cosphi * coslam;
|
|
2236
|
if (y <= Proj4js.common.EPSLN) {
|
|
2237
|
F_ERROR;
|
|
2238
|
}
|
|
2239
|
y = this.akm1 / y;
|
|
2240
|
x = y * cosphi * sinlam;
|
|
2241
|
y *= sinphi;
|
|
2242
|
break;
|
|
2243
|
case this.OBLIQ:
|
|
2244
|
y = 1. + this.sinph0 * sinphi + this.cosph0 * cosphi * coslam;
|
|
2245
|
if (y <= Proj4js.common.EPSLN) {
|
|
2246
|
F_ERROR;
|
|
2247
|
}
|
|
2248
|
y = this.akm1 / y;
|
|
2249
|
x = y * cosphi * sinlam;
|
|
2250
|
y *= this.cosph0 * sinphi - this.sinph0 * cosphi * coslam;
|
|
2251
|
break;
|
|
2252
|
case this.N_POLE:
|
|
2253
|
coslam = -coslam;
|
|
2254
|
lat = -lat;
|
|
2255
|
//Note: no break here so it conitnues through S_POLE
|
|
2256
|
case this.S_POLE:
|
|
2257
|
if (Math.abs(lat - Proj4js.common.HALF_PI) < this.TOL) {
|
|
2258
|
F_ERROR;
|
|
2259
|
}
|
|
2260
|
y = this.akm1 * Math.tan(Proj4js.common.FORTPI + .5 * lat)
|
|
2261
|
x = sinlam * y;
|
|
2262
|
y *= coslam;
|
|
2263
|
break;
|
|
2264
|
}
|
|
2265
|
} else {
|
|
2266
|
coslam = Math.cos(lon);
|
|
2267
|
sinlam = Math.sin(lon);
|
|
2268
|
sinphi = Math.sin(lat);
|
|
2269
|
if (this.mode == this.OBLIQ || this.mode == this.EQUIT) {
|
|
2270
|
X = 2. * Math.atan(this.ssfn_(lat, sinphi, this.e));
|
|
2271
|
sinX = Math.sin(X - Proj4js.common.HALF_PI);
|
|
2272
|
cosX = Math.cos(X);
|
|
2273
|
}
|
|
2274
|
switch (this.mode) {
|
|
2275
|
case this.OBLIQ:
|
|
2276
|
A = this.akm1 / (this.cosX1 * (1. + this.sinX1 * sinX + this.cosX1 * cosX * coslam));
|
|
2277
|
y = A * (this.cosX1 * sinX - this.sinX1 * cosX * coslam);
|
|
2278
|
x = A * cosX;
|
|
2279
|
break;
|
|
2280
|
case this.EQUIT:
|
|
2281
|
A = 2. * this.akm1 / (1. + cosX * coslam);
|
|
2282
|
y = A * sinX;
|
|
2283
|
x = A * cosX;
|
|
2284
|
break;
|
|
2285
|
case this.S_POLE:
|
|
2286
|
lat = -lat;
|
|
2287
|
coslam = - coslam;
|
|
2288
|
sinphi = -sinphi;
|
|
2289
|
case this.N_POLE:
|
|
2290
|
x = this.akm1 * Proj4js.common.tsfnz(this.e, lat, sinphi);
|
|
2291
|
y = - x * coslam;
|
|
2292
|
break;
|
|
2293
|
}
|
|
2294
|
x = x * sinlam;
|
|
2295
|
}
|
|
2296
|
p.x = x*this.a + this.x0;
|
|
2297
|
p.y = y*this.a + this.y0;
|
|
2298
|
return p;
|
|
2299
|
},
|
|
2300
|
|
|
2301
|
|
|
2302
|
//* Stereographic inverse equations--mapping x,y to lat/long
|
|
2303
|
inverse: function(p) {
|
|
2304
|
var x = (p.x - this.x0)/this.a; /* descale and de-offset */
|
|
2305
|
var y = (p.y - this.y0)/this.a;
|
|
2306
|
var lon, lat
|
|
2307
|
|
|
2308
|
var cosphi, sinphi, tp=0.0, phi_l=0.0, rho, halfe=0.0, pi2=0.0;
|
|
2309
|
var i;
|
|
2310
|
|
|
2311
|
if (this.sphere) {
|
|
2312
|
var c, rh, sinc, cosc;
|
|
2313
|
|
|
2314
|
rh = Math.sqrt(x*x + y*y);
|
|
2315
|
c = 2. * Math.atan(rh / this.akm1);
|
|
2316
|
sinc = Math.sin(c);
|
|
2317
|
cosc = Math.cos(c);
|
|
2318
|
lon = 0.;
|
|
2319
|
switch (this.mode) {
|
|
2320
|
case this.EQUIT:
|
|
2321
|
if (Math.abs(rh) <= Proj4js.common.EPSLN) {
|
|
2322
|
lat = 0.;
|
|
2323
|
} else {
|
|
2324
|
lat = Math.asin(y * sinc / rh);
|
|
2325
|
}
|
|
2326
|
if (cosc != 0. || x != 0.) lon = Math.atan2(x * sinc, cosc * rh);
|
|
2327
|
break;
|
|
2328
|
case this.OBLIQ:
|
|
2329
|
if (Math.abs(rh) <= Proj4js.common.EPSLN) {
|
|
2330
|
lat = this.phi0;
|
|
2331
|
} else {
|
|
2332
|
lat = Math.asin(cosc * sinph0 + y * sinc * cosph0 / rh);
|
|
2333
|
}
|
|
2334
|
c = cosc - sinph0 * Math.sin(lat);
|
|
2335
|
if (c != 0. || x != 0.) {
|
|
2336
|
lon = Math.atan2(x * sinc * cosph0, c * rh);
|
|
2337
|
}
|
|
2338
|
break;
|
|
2339
|
case this.N_POLE:
|
|
2340
|
y = -y;
|
|
2341
|
case this.S_POLE:
|
|
2342
|
if (Math.abs(rh) <= Proj4js.common.EPSLN) {
|
|
2343
|
lat = this.phi0;
|
|
2344
|
} else {
|
|
2345
|
lat = Math.asin(this.mode == this.S_POLE ? -cosc : cosc);
|
|
2346
|
}
|
|
2347
|
lon = (x == 0. && y == 0.) ? 0. : Math.atan2(x, y);
|
|
2348
|
break;
|
|
2349
|
}
|
|
2350
|
} else {
|
|
2351
|
rho = Math.sqrt(x*x + y*y);
|
|
2352
|
switch (this.mode) {
|
|
2353
|
case this.OBLIQ:
|
|
2354
|
case this.EQUIT:
|
|
2355
|
tp = 2. * Math.atan2(rho * this.cosX1 , this.akm1);
|
|
2356
|
cosphi = Math.cos(tp);
|
|
2357
|
sinphi = Math.sin(tp);
|
|
2358
|
if( rho == 0.0 ) {
|
|
2359
|
phi_l = Math.asin(cosphi * this.sinX1);
|
|
2360
|
} else {
|
|
2361
|
phi_l = Math.asin(cosphi * this.sinX1 + (y * sinphi * this.cosX1 / rho));
|
|
2362
|
}
|
|
2363
|
|
|
2364
|
tp = Math.tan(.5 * (Proj4js.common.HALF_PI + phi_l));
|
|
2365
|
x *= sinphi;
|
|
2366
|
y = rho * this.cosX1 * cosphi - y * this.sinX1* sinphi;
|
|
2367
|
pi2 = Proj4js.common.HALF_PI;
|
|
2368
|
halfe = .5 * this.e;
|
|
2369
|
break;
|
|
2370
|
case this.N_POLE:
|
|
2371
|
y = -y;
|
|
2372
|
case this.S_POLE:
|
|
2373
|
tp = - rho / this.akm1
|
|
2374
|
phi_l = Proj4js.common.HALF_PI - 2. * Math.atan(tp);
|
|
2375
|
pi2 = -Proj4js.common.HALF_PI;
|
|
2376
|
halfe = -.5 * this.e;
|
|
2377
|
break;
|
|
2378
|
}
|
|
2379
|
for (i = this.NITER; i--; phi_l = lat) { //check this
|
|
2380
|
sinphi = this.e * Math.sin(phi_l);
|
|
2381
|
lat = 2. * Math.atan(tp * Math.pow((1.+sinphi)/(1.-sinphi), halfe)) - pi2;
|
|
2382
|
if (Math.abs(phi_l - lat) < this.CONV) {
|
|
2383
|
if (this.mode == this.S_POLE) lat = -lat;
|
|
2384
|
lon = (x == 0. && y == 0.) ? 0. : Math.atan2(x, y);
|
|
2385
|
p.x = lon;
|
|
2386
|
p.y = lat
|
|
2387
|
return p;
|
|
2388
|
}
|
|
2389
|
}
|
|
2390
|
}
|
|
2391
|
}
|
|
2392
|
};
|
|
2393
|
/* ======================================================================
|
|
2394
|
projCode/mill.js
|
|
2395
|
====================================================================== */
|
|
2396
|
|
|
2397
|
/*******************************************************************************
|
|
2398
|
NAME MILLER CYLINDRICAL
|
|
2399
|
|
|
2400
|
PURPOSE: Transforms input longitude and latitude to Easting and
|
|
2401
|
Northing for the Miller Cylindrical projection. The
|
|
2402
|
longitude and latitude must be in radians. The Easting
|
|
2403
|
and Northing values will be returned in meters.
|
|
2404
|
|
|
2405
|
PROGRAMMER DATE
|
|
2406
|
---------- ----
|
|
2407
|
T. Mittan March, 1993
|
|
2408
|
|
|
2409
|
This function was adapted from the Lambert Azimuthal Equal Area projection
|
|
2410
|
code (FORTRAN) in the General Cartographic Transformation Package software
|
|
2411
|
which is available from the U.S. Geological Survey National Mapping Division.
|
|
2412
|
|
|
2413
|
ALGORITHM REFERENCES
|
|
2414
|
|
|
2415
|
1. "New Equal-Area Map Projections for Noncircular Regions", John P. Snyder,
|
|
2416
|
The American Cartographer, Vol 15, No. 4, October 1988, pp. 341-355.
|
|
2417
|
|
|
2418
|
2. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
|
|
2419
|
Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
|
|
2420
|
State Government Printing Office, Washington D.C., 1987.
|
|
2421
|
|
|
2422
|
3. "Software Documentation for GCTP General Cartographic Transformation
|
|
2423
|
Package", U.S. Geological Survey National Mapping Division, May 1982.
|
|
2424
|
*******************************************************************************/
|
|
2425
|
|
|
2426
|
Proj4js.Proj.mill = {
|
|
2427
|
|
|
2428
|
/* Initialize the Miller Cylindrical projection
|
|
2429
|
-------------------------------------------*/
|
|
2430
|
init: function() {
|
|
2431
|
//no-op
|
|
2432
|
},
|
|
2433
|
|
|
2434
|
|
|
2435
|
/* Miller Cylindrical forward equations--mapping lat,long to x,y
|
|
2436
|
------------------------------------------------------------*/
|
|
2437
|
forward: function(p) {
|
|
2438
|
var lon=p.x;
|
|
2439
|
var lat=p.y;
|
|
2440
|
/* Forward equations
|
|
2441
|
-----------------*/
|
|
2442
|
dlon = Proj4js.common.adjust_lon(lon -this.long0);
|
|
2443
|
var x = this.x0 +this.R * dlon;
|
|
2444
|
var y = this.y0 + this.R *Math.log(Math.tan((Proj4js.common.PI / 4.0) + (lat / 2.5))) * 1.25;
|
|
2445
|
|
|
2446
|
p.x=x;
|
|
2447
|
p.y=y;
|
|
2448
|
return p;
|
|
2449
|
},//millFwd()
|
|
2450
|
|
|
2451
|
/* Miller Cylindrical inverse equations--mapping x,y to lat/long
|
|
2452
|
------------------------------------------------------------*/
|
|
2453
|
inverse: function(p) {
|
|
2454
|
p. x -= this.x0;
|
|
2455
|
p. y -= this.y0;
|
|
2456
|
|
|
2457
|
var lon = Proj4js.common.adjust_lon(this.long0 + p.x /this.R);
|
|
2458
|
var lat = 2.5 * (Math.atan(Math.exp(p.y/ this.R / 1.25)) - Proj4js.common.PI / 4.0);
|
|
2459
|
|
|
2460
|
p.x=lon;
|
|
2461
|
p.y=lat;
|
|
2462
|
return p;
|
|
2463
|
}//millInv()
|
|
2464
|
};
|
|
2465
|
/* ======================================================================
|
|
2466
|
projCode/sinu.js
|
|
2467
|
====================================================================== */
|
|
2468
|
|
|
2469
|
/*******************************************************************************
|
|
2470
|
NAME SINUSOIDAL
|
|
2471
|
|
|
2472
|
PURPOSE: Transforms input longitude and latitude to Easting and
|
|
2473
|
Northing for the Sinusoidal projection. The
|
|
2474
|
longitude and latitude must be in radians. The Easting
|
|
2475
|
and Northing values will be returned in meters.
|
|
2476
|
|
|
2477
|
PROGRAMMER DATE
|
|
2478
|
---------- ----
|
|
2479
|
D. Steinwand, EROS May, 1991
|
|
2480
|
|
|
2481
|
This function was adapted from the Sinusoidal projection code (FORTRAN) in the
|
|
2482
|
General Cartographic Transformation Package software which is available from
|
|
2483
|
the U.S. Geological Survey National Mapping Division.
|
|
2484
|
|
|
2485
|
ALGORITHM REFERENCES
|
|
2486
|
|
|
2487
|
1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
|
|
2488
|
Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
|
|
2489
|
State Government Printing Office, Washington D.C., 1987.
|
|
2490
|
|
|
2491
|
2. "Software Documentation for GCTP General Cartographic Transformation
|
|
2492
|
Package", U.S. Geological Survey National Mapping Division, May 1982.
|
|
2493
|
*******************************************************************************/
|
|
2494
|
|
|
2495
|
Proj4js.Proj.sinu = {
|
|
2496
|
|
|
2497
|
/* Initialize the Sinusoidal projection
|
|
2498
|
------------------------------------*/
|
|
2499
|
init: function() {
|
|
2500
|
/* Place parameters in static storage for common use
|
|
2501
|
-------------------------------------------------*/
|
|
2502
|
this.R = 6370997.0; //Radius of earth
|
|
2503
|
},
|
|
2504
|
|
|
2505
|
/* Sinusoidal forward equations--mapping lat,long to x,y
|
|
2506
|
-----------------------------------------------------*/
|
|
2507
|
forward: function(p) {
|
|
2508
|
var x,y,delta_lon;
|
|
2509
|
var lon=p.x;
|
|
2510
|
var lat=p.y;
|
|
2511
|
/* Forward equations
|
|
2512
|
-----------------*/
|
|
2513
|
delta_lon = Proj4js.common.adjust_lon(lon - this.long0);
|
|
2514
|
x = this.R * delta_lon * Math.cos(lat) + this.x0;
|
|
2515
|
y = this.R * lat + this.y0;
|
|
2516
|
|
|
2517
|
p.x=x;
|
|
2518
|
p.y=y;
|
|
2519
|
return p;
|
|
2520
|
},
|
|
2521
|
|
|
2522
|
inverse: function(p) {
|
|
2523
|
var lat,temp,lon;
|
|
2524
|
|
|
2525
|
/* Inverse equations
|
|
2526
|
-----------------*/
|
|
2527
|
p.x -= this.x0;
|
|
2528
|
p.y -= this.y0;
|
|
2529
|
lat = p.y / this.R;
|
|
2530
|
if (Math.abs(lat) > Proj4js.common.HALF_PI) {
|
|
2531
|
Proj4js.reportError("sinu:Inv:DataError");
|
|
2532
|
}
|
|
2533
|
temp = Math.abs(lat) - Proj4js.common.HALF_PI;
|
|
2534
|
if (Math.abs(temp) > Proj4js.common.EPSLN) {
|
|
2535
|
temp = this.long0+ p.x / (this.R *Math.cos(lat));
|
|
2536
|
lon = Proj4js.common.adjust_lon(temp);
|
|
2537
|
} else {
|
|
2538
|
lon = this.long0;
|
|
2539
|
}
|
|
2540
|
|
|
2541
|
p.x=lon;
|
|
2542
|
p.y=lat;
|
|
2543
|
return p;
|
|
2544
|
}
|
|
2545
|
};
|
|
2546
|
|
|
2547
|
|
|
2548
|
/* ======================================================================
|
|
2549
|
projCode/geocent.js
|
|
2550
|
====================================================================== */
|
|
2551
|
|
|
2552
|
/*
|
|
2553
|
Author: Richard Greenwood rich@greenwoodmap.com
|
|
2554
|
License: LGPL as per: http://www.gnu.org/copyleft/lesser.html
|
|
2555
|
*/
|
|
2556
|
|
|
2557
|
/**
|
|
2558
|
* convert between geodetic coordinates (longitude, latitude, height)
|
|
2559
|
* and gecentric coordinates (X, Y, Z)
|
|
2560
|
* ported from Proj 4.9.9 geocent.c
|
|
2561
|
*/
|
|
2562
|
|
|
2563
|
|
|
2564
|
// following constants #define'd in geocent.h
|
|
2565
|
// var GEOCENT_NO_ERROR = 0x0000;
|
|
2566
|
var GEOCENT_LAT_ERROR = 0x0001;
|
|
2567
|
// var GEOCENT_LON_ERROR = 0x0002;
|
|
2568
|
// var cs.a_ERROR = 0x0004;
|
|
2569
|
// var cs.b_ERROR = 0x0008;
|
|
2570
|
// var cs.a_LESS_B_ERROR = 0x0010;
|
|
2571
|
|
|
2572
|
// following constants from geocent.c
|
|
2573
|
var COS_67P5 = 0.38268343236508977; /* cosine of 67.5 degrees */
|
|
2574
|
var AD_C = 1.0026000; /* Toms region 1 constant */
|
|
2575
|
|
|
2576
|
function cs_geodetic_to_geocentric (cs, p) {
|
|
2577
|
|
|
2578
|
/*
|
|
2579
|
* The function Convert_Geodetic_To_Geocentric converts geodetic coordinates
|
|
2580
|
* (latitude, longitude, and height) to geocentric coordinates (X, Y, Z),
|
|
2581
|
* according to the current ellipsoid parameters.
|
|
2582
|
*
|
|
2583
|
* Latitude : Geodetic latitude in radians (input)
|
|
2584
|
* Longitude : Geodetic longitude in radians (input)
|
|
2585
|
* Height : Geodetic height, in meters (input)
|
|
2586
|
* X : Calculated Geocentric X coordinate, in meters (output)
|
|
2587
|
* Y : Calculated Geocentric Y coordinate, in meters (output)
|
|
2588
|
* Z : Calculated Geocentric Z coordinate, in meters (output)
|
|
2589
|
*
|
|
2590
|
*/
|
|
2591
|
|
|
2592
|
var Longitude = p.x;
|
|
2593
|
var Latitude = p.y;
|
|
2594
|
var Height = p.z;
|
|
2595
|
var X; // output
|
|
2596
|
var Y;
|
|
2597
|
var Z;
|
|
2598
|
|
|
2599
|
var Error_Code=0; // GEOCENT_NO_ERROR;
|
|
2600
|
var Rn; /* Earth radius at location */
|
|
2601
|
var Sin_Lat; /* Math.sin(Latitude) */
|
|
2602
|
var Sin2_Lat; /* Square of Math.sin(Latitude) */
|
|
2603
|
var Cos_Lat; /* Math.cos(Latitude) */
|
|
2604
|
|
|
2605
|
/*
|
|
2606
|
** Don't blow up if Latitude is just a little out of the value
|
|
2607
|
** range as it may just be a rounding issue. Also removed longitude
|
|
2608
|
** test, it should be wrapped by Math.cos() and Math.sin(). NFW for PROJ.4, Sep/2001.
|
|
2609
|
*/
|
|
2610
|
if( Latitude < -HALF_PI && Latitude > -1.001 * HALF_PI )
|
|
2611
|
Latitude = -HALF_PI;
|
|
2612
|
else if( Latitude > HALF_PI && Latitude < 1.001 * HALF_PI )
|
|
2613
|
Latitude = HALF_PI;
|
|
2614
|
else if ((Latitude < -HALF_PI) || (Latitude > HALF_PI))
|
|
2615
|
{ /* Latitude out of range */
|
|
2616
|
Error_Code |= GEOCENT_LAT_ERROR;
|
|
2617
|
}
|
|
2618
|
|
|
2619
|
if (!Error_Code)
|
|
2620
|
{ /* no errors */
|
|
2621
|
if (Longitude > PI)
|
|
2622
|
Longitude -= (2*PI);
|
|
2623
|
Sin_Lat = Math.sin(Latitude);
|
|
2624
|
Cos_Lat = Math.cos(Latitude);
|
|
2625
|
Sin2_Lat = Sin_Lat * Sin_Lat;
|
|
2626
|
Rn = cs.a / (Math.sqrt(1.0e0 - cs.es * Sin2_Lat));
|
|
2627
|
X = (Rn + Height) * Cos_Lat * Math.cos(Longitude);
|
|
2628
|
Y = (Rn + Height) * Cos_Lat * Math.sin(Longitude);
|
|
2629
|
Z = ((Rn * (1 - cs.es)) + Height) * Sin_Lat;
|
|
2630
|
|
|
2631
|
}
|
|
2632
|
|
|
2633
|
p.x = X;
|
|
2634
|
p.y = Y;
|
|
2635
|
p.z = Z;
|
|
2636
|
return Error_Code;
|
|
2637
|
} // cs_geodetic_to_geocentric()
|
|
2638
|
|
|
2639
|
|
|
2640
|
/** Convert_Geocentric_To_Geodetic
|
|
2641
|
* The method used here is derived from 'An Improved Algorithm for
|
|
2642
|
* Geocentric to Geodetic Coordinate Conversion', by Ralph Toms, Feb 1996
|
|
2643
|
*/
|
|
2644
|
|
|
2645
|
function cs_geocentric_to_geodetic (cs, p) {
|
|
2646
|
|
|
2647
|
var X =p.x;
|
|
2648
|
var Y = p.y;
|
|
2649
|
var Z = p.z;
|
|
2650
|
var Longitude;
|
|
2651
|
var Latitude;
|
|
2652
|
var Height;
|
|
2653
|
|
|
2654
|
var W; /* distance from Z axis */
|
|
2655
|
var W2; /* square of distance from Z axis */
|
|
2656
|
var T0; /* initial estimate of vertical component */
|
|
2657
|
var T1; /* corrected estimate of vertical component */
|
|
2658
|
var S0; /* initial estimate of horizontal component */
|
|
2659
|
var S1; /* corrected estimate of horizontal component */
|
|
2660
|
var Sin_B0; /* Math.sin(B0), B0 is estimate of Bowring aux variable */
|
|
2661
|
var Sin3_B0; /* cube of Math.sin(B0) */
|
|
2662
|
var Cos_B0; /* Math.cos(B0) */
|
|
2663
|
var Sin_p1; /* Math.sin(phi1), phi1 is estimated latitude */
|
|
2664
|
var Cos_p1; /* Math.cos(phi1) */
|
|
2665
|
var Rn; /* Earth radius at location */
|
|
2666
|
var Sum; /* numerator of Math.cos(phi1) */
|
|
2667
|
var At_Pole; /* indicates location is in polar region */
|
|
2668
|
|
|
2669
|
X = parseFloat(X); // cast from string to float
|
|
2670
|
Y = parseFloat(Y);
|
|
2671
|
Z = parseFloat(Z);
|
|
2672
|
|
|
2673
|
At_Pole = false;
|
|
2674
|
if (X != 0.0)
|
|
2675
|
{
|
|
2676
|
Longitude = Math.atan2(Y,X);
|
|
2677
|
}
|
|
2678
|
else
|
|
2679
|
{
|
|
2680
|
if (Y > 0)
|
|
2681
|
{
|
|
2682
|
Longitude = HALF_PI;
|
|
2683
|
}
|
|
2684
|
else if (Y < 0)
|
|
2685
|
{
|
|
2686
|
Longitude = -HALF_PI;
|
|
2687
|
}
|
|
2688
|
else
|
|
2689
|
{
|
|
2690
|
At_Pole = true;
|
|
2691
|
Longitude = 0.0;
|
|
2692
|
if (Z > 0.0)
|
|
2693
|
{ /* north pole */
|
|
2694
|
Latitude = HALF_PI;
|
|
2695
|
}
|
|
2696
|
else if (Z < 0.0)
|
|
2697
|
{ /* south pole */
|
|
2698
|
Latitude = -HALF_PI;
|
|
2699
|
}
|
|
2700
|
else
|
|
2701
|
{ /* center of earth */
|
|
2702
|
Latitude = HALF_PI;
|
|
2703
|
Height = -cs.b;
|
|
2704
|
return;
|
|
2705
|
}
|
|
2706
|
}
|
|
2707
|
}
|
|
2708
|
W2 = X*X + Y*Y;
|
|
2709
|
W = Math.sqrt(W2);
|
|
2710
|
T0 = Z * AD_C;
|
|
2711
|
S0 = Math.sqrt(T0 * T0 + W2);
|
|
2712
|
Sin_B0 = T0 / S0;
|
|
2713
|
Cos_B0 = W / S0;
|
|
2714
|
Sin3_B0 = Sin_B0 * Sin_B0 * Sin_B0;
|
|
2715
|
T1 = Z + cs.b * cs.ep2 * Sin3_B0;
|
|
2716
|
Sum = W - cs.a * cs.es * Cos_B0 * Cos_B0 * Cos_B0;
|
|
2717
|
S1 = Math.sqrt(T1*T1 + Sum * Sum);
|
|
2718
|
Sin_p1 = T1 / S1;
|
|
2719
|
Cos_p1 = Sum / S1;
|
|
2720
|
Rn = cs.a / Math.sqrt(1.0 - cs.es * Sin_p1 * Sin_p1);
|
|
2721
|
if (Cos_p1 >= COS_67P5)
|
|
2722
|
{
|
|
2723
|
Height = W / Cos_p1 - Rn;
|
|
2724
|
}
|
|
2725
|
else if (Cos_p1 <= -COS_67P5)
|
|
2726
|
{
|
|
2727
|
Height = W / -Cos_p1 - Rn;
|
|
2728
|
}
|
|
2729
|
else
|
|
2730
|
{
|
|
2731
|
Height = Z / Sin_p1 + Rn * (cs.es - 1.0);
|
|
2732
|
}
|
|
2733
|
if (At_Pole == false)
|
|
2734
|
{
|
|
2735
|
Latitude = Math.atan(Sin_p1 / Cos_p1);
|
|
2736
|
}
|
|
2737
|
|
|
2738
|
p.x = Longitude;
|
|
2739
|
p.y =Latitude;
|
|
2740
|
p.z = Height;
|
|
2741
|
return 0;
|
|
2742
|
} // cs_geocentric_to_geodetic()
|
|
2743
|
|
|
2744
|
|
|
2745
|
|
|
2746
|
/****************************************************************/
|
|
2747
|
// pj_geocentic_to_wgs84(defn, p )
|
|
2748
|
// defn = coordinate system definition,
|
|
2749
|
// p = point to transform in geocentric coordinates (x,y,z)
|
|
2750
|
function cs_geocentric_to_wgs84( defn, p ) {
|
|
2751
|
|
|
2752
|
if( defn.datum_type == PJD_3PARAM )
|
|
2753
|
{
|
|
2754
|
// if( x[io] == HUGE_VAL )
|
|
2755
|
// continue;
|
|
2756
|
p.x += defn.datum_params[0];
|
|
2757
|
p.y += defn.datum_params[1];
|
|
2758
|
p.z += defn.datum_params[2];
|
|
2759
|
|
|
2760
|
}
|
|
2761
|
else // if( defn.datum_type == PJD_7PARAM )
|
|
2762
|
{
|
|
2763
|
var Dx_BF =defn.datum_params[0];
|
|
2764
|
var Dy_BF =defn.datum_params[1];
|
|
2765
|
var Dz_BF =defn.datum_params[2];
|
|
2766
|
var Rx_BF =defn.datum_params[3];
|
|
2767
|
var Ry_BF =defn.datum_params[4];
|
|
2768
|
var Rz_BF =defn.datum_params[5];
|
|
2769
|
var M_BF =defn.datum_params[6];
|
|
2770
|
// if( x[io] == HUGE_VAL )
|
|
2771
|
// continue;
|
|
2772
|
var x_out = M_BF*( p.x - Rz_BF*p.y + Ry_BF*p.z) + Dx_BF;
|
|
2773
|
var y_out = M_BF*( Rz_BF*p.x + p.y - Rx_BF*p.z) + Dy_BF;
|
|
2774
|
var z_out = M_BF*(-Ry_BF*p.x + Rx_BF*p.y + p.z) + Dz_BF;
|
|
2775
|
p.x = x_out;
|
|
2776
|
p.y = y_out;
|
|
2777
|
p.z = z_out;
|
|
2778
|
}
|
|
2779
|
} // cs_geocentric_to_wgs84
|
|
2780
|
|
|
2781
|
/****************************************************************/
|
|
2782
|
// pj_geocentic_from_wgs84()
|
|
2783
|
// coordinate system definition,
|
|
2784
|
// point to transform in geocentric coordinates (x,y,z)
|
|
2785
|
function cs_geocentric_from_wgs84( defn, p ) {
|
|
2786
|
|
|
2787
|
if( defn.datum_type == PJD_3PARAM )
|
|
2788
|
{
|
|
2789
|
//if( x[io] == HUGE_VAL )
|
|
2790
|
// continue;
|
|
2791
|
p.x -= defn.datum_params[0];
|
|
2792
|
p.y -= defn.datum_params[1];
|
|
2793
|
p.z -= defn.datum_params[2];
|
|
2794
|
|
|
2795
|
}
|
|
2796
|
else // if( defn.datum_type == PJD_7PARAM )
|
|
2797
|
{
|
|
2798
|
var Dx_BF =defn.datum_params[0];
|
|
2799
|
var Dy_BF =defn.datum_params[1];
|
|
2800
|
var Dz_BF =defn.datum_params[2];
|
|
2801
|
var Rx_BF =defn.datum_params[3];
|
|
2802
|
var Ry_BF =defn.datum_params[4];
|
|
2803
|
var Rz_BF =defn.datum_params[5];
|
|
2804
|
var M_BF =defn.datum_params[6];
|
|
2805
|
var x_tmp = (p.x - Dx_BF) / M_BF;
|
|
2806
|
var y_tmp = (p.y - Dy_BF) / M_BF;
|
|
2807
|
var z_tmp = (p.z - Dz_BF) / M_BF;
|
|
2808
|
//if( x[io] == HUGE_VAL )
|
|
2809
|
// continue;
|
|
2810
|
|
|
2811
|
p.x = x_tmp + Rz_BF*y_tmp - Ry_BF*z_tmp;
|
|
2812
|
p.y = -Rz_BF*x_tmp + y_tmp + Rx_BF*z_tmp;
|
|
2813
|
p.z = Ry_BF*x_tmp - Rx_BF*y_tmp + z_tmp;
|
|
2814
|
}
|
|
2815
|
} //cs_geocentric_from_wgs84()
|
|
2816
|
/* ======================================================================
|
|
2817
|
projCode/vandg.js
|
|
2818
|
====================================================================== */
|
|
2819
|
|
|
2820
|
/*******************************************************************************
|
|
2821
|
NAME VAN DER GRINTEN
|
|
2822
|
|
|
2823
|
PURPOSE: Transforms input Easting and Northing to longitude and
|
|
2824
|
latitude for the Van der Grinten projection. The
|
|
2825
|
Easting and Northing must be in meters. The longitude
|
|
2826
|
and latitude values will be returned in radians.
|
|
2827
|
|
|
2828
|
PROGRAMMER DATE
|
|
2829
|
---------- ----
|
|
2830
|
T. Mittan March, 1993
|
|
2831
|
|
|
2832
|
This function was adapted from the Van Der Grinten projection code
|
|
2833
|
(FORTRAN) in the General Cartographic Transformation Package software
|
|
2834
|
which is available from the U.S. Geological Survey National Mapping Division.
|
|
2835
|
|
|
2836
|
ALGORITHM REFERENCES
|
|
2837
|
|
|
2838
|
1. "New Equal-Area Map Projections for Noncircular Regions", John P. Snyder,
|
|
2839
|
The American Cartographer, Vol 15, No. 4, October 1988, pp. 341-355.
|
|
2840
|
|
|
2841
|
2. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
|
|
2842
|
Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
|
|
2843
|
State Government Printing Office, Washington D.C., 1987.
|
|
2844
|
|
|
2845
|
3. "Software Documentation for GCTP General Cartographic Transformation
|
|
2846
|
Package", U.S. Geological Survey National Mapping Division, May 1982.
|
|
2847
|
*******************************************************************************/
|
|
2848
|
|
|
2849
|
Proj4js.Proj.vandg = {
|
|
2850
|
|
|
2851
|
/* Initialize the Van Der Grinten projection
|
|
2852
|
----------------------------------------*/
|
|
2853
|
init: function() {
|
|
2854
|
this.R = 6370997.0; //Radius of earth
|
|
2855
|
},
|
|
2856
|
|
|
2857
|
forward: function(p) {
|
|
2858
|
|
|
2859
|
var lon=p.x;
|
|
2860
|
var lat=p.y;
|
|
2861
|
|
|
2862
|
/* Forward equations
|
|
2863
|
-----------------*/
|
|
2864
|
var dlon = Proj4js.common.adjust_lon(lon - this.long0);
|
|
2865
|
var x,y;
|
|
2866
|
|
|
2867
|
if (Math.abs(lat) <= Proj4js.common.EPSLN) {
|
|
2868
|
x = this.x0 + this.R * dlon;
|
|
2869
|
y = this.y0;
|
|
2870
|
}
|
|
2871
|
var theta = Proj4js.common.asinz(2.0 * Math.abs(lat / Proj4js.common.PI));
|
|
2872
|
if ((Math.abs(dlon) <= Proj4js.common.EPSLN) || (Math.abs(Math.abs(lat) - Proj4js.common.HALF_PI) <= Proj4js.common.EPSLN)) {
|
|
2873
|
x = this.x0;
|
|
2874
|
if (lat >= 0) {
|
|
2875
|
y = this.y0 + Proj4js.common.PI * this.R * Math.tan(.5 * theta);
|
|
2876
|
} else {
|
|
2877
|
y = this.y0 + Proj4js.common.PI * this.R * - Math.tan(.5 * theta);
|
|
2878
|
}
|
|
2879
|
// return(OK);
|
|
2880
|
}
|
|
2881
|
var al = .5 * Math.abs((Proj4js.common.PI / dlon) - (dlon / Proj4js.common.PI));
|
|
2882
|
var asq = al * al;
|
|
2883
|
var sinth = Math.sin(theta);
|
|
2884
|
var costh = Math.cos(theta);
|
|
2885
|
|
|
2886
|
var g = costh / (sinth + costh - 1.0);
|
|
2887
|
var gsq = g * g;
|
|
2888
|
var m = g * (2.0 / sinth - 1.0);
|
|
2889
|
var msq = m * m;
|
|
2890
|
var con = Proj4js.common.PI * this.R * (al * (g - msq) + Math.sqrt(asq * (g - msq) * (g - msq) - (msq + asq) * (gsq - msq))) / (msq + asq);
|
|
2891
|
if (dlon < 0) {
|
|
2892
|
con = -con;
|
|
2893
|
}
|
|
2894
|
x = this.x0 + con;
|
|
2895
|
con = Math.abs(con / (Proj4js.common.PI * this.R));
|
|
2896
|
if (lat >= 0) {
|
|
2897
|
y = this.y0 + Proj4js.common.PI * this.R * Math.sqrt(1.0 - con * con - 2.0 * al * con);
|
|
2898
|
} else {
|
|
2899
|
y = this.y0 - Proj4js.common.PI * this.R * Math.sqrt(1.0 - con * con - 2.0 * al * con);
|
|
2900
|
}
|
|
2901
|
p.x = x;
|
|
2902
|
p.y = y;
|
|
2903
|
return p;
|
|
2904
|
},
|
|
2905
|
|
|
2906
|
/* Van Der Grinten inverse equations--mapping x,y to lat/long
|
|
2907
|
---------------------------------------------------------*/
|
|
2908
|
inverse: function(p) {
|
|
2909
|
var dlon;
|
|
2910
|
var xx,yy,xys,c1,c2,c3;
|
|
2911
|
var al,asq;
|
|
2912
|
var a1;
|
|
2913
|
var m1;
|
|
2914
|
var con;
|
|
2915
|
var th1;
|
|
2916
|
var d;
|
|
2917
|
|
|
2918
|
/* inverse equations
|
|
2919
|
-----------------*/
|
|
2920
|
p.x -= this.x0;
|
|
2921
|
p.y -= this.y0;
|
|
2922
|
con = Proj4js.common.PI * this.R;
|
|
2923
|
xx = p.x / con;
|
|
2924
|
yy =p.y / con;
|
|
2925
|
xys = xx * xx + yy * yy;
|
|
2926
|
c1 = -Math.abs(yy) * (1.0 + xys);
|
|
2927
|
c2 = c1 - 2.0 * yy * yy + xx * xx;
|
|
2928
|
c3 = -2.0 * c1 + 1.0 + 2.0 * yy * yy + xys * xys;
|
|
2929
|
d = yy * yy / c3 + (2.0 * c2 * c2 * c2 / c3 / c3 / c3 - 9.0 * c1 * c2 / c3 /c3) / 27.0;
|
|
2930
|
a1 = (c1 - c2 * c2 / 3.0 / c3) / c3;
|
|
2931
|
m1 = 2.0 * Math.sqrt( -a1 / 3.0);
|
|
2932
|
con = ((3.0 * d) / a1) / m1;
|
|
2933
|
if (Math.abs(con) > 1.0) {
|
|
2934
|
if (con >= 0.0) {
|
|
2935
|
con = 1.0;
|
|
2936
|
} else {
|
|
2937
|
con = -1.0;
|
|
2938
|
}
|
|
2939
|
}
|
|
2940
|
th1 = Math.acos(con) / 3.0;
|
|
2941
|
if (p.y >= 0) {
|
|
2942
|
lat = (-m1 *Math.cos(th1 + Proj4js.common.PI / 3.0) - c2 / 3.0 / c3) * Proj4js.common.PI;
|
|
2943
|
} else {
|
|
2944
|
lat = -(-m1 * Math.cos(th1 + PI / 3.0) - c2 / 3.0 / c3) * Proj4js.common.PI;
|
|
2945
|
}
|
|
2946
|
|
|
2947
|
if (Math.abs(xx) < Proj4js.common.EPSLN) {
|
|
2948
|
lon = this.long0;
|
|
2949
|
}
|
|
2950
|
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);
|
|
2951
|
|
|
2952
|
p.x=lon;
|
|
2953
|
p.y=lat;
|
|
2954
|
return p;
|
|
2955
|
}
|
|
2956
|
};
|
|
2957
|
/* ======================================================================
|
|
2958
|
projCode/gauss.js
|
|
2959
|
====================================================================== */
|
|
2960
|
|
|
2961
|
|
|
2962
|
Proj4js.Proj.gauss = {
|
|
2963
|
|
|
2964
|
init : function() {
|
|
2965
|
sphi = Math.sin(this.lat0);
|
|
2966
|
cphi = Math.cos(this.lat0);
|
|
2967
|
cphi *= cphi;
|
|
2968
|
this.rc = Math.sqrt(1.0 - this.es) / (1.0 - this.es * sphi * sphi);
|
|
2969
|
this.C = Math.sqrt(1.0 + this.es * cphi * cphi / (1.0 - this.es));
|
|
2970
|
this.phic0 = Math.asin(sphi / this.C);
|
|
2971
|
this.ratexp = 0.5 * this.C * this.e;
|
|
2972
|
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));
|
|
2973
|
},
|
|
2974
|
|
|
2975
|
forward : function(p) {
|
|
2976
|
var lon = p.x;
|
|
2977
|
var lat = p.y;
|
|
2978
|
|
|
2979
|
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;
|
|
2980
|
p.x = this.C * lon;
|
|
2981
|
return p;
|
|
2982
|
},
|
|
2983
|
|
|
2984
|
inverse : function(p) {
|
|
2985
|
var DEL_TOL = 1e-14;
|
|
2986
|
var lon = p.x / this.C;
|
|
2987
|
var lat = p.y;
|
|
2988
|
num = Math.pow(Math.tan(0.5 * lat + Proj4js.common.FORTPI)/this.K, 1./this.C);
|
|
2989
|
for (var i = Proj4js.common.MAX_ITER; i>0; --i) {
|
|
2990
|
lat = 2.0 * Math.atan(num * Proj4js.common.srat(this.e * Math.sin(p.y), -0.5 * this.e)) - Proj4js.common.HALF_PI;
|
|
2991
|
if (Math.abs(lat - p.y) < DEL_TOL) break;
|
|
2992
|
p.y = lat;
|
|
2993
|
}
|
|
2994
|
/* convergence failed */
|
|
2995
|
if (!i) {
|
|
2996
|
Proj4js.reportError("gauss:inverse:convergence failed");
|
|
2997
|
return null;
|
|
2998
|
}
|
|
2999
|
p.x = lon;
|
|
3000
|
p.y = lat;
|
|
3001
|
return p;
|
|
3002
|
}
|
|
3003
|
};
|
|
3004
|
|
|
3005
|
/* ======================================================================
|
|
3006
|
projCode/omerc.js
|
|
3007
|
====================================================================== */
|
|
3008
|
|
|
3009
|
/*******************************************************************************
|
|
3010
|
NAME OBLIQUE MERCATOR (HOTINE)
|
|
3011
|
|
|
3012
|
PURPOSE: Transforms input longitude and latitude to Easting and
|
|
3013
|
Northing for the Oblique Mercator projection. The
|
|
3014
|
longitude and latitude must be in radians. The Easting
|
|
3015
|
and Northing values will be returned in meters.
|
|
3016
|
|
|
3017
|
PROGRAMMER DATE
|
|
3018
|
---------- ----
|
|
3019
|
T. Mittan Mar, 1993
|
|
3020
|
|
|
3021
|
ALGORITHM REFERENCES
|
|
3022
|
|
|
3023
|
1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
|
|
3024
|
Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
|
|
3025
|
State Government Printing Office, Washington D.C., 1987.
|
|
3026
|
|
|
3027
|
2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections",
|
|
3028
|
U.S. Geological Survey Professional Paper 1453 , United State Government
|
|
3029
|
Printing Office, Washington D.C., 1989.
|
|
3030
|
*******************************************************************************/
|
|
3031
|
|
|
3032
|
Proj4js.Proj.omerc = {
|
|
3033
|
|
|
3034
|
/* Initialize the Oblique Mercator projection
|
|
3035
|
------------------------------------------*/
|
|
3036
|
init: function() {
|
|
3037
|
if (!this.mode) this.mode=0;
|
|
3038
|
if (!this.lon1) {this.lon1=0;this.mode=1;}
|
|
3039
|
if (!this.lon2) this.lon2=0;
|
|
3040
|
if (!this.lat2) this.lat2=0;
|
|
3041
|
|
|
3042
|
/* Place parameters in static storage for common use
|
|
3043
|
-------------------------------------------------*/
|
|
3044
|
var temp = this.b/ this.a;
|
|
3045
|
var es = 1.0 - Math.pow(temp,2);
|
|
3046
|
var e = Math.sqrt(es);
|
|
3047
|
|
|
3048
|
this.sin_p20=Math.sin(this.lat0);
|
|
3049
|
this.cos_p20=Math.cos(this.lat0);
|
|
3050
|
|
|
3051
|
this.con = 1.0 - this.es * this.sin_p20 * this.sin_p20;
|
|
3052
|
this.com = Math.sqrt(1.0 - es);
|
|
3053
|
this.bl = Math.sqrt(1.0 + this.es * Math.pow(this.cos_p20,4.0)/(1.0 - es));
|
|
3054
|
this.al = this.a * this.bl * this.k0 * this.com / this.con;
|
|
3055
|
if (Math.abs(this.lat0) < Proj4js.common.EPSLN) {
|
|
3056
|
this.ts = 1.0;
|
|
3057
|
this.d = 1.0;
|
|
3058
|
this.el = 1.0;
|
|
3059
|
} else {
|
|
3060
|
this.ts = Proj4js.common.tsfnz(this.e,this.lat0,this.sin_p20);
|
|
3061
|
this.con = Math.sqrt(this.con);
|
|
3062
|
this.d = this.bl * this.com / (this.cos_p20 * this.con);
|
|
3063
|
if ((this.d * this.d - 1.0) > 0.0) {
|
|
3064
|
if (this.lat0 >= 0.0) {
|
|
3065
|
this.f = this.d + Math.sqrt(this.d * this.d - 1.0);
|
|
3066
|
} else {
|
|
3067
|
this.f = this.d - Math.sqrt(this.d * this.d - 1.0);
|
|
3068
|
}
|
|
3069
|
} else {
|
|
3070
|
this.f = this.d;
|
|
3071
|
}
|
|
3072
|
this.el = this.f * Math.pow(this.ts,this.bl);
|
|
3073
|
}
|
|
3074
|
|
|
3075
|
//this.longc=52.60353916666667;
|
|
3076
|
|
|
3077
|
if (this.mode != 0) {
|
|
3078
|
this.g = .5 * (this.f - 1.0/this.f);
|
|
3079
|
this.gama = Proj4js.common.asinz(Math.sin(this.alpha) / this.d);
|
|
3080
|
this.longc= this.longc - Proj4js.common.asinz(this.g * Math.tan(this.gama))/this.bl;
|
|
3081
|
|
|
3082
|
/* Report parameters common to format B
|
|
3083
|
-------------------------------------*/
|
|
3084
|
//genrpt(azimuth * R2D,"Azimuth of Central Line: ");
|
|
3085
|
//cenlon(lon_origin);
|
|
3086
|
// cenlat(lat_origin);
|
|
3087
|
|
|
3088
|
this.con = Math.abs(this.lat0);
|
|
3089
|
if ((this.con > Proj4js.common.EPSLN) && (Math.abs(this.con - Proj4js.common.HALF_PI) > Proj4js.common.EPSLN)) {
|
|
3090
|
this.singam=Math.sin(this.gama);
|
|
3091
|
this.cosgam=Math.cos(this.gama);
|
|
3092
|
|
|
3093
|
this.sinaz=Math.sin(this.alpha);
|
|
3094
|
this.cosaz=Math.cos(this.alpha);
|
|
3095
|
|
|
3096
|
if (this.lat0>= 0) {
|
|
3097
|
this.u = (this.al / this.bl) * Math.atan(Math.sqrt(this.d*this.d - 1.0)/this.cosaz);
|
|
3098
|
} else {
|
|
3099
|
this.u = -(this.al / this.bl) *Math.atan(Math.sqrt(this.d*this.d - 1.0)/this.cosaz);
|
|
3100
|
}
|
|
3101
|
} else {
|
|
3102
|
Proj4js.reportError("omerc:Init:DataError");
|
|
3103
|
}
|
|
3104
|
} else {
|
|
3105
|
this.sinphi =Math. sin(this.at1);
|
|
3106
|
this.ts1 = Proj4js.common.tsfnz(this.e,this.lat1,this.sinphi);
|
|
3107
|
this.sinphi = Math.sin(this.lat2);
|
|
3108
|
this.ts2 = Proj4js.common.tsfnz(this.e,this.lat2,this.sinphi);
|
|
3109
|
this.h = Math.pow(this.ts1,this.bl);
|
|
3110
|
this.l = Math.pow(this.ts2,this.bl);
|
|
3111
|
this.f = this.el/this.h;
|
|
3112
|
this.g = .5 * (this.f - 1.0/this.f);
|
|
3113
|
this.j = (this.el * this.el - this.l * this.h)/(this.el * this.el + this.l * this.h);
|
|
3114
|
this.p = (this.l - this.h) / (this.l + this.h);
|
|
3115
|
this.dlon = this.lon1 - this.lon2;
|
|
3116
|
if (this.dlon < -Proj4js.common.PI) this.lon2 = this.lon2 - 2.0 * Proj4js.common.PI;
|
|
3117
|
if (this.dlon > Proj4js.common.PI) this.lon2 = this.lon2 + 2.0 * Proj4js.common.PI;
|
|
3118
|
this.dlon = this.lon1 - this.lon2;
|
|
3119
|
this.longc = .5 * (this.lon1 + this.lon2) -Math.atan(this.j * Math.tan(.5 * this.bl * this.dlon)/this.p)/this.bl;
|
|
3120
|
this.dlon = Proj4js.common.adjust_lon(this.lon1 - this.longc);
|
|
3121
|
this.gama = Math.atan(Math.sin(this.bl * this.dlon)/this.g);
|
|
3122
|
this.alpha = Proj4js.common.asinz(this.d * Math.sin(this.gama));
|
|
3123
|
|
|
3124
|
/* Report parameters common to format A
|
|
3125
|
-------------------------------------*/
|
|
3126
|
|
|
3127
|
if (Math.abs(this.lat1 - this.lat2) <= Proj4js.common.EPSLN) {
|
|
3128
|
Proj4js.reportError("omercInitDataError");
|
|
3129
|
//return(202);
|
|
3130
|
} else {
|
|
3131
|
this.con = Math.abs(this.lat1);
|
|
3132
|
}
|
|
3133
|
if ((this.con <= Proj4js.common.EPSLN) || (Math.abs(this.con - HALF_PI) <= Proj4js.common.EPSLN)) {
|
|
3134
|
Proj4js.reportError("omercInitDataError");
|
|
3135
|
//return(202);
|
|
3136
|
} else {
|
|
3137
|
if (Math.abs(Math.abs(this.lat0) - Proj4js.common.HALF_PI) <= Proj4js.common.EPSLN) {
|
|
3138
|
Proj4js.reportError("omercInitDataError");
|
|
3139
|
//return(202);
|
|
3140
|
}
|
|
3141
|
}
|
|
3142
|
|
|
3143
|
this.singam=Math.sin(this.gam);
|
|
3144
|
this.cosgam=Math.cos(this.gam);
|
|
3145
|
|
|
3146
|
this.sinaz=Math.sin(this.alpha);
|
|
3147
|
this.cosaz=Math.cos(this.alpha);
|
|
3148
|
|
|
3149
|
|
|
3150
|
if (this.lat0 >= 0) {
|
|
3151
|
this.u = (this.al/this.bl) * Math.atan(Math.sqrt(this.d * this.d - 1.0)/this.cosaz);
|
|
3152
|
} else {
|
|
3153
|
this.u = -(this.al/this.bl) * Math.atan(Math.sqrt(this.d * this.d - 1.0)/this.cosaz);
|
|
3154
|
}
|
|
3155
|
}
|
|
3156
|
},
|
|
3157
|
|
|
3158
|
|
|
3159
|
/* Oblique Mercator forward equations--mapping lat,long to x,y
|
|
3160
|
----------------------------------------------------------*/
|
|
3161
|
forward: function(p) {
|
|
3162
|
var theta; /* angle */
|
|
3163
|
var sin_phi, cos_phi;/* sin and cos value */
|
|
3164
|
var b; /* temporary values */
|
|
3165
|
var c, t, tq; /* temporary values */
|
|
3166
|
var con, n, ml; /* cone constant, small m */
|
|
3167
|
var q,us,vl;
|
|
3168
|
var ul,vs;
|
|
3169
|
var s;
|
|
3170
|
var dlon;
|
|
3171
|
var ts1;
|
|
3172
|
|
|
3173
|
var lon=p.x;
|
|
3174
|
var lat=p.y;
|
|
3175
|
/* Forward equations
|
|
3176
|
-----------------*/
|
|
3177
|
sin_phi = Math.sin(lat);
|
|
3178
|
dlon = Proj4js.common.adjust_lon(lon - this.longc);
|
|
3179
|
vl = Math.sin(this.bl * dlon);
|
|
3180
|
if (Math.abs(Math.abs(lat) - Proj4js.common.HALF_PI) > Proj4js.common.EPSLN) {
|
|
3181
|
ts1 = Proj4js.common.tsfnz(this.e,lat,sin_phi);
|
|
3182
|
q = this.el / (Math.pow(ts1,this.bl));
|
|
3183
|
s = .5 * (q - 1.0 / q);
|
|
3184
|
t = .5 * (q + 1.0/ q);
|
|
3185
|
ul = (s * this.singam - vl * this.cosgam) / t;
|
|
3186
|
con = Math.cos(this.bl * dlon);
|
|
3187
|
if (Math.abs(con) < .0000001) {
|
|
3188
|
us = this.al * this.bl * dlon;
|
|
3189
|
} else {
|
|
3190
|
us = this.al * Math.atan((s * this.cosgam + vl * this.singam) / con)/this.bl;
|
|
3191
|
if (con < 0) us = us + Proj4js.common.PI * this.al / this.bl;
|
|
3192
|
}
|
|
3193
|
} else {
|
|
3194
|
if (lat >= 0) {
|
|
3195
|
ul = this.singam;
|
|
3196
|
} else {
|
|
3197
|
ul = -this.singam;
|
|
3198
|
}
|
|
3199
|
us = this.al * lat / this.bl;
|
|
3200
|
}
|
|
3201
|
if (Math.abs(Math.abs(ul) - 1.0) <= Proj4js.common.EPSLN) {
|
|
3202
|
//alert("Point projects into infinity","omer-for");
|
|
3203
|
Proj4js.reportError("omercFwdInfinity");
|
|
3204
|
//return(205);
|
|
3205
|
}
|
|
3206
|
vs = .5 * this.al * Math.log((1.0 - ul)/(1.0 + ul)) / this.bl;
|
|
3207
|
us = us - this.u;
|
|
3208
|
var x = this.x0 + vs * this.cosaz + us * this.sinaz;
|
|
3209
|
var y = this.y0 + us * this.cosaz - vs * this.sinaz;
|
|
3210
|
|
|
3211
|
p.x=x;
|
|
3212
|
p.y=y;
|
|
3213
|
return p;
|
|
3214
|
},
|
|
3215
|
|
|
3216
|
inverse: function(p) {
|
|
3217
|
var delta_lon; /* Delta longitude (Given longitude - center */
|
|
3218
|
var theta; /* angle */
|
|
3219
|
var delta_theta; /* adjusted longitude */
|
|
3220
|
var sin_phi, cos_phi;/* sin and cos value */
|
|
3221
|
var b; /* temporary values */
|
|
3222
|
var c, t, tq; /* temporary values */
|
|
3223
|
var con, n, ml; /* cone constant, small m */
|
|
3224
|
var vs,us,q,s,ts1;
|
|
3225
|
var vl,ul,bs;
|
|
3226
|
var dlon;
|
|
3227
|
var flag;
|
|
3228
|
|
|
3229
|
/* Inverse equations
|
|
3230
|
-----------------*/
|
|
3231
|
p.x -= this.x0;
|
|
3232
|
p.y -= this.y0;
|
|
3233
|
flag = 0;
|
|
3234
|
vs = p.x * this.cosaz - p.y * this.sinaz;
|
|
3235
|
us = p.y * this.cosaz + p.x * this.sinaz;
|
|
3236
|
us = us + this.u;
|
|
3237
|
q = Math.exp(-this.bl * vs / this.al);
|
|
3238
|
s = .5 * (q - 1.0/q);
|
|
3239
|
t = .5 * (q + 1.0/q);
|
|
3240
|
vl = Math.sin(this.bl * us / this.al);
|
|
3241
|
ul = (vl * this.cosgam + s * this.singam)/t;
|
|
3242
|
if (Math.abs(Math.abs(ul) - 1.0) <= Proj4js.common.EPSLN)
|
|
3243
|
{
|
|
3244
|
lon = this.longc;
|
|
3245
|
if (ul >= 0.0) {
|
|
3246
|
lat = Proj4js.common.HALF_PI;
|
|
3247
|
} else {
|
|
3248
|
lat = -Proj4js.common.HALF_PI;
|
|
3249
|
}
|
|
3250
|
} else {
|
|
3251
|
con = 1.0 / this.bl;
|
|
3252
|
ts1 =Math.pow((this.el / Math.sqrt((1.0 + ul) / (1.0 - ul))),con);
|
|
3253
|
lat = Proj4js.common.phi2z(this.e,ts1);
|
|
3254
|
//if (flag != 0)
|
|
3255
|
//return(flag);
|
|
3256
|
//~ con = Math.cos(this.bl * us /al);
|
|
3257
|
theta = this.longc - Math.atan2((s * this.cosgam - vl * this.singam) , con)/this.bl;
|
|
3258
|
lon = Proj4js.common.adjust_lon(theta);
|
|
3259
|
}
|
|
3260
|
p.x=lon;
|
|
3261
|
p.y=lat;
|
|
3262
|
return p;
|
|
3263
|
}
|
|
3264
|
};
|
|
3265
|
/* ======================================================================
|
|
3266
|
projCode/lcc.js
|
|
3267
|
====================================================================== */
|
|
3268
|
|
|
3269
|
/*******************************************************************************
|
|
3270
|
NAME LAMBERT CONFORMAL CONIC
|
|
3271
|
|
|
3272
|
PURPOSE: Transforms input longitude and latitude to Easting and
|
|
3273
|
Northing for the Lambert Conformal Conic projection. The
|
|
3274
|
longitude and latitude must be in radians. The Easting
|
|
3275
|
and Northing values will be returned in meters.
|
|
3276
|
|
|
3277
|
|
|
3278
|
ALGORITHM REFERENCES
|
|
3279
|
|
|
3280
|
1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
|
|
3281
|
Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
|
|
3282
|
State Government Printing Office, Washington D.C., 1987.
|
|
3283
|
|
|
3284
|
2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections",
|
|
3285
|
U.S. Geological Survey Professional Paper 1453 , United State Government
|
|
3286
|
*******************************************************************************/
|
|
3287
|
|
|
3288
|
|
|
3289
|
//<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
|
|
3290
|
|
|
3291
|
// Initialize the Lambert Conformal conic projection
|
|
3292
|
// -----------------------------------------------------------------
|
|
3293
|
|
|
3294
|
//Proj4js.Proj.lcc = Class.create();
|
|
3295
|
Proj4js.Proj.lcc = {
|
|
3296
|
init : function() {
|
|
3297
|
|
|
3298
|
// array of: r_maj,r_min,lat1,lat2,c_lon,c_lat,false_east,false_north
|
|
3299
|
//double c_lat; /* center latitude */
|
|
3300
|
//double c_lon; /* center longitude */
|
|
3301
|
//double lat1; /* first standard parallel */
|
|
3302
|
//double lat2; /* second standard parallel */
|
|
3303
|
//double r_maj; /* major axis */
|
|
3304
|
//double r_min; /* minor axis */
|
|
3305
|
//double false_east; /* x offset in meters */
|
|
3306
|
//double false_north; /* y offset in meters */
|
|
3307
|
|
|
3308
|
if (!this.lat2){this.lat2=this.lat0;}//if lat2 is not defined
|
|
3309
|
if (!this.k0) this.k0 = 1.0;
|
|
3310
|
|
|
3311
|
// Standard Parallels cannot be equal and on opposite sides of the equator
|
|
3312
|
if (Math.abs(this.lat1+this.lat2) < Proj4js.common.EPSLN) {
|
|
3313
|
Proj4js.reportError("lcc:init: Equal Latitudes");
|
|
3314
|
return;
|
|
3315
|
}
|
|
3316
|
|
|
3317
|
var temp = this.b / this.a;
|
|
3318
|
this.e = Math.sqrt(1.0 - temp*temp);
|
|
3319
|
|
|
3320
|
var sin1 = Math.sin(this.lat1);
|
|
3321
|
var cos1 = Math.cos(this.lat1);
|
|
3322
|
var ms1 = Proj4js.common.msfnz(this.e, sin1, cos1);
|
|
3323
|
var ts1 = Proj4js.common.tsfnz(this.e, this.lat1, sin1);
|
|
3324
|
|
|
3325
|
var sin2 = Math.sin(this.lat2);
|
|
3326
|
var cos2 = Math.cos(this.lat2);
|
|
3327
|
var ms2 = Proj4js.common.msfnz(this.e, sin2, cos2);
|
|
3328
|
var ts2 = Proj4js.common.tsfnz(this.e, this.lat2, sin2);
|
|
3329
|
|
|
3330
|
var ts0 = Proj4js.common.tsfnz(this.e, this.lat0, Math.sin(this.lat0));
|
|
3331
|
|
|
3332
|
if (Math.abs(this.lat1 - this.lat2) > Proj4js.common.EPSLN) {
|
|
3333
|
this.ns = Math.log(ms1/ms2)/Math.log(ts1/ts2);
|
|
3334
|
} else {
|
|
3335
|
this.ns = sin1;
|
|
3336
|
}
|
|
3337
|
this.f0 = ms1 / (this.ns * Math.pow(ts1, this.ns));
|
|
3338
|
this.rh = this.a * this.f0 * Math.pow(ts0, this.ns);
|
|
3339
|
if (!this.title) this.title = "Lambert Conformal Conic";
|
|
3340
|
},
|
|
3341
|
|
|
3342
|
|
|
3343
|
// Lambert Conformal conic forward equations--mapping lat,long to x,y
|
|
3344
|
// -----------------------------------------------------------------
|
|
3345
|
forward : function(p) {
|
|
3346
|
|
|
3347
|
var lon = p.x;
|
|
3348
|
var lat = p.y;
|
|
3349
|
|
|
3350
|
// convert to radians
|
|
3351
|
if ( lat <= 90.0 && lat >= -90.0 && lon <= 180.0 && lon >= -180.0) {
|
|
3352
|
//lon = lon * Proj4js.common.D2R;
|
|
3353
|
//lat = lat * Proj4js.common.D2R;
|
|
3354
|
} else {
|
|
3355
|
Proj4js.reportError("lcc:forward: llInputOutOfRange: "+ lon +" : " + lat);
|
|
3356
|
return null;
|
|
3357
|
}
|
|
3358
|
|
|
3359
|
var con = Math.abs( Math.abs(lat) - Proj4js.common.HALF_PI);
|
|
3360
|
var ts;
|
|
3361
|
if (con > Proj4js.common.EPSLN) {
|
|
3362
|
ts = Proj4js.common.tsfnz(this.e, lat, Math.sin(lat) );
|
|
3363
|
rh1 = this.a * this.f0 * Math.pow(ts, this.ns);
|
|
3364
|
} else {
|
|
3365
|
con = lat * this.ns;
|
|
3366
|
if (con <= 0) {
|
|
3367
|
Proj4js.reportError("lcc:forward: No Projection");
|
|
3368
|
return null;
|
|
3369
|
}
|
|
3370
|
rh1 = 0;
|
|
3371
|
}
|
|
3372
|
var theta = this.ns * Proj4js.common.adjust_lon(lon - this.long0);
|
|
3373
|
p.x = this.k0 * (rh1 * Math.sin(theta)) + this.x0;
|
|
3374
|
p.y = this.k0 * (this.rh - rh1 * Math.cos(theta)) + this.y0;
|
|
3375
|
|
|
3376
|
return p;
|
|
3377
|
},
|
|
3378
|
|
|
3379
|
// Lambert Conformal Conic inverse equations--mapping x,y to lat/long
|
|
3380
|
// -----------------------------------------------------------------
|
|
3381
|
inverse : function(p) {
|
|
3382
|
|
|
3383
|
var rh1, con, ts;
|
|
3384
|
var lat, lon;
|
|
3385
|
x = (p.x - this.x0)/this.k0;
|
|
3386
|
y = (this.rh - (p.y - this.y0)/this.k0);
|
|
3387
|
if (this.ns > 0) {
|
|
3388
|
rh1 = Math.sqrt (x * x + y * y);
|
|
3389
|
con = 1.0;
|
|
3390
|
} else {
|
|
3391
|
rh1 = -Math.sqrt (x * x + y * y);
|
|
3392
|
con = -1.0;
|
|
3393
|
}
|
|
3394
|
var theta = 0.0;
|
|
3395
|
if (rh1 != 0) {
|
|
3396
|
theta = Math.atan2((con * x),(con * y));
|
|
3397
|
}
|
|
3398
|
if ((rh1 != 0) || (this.ns > 0.0)) {
|
|
3399
|
con = 1.0/this.ns;
|
|
3400
|
ts = Math.pow((rh1/(this.a * this.f0)), con);
|
|
3401
|
lat = Proj4js.common.phi2z(this.e, ts);
|
|
3402
|
if (lat == -9999) return null;
|
|
3403
|
} else {
|
|
3404
|
lat = -Proj4js.common.HALF_PI;
|
|
3405
|
}
|
|
3406
|
lon = Proj4js.common.adjust_lon(theta/this.ns + this.long0);
|
|
3407
|
|
|
3408
|
p.x = lon;
|
|
3409
|
p.y = lat;
|
|
3410
|
return p;
|
|
3411
|
}
|
|
3412
|
};
|
|
3413
|
|
|
3414
|
|
|
3415
|
|
|
3416
|
|
|
3417
|
/* ======================================================================
|
|
3418
|
projCode/laea.js
|
|
3419
|
====================================================================== */
|
|
3420
|
|
|
3421
|
/*******************************************************************************
|
|
3422
|
NAME LAMBERT AZIMUTHAL EQUAL-AREA
|
|
3423
|
|
|
3424
|
PURPOSE: Transforms input longitude and latitude to Easting and
|
|
3425
|
Northing for the Lambert Azimuthal Equal-Area projection. The
|
|
3426
|
longitude and latitude must be in radians. The Easting
|
|
3427
|
and Northing values will be returned in meters.
|
|
3428
|
|
|
3429
|
PROGRAMMER DATE
|
|
3430
|
---------- ----
|
|
3431
|
D. Steinwand, EROS March, 1991
|
|
3432
|
|
|
3433
|
This function was adapted from the Lambert Azimuthal Equal Area projection
|
|
3434
|
code (FORTRAN) in the General Cartographic Transformation Package software
|
|
3435
|
which is available from the U.S. Geological Survey National Mapping Division.
|
|
3436
|
|
|
3437
|
ALGORITHM REFERENCES
|
|
3438
|
|
|
3439
|
1. "New Equal-Area Map Projections for Noncircular Regions", John P. Snyder,
|
|
3440
|
The American Cartographer, Vol 15, No. 4, October 1988, pp. 341-355.
|
|
3441
|
|
|
3442
|
2. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
|
|
3443
|
Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
|
|
3444
|
State Government Printing Office, Washington D.C., 1987.
|
|
3445
|
|
|
3446
|
3. "Software Documentation for GCTP General Cartographic Transformation
|
|
3447
|
Package", U.S. Geological Survey National Mapping Division, May 1982.
|
|
3448
|
*******************************************************************************/
|
|
3449
|
|
|
3450
|
Proj4js.Proj.laea = {
|
|
3451
|
|
|
3452
|
|
|
3453
|
/* Initialize the Lambert Azimuthal Equal Area projection
|
|
3454
|
------------------------------------------------------*/
|
|
3455
|
init: function() {
|
|
3456
|
this.sin_lat_o=Math.sin(this.lat0);
|
|
3457
|
this.cos_lat_o=Math.cos(this.lat0);
|
|
3458
|
},
|
|
3459
|
|
|
3460
|
/* Lambert Azimuthal Equal Area forward equations--mapping lat,long to x,y
|
|
3461
|
-----------------------------------------------------------------------*/
|
|
3462
|
forward: function(p) {
|
|
3463
|
|
|
3464
|
/* Forward equations
|
|
3465
|
-----------------*/
|
|
3466
|
var lon=p.x;
|
|
3467
|
var lat=p.y;
|
|
3468
|
var delta_lon = Proj4js.common.adjust_lon(lon - this.long0);
|
|
3469
|
|
|
3470
|
//v 1.0
|
|
3471
|
var sin_lat=Math.sin(lat);
|
|
3472
|
var cos_lat=Math.cos(lat);
|
|
3473
|
|
|
3474
|
var sin_delta_lon=Math.sin(delta_lon);
|
|
3475
|
var cos_delta_lon=Math.cos(delta_lon);
|
|
3476
|
|
|
3477
|
var g =this.sin_lat_o * sin_lat +this.cos_lat_o * cos_lat * cos_delta_lon;
|
|
3478
|
if (g == -1.0) {
|
|
3479
|
Proj4js.reportError("laea:fwd:Point projects to a circle of radius "+ 2.0 * R);
|
|
3480
|
return null;
|
|
3481
|
}
|
|
3482
|
var ksp = this.a * Math.sqrt(2.0 / (1.0 + g));
|
|
3483
|
var x = ksp * cos_lat * sin_delta_lon + this.x0;
|
|
3484
|
var y = ksp * (this.cos_lat_o * sin_lat - this.sin_lat_o * cos_lat * cos_delta_lon) + this.x0;
|
|
3485
|
p.x = x;
|
|
3486
|
p.y = y
|
|
3487
|
return p;
|
|
3488
|
},//lamazFwd()
|
|
3489
|
|
|
3490
|
/* Inverse equations
|
|
3491
|
-----------------*/
|
|
3492
|
inverse: function(p) {
|
|
3493
|
p.x -= this.x0;
|
|
3494
|
p.y -= this.y0;
|
|
3495
|
|
|
3496
|
var Rh = Math.sqrt(p.x *p.x +p.y * p.y);
|
|
3497
|
var temp = Rh / (2.0 * this.a);
|
|
3498
|
|
|
3499
|
if (temp > 1) {
|
|
3500
|
Proj4js.reportError("laea:Inv:DataError");
|
|
3501
|
return null;
|
|
3502
|
}
|
|
3503
|
|
|
3504
|
var z = 2.0 * Proj4js.common.asinz(temp);
|
|
3505
|
var sin_z=Math.sin(z);
|
|
3506
|
var cos_z=Math.cos(z);
|
|
3507
|
|
|
3508
|
var lon =this.long0;
|
|
3509
|
if (Math.abs(Rh) > Proj4js.common.EPSLN) {
|
|
3510
|
var lat = Proj4js.common.asinz(this.sin_lat_o * cos_z +this. cos_lat_o * sin_z *p.y / Rh);
|
|
3511
|
var temp =Math.abs(this.lat0) - Proj4js.common.HALF_PI;
|
|
3512
|
if (Math.abs(temp) > Proj4js.common.EPSLN) {
|
|
3513
|
temp = cos_z -this.sin_lat_o * Math.sin(lat);
|
|
3514
|
if(temp!=0.0) lon=Proj4js.common.adjust_lon(this.long0+Math.atan2(p.x*sin_z*this.cos_lat_o,temp*Rh));
|
|
3515
|
} else if (this.lat0 < 0.0) {
|
|
3516
|
lon = Proj4js.common.adjust_lon(this.long0 - Math.atan2(-p.x,p.y));
|
|
3517
|
} else {
|
|
3518
|
lon = Proj4js.common.adjust_lon(this.long0 + Math.atan2(p.x, -p.y));
|
|
3519
|
}
|
|
3520
|
} else {
|
|
3521
|
lat = this.lat0;
|
|
3522
|
}
|
|
3523
|
//return(OK);
|
|
3524
|
p.x = lon;
|
|
3525
|
p.y = lat;
|
|
3526
|
return p;
|
|
3527
|
}//lamazInv()
|
|
3528
|
};
|
|
3529
|
|
|
3530
|
|
|
3531
|
|
|
3532
|
/* ======================================================================
|
|
3533
|
projCode/aeqd.js
|
|
3534
|
====================================================================== */
|
|
3535
|
|
|
3536
|
Proj4js.Proj.aeqd = {
|
|
3537
|
|
|
3538
|
init : function() {
|
|
3539
|
this.sin_p12=Math.sin(this.lat0)
|
|
3540
|
this.cos_p12=Math.cos(this.lat0)
|
|
3541
|
},
|
|
3542
|
|
|
3543
|
forward: function(p) {
|
|
3544
|
var lon=p.x;
|
|
3545
|
var lat=p.y;
|
|
3546
|
var ksp;
|
|
3547
|
|
|
3548
|
var sinphi=Math.sin(p.y);
|
|
3549
|
var cosphi=Math.cos(p.y);
|
|
3550
|
var dlon = Proj4js.common.adjust_lon(lon - this.long0);
|
|
3551
|
var coslon = Math.cos(dlon);
|
|
3552
|
var g = this.sin_p12 * sinphi + this.cos_p12 * cosphi * coslon;
|
|
3553
|
if (Math.abs(Math.abs(g) - 1.0) < Proj4js.common.EPSLN) {
|
|
3554
|
ksp = 1.0;
|
|
3555
|
if (g < 0.0) {
|
|
3556
|
Proj4js.reportError("aeqd:Fwd:PointError");
|
|
3557
|
return;
|
|
3558
|
}
|
|
3559
|
} else {
|
|
3560
|
var z = Math.acos(g);
|
|
3561
|
ksp = z/Math.sin(z);
|
|
3562
|
}
|
|
3563
|
p.x = this.x0 + this.a * ksp * cosphi * Math.sin(dlon);
|
|
3564
|
p.y = this.y0 + this.a * ksp * (this.cos_p12 * sinphi - this.sin_p12 * cosphi * coslon);
|
|
3565
|
return p;
|
|
3566
|
},
|
|
3567
|
|
|
3568
|
inverse: function(p){
|
|
3569
|
p.x -= this.x0;
|
|
3570
|
p.y -= this.y0;
|
|
3571
|
|
|
3572
|
var rh = Math.sqrt(p.x * p.x + p.y *p.y);
|
|
3573
|
if (rh > (2.0 * Proj4js.common.HALF_PI * this.a)) {
|
|
3574
|
Proj4js.reportError("aeqdInvDataError");
|
|
3575
|
return;
|
|
3576
|
}
|
|
3577
|
var z = rh / this.a;
|
|
3578
|
|
|
3579
|
var sinz=Math.sin(z)
|
|
3580
|
var cosz=Math.cos(z)
|
|
3581
|
|
|
3582
|
var lon = this.long0;
|
|
3583
|
var lat;
|
|
3584
|
if (Math.abs(rh) <= Proj4js.common.EPSLN) {
|
|
3585
|
lat = this.lat0;
|
|
3586
|
} else {
|
|
3587
|
lat = Proj4js.common.asinz(cosz * this.sin_p12 + (p.y * sinz * this.cos_p12) / rh);
|
|
3588
|
var con = Math.abs(this.lat0) - Proj4js.common.HALF_PI;
|
|
3589
|
if (Math.abs(con) <= Proj4js.common.EPSLN) {
|
|
3590
|
if (lat0 >= 0.0) {
|
|
3591
|
lon = Proj4js.common.adjust_lon(this.long0 + Math.atan2(p.x , -p.y));
|
|
3592
|
} else {
|
|
3593
|
lon = Proj4js.common.adjust_lon(this.long0 - Math.atan2(-p.x , p.y));
|
|
3594
|
}
|
|
3595
|
} else {
|
|
3596
|
con = cosz - this.sin_p12 * Math.sin(lat);
|
|
3597
|
if ((Math.abs(con) < Proj4js.common.EPSLN) && (Math.abs(p.x) < Proj4js.common.EPSLN)) {
|
|
3598
|
//no-op, just keep the lon value as is
|
|
3599
|
} else {
|
|
3600
|
var temp = Math.atan2((p.x * sinz * this.cos_p12), (con * rh));
|
|
3601
|
lon = Proj4js.common.adjust_lon(this.long0 + Math.atan2((p.x * sinz * this.cos_p12), (con * rh)));
|
|
3602
|
}
|
|
3603
|
}
|
|
3604
|
}
|
|
3605
|
|
|
3606
|
p.x = lon;
|
|
3607
|
p.y = lat;
|
|
3608
|
return p;
|
|
3609
|
}
|
|
3610
|
};
|
|
3611
|
/* ======================================================================
|
|
3612
|
projCode/moll.js
|
|
3613
|
====================================================================== */
|
|
3614
|
|
|
3615
|
/*******************************************************************************
|
|
3616
|
NAME MOLLWEIDE
|
|
3617
|
|
|
3618
|
PURPOSE: Transforms input longitude and latitude to Easting and
|
|
3619
|
Northing for the MOllweide projection. The
|
|
3620
|
longitude and latitude must be in radians. The Easting
|
|
3621
|
and Northing values will be returned in meters.
|
|
3622
|
|
|
3623
|
PROGRAMMER DATE
|
|
3624
|
---------- ----
|
|
3625
|
D. Steinwand, EROS May, 1991; Updated Sept, 1992; Updated Feb, 1993
|
|
3626
|
S. Nelson, EDC Jun, 2993; Made corrections in precision and
|
|
3627
|
number of iterations.
|
|
3628
|
|
|
3629
|
ALGORITHM REFERENCES
|
|
3630
|
|
|
3631
|
1. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections",
|
|
3632
|
U.S. Geological Survey Professional Paper 1453 , United State Government
|
|
3633
|
Printing Office, Washington D.C., 1989.
|
|
3634
|
|
|
3635
|
2. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
|
|
3636
|
Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
|
|
3637
|
State Government Printing Office, Washington D.C., 1987.
|
|
3638
|
*******************************************************************************/
|
|
3639
|
|
|
3640
|
Proj4js.Proj.moll = {
|
|
3641
|
|
|
3642
|
/* Initialize the Mollweide projection
|
|
3643
|
------------------------------------*/
|
|
3644
|
init: function(){
|
|
3645
|
//no-op
|
|
3646
|
},
|
|
3647
|
|
|
3648
|
/* Mollweide forward equations--mapping lat,long to x,y
|
|
3649
|
----------------------------------------------------*/
|
|
3650
|
forward: function(p) {
|
|
3651
|
|
|
3652
|
/* Forward equations
|
|
3653
|
-----------------*/
|
|
3654
|
var lon=p.x;
|
|
3655
|
var lat=p.y;
|
|
3656
|
|
|
3657
|
var delta_lon = Proj4js.common.adjust_lon(lon - this.long0);
|
|
3658
|
var theta = lat;
|
|
3659
|
var con = Proj4js.common.PI * Math.sin(lat);
|
|
3660
|
|
|
3661
|
/* Iterate using the Newton-Raphson method to find theta
|
|
3662
|
-----------------------------------------------------*/
|
|
3663
|
for (var i=0;;i++) {
|
|
3664
|
var delta_theta = -(theta + Math.sin(theta) - con)/ (1.0 + Math.cos(theta));
|
|
3665
|
theta += delta_theta;
|
|
3666
|
if (Math.abs(delta_theta) < Proj4js.common.EPSLN) break;
|
|
3667
|
if (i >= 50) {
|
|
3668
|
Proj4js.reportError("moll:Fwd:IterationError");
|
|
3669
|
//return(241);
|
|
3670
|
}
|
|
3671
|
}
|
|
3672
|
theta /= 2.0;
|
|
3673
|
|
|
3674
|
/* If the latitude is 90 deg, force the x coordinate to be "0 + false easting"
|
|
3675
|
this is done here because of precision problems with "cos(theta)"
|
|
3676
|
--------------------------------------------------------------------------*/
|
|
3677
|
if (Proj4js.common.PI/2 - Math.abs(lat) < Proj4js.common.EPSLN) delta_lon =0;
|
|
3678
|
var x = 0.900316316158 * this.R * delta_lon * Math.cos(theta) + this.x0;
|
|
3679
|
var y = 1.4142135623731 * this.R * Math.sin(theta) + this.y0;
|
|
3680
|
|
|
3681
|
p.x=x;
|
|
3682
|
p.y=y;
|
|
3683
|
return p;
|
|
3684
|
},
|
|
3685
|
|
|
3686
|
inverse: function(p){
|
|
3687
|
var theta;
|
|
3688
|
var arg;
|
|
3689
|
|
|
3690
|
/* Inverse equations
|
|
3691
|
-----------------*/
|
|
3692
|
p.x-= this.x0;
|
|
3693
|
//~ p.y -= this.y0;
|
|
3694
|
var arg = p.y / (1.4142135623731 * this.R);
|
|
3695
|
|
|
3696
|
/* Because of division by zero problems, 'arg' can not be 1.0. Therefore
|
|
3697
|
a number very close to one is used instead.
|
|
3698
|
-------------------------------------------------------------------*/
|
|
3699
|
if(Math.abs(arg) > 0.999999999999) arg=0.999999999999;
|
|
3700
|
var theta =Math.asin(arg);
|
|
3701
|
var lon = Proj4js.common.adjust_lon(this.long0 + (p.x / (0.900316316158 * this.R * Math.cos(theta))));
|
|
3702
|
if(lon < (-Proj4js.common.PI)) lon= -Proj4js.common.PI;
|
|
3703
|
if(lon > Proj4js.common.PI) lon= Proj4js.common.PI;
|
|
3704
|
arg = (2.0 * theta + Math.sin(2.0 * theta)) / Proj4js.common.PI;
|
|
3705
|
if(Math.abs(arg) > 1.0)arg=1.0;
|
|
3706
|
var lat = Math.asin(arg);
|
|
3707
|
//return(OK);
|
|
3708
|
|
|
3709
|
p.x=lon;
|
|
3710
|
p.y=lat;
|
|
3711
|
return p;
|
|
3712
|
}
|
|
3713
|
};
|
|
3714
|
|