#!/usr/bin/perl -W # use Image::Magick; use POSIX; # Version 1.0 21/12/2006 # Version 1.1 14/3/2007 web form use CGI; ##vars # my $maj1830=6377563.396; # Airy1830 my $min1830=6356256.909; my $maj84=6378137.0; # GRS80 aka WGS84 my $min84=6356752.3141; my $E1830=0.00667054007414908; # ref ellipsoid 1830 my $E84=0.00669438003551284; # ref ellipsoid 84 my $pi=3.14159265358979; my @ns=("S","N"); my @ew=("W","E"); my $mapgen="lmap.htm"; my $mapjoin="lmap.bat"; my $mapinfo="lmap.inf"; my $leadin="montage -geometry +0+0 -tile "; # 13x8 NT07NE08.gif my $leadout=" lmap.gif"; my ($topleft,$g); my ($East,$North); my ($Lat,$Long); my ($UTMEasting, $UTMNorthing, $UTMZone, $longitudeZone); my $gifs=""; my $mapname="Map Name"; # ##eovars my $ostl='NN200300'; my $osbr='NT500300'; my $errortext=''; ##check incoming, make ostl and osbr be what came in subject to error checks $query = new CGI; $ostlin = $query->param('ostlin'); if ($ostlin) { if ($ostlin =~/[A-Z]{2}\d{6}/) { $ostl=$ostlin; } else { $errortext.="Expecting OSGB top left in form NT230450!
"; } } $osbrin = $query->param('osbrin'); if ($osbrin) { if ($osbrin =~/[A-Z]{2}\d{6}/) { $osbr=$osbrin; } else { $errortext.="Expecting OSGB bottom right in form NT450230!
"; } } my $header=" OSGB to Lat/Long for UI-View32 Maps
Back Home -> Back to Project Index -> Back to APRS Index http://graham.auld.me.uk

Convert OSGB to Decimal Lat & Lon

\n"; my $footer="
Author: osgb2ll.pl v1.1
Dated: 14th March 2007\n\n"; print "Content-type: text/html\nCache-Control: max-age=1\n\n"; print "\n"; print $header; if ($errortext ne '') { print "Error!
$errortext
\n"; } print "
\n"; print "
Input OS Grid Corners Output for UI-View32 inf file
"; #input side print "

OSGB Top Left
OSGB Bottom Right
"; print "		"; #output side print "
"; print "mapname.inf:
\n"; print ""; print "
"; print "
Want to find the OSGB Grid reference of your map image?

Then go to http://www.streetmap.co.uk, find the map area of interest which matches the extent of your own image, clicking on the required position, which will then be marked by an orange coloured arrow. Then look for \"Click here to convert/measure coordinates\" The link will take you to a Grid Conversion Results page. The 'LR' row in the table will provide you with the OS Grid reference of your selected position. Enter the OSGB coordinates into the the appropriate corner of the map above. After pressing the convert button, a fully formatted UI-View32 inf file will be displayed on the right hand side. Just copy and paste this into a .inf file in the UI-View32 map area, along with a map image file of the same name."; print "

Note: You are fitting a 2-dimentional map to the 3-dimentional curvature of the earth. On small maps there is little problem but on larger maps, the deviation becomes more obvious. Thus while the identified corners and its major diagonal will be correct, the more distant areas accumulate a greater degree of error. Adjustments may be required to reduce the errors in the primary area of interest."; print $footer; ###theend### exit(); ############below here is conversion stuff################# # call like so to generate some text as should appear in .inf file... # genmap ("NN200300","NT500300"); # # my $pos="NT090790"; #$pos="NT500300"; #$pos="NT080770"; #getOSRefFromSixFigureReference($pos); #print "Ref $East, $North\n"; #$ngr=toSixFigureString($East,$North); #print "$pos = $ngr\n"; #exit(0); #testpos("NT080770"); # eastings 308000, 677000 - 55.58.632N 3.26.467W for APRS mapping #exit(0); #genmap ("NN200300","NT500300"); # # Given the 6-figure NGR reference of the top left and bottom right squares, it will generate a 'lmap.htm' page # which allows all the 'Streetmap.co.uk' 1k x 1k images to be displayed in a large table. # It will then generate a 'lmap.bat' file which can subsequentially run to stitch all the 1k x 1k images together # into one large 'lmap.gif' but all the required images must be available in the local folder first, hence check the HTM display. # A 'lmap.inf' file will be generated which, after name editing, can be used with UI-View32. ## sub genmap { my($fromref,$toref)=@_; getOSRefFromSixFigureReference($toref); # bottom right of area my $gx1=$East; my $gy1=$North; getOSRefFromSixFigureReference($fromref); # top left of area my $gx=$East; my $gy=$North; # print "$fromref to $toref - ($gx,$gy) to ($gx1,$gy1)\n"; # #conversion from... # open (H, ">$mapgen") or die "Can't open $mapgen $!"; # gen web file to check map # print H "\n\n"; # my $gifs=''; # for ($j=$gy;$j>=$gy1;$j-=10000) { # top left down # print H ""; # for ($i=$gx;$i<=$gx1;$i+=10000) { # top left right #print "$i,$j "; # $ngr=toSixFigureString($i,$j); #my $s1=substr($ngr, 0, 2); # my $s2=substr($ngr, 2, 1); # my $s2a=$s2; # my $s3=substr($ngr, 5, 1); # my $smngr=$s1.$s2.$s3; # print "$smngr "; # print H "\n"; # } # print H "
"; # $gifs.=' ..\\squares\\'.$smngr.'.gif'; # } # print "\n"; # print H "
\n\n"; # close(H); # my $hor=($gx1-$gx)/10000+1; # my $ver=($gy-$gy1)/10000+1; # open (B, ">$mapjoin") or die "Can't open $mapjoin $!"; # batch file to gen full map # print B $leadin.$hor.'x'.$ver.$gifs.$leadout."\n"; # close(B); #my $m1=OSreftoDEClatlong($fromref,0,10000-($ver*0.5)); # middle left #my $m2=OSreftoDEClatlong($fromref,($hor*0.5),10000); # middle top #my $m3=OSreftoDEClatlong($toref,10000,($ver*0.5)); # middle right #my $m4=OSreftoDEClatlong($toref,10000-($hor*0.5),0); # middle bottom # eastings are ~100m out $s4=OSreftoDEClatlong($fromref,0,10000); # move to top left of first square toOSRef($Lat,$Long); # Convert to Lat/Long $East=int($East+0.5); $North=int($North+0.5); $s6=toSixFigureString($East,$North); #print "$s6 ($East,$North)"; $s5=OSreftoDEClatlong($toref,10000,0); # move bottom right of last square toOSRef($Lat,$Long); # Convert to Lat/Long $East=int($East+0.5); $North=int($North+0.5); $s6=toSixFigureString(($East),$North); #print " - $s6 ($East,$North)\n"; #print "$s4 - $s5\n"; #open (R, ">$mapinfo") or die "Can't open $mapinfo $!"; # APRS map info # print R "$s4\n$s5\nNew LMap\n"; print "$s4\n$s5\n$mapname\n"; #close(R); # print $footer; } # # Position check routine # sub testpos { my ($pos)=@_; $pos=OSreftoDEClatlong($pos,0,0); printf "$pos for APRS mapping\n"; } # ------------------------------ sub deg2rad { my ($deg)=@_; return $deg*$pi/180 } sub rad2deg { my ($rad)=@_; return $rad*180/$pi } sub sinSquared { my ($x)=@_; return sin($x)*sin($x) } sub cosSquared { my ($x)=@_; return cos($x)*cos($x) } sub tanSquared { my ($x)=@_; return tan($x)*tan($x) } sub sec { my ($x)=@_; return 1.0 / cos($x) } # # numeric latitide, longditude to decimal minutes for APRS # sub OSreftoDEClatlong { my ($pos,$xdelta,$ydelta)=@_; getOSRefFromSixFigureReference($pos); # Alpha-numeric pos to fully numeric $East+=$xdelta; $North+=$ydelta; gridtoLatLng($East,$North); # Convert to Lat/Long #OSGB36ToWGS84($Lat,$Long); # Geodetic Transform GB National grid to World Geodetic system my $NS='N'; my $EW; my $lo2; my $la1=int($Lat); my $la2=($Lat-$la1)*60; if ($Long<0.0) { $EW='W'; $lo2=-$Long; } else { $EW='E'; $lo2=$Long; } my $lo1=int($lo2); $lo2=($lo2-$lo1)*60; $lo1= sprintf "%d.%.3f%s,%d.%.3f%s",$la1,$la2,$NS,$lo1,$lo2,$EW; return $lo1; } # Convert this Latitude Longditude from OSGB36 datum to WGS84 datum. # input numeric (latitude, longditude) # result to global $Lat, $Long # sub OSGB36ToWGS84 { my ($lat,$lng)=@_; my $a = $maj1830; my $b = $min1830; my $eSquared = (($a * $a) - ($b * $b)) / ($a * $a); $phi = deg2rad($lat); $lambda = deg2rad($lng); $v = $a / (sqrt(1 - $eSquared * sinSquared($phi))); $H = 0; # height my $x = ($v + $H) * cos($phi) * cos($lambda); my $y = ($v + $H) * cos($phi) * sin($lambda); my $z = ((1 - $eSquared) * $v + $H) * sin($phi); my $tx = 446.448; # my $ty = -124.157; my $ty = -125.157; my $tz = 542.060; my $s = -0.0000204894; my $rx = deg2rad( 0.00004172222); my $ry = deg2rad( 0.00006861111); my $rz = deg2rad( 0.00023391666); my $xB = $tx + ($x * (1 + $s)) + (-$rx * $y) + ($ry * $z); my $yB = $ty + ($rz * $x) + ($y * (1 + $s)) + (-$rx * $z); my $zB = $tz + (-$ry * $x) + ($rx * $y) + ($z * (1 + $s)); $a = $maj84; $b = $min84; $eSquared = (($a * $a) - ($b * $b)) / ($a * $a); my $lambdaB = rad2deg(atan($yB / $xB)); my $p = sqrt(($xB * $xB) + ($yB * $yB)); my $phiN = atan($zB / ($p * (1 - $eSquared))); for ($i = 1; $i < 10; $i++) { my $v = $a / (sqrt(1 - $eSquared * sinSquared($phiN))); $phiN1 = atan(($zB + ($eSquared * $v * sin($phiN))) / $p); $phiN = $phiN1; } my $phiB = rad2deg($phiN); $Lat = $phiB; $Long = $lambdaB; } # # Convert this grid reference into a latitude and longitude # input numeric (Eastings, Northings) # result to global $Lat, $Long # sub gridtoLatLng { my ($E,$N)=@_; my $OSGB_F0 = 0.9996012717; # National grid my $N0 = -100000.0; my $E0 = 400000.0; my $phi0 = deg2rad(49.0); my $lambda0 = deg2rad(-2.0); my $a = $maj1830; my $b = $min1830; my $eSquared = (($a * $a) - ($b * $b)) / ($a * $a); my $phi = 0.0; my $lambda = 0.0; my $n = ($a - $b) / ($a + $b); my $M = 0.0; my $phiPrime = (($N - $N0) / ($a * $OSGB_F0)) + $phi0; do { $M = ($b * $OSGB_F0) * (((1 + $n + ((5.0 / 4.0) * $n * $n) + ((5.0 / 4.0) * $n * $n * $n)) * ($phiPrime - $phi0)) - (((3 * $n) + (3 * $n * $n) + ((21.0 / 8.0) * $n * $n * $n)) * sin($phiPrime - $phi0) * cos($phiPrime + $phi0)) + ((((15.0 / 8.0) * $n * $n) + ((15.0 / 8.0) * $n * $n * $n)) * sin(2.0 * ($phiPrime - $phi0)) * cos(2.0 * ($phiPrime + $phi0))) - (((35.0 / 24.0) * $n * $n * $n) * sin(3.0 * ($phiPrime - $phi0)) * cos(3.0 * ($phiPrime + $phi0)))); $phiPrime += ($N - $N0 - $M) / ($a * $OSGB_F0); } while (($N - $N0 - $M) >= 0.001); my $v = $a * $OSGB_F0 * pow(1.0 - $eSquared * sinSquared($phiPrime), -0.5); my $rho = $a * $OSGB_F0 * (1.0 - $eSquared) * pow(1.0 - $eSquared * sinSquared($phiPrime), -1.5); my $etaSquared = ($v / $rho) - 1.0; my $VII = tan($phiPrime) / (2 * $rho * $v); my $VIII = (tan($phiPrime) / (24.0 * $rho * pow($v, 3.0))) * (5.0 + (3.0 * tanSquared($phiPrime)) + $etaSquared - (9.0 * tanSquared($phiPrime) * $etaSquared)); my $IX = (tan($phiPrime) / (720.0 * $rho * pow($v, 5.0))) * (61.0 + (90.0 * tanSquared($phiPrime)) + (45.0 * tanSquared($phiPrime) * tanSquared($phiPrime))); my $X = sec($phiPrime) / $v; my $XI = (sec($phiPrime) / (6.0 * $v * $v * $v)) * (($v / $rho) + (2 * tanSquared($phiPrime))); my $XII = (sec($phiPrime) / (120.0 * pow($v, 5.0))) * (5.0 + (28.0 * tanSquared($phiPrime)) + (24.0 * tanSquared($phiPrime) * tanSquared($phiPrime))); my $XIIA = (sec($phiPrime) / (5040.0 * pow($v, 7.0))) * (61.0 + (662.0 * tanSquared($phiPrime)) + (1320.0 * tanSquared($phiPrime) * tanSquared($phiPrime)) + (720.0 * tanSquared($phiPrime) * tanSquared($phiPrime) * tanSquared($phiPrime))); $phi = $phiPrime - ($VII * pow($E - $E0, 2.0)) + ($VIII * pow($E - $E0, 4.0)) - ($IX * pow($E - $E0, 6.0)); $lambda = $lambda0 + ($X * ($E - $E0)) - ($XI * pow($E - $E0, 3.0)) + ($XII * pow($E - $E0, 5.0)) - ($XIIA * pow($E - $E0, 7.0)); $Lat=rad2deg($phi); # global lat return $Long=rad2deg($lambda); # global long return } # # Convert this UTM reference to a latitude and longitude # input numeric (Eastings, Northings) # result to global $Lat, $Long # sub UTMtoLatLng { my ($x,$y,$zoneNumber,$zoneLetter)=@_; $x-=500000.0; my $UTM_F0 = 0.9996; my $a=$maj84; my $b=$min84; my $eSquared = (($a * $a) - ($b * $b)) / ($a * $a); my $ePrimeSquared = $eSquared / (1.0 - $eSquared); my $e1 = (1 - sqrt(1 - $eSquared)) / (1 + sqrt(1 - $eSquared)); #$x = $this->easting - 500000.0;; #$y = $this->northing; #$zoneNumber = $this->lngZone; #$zoneLetter = $this->latZone; my $longitudeOrigin = ($zoneNumber - 1.0) * 6.0 - 180.0 + 3.0; # Correct y for southern hemisphere if ((ord($zoneLetter) - ord("N")) < 0) { $y -= 10000000.0; } my $m = $y / $UTM_F0; my $mu = $m / ($a * (1.0 - $eSquared / 4.0 - 3.0 * $eSquared * $eSquared / 64.0 - 5.0 * pow($eSquared, 3.0) / 256.0)); my $phi1Rad = $mu + (3.0 * $e1 / 2.0 - 27.0 * pow($e1, 3.0) / 32.0) * sin(2.0 * $mu) + (21.0 * $e1 * $e1 / 16.0 - 55.0 * pow($e1, 4.0) / 32.0) * sin(4.0 * $mu) + (151.0 * pow($e1, 3.0) / 96.0) * sin(6.0 * $mu); my $n = $a / sqrt(1.0 - $eSquared * sin($phi1Rad) * sin($phi1Rad)); my $t = tan($phi1Rad) * tan($phi1Rad); my $c = $ePrimeSquared * cos($phi1Rad) * cos($phi1Rad); my $r = $a * (1.0 - $eSquared) / pow(1.0 - $eSquared * sin($phi1Rad) * sin($phi1Rad), 1.5); my $d = $x / ($n * $UTM_F0); my $latitude = ($phi1Rad - ($n * tan($phi1Rad) / $r) * ($d * $d / 2.0 - (5.0 + (3.0 * $t) + (10.0 * $c) - (4.0 * $c * $c) - (9.0 * $ePrimeSquared)) * pow($d, 4.0) / 24.0 + (61.0 + (90.0 * $t) + (298.0 * $c) + (45.0 * $t * $t) - (252.0 * $ePrimeSquared) - (3.0 * $c * $c)) * pow($d, 6.0) / 720.0)) * (180.0/$pi); my $longitude = $longitudeOrigin + (($d - (1.0 + 2.0 * $t + $c) * pow($d, 3.0)/6.0 + (5.0 - (2.0 * $c) + (28.0 * $t) - (3.0 * $c * $c) + (8.0 * $ePrimeSquared) + (24.0 * $t * $t)) * pow($d, 5.0) / 120.0) / cos($phi1Rad)) * (180.0/$pi); $Lat=$latitude; # global Lat return $Long=$longitude; # global Long return } # # Take a string formatted as a six-figure OS grid reference (e.g. # "TG514131") and return a reference to an OSRef object that represents # that grid reference. The first character must be H, N, S, O or T. # The second character can be any uppercase character from A through Z # excluding I. # input map reference string (ref) # result to global $East, $North # sub getOSRefFromSixFigureReference { my ($ref)=@_; my $char1 = substr($ref, 0, 1); my $char2 = substr($ref, 1, 1); my $east = int(substr($ref, 2, 3)) * 100; my $north = int(substr($ref, 5, 3)) * 100; if ($char1 eq 'H') { $north += 1000000; } elsif ($char1 eq 'N') { $north += 500000; } elsif ($char1 eq 'O') { $north += 500000; $east += 500000; } elsif ($char1 eq 'T') { $east += 500000; } my $char2ord = ord($char2); if ($char2ord > 73) {$char2ord--;} # Adjust for no I my $nx = (($char2ord - 65) % 5) * 100000; my $ny = (4 - floor(($char2ord - 65) / 5)) * 100000; $East=$east + $nx; # global OSRef return $North=$north + $ny; # global OSRef return } # # Convert this grid reference into a string using a standard six-figure # grid reference including the two-character designation for the 100km sub toSixFigureString { my ($easting,$northing)=@_; my $hundredkmE=floor($easting / 100000); my $hundredkmN=floor($northing / 100000); my $firstLetter=""; if ($hundredkmN < 5) { if ($hundredkmE < 5) { $firstLetter="S"; } else { $firstLetter="T"; } } elsif ($hundredkmN < 10) { if ($hundredkmE < 5) { $firstLetter="N"; } else { $firstLetter="O"; } } else { $firstLetter="H"; } my $secondLetter=""; my $index=65 + ((4 - ($hundredkmN % 5)) * 5) + ($hundredkmE % 5); my $ti=$index; if ($index >=73) {$index++}; $secondLetter=chr($index); my $e=floor(($easting - (100000 * $hundredkmE)) / 100); my $n=floor(($northing - (100000 * $hundredkmN)) / 100); my $es=$e; if ($e < 100) {$es="0$es"}; if ($e < 10) {$es="0$es"}; my $ns=$n; if ($n < 100) {$ns="0$ns"}; if ($n < 10) {$ns="0$ns"}; return $firstLetter . $secondLetter . $es . $ns; } # # Work out the UTM latitude zone from the latitude # input latitude (latitude) # return from subroutine zone letter # sub getUTMLatitudeZoneLetter { my ($latitude)=@_; if ((84 >= $latitude) && ($latitude >= 72)) {return "X"} elsif (( 72 > $latitude) && ($latitude >= 64)) {return "W"} elsif (( 64 > $latitude) && ($latitude >= 56)) {return "V"} elsif (( 56 > $latitude) && ($latitude >= 48)) {return "U"} elsif (( 48 > $latitude) && ($latitude >= 40)) {return "T"} elsif (( 40 > $latitude) && ($latitude >= 32)) {return "S"} elsif (( 32 > $latitude) && ($latitude >= 24)) {return "R"} elsif (( 24 > $latitude) && ($latitude >= 16)) {return "Q"} elsif (( 16 > $latitude) && ($latitude >= 8)) {return "P"} elsif (( 8 > $latitude) && ($latitude >= 0)) {return "N"} elsif (( 0 > $latitude) && ($latitude >= -8)) {return "M"} elsif (( -8 > $latitude) && ($latitude >= -16)) {return "L"} elsif ((-16 > $latitude) && ($latitude >= -24)) {return "K"} elsif ((-24 > $latitude) && ($latitude >= -32)) {return "J"} elsif ((-32 > $latitude) && ($latitude >= -40)) {return "H"} elsif ((-40 > $latitude) && ($latitude >= -48)) {return "G"} elsif ((-48 > $latitude) && ($latitude >= -56)) {return "F"} elsif ((-56 > $latitude) && ($latitude >= -64)) {return "E"} elsif ((-64 > $latitude) && ($latitude >= -72)) {return "D"} elsif ((-72 > $latitude) && ($latitude >= -80)) {return "C"} else {return 'Z'} } # # Convert this LatLng object into an OSGB grid reference. Note that this # function does not take into account the bounds of the OSGB grid - # beyond the bounds of the OSGB grid, the resulting OSRef object has no # meaning # sub toOSRef { my ($lat,$long)=@_; my $OSGB_F0=0.9996012717; my $N0=-100000.0; my $E0=400000.0; my $phi0=deg2rad(49.0); my $lambda0=deg2rad(-2.0); my $a=$maj1830; my $b=$min1830; my $eSquared=$E1830; my $phi=deg2rad($lat); my $lambda=deg2rad($long); my $E=0.0; my $N=0.0; my $n=($a - $b) / ($a + $b); my $v=$a * $OSGB_F0 * pow(1.0 - $eSquared * sinSquared($phi), -0.5); my $rho=$a * $OSGB_F0 * (1.0 - $eSquared) * pow(1.0 - $eSquared * sinSquared($phi), -1.5); my $etaSquared=($v / $rho) - 1.0; my $M=($b * $OSGB_F0) * (((1 + $n + ((5.0 / 4.0) * $n * $n) + ((5.0 / 4.0) * $n * $n * $n)) * ($phi - $phi0)) - (((3 * $n) + (3 * $n * $n) + ((21.0 / 8.0) * $n * $n * $n)) * sin($phi - $phi0) * cos($phi + $phi0)) + ((((15.0 / 8.0) * $n * $n) + ((15.0 / 8.0) * $n * $n * $n)) * sin(2.0 * ($phi - $phi0)) * cos(2.0 * ($phi + $phi0))) - (((35.0 / 24.0) * $n * $n * $n) * sin(3.0 * ($phi - $phi0)) * cos(3.0 * ($phi + $phi0)))); my $I=$M + $N0; my $II=($v / 2.0) * sin($phi) * cos($phi); my $III=($v / 24.0) * sin($phi) * pow(cos($phi), 3.0) * (5.0 - tanSquared($phi) + (9.0 * $etaSquared)); my $IIIA=($v / 720.0) * sin($phi) * pow(cos($phi), 5.0) * (61.0 - (58.0 * tanSquared($phi)) + pow(tan($phi), 4.0)); my $IV=$v * cos($phi); my $V=($v / 6.0) * pow(cos($phi), 3.0) * (($v / $rho) - tanSquared($phi)); my $VI=($v / 120.0) * pow(cos($phi), 5.0) * (5.0 - (18.0 * tanSquared($phi)) + (pow(tan($phi), 4.0)) + (14 * $etaSquared) - (58 * tanSquared($phi) * $etaSquared)); $North=$I + ($II * pow($lambda - $lambda0, 2.0)) + ($III * pow($lambda - $lambda0, 4.0)) + ($IIIA * pow($lambda - $lambda0, 6.0)); $East=$E0 + ($IV * ($lambda - $lambda0)) + ($V * pow($lambda - $lambda0, 3.0)) + ($VI * pow($lambda - $lambda0, 5.0)); # Global North, East return } # Convert a latitude and longitude to an UTM reference # # sub toUTMRef { my ($latitude,$longitude)=@_; my $UTM_F0 =0.9996; my $a=$maj84; my $eSquared=$E84; my $latitudeRad=$latitude * ($pi / 180.0); my $longitudeRad=$longitude * ($pi / 180.0); $longitudeZone=int(($longitude + 180.0) / 6.0) + 1; # Special zone for Norway if ($latitude>=56.0 && $latitude<64.0 && $longitude>=3.0 && $longitude<12.0) { $longitudeZone=32; } # Special zones for Svalbard if ($latitude>=72.0 && $latitude<84.0) { if ($longitude>=0.0 && $longitude<9.0) { $longitudeZone=31; } elsif ($longitude>=9.0 && $longitude<21.0) { $longitudeZone=33; } elsif ($longitude>=21.0 && $longitude<33.0) { $longitudeZone=35; } elsif ($longitude>=33.0 && $longitude<42.0) { $longitudeZone=37; } } my $longitudeOrigin=($longitudeZone - 1) * 6 - 180 + 3; my $longitudeOriginRad=$longitudeOrigin * (pi() / 180.0); $UTMZone=getUTMLatitudeZoneLetter($latitude); my $ePrimeSquared=($eSquared) / (1 - $eSquared); my $n=$a / sqrt(1 - $eSquared * sin($latitudeRad) * sin($latitudeRad)); my $t=tan($latitudeRad) * tan($latitudeRad); my $c=$ePrimeSquared * cos($latitudeRad) * cos($latitudeRad); my $A=cos($latitudeRad) * ($longitudeRad - $longitudeOriginRad); my $M=$a * ((1 - $eSquared / 4 - 3 * $eSquared * $eSquared / 64 - 5 * $eSquared * $eSquared * $eSquared / 256) * $latitudeRad - (3 * $eSquared / 8 + 3 * $eSquared * $eSquared / 32 + 45 * $eSquared * $eSquared * $eSquared / 1024) * sin(2 * $latitudeRad) + (15 * $eSquared * $eSquared / 256 + 45 * $eSquared * $eSquared * $eSquared / 1024) * sin(4 * $latitudeRad) - (35 * $eSquared * $eSquared * $eSquared / 3072) * sin(6 * $latitudeRad)); $UTMEasting=$UTM_F0 * $n * ($A + (1 - $t + $c) * pow($A, 3.0) / 6 + (5 - 18 * $t + $t * $t + 72 * $c - 58 * $ePrimeSquared) * pow($A, 5.0) / 120) + 500000.0; $UTMNorthing=$UTM_F0 * ($M + $n * tan($latitudeRad) * ($A * $A / 2 + (5 - $t + (9 * $c) + (4 * $c * $c)) * pow($A, 4.0) / 24 + (61 - (58 * $t) + ($t * $t) + (600 * $c) - (330 * $ePrimeSquared)) * pow($A, 6.0) / 720)); # Adjust for the southern hemisphere if ($latitude < 0) { $UTMNorthing +=10000000.0; } # ($UTMEasting, $UTMNorthing, $UTMZone, $longitudeZone); # Global values returned. }