| 1001 |
delphine |
1 |
<?php
|
|
|
2 |
|
|
|
3 |
/*******************************************************************************
|
|
|
4 |
NAME CASSINI
|
|
|
5 |
|
|
|
6 |
PURPOSE: Transforms input longitude and latitude to Easting and
|
|
|
7 |
Northing for the Cassini projection. The
|
|
|
8 |
longitude and latitude must be in radians. The Easting
|
|
|
9 |
and Northing values will be returned in meters.
|
|
|
10 |
Ported from PROJ.4.
|
|
|
11 |
|
|
|
12 |
|
|
|
13 |
ALGORITHM REFERENCES
|
|
|
14 |
|
|
|
15 |
1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
|
|
|
16 |
Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
|
|
|
17 |
State Government Printing Office, Washington D.C., 1987.
|
|
|
18 |
|
|
|
19 |
2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections",
|
|
|
20 |
U.S. Geological Survey Professional Paper 1453 , United State Government
|
|
|
21 |
****************************************************************************** */
|
|
|
22 |
|
|
|
23 |
/**
|
|
|
24 |
* Author : Julien Moquet
|
|
|
25 |
*
|
|
|
26 |
* Inspired by Proj4php from Mike Adair madairATdmsolutions.ca
|
|
|
27 |
* and Richard Greenwood rich@greenwoodma$p->com
|
|
|
28 |
* License: LGPL as per: http://www.gnu.org/copyleft/lesser.html
|
|
|
29 |
*/
|
|
|
30 |
//Proj4php.defs["EPSG:28191"] = "+proj=cass +lat_0=31.73409694444445 +lon_0=35.21208055555556 +x_0=170251.555 +y_0=126867.909 +a=6378300.789 +b=6356566.435 +towgs84=-275.722,94.7824,340.894,-8.001,-4.42,-11.821,1 +units=m +no_defs";
|
|
|
31 |
// Initialize the Cassini projection
|
|
|
32 |
// -----------------------------------------------------------------
|
|
|
33 |
|
|
|
34 |
class Proj4phpProjCass {
|
|
|
35 |
|
|
|
36 |
public function init() {
|
|
|
37 |
if( !$this->sphere ) {
|
|
|
38 |
$this->en = Proj4php::$common->pj_enfn( $this->es );
|
|
|
39 |
$this->m0 = Proj4php::$common->pj_mlfn( $this->lat0, sin( $this->lat0 ), cos( $this->lat0 ), $this->en );
|
|
|
40 |
}
|
|
|
41 |
}
|
|
|
42 |
|
|
|
43 |
protected $C1 = .16666666666666666666;
|
|
|
44 |
protected $C2 = .00833333333333333333;
|
|
|
45 |
protected $C3 = .04166666666666666666;
|
|
|
46 |
protected $C4 = .33333333333333333333;
|
|
|
47 |
protected $C5 = .06666666666666666666;
|
|
|
48 |
|
|
|
49 |
/* Cassini forward equations--mapping lat,long to x,y
|
|
|
50 |
----------------------------------------------------------------------- */
|
|
|
51 |
public function forward( $p ) {
|
|
|
52 |
|
|
|
53 |
/* Forward equations
|
|
|
54 |
----------------- */
|
|
|
55 |
#$x;
|
|
|
56 |
#$y;
|
|
|
57 |
$lam = $p->x;
|
|
|
58 |
$phi = $p->y;
|
|
|
59 |
$lam = Proj4php::$common->adjust_lon( $lam - $this->long0 );
|
|
|
60 |
|
|
|
61 |
if( $this->sphere ) {
|
|
|
62 |
$x = asin( cos( $phi ) * sin( $lam ) );
|
|
|
63 |
$y = atan2( tan( $phi ), cos( $lam ) ) - $this->phi0;
|
|
|
64 |
} else {
|
|
|
65 |
//ellipsoid
|
|
|
66 |
$this->n = sin( $phi );
|
|
|
67 |
$this->c = cos( $phi );
|
|
|
68 |
$y = $this->pj_mlfn( $phi, $this->n, $this->c, $this->en );
|
|
|
69 |
$this->n = 1. / sqrt( 1. - $this->es * $this->n * $this->n );
|
|
|
70 |
$this->tn = tan( $phi );
|
|
|
71 |
$this->t = $this->tn * $this->tn;
|
|
|
72 |
$this->a1 = $lam * $this->c;
|
|
|
73 |
$this->c *= $this->es * $this->c / (1 - $this->es);
|
|
|
74 |
$this->a2 = $this->a1 * $this->a1;
|
|
|
75 |
$x = $this->n * $this->a1 * (1. - $this->a2 * $this->t * ($this->C1 - (8. - $this->t + 8. * $this->c) * $this->a2 * $this->C2));
|
|
|
76 |
$y -= $this->m0 - $this->n * $this->tn * $this->a2 * (.5 + (5. - $this->t + 6. * $this->c) * $this->a2 * $this->C3);
|
|
|
77 |
}
|
|
|
78 |
|
|
|
79 |
$p->x = $this->a * $x + $this->x0;
|
|
|
80 |
$p->y = $this->a * $y + $this->y0;
|
|
|
81 |
|
|
|
82 |
return $p;
|
|
|
83 |
}
|
|
|
84 |
|
|
|
85 |
/* Inverse equations
|
|
|
86 |
----------------- */
|
|
|
87 |
public function inverse( $p ) {
|
|
|
88 |
$p->x -= $this->x0;
|
|
|
89 |
$p->y -= $this->y0;
|
|
|
90 |
$x = $p->x / $this->a;
|
|
|
91 |
$y = $p->y / $this->a;
|
|
|
92 |
|
|
|
93 |
if( $this->sphere ) {
|
|
|
94 |
$this->dd = $y + $this->lat0;
|
|
|
95 |
$phi = asin( sin( $this->dd ) * cos( $x ) );
|
|
|
96 |
$lam = atan2( tan( $x ), cos( $this->dd ) );
|
|
|
97 |
} else {
|
|
|
98 |
/* ellipsoid */
|
|
|
99 |
$ph1 = Proj4php::$common->pj_inv_mlfn( $this->m0 + $y, $this->es, $this->en );
|
|
|
100 |
$this->tn = tan( $ph1 );
|
|
|
101 |
$this->t = $this->tn * $this->tn;
|
|
|
102 |
$this->n = sin( $ph1 );
|
|
|
103 |
$this->r = 1. / (1. - $this->es * $this->n * $this->n);
|
|
|
104 |
$this->n = sqrt( $this->r );
|
|
|
105 |
$this->r *= (1. - $this->es) * $this->n;
|
|
|
106 |
$this->dd = $x / $this->n;
|
|
|
107 |
$this->d2 = $this->dd * $this->dd;
|
|
|
108 |
$phi = $ph1 - ($this->n * $this->tn / $this->r) * $this->d2 * (.5 - (1. + 3. * $this->t) * $this->d2 * $this->C3);
|
|
|
109 |
$lam = $this->dd * (1. + $this->t * $this->d2 * (-$this->C4 + (1. + 3. * $this->t) * $this->d2 * $this->C5)) / cos( $ph1 );
|
|
|
110 |
}
|
|
|
111 |
$p->x = Proj4php::$common->adjust_lon( $this->long0 + $lam );
|
|
|
112 |
$p->y = $phi;
|
|
|
113 |
|
|
|
114 |
return $p;
|
|
|
115 |
}
|
|
|
116 |
}
|
|
|
117 |
|
|
|
118 |
Proj4php::$proj['cass'] = new Proj4phpProjCass();
|