Code Red > AutoLISP (Vanilla / Visual)
C++ to lisp
irneb:
First you need to know those functions. E.g. C#'s Math.Exp raises the math constant e to the specified power: http://msdn.microsoft.com/en-us/library/system.math.exp(v=vs.110).aspx
The same function in AutoLisp is simply exp: http://help.autodesk.com/view/ACD/2015/ENU/?guid=GUID-FD0C918B-A162-4939-9F4E-FFE8863C3FC8
Same goes for:
* cos - http://help.autodesk.com/view/ACD/2015/ENU/?guid=GUID-5F39D690-A986-498D-929A-56BC891E08CD
* atan - http://help.autodesk.com/view/ACD/2015/ENU/?guid=GUID-1A40A16E-8ABE-437D-9888-2F43CEA0B5CE
Unfortunately AutoLisp doesn't come standard with the tan function. But it's reasonably easy to make one yourself. Basically tan = sin/cos. So
--- Code - Auto/Visual Lisp: ---(defun tan (x) (/ (sin x) (cos x)))
Though that means you run the risk of getting a divide-by-zero error on angles like 90 degrees. So a safer way might be:
--- Code - Auto/Visual Lisp: ---(defun tan (x / c) (if (= (setq c (cos x)) 0.0) (* (if (< x 0) -1 1) 1.7976931348623158e+308) ;Closest we can get to infinity using a double floating point (/ (sin x) c)))
Then those 3 lines are:
--- Code - Auto/Visual Lisp: ---(setq senoheps (/ (- (exp eps) (exp (- eps))) 2.0))(setq delt (atan (/ senoheps (cos nab))))(setq tao (atan (* (cos delt) (tab nab))))
pedroantonio123:
Hello, I have now assembled everything.
It also works without error but the result is not correct.
(ToLatLon 451602.19 4519076.99 "33N")
should Latitud: 40.821287 Longitud; 14.426080
but with (ToLatLon) Result is Latitud: 40.5467 Longitud: 14.4284
Can someone please check the code C#. Whether this outputs the correct result :evil: :-P :oops: :cry: :-D
Thanks
Lisp Code
--- Code: ---(defun ToLatLon (utmX utmY utmZone)
(defun tan (x / c)
(if (= (setq c (cos x)) 0.0)
(* (if (< x 0) -1 1) 1.7976931348623158e+308) ;Closest we can get to infinity using a double floating point
(/ (sin x) c)
)
)
(setq latitude 0.0)
(setq longitude 0.0)
;;bool isNorthHemisphere = utmZone.Last() >= 'N';
(setq isNorthHemisphere (wcmatch utmZone "*[Nn]"))
;;var diflat = -0.00066286966871111111111111111111111111;
(setq diflat -0.00066286966871111111111111111111111111)
;;var diflon = -0.0003868060578;
(setq diflon -0.0003868060578)
;;var zone = int.Parse(utmZone.Remove(utmZone.Length - 1));
(setq zone (atoi (substr utmZone 1 (1- (strlen utmZone)))))
;;var c_sa = 6378137.000000;
(setq c_sa 6378137.0)
;;var c_sb = 6356752.314245;
(setq c_sb 6356752.314245)
;;var e2 = Math.Pow((Math.Pow(c_sa, 2) - Math.Pow(c_sb, 2)), 0.5) / c_sb;
(setq e2 (/ (expt (- (expt c_sa 2.0) (expt c_sb 2.0)) 0.5) c_sb))
;;var e2cuadrada = Math.Pow(e2, 2);
(setq e2cuadrada (expt e2 2.0))
;;var c = Math.Pow(c_sa, 2) / c_sb;
(setq c (/ (expt c_sa 2.0) c_sb))
;;var x = utmX - 500000;
(setq x (- utmX 500000))
;;var y = isNorthHemisphere ? utmY : utmY - 10000000;
(setq y (if isNorthHemisphere utmY (- utmY 10000000)))
;;var s = ((zone * 6.0) - 183.0);
(setq s (- (* zone 6.0) 183.0))
;;var lat = y / (6366197.724 * 0.9996); // Change c_sa for 6366197.724
(setq lat (/ y (* 6366197.724 0.9996)))
;;var v = (c / Math.Pow(1 + (e2cuadrada * Math.Pow(Math.Cos(lat), 2)), 0.5)) * 0.9996;
(setq v (* (/ c (expt (+ 1 (* e2cuadrada (expt (cos lat) 2) )) 0.5) ) 0.9996))
;;var a = x / v;
(setq a (/ x v))
;;var a1 = Math.Sin(2 * lat);
(setq a1 (sin (* 2.0 lat)))
;;var a2 = a1 * Math.Pow((Math.Cos(lat)), 2);
(setq a2 (* a1 (expt (cos lat) 2.0)))
;;var j2 = lat + (a1 / 2.0);
(setq j2 (+ lat (/ a1 2.0)))
;;var j4 = ((3 * j2) + a2) / 4.0;
(setq j4 (/ (+ (* 3.0 j2) a2) 4.0))
;;var j6 = (5 * j4 + a2 * Math.Pow((Math.Cos(lat)), 2)) / 3.0;
(setq j6 (/ (+ (* 5.0 j4) (* a2 (expt (cos lat) 2.0))) 3.0))
;;var alfa = (3.0 / 4.0) * e2cuadrada;
(setq alfa (* (/ 3.0 4.0) e2cuadrada))
;;var beta = (5.0 / 3.0) * Math.Pow(alfa, 2);
(setq beta (* (/ 5.0 3.0) (expt alfa 2.0)))
;;var gama = (35.0 / 27.0) * Math.Pow(alfa, 3);
(setq gama (* (/ 35.0 27.0) (expt alfa 3.0)))
;;var bm = 0.9996 * c * (lat - alfa * j2 + beta * j4 - gama * j6);
(setq bm (* 0.9996 c (- lat (+ (* alfa j2) (* beta j4)) (* gama j6)))) ;;;not shore
;;var b = (y - bm) / v;
(setq b (/ (- y bm) v))
;;var epsi = ((e2cuadrada * Math.Pow(a, 2)) / 2.0) * Math.Pow((Math.Cos(lat)), 2);
(setq epsi (* (/ (* e2cuadrada (expt a 2.0)) 2.0) (expt (cos lat) 2.0)))
;;var eps = a * (1 - (epsi / 3.0));
(setq eps (* a (- 1 (/ epsi 3.0))))
;;var nab = (b * (1 - epsi)) + lat;
(setq nab (+ (* b (1- epsi)) lat))
;;var senoheps = (Math.Exp(eps) - Math.Exp(-eps)) / 2.0;
(setq senoheps (/ (- (exp eps) (exp (- eps))) 2.0))
;;var delt = Math.Atan(senoheps / (Math.Cos(nab)));
(setq delt (atan (/ senoheps (cos nab))))
;;var tao = Math.Atan(Math.Cos(delt) * Math.Tan(nab));
(setq tao (atan (* (cos delt) (tan nab))))
;;longitude = (delt / Math.PI) * 180 + s;
(setq longitude (+ (* (/ delt pi) 180.0) s))
;;latitude = (((lat + (1 + e2cuadrada * Math.Pow(Math.Cos(lat), 2) - (3.0 / 2.0) * e2cuadrada * Math.Sin(lat) * Math.Cos(lat) * (tao - lat)) * (tao - lat))) / Math.PI) * 180; // era incorrecto el calculo
(setq latitude (* (/ (+ lat (* (- (+ 1 (* e2cuadrada (expt (cos lat) 2.0))) (* (/ 3.0 2.0) e2cuadrada (sin lat) (cos lat) (- tao lat))) (- tao lat))) pi) 180.0))
(princ "\nLatitud: ")
(princ latitude)
(princ "\nLongitud: ")
(princ longitude)
(princ)
)
--- End code ---
irneb:
Nope, it seems as if the C# code is running it properly. E.g. running this in MonoDevelop in Linux:
--- Code - C#: ---using System; namespace ToLatLon_Test{ class MainClass { public static void Main(string[] args) { Console.WriteLine("UTM to LatLon test"); string input; double utmX, utmY; string Zone; try { Console.Write("Enter utmX: "); input = Console.ReadLine(); utmX = double.Parse(input); Console.Write("Enter utmY: "); input = Console.ReadLine(); utmY = double.Parse(input); Console.Write("Enter UTM Zone: "); Zone = Console.ReadLine(); Console.WriteLine(); ToLatLon(utmX, utmY, Zone); } catch(Exception e) { Console.Write("Error: "); Console.WriteLine(e.Message); } } public static void ToLatLon(double utmX, double utmY, string utmZone) { double latitude = 0; double longitude = 0; // Changed from utmZone.Last as the Linq version is much slower than simply indexing // the last position in the string's array of characters bool isNorthHemisphere = utmZone[utmZone.Length-1] >= 'N'; var diflat = -0.00066286966871111111111111111111111111; var diflon = -0.0003868060578; var zone = int.Parse(utmZone.Remove(utmZone.Length - 1)); var c_sa = 6378137.000000; var c_sb = 6356752.314245; var e2 = Math.Pow((Math.Pow(c_sa, 2) - Math.Pow(c_sb, 2)), 0.5) / c_sb; var e2cuadrada = Math.Pow(e2, 2); var c = Math.Pow(c_sa, 2) / c_sb; var x = utmX - 500000; var y = isNorthHemisphere ? utmY : utmY - 10000000; var s = ((zone * 6.0) - 183.0); var lat = y / (6366197.724 * 0.9996); // Change c_sa for 6366197.724 var v = (c / Math.Pow(1 + (e2cuadrada * Math.Pow(Math.Cos(lat), 2)), 0.5)) * 0.9996; var a = x / v; var a1 = Math.Sin(2 * lat); var a2 = a1 * Math.Pow((Math.Cos(lat)), 2); var j2 = lat + (a1 / 2.0); var j4 = ((3 * j2) + a2) / 4.0; var j6 = (5 * j4 + a2 * Math.Pow((Math.Cos(lat)), 2)) / 3.0; // saque a2 de multiplicar por el coseno de lat y elevar al cuadrado var alfa = (3.0 / 4.0) * e2cuadrada; var beta = (5.0 / 3.0) * Math.Pow(alfa, 2); var gama = (35.0 / 27.0) * Math.Pow(alfa, 3); var bm = 0.9996 * c * (lat - alfa * j2 + beta * j4 - gama * j6); var b = (y - bm) / v; var epsi = ((e2cuadrada * Math.Pow(a, 2)) / 2.0) * Math.Pow((Math.Cos(lat)), 2); var eps = a * (1 - (epsi / 3.0)); var nab = (b * (1 - epsi)) + lat; var senoheps = (Math.Exp(eps) - Math.Exp(-eps)) / 2.0; var delt = Math.Atan(senoheps / (Math.Cos(nab))); var tao = Math.Atan(Math.Cos(delt) * Math.Tan(nab)); longitude = (delt / Math.PI) * 180 + s; latitude = (((lat + (1 + e2cuadrada * Math.Pow(Math.Cos(lat), 2) - (3.0 / 2.0) * e2cuadrada * Math.Sin(lat) * Math.Cos(lat) * (tao - lat)) * (tao - lat))) / Math.PI) * 180; // era incorrecto el calculo Console.WriteLine("Latitud: " + latitude.ToString() + "\nLongitud: " + longitude.ToString()); } }}Produces the following:
--- Code: ---UTM to LatLon test
Enter utmX: 451602.19
Enter utmY: 4519076.99
Enter UTM Zone: 33N
Latitud: 40.8212871342299
Longitud: 14.4260799700565
Press any key to continue...
--- End code ---
irneb:
Looking through your Lisp code ... this might be the culprit:
--- Code - Auto/Visual Lisp: ---;;var bm = 0.9996 * c * (lat - alfa * j2 + beta * j4 - gama * j6);(setq bm (* 0.9996 c (- lat (+ (* alfa j2) (* beta j4)) (* gama j6)))) ;;;not shore;; Should be something like this (setq bm (* 0.9996 c (- (+ (- lat (* alpha j2)) (* beta j4)) (* gama j6)))))
pedroantonio124:
Hello irneb,
Thanks for your help :-)
Unfortunately, I still receive the wrong result with the lisp code.
Maybe I can find the mistake if I know the result of every line.
Could you let me the C # run through over again code? Please
With a Console.WriteLine for every line.
PedroAntonio :evil:
--- Code: --- public static void ToLatLon(double utmX, double utmY, string utmZone)
{
double latitude = 0;
double longitude = 0;
bool isNorthHemisphere = utmZone.Last() >= 'N';
var diflat = -0.00066286966871111111111111111111111111;
var diflon = -0.0003868060578;
var zone = int.Parse(utmZone.Remove(utmZone.Length - 1));
var c_sa = 6378137.000000;
var c_sb = 6356752.314245;
var e2 = Math.Pow((Math.Pow(c_sa, 2) - Math.Pow(c_sb, 2)), 0.5) / c_sb;
Console.WriteLine ("\n" e2.ToString() )
var e2cuadrada = Math.Pow(e2, 2);
Console.WriteLine ("\n" e2cuadrada.ToString() )
var c = Math.Pow(c_sa, 2) / c_sb;
Console.WriteLine ("\n" c.ToString() )
var x = utmX - 500000;
Console.WriteLine ("\n" x.ToString() )
var y = isNorthHemisphere ? utmY : utmY - 10000000;
Console.WriteLine ("\n" y.ToString() )
var s = ((zone * 6.0) - 183.0);
Console.WriteLine ("\n" s.ToString() )
var lat = y / (6366197.724 * 0.9996); // Change c_sa for 6366197.724
Console.WriteLine ("\n" lat.ToString() )
var v = (c / Math.Pow(1 + (e2cuadrada * Math.Pow(Math.Cos(lat), 2)), 0.5)) * 0.9996;
Console.WriteLine ("\n" v.ToString() )
var a = x / v;
Console.WriteLine ("\n" a.ToString() )
var a1 = Math.Sin(2 * lat);
Console.WriteLine ("\n" a1.ToString() )
var a2 = a1 * Math.Pow((Math.Cos(lat)), 2);
Console.WriteLine ("\n" a2.ToString() )
var j2 = lat + (a1 / 2.0);
Console.WriteLine ("\n" j2.ToString() )
var j4 = ((3 * j2) + a2) / 4.0;
Console.WriteLine ("\n" j4.ToString() )
var j6 = (5 * j4 + a2 * Math.Pow((Math.Cos(lat)), 2)) / 3.0; // saque a2 de multiplicar por el coseno de lat y elevar al cuadrado
Console.WriteLine ("\n" j6.ToString() )
var alfa = (3.0 / 4.0) * e2cuadrada;
Console.WriteLine ("\n" alfa.ToString() )
var beta = (5.0 / 3.0) * Math.Pow(alfa, 2);
Console.WriteLine ("\n" beta.ToString() )
var gama = (35.0 / 27.0) * Math.Pow(alfa, 3);
Console.WriteLine ("\n" gama.ToString() )
var bm = 0.9996 * c * (lat - alfa * j2 + beta * j4 - gama * j6);
Console.WriteLine ("\n" bm.ToString() )
var b = (y - bm) / v;
Console.WriteLine ("\n" b.ToString() )
var epsi = ((e2cuadrada * Math.Pow(a, 2)) / 2.0) * Math.Pow((Math.Cos(lat)), 2);
Console.WriteLine ("\n" epsi.ToString() )
var eps = a * (1 - (epsi / 3.0));
Console.WriteLine ("\n" eps.ToString() )
var nab = (b * (1 - epsi)) + lat;
Console.WriteLine ("\n" nab.ToString() )
var senoheps = (Math.Exp(eps) - Math.Exp(-eps)) / 2.0;
Console.WriteLine ("\n" senoheps.ToString() )
var delt = Math.Atan(senoheps / (Math.Cos(nab)));
Console.WriteLine ("\n" delt.ToString() )
var tao = Math.Atan(Math.Cos(delt) * Math.Tan(nab));
Console.WriteLine ("\n" tao.ToString() )
longitude = (delt / Math.PI) * 180 + s;
latitude = (((lat + (1 + e2cuadrada * Math.Pow(Math.Cos(lat), 2) - (3.0 / 2.0) * e2cuadrada * Math.Sin(lat) * Math.Cos(lat) * (tao - lat)) * (tao - lat))) / Math.PI) * 180; // era incorrecto el calculo
Console.WriteLine("Latitud: " + latitude.ToString() + "\nLongitud: " + longitude.ToString());
}
--- End code ---
Navigation
[0] Message Index
[#] Next page
[*] Previous page
Go to full version