| 1001 |
delphine |
1 |
<?php
|
|
|
2 |
/**
|
|
|
3 |
* Author : Julien Moquet
|
|
|
4 |
*
|
|
|
5 |
* Inspired by Proj4php from Mike Adair madairATdmsolutions.ca
|
|
|
6 |
* and Richard Greenwood rich@greenwoodma$p->com
|
|
|
7 |
* License: LGPL as per: http://www.gnu.org/copyleft/lesser.html
|
|
|
8 |
*/
|
|
|
9 |
/* * *****************************************************************************
|
|
|
10 |
NAME OBLIQUE MERCATOR (HOTINE)
|
|
|
11 |
|
|
|
12 |
PURPOSE: Transforms input longitude and latitude to Easting and
|
|
|
13 |
Northing for the Oblique Mercator projection. The
|
|
|
14 |
longitude and latitude must be in radians. The Easting
|
|
|
15 |
and Northing values will be returned in meters.
|
|
|
16 |
|
|
|
17 |
PROGRAMMER DATE
|
|
|
18 |
---------- ----
|
|
|
19 |
T. Mittan Mar, 1993
|
|
|
20 |
|
|
|
21 |
ALGORITHM REFERENCES
|
|
|
22 |
|
|
|
23 |
1. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
|
|
|
24 |
Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
|
|
|
25 |
State Government Printing Office, Washington D.C., 1987.
|
|
|
26 |
|
|
|
27 |
2. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections",
|
|
|
28 |
U.S. Geological Survey Professional Paper 1453 , United State Government
|
|
|
29 |
Printing Office, Washington D.C., 1989.
|
|
|
30 |
* ***************************************************************************** */
|
|
|
31 |
|
|
|
32 |
class Proj4phpProjOmerc {
|
|
|
33 |
/* Initialize the Oblique Mercator projection
|
|
|
34 |
------------------------------------------ */
|
|
|
35 |
|
|
|
36 |
public function init() {
|
|
|
37 |
if( !$this->mode )
|
|
|
38 |
$this->mode = 0;
|
|
|
39 |
if( !$this->lon1 ) {
|
|
|
40 |
$this->lon1 = 0;
|
|
|
41 |
$this->mode = 1;
|
|
|
42 |
}
|
|
|
43 |
if( !$this->lon2 )
|
|
|
44 |
$this->lon2 = 0;
|
|
|
45 |
if( !$this->lat2 )
|
|
|
46 |
$this->lat2 = 0;
|
|
|
47 |
|
|
|
48 |
/* Place parameters in static storage for common use
|
|
|
49 |
------------------------------------------------- */
|
|
|
50 |
$temp = $this->b / $this->a;
|
|
|
51 |
$es = 1.0 - pow( $temp, 2 );
|
|
|
52 |
$e = sqrt( $es );
|
|
|
53 |
|
|
|
54 |
$this->sin_p20 = sin( $this->lat0 );
|
|
|
55 |
$this->cos_p20 = cos( $this->lat0 );
|
|
|
56 |
|
|
|
57 |
$this->con = 1.0 - $this->es * $this->sin_p20 * $this->sin_p20;
|
|
|
58 |
$this->com = sqrt( 1.0 - $es );
|
|
|
59 |
$this->bl = sqrt( 1.0 + $this->es * pow( $this->cos_p20, 4.0 ) / (1.0 - $es) );
|
|
|
60 |
$this->al = $this->a * $this->bl * $this->k0 * $this->com / $this->con;
|
|
|
61 |
if( abs( $this->lat0 ) < Proj4php::$common->EPSLN ) {
|
|
|
62 |
$this->ts = 1.0;
|
|
|
63 |
$this->d = 1.0;
|
|
|
64 |
$this->el = 1.0;
|
|
|
65 |
} else {
|
|
|
66 |
$this->ts = Proj4php::$common->tsfnz( $this->e, $this->lat0, $this->sin_p20 );
|
|
|
67 |
$this->con = sqrt( $this->con );
|
|
|
68 |
$this->d = $this->bl * $this->com / ($this->cos_p20 * $this->con);
|
|
|
69 |
if( ($this->d * $this->d - 1.0) > 0.0 ) {
|
|
|
70 |
if( $this->lat0 >= 0.0 ) {
|
|
|
71 |
$this->f = $this->d + sqrt( $this->d * $this->d - 1.0 );
|
|
|
72 |
} else {
|
|
|
73 |
$this->f = $this->d - sqrt( $this->d * $this->d - 1.0 );
|
|
|
74 |
}
|
|
|
75 |
} else {
|
|
|
76 |
$this->f = $this->d;
|
|
|
77 |
}
|
|
|
78 |
$this->el = $this->f * pow( $this->ts, $this->bl );
|
|
|
79 |
}
|
|
|
80 |
|
|
|
81 |
//$this->longc=52.60353916666667;
|
|
|
82 |
|
|
|
83 |
if( $this->mode != 0 ) {
|
|
|
84 |
$this->g = .5 * ($this->f - 1.0 / $this->f);
|
|
|
85 |
$this->gama = Proj4php::$common->asinz( sin( $this->alpha ) / $this->d );
|
|
|
86 |
$this->longc = $this->longc - Proj4php::$common->asinz( $this->g * tan( $this->gama ) ) / $this->bl;
|
|
|
87 |
|
|
|
88 |
/* Report parameters common to format B
|
|
|
89 |
------------------------------------- */
|
|
|
90 |
//genrpt(azimuth * R2D,"Azimuth of Central Line: ");
|
|
|
91 |
//cenlon(lon_origin);
|
|
|
92 |
// cenlat(lat_origin);
|
|
|
93 |
|
|
|
94 |
$this->con = abs( $this->lat0 );
|
|
|
95 |
if( ($this->con > Proj4php::$common->EPSLN) && (abs( $this->con - Proj4php::$common->HALF_PI ) > Proj4php::$common->EPSLN) ) {
|
|
|
96 |
$this->singam = sin( $this->gama );
|
|
|
97 |
$this->cosgam = cos( $this->gama );
|
|
|
98 |
|
|
|
99 |
$this->sinaz = sin( $this->alpha );
|
|
|
100 |
$this->cosaz = cos( $this->alpha );
|
|
|
101 |
|
|
|
102 |
if( $this->lat0 >= 0 ) {
|
|
|
103 |
$this->u = ($this->al / $this->bl) * atan( sqrt( $this->d * $this->d - 1.0 ) / $this->cosaz );
|
|
|
104 |
} else {
|
|
|
105 |
$this->u = -($this->al / $this->bl) * atan( sqrt( $this->d * $this->d - 1.0 ) / $this->cosaz );
|
|
|
106 |
}
|
|
|
107 |
} else {
|
|
|
108 |
Proj4php::reportError( "omerc:Init:DataError" );
|
|
|
109 |
}
|
|
|
110 |
} else {
|
|
|
111 |
$this->sinphi = sin( $this->at1 );
|
|
|
112 |
$this->ts1 = Proj4php::$common->tsfnz( $this->e, $this->lat1, $this->sinphi );
|
|
|
113 |
$this->sinphi = sin( $this->lat2 );
|
|
|
114 |
$this->ts2 = Proj4php::$common->tsfnz( $this->e, $this->lat2, $this->sinphi );
|
|
|
115 |
$this->h = pow( $this->ts1, $this->bl );
|
|
|
116 |
$this->l = pow( $this->ts2, $this->bl );
|
|
|
117 |
$this->f = $this->el / $this->h;
|
|
|
118 |
$this->g = .5 * ($this->f - 1.0 / $this->f);
|
|
|
119 |
$this->j = ($this->el * $this->el - $this->l * $this->h) / ($this->el * $this->el + $this->l * $this->h);
|
|
|
120 |
$this->p = ($this->l - $this->h) / ($this->l + $this->h);
|
|
|
121 |
$this->dlon = $this->lon1 - $this->lon2;
|
|
|
122 |
if( $this->dlon < -Proj4php::$common->PI )
|
|
|
123 |
$this->lon2 = $this->lon2 - 2.0 * Proj4php::$common->PI;
|
|
|
124 |
if( $this->dlon > Proj4php::$common->PI )
|
|
|
125 |
$this->lon2 = $this->lon2 + 2.0 * Proj4php::$common->PI;
|
|
|
126 |
$this->dlon = $this->lon1 - $this->lon2;
|
|
|
127 |
$this->longc = .5 * ($this->lon1 + $this->lon2) - atan( $this->j * tan( .5 * $this->bl * $this->dlon ) / $this->p ) / $this->bl;
|
|
|
128 |
$this->dlon = Proj4php::$common->adjust_lon( $this->lon1 - $this->longc );
|
|
|
129 |
$this->gama = atan( sin( $this->bl * $this->dlon ) / $this->g );
|
|
|
130 |
$this->alpha = Proj4php::$common->asinz( $this->d * sin( $this->gama ) );
|
|
|
131 |
|
|
|
132 |
/* Report parameters common to format A
|
|
|
133 |
------------------------------------- */
|
|
|
134 |
if( abs( $this->lat1 - $this->lat2 ) <= Proj4php::$common->EPSLN ) {
|
|
|
135 |
Proj4php::reportError( "omercInitDataError" );
|
|
|
136 |
//return(202);
|
|
|
137 |
} else {
|
|
|
138 |
$this->con = abs( $this->lat1 );
|
|
|
139 |
}
|
|
|
140 |
if( ($this->con <= Proj4php::$common->EPSLN) || (abs( $this->con - Proj4php::$common->HALF_PI ) <= Proj4php::$common->EPSLN) ) {
|
|
|
141 |
Proj4php::reportError( "omercInitDataError" );
|
|
|
142 |
//return(202);
|
|
|
143 |
} else {
|
|
|
144 |
if( abs( abs( $this->lat0 ) - Proj4php::$common->HALF_PI ) <= Proj4php::$common->EPSLN ) {
|
|
|
145 |
Proj4php::reportError( "omercInitDataError" );
|
|
|
146 |
//return(202);
|
|
|
147 |
}
|
|
|
148 |
}
|
|
|
149 |
|
|
|
150 |
$this->singam = sin( $this->gam );
|
|
|
151 |
$this->cosgam = cos( $this->gam );
|
|
|
152 |
|
|
|
153 |
$this->sinaz = sin( $this->alpha );
|
|
|
154 |
$this->cosaz = cos( $this->alpha );
|
|
|
155 |
|
|
|
156 |
|
|
|
157 |
if( $this->lat0 >= 0 ) {
|
|
|
158 |
$this->u = ($this->al / $this->bl) * atan( sqrt( $this->d * $this->d - 1.0 ) / $this->cosaz );
|
|
|
159 |
} else {
|
|
|
160 |
$this->u = -($this->al / $this->bl) * atan( sqrt( $this->d * $this->d - 1.0 ) / $this->cosaz );
|
|
|
161 |
}
|
|
|
162 |
}
|
|
|
163 |
}
|
|
|
164 |
|
|
|
165 |
/* Oblique Mercator forward equations--mapping lat,long to x,y
|
|
|
166 |
---------------------------------------------------------- */
|
|
|
167 |
public function forward( $p ) {
|
|
|
168 |
|
|
|
169 |
/*
|
|
|
170 |
$theta; // angle
|
|
|
171 |
$sin_phi;
|
|
|
172 |
$cos_phi; // sin and cos value
|
|
|
173 |
$b; // temporary values
|
|
|
174 |
$c;
|
|
|
175 |
$t;
|
|
|
176 |
$tq; // temporary values
|
|
|
177 |
$con;
|
|
|
178 |
$n;
|
|
|
179 |
$ml; // cone constant, small m
|
|
|
180 |
$q;
|
|
|
181 |
$us;
|
|
|
182 |
$vl;
|
|
|
183 |
$ul;
|
|
|
184 |
$vs;
|
|
|
185 |
$s;
|
|
|
186 |
$dlon;
|
|
|
187 |
$ts1;
|
|
|
188 |
*/
|
|
|
189 |
|
|
|
190 |
$lon = $p->x;
|
|
|
191 |
$lat = $p->y;
|
|
|
192 |
|
|
|
193 |
/* Forward equations
|
|
|
194 |
----------------- */
|
|
|
195 |
$sin_phi = sin( $lat );
|
|
|
196 |
$dlon = Proj4php::$common->adjust_lon( $lon - $this->longc );
|
|
|
197 |
$vl = sin( $this->bl * $dlon );
|
|
|
198 |
if( abs( abs( $lat ) - Proj4php::$common->HALF_PI ) > Proj4php::$common->EPSLN ) {
|
|
|
199 |
$ts1 = Proj4php::$common->tsfnz( $this->e, $lat, $sin_phi );
|
|
|
200 |
$q = $this->el / (pow( $ts1, $this->bl ));
|
|
|
201 |
$s = .5 * ($q - 1.0 / $q);
|
|
|
202 |
$t = .5 * ($q + 1.0 / $q);
|
|
|
203 |
$ul = ($s * $this->singam - $vl * $this->cosgam) / $t;
|
|
|
204 |
$con = cos( $this->bl * $dlon );
|
|
|
205 |
if( abs( con ) < .0000001 ) {
|
|
|
206 |
$us = $this->al * $this->bl * $dlon;
|
|
|
207 |
} else {
|
|
|
208 |
$us = $this->al * atan( ($s * $this->cosgam + $vl * $this->singam) / $con ) / $this->bl;
|
|
|
209 |
if( $con < 0 )
|
|
|
210 |
$us = $us + Proj4php::$common->PI * $this->al / $this->bl;
|
|
|
211 |
}
|
|
|
212 |
} else {
|
|
|
213 |
if( $lat >= 0 ) {
|
|
|
214 |
$ul = $this->singam;
|
|
|
215 |
} else {
|
|
|
216 |
$ul = -$this->singam;
|
|
|
217 |
}
|
|
|
218 |
$us = $this->al * $lat / $this->bl;
|
|
|
219 |
}
|
|
|
220 |
if( abs( abs( $ul ) - 1.0 ) <= Proj4php::$common->EPSLN ) {
|
|
|
221 |
//alert("Point projects into infinity","omer-for");
|
|
|
222 |
Proj4php::reportError( "omercFwdInfinity" );
|
|
|
223 |
//return(205);
|
|
|
224 |
}
|
|
|
225 |
$vs = .5 * $this->al * log( (1.0 - $ul) / (1.0 + $ul) ) / $this->bl;
|
|
|
226 |
$us = $us - $this->u;
|
|
|
227 |
$p->x = $this->x0 + $vs * $this->cosaz + $us * $this->sinaz;
|
|
|
228 |
$p->y = $this->y0 + $us * $this->cosaz - $vs * $this->sinaz;
|
|
|
229 |
|
|
|
230 |
return $p;
|
|
|
231 |
}
|
|
|
232 |
|
|
|
233 |
/**
|
|
|
234 |
*
|
|
|
235 |
* @param type $p
|
|
|
236 |
* @return type
|
|
|
237 |
*/
|
|
|
238 |
public function inverse( $p ) {
|
|
|
239 |
/*
|
|
|
240 |
$delta_lon; /* Delta longitude (Given longitude - center
|
|
|
241 |
$theta; /* angle
|
|
|
242 |
$delta_theta; /* adjusted longitude
|
|
|
243 |
$sin_phi;
|
|
|
244 |
$cos_phi; /* sin and cos value
|
|
|
245 |
$b; /* temporary values
|
|
|
246 |
$c;
|
|
|
247 |
$t;
|
|
|
248 |
$tq; /* temporary values
|
|
|
249 |
$con;
|
|
|
250 |
$n;
|
|
|
251 |
$ml; /* cone constant, small m
|
|
|
252 |
$vs;
|
|
|
253 |
$us;
|
|
|
254 |
$q;
|
|
|
255 |
$s;
|
|
|
256 |
$ts1;
|
|
|
257 |
$vl;
|
|
|
258 |
$ul;
|
|
|
259 |
$bs;
|
|
|
260 |
$dlon;
|
|
|
261 |
$flag;
|
|
|
262 |
*/
|
|
|
263 |
|
|
|
264 |
/* Inverse equations
|
|
|
265 |
----------------- */
|
|
|
266 |
$p->x -= $this->x0;
|
|
|
267 |
$p->y -= $this->y0;
|
|
|
268 |
#$flag = 0;
|
|
|
269 |
$vs = $p->x * $this->cosaz - $p->y * $this->sinaz;
|
|
|
270 |
$us = $p->y * $this->cosaz + $p->x * $this->sinaz;
|
|
|
271 |
$us = $us + $this->u;
|
|
|
272 |
$q = exp( -$this->bl * $vs / $this->al );
|
|
|
273 |
$s = .5 * ($q - 1.0 / $q);
|
|
|
274 |
$t = .5 * ($q + 1.0 / $q);
|
|
|
275 |
$vl = sin( $this->bl * $us / $this->al );
|
|
|
276 |
$ul = ($vl * $this->cosgam + $s * $this->singam) / $t;
|
|
|
277 |
if( abs( abs( $ul ) - 1.0 ) <= Proj4php::$common->EPSLN ) {
|
|
|
278 |
$lon = $this->longc;
|
|
|
279 |
if( ul >= 0.0 ) {
|
|
|
280 |
$lat = Proj4php::$common->HALF_PI;
|
|
|
281 |
} else {
|
|
|
282 |
$lat = -Proj4php::$common->HALF_PI;
|
|
|
283 |
}
|
|
|
284 |
} else {
|
|
|
285 |
$con = 1.0 / $this->bl;
|
|
|
286 |
$ts1 = pow( ($this->el / sqrt( (1.0 + $ul) / (1.0 - $ul) ) ), $con );
|
|
|
287 |
$lat = Proj4php::$common->phi2z( $this->e, $ts1 );
|
|
|
288 |
//if ($flag != 0)
|
|
|
289 |
//return($flag);
|
|
|
290 |
//~ con = cos($this->bl * us /al);
|
|
|
291 |
$theta = $this->longc - atan2( ($s * $this->cosgam - $vl * $this->singam ), $con ) / $this->bl;
|
|
|
292 |
$lon = Proj4php::$common->adjust_lon( $theta );
|
|
|
293 |
}
|
|
|
294 |
$p->x = $lon;
|
|
|
295 |
$p->y = $lat;
|
|
|
296 |
return $p;
|
|
|
297 |
}
|
|
|
298 |
|
|
|
299 |
}
|
|
|
300 |
|
|
|
301 |
Proj4php::$proj['omerc'] = new Proj4phpProjOmerc();
|