Skip to content

Instantly share code, notes, and snippets.

@martinsbalodis
Created February 2, 2012 19:51
Show Gist options
  • Save martinsbalodis/1725377 to your computer and use it in GitHub Desktop.
Save martinsbalodis/1725377 to your computer and use it in GitHub Desktop.
public function LKS2LatLong ($east, $north)
{
$K0 = 0.9996;
$FALSE_E = 500000.0;
$FALSE_N = 0;
$f = 1/298.257223563;//298.257222101;//298.257223563
$a = 6378137.0;
$e = sqrt(2.0 * $f - $f * $f);
$e_dash = sqrt($e*$e + pow($e,4)/(1.0 - $e * $e));
$north_d = $north - $FALSE_N;
$east_d = $east - $FALSE_E;
$m = $north_d/$K0;
$A0 = 1.0-(($e*$e)/4.0)-((3.0*pow($e,4))/64.0)-((5.0*pow($e,6))/256.0);
$A2 = (3.0/8.0)*($e*$e+Pow($e,4)/4.0+(15.0*Pow($e,6))/128.0);
$A4 = (15.0/256.0)*(Pow($e,4)+(3.0*Pow($e,6))/4.0);
$A6 = (35.0*Pow($e,6))/3072.0;
$clat = $m/($a*$A0);
$x1 = $a*$A0*$clat-$a*$A2*sin(2.0*$clat)+$a*$A4*sin(4.0*$clat)-$a*$A6*sin(6.0*$clat);
//echo abs($m - $x);
while (abs($m - $x1)>0.0001) {
#echo $i++;
$clat = $clat+($m-$x1)/30.81851/3600*pi()/180.0;
$x1 = $a*$A0*$clat-$a*$A2*sin(2.0*$clat)+$a*$A4*sin(4.0*$clat)-$a*$A6*sin(6.0*$clat);
}
$cent_mer = 24*pi()/180.0;
$c=$a/sqrt(1.0-$e*$e);
$cos_clat=cos($clat);
$V=sqrt(1.0+Pow($e_dash,2)*Pow($cos_clat,2));
$p=$c/Pow($V,3);
$v1=$c/$V;
$t=tan($clat);
$vp=$v1/$p;
$Result = array();
$A0 = $clat-$t/($K0*$p)*pow($east_d,2)/(2.0*$K0*$v1);
$A2 = $t/($K0*$p)*pow($east_d,4)/(24.0*pow($K0,3)*pow($v1,3))*(-4.0*pow($vp,3)+9.0*$vp*(1-$t*$t)+12.0*$t*$t);
$A4 = $t/($K0*$p)*pow($east_d,6)/(720.0*pow($K0,5)*pow($v1,5))*
(8.0*pow($vp,4)*(11.0-24.0*$t*$t)-12.0*pow($vp,3)*(21.0-71.0*$t*$t)+
15.0*$vp*$vp*(15.0-98.0*$t*$t+15.0*pow($t,4))+
180.0*$vp*(5.0*$t*$t-3.0*pow($t,4))+360.0*pow($t,4));
$A6 = $t/($K0*$p)*pow($east_d,8)/(40320.0*pow($K0,7)*pow($v1,7))*(1385.0+
3633.0*$t*$t+4095.0*pow($t,4)+1575.0*pow($t,6));
$Result['x'] = ($A0+$A2+$A4+$A6)*180.0/pi();
$x1 = 1.0/$cos_clat;
$A0 = $x1*$east_d/($K0*$v1);
$A2 = $x1*pow($east_d,3)/(6.0*pow($K0,3)*pow($v1,3))*($vp+2.0*$t*$t);
$A4 = $x1*pow($east_d,5)/(120.0*pow($K0,5)*pow($v1,5))*(-4.0*pow($vp,3)*(1.0-6.0*$t*$t)+$vp*$vp*(9.0-68.0*$t*$t)+
72.0*$vp*$t*$t+24.0*pow($t,4));
$A6 = $x1*pow($east_d,7)/(5040.0*pow($K0,7)*pow($v1,7))*(61.0+662.0*$t*$t+1320.0*pow($t,4)+720.0*pow($t,6));
$Result['y'] = ($cent_mer+$A0-$A2+$A4-$A6)*180.0/pi();
return $Result;
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment