Lisp在GP之间绘制一条线

Cristi_an_24 · 发表于 2022-7-6 07:41:08

你好

我试着写一个lisp在两个GPs点之间画一条线，它们之间的距离是用havesine公式计算的，但我得到了一个错误[；错误：错误的参数类型：consp“56.06613,3.4514”]，我认为它由于某种原因无法识别坐标，有人能帮我吗。
我开始为自己写这个lisp是为了好玩，因为我有一部智能手机，可以从中提取GPS坐标，我想看看是否可以更准确地测量花园的大小
我不得不承认我是一个Lisp程序写作的乞丐。

顺致敬意，
克里斯蒂安

copy

 
;draw real distances from GPS coordinates
;load points from csv files
(Defun C:gpsd ( / dataL csv_name ofil rd-pos n Pi1 P1 cntr P2n ang ED lat1 lat2 Dlat Dlong a c RD ent entdata)
;*******************load files***************
(setq dataL '());define an empty list
(setq csv_name "");define an empty name for file
(setq csv_name (getfiled "Select Point List .csv file" "" "csv" 2));get the file name 
(setq ofil (open CSV_name "r")); define "ofile" as open the file
(while (setq rd-pos (read-line ofil));read each line from file
(setq dataL (cons rd-pos dataL));add each row to list
);close while
(close ofil);close file
(setq dataL (reverse dataL));put list in right order
(setq n (length dataL));get the number of points
(setq Pi1 (nth 0 dataL));get initial GPS point
(setq P1 '(0 0));set initial Real point to "0,0"
(setq xi1 (car Pi1));get initial x from first point
(setq yi1 (cadr Pi1));get initial y from first point
(setq cntr 1)
(while (<= cntr n)
(setq P2n (nth cntr dataL));choose second GPS point by cntr
(setq Pxn (car P2n));get x from second point
(setq Pyn (cadr P2n));get y from second point
(setq ang (angle Pi1 P2n));get angle between points
(setq ang (* 180.0 (/ ang pi)));convert angles from radians to degrees
;calculate real distance
(setq ED 6371000);earth diameter in metres
(setq lat1 (* pi (/ xi1 180.0)));lat in rad
(setq lat2 (* pi (/ Pxn 180.0)))
(setq Dlat (* pi(/ (- Pxn xi1) 180.0)))
(setq Dlong (* pi(/ (- Pyn yi1) 180.0))) 
(setq a (+ (sin(/ 2 Dlat)^2) (* (* (cos lat1) (cos lat2)) (sin (/ 2 Dlong)^2))))
(setq c (* 2 (atan (sqrt a) (sqrt (- 1 a)))));calculate earth angle
(setq RD (* ED c));calculate real distance between GPs Points
(setq RD (* RD 1000));convert distance from metres to mm
(setq P2 (polar P1 ang RD));calculate Real P2 by distance and angle
(command ".line" P1 P2 "");draw a line as real distance between GPS points
(setq ent (entlast));select last entoty created
(setq entdata (entget ent));extract info from it
(setq P1 (cdr (assoc 11 entdata)));reset Real P1 as (X2,Y2) from last line
(setq Pi1 P2n);reset GPS Pi1 as next point
(setq xi1 (car Pi1));recalculate x from new point
(setq yi1 (cadr Pi1));recalculate y from new point
;reset counter
(setq cntr (+ 1 cntr)) 
);close while loop
);close lisp

SLW210 · 发表于 2022-7-6 07:45:05

谢谢，我明天会试试Marko的版本，我想我是从csv中读取每一行，作为每个点的坐标列表，显然我错了。

当做
克里斯蒂安

marko_ribar · 发表于 2022-7-6 07:48:03

我再次需要你的帮助，这次我相信角度是负值，我的程序停止了，我试图用“abs”强迫我的计算得到正值，但仍然没有结果
我也改变了我的代码，我现在计算相对于第一个点的距离，得到真正的x和y，我把它们存储到一个列表中，最后我想画一条线来连接所有的点。

copy


;draw real distances from GPS coordinates
;load points from csv files
(Defun C:gpsd ( / Xcoord Ycoord dataL csv_name ofil rd-pos n Pi1 P1 cntr P2n ang ED lat1 lat2 Dlat Dlong a c RD ent entdata)
;*******************load files***************
(setq dataL '());define an empty list
(setq csv_name "");define an empty name for file
(setq csv_name (getfiled "Select Point List .csv file" "" "csv" 2));get the file name 
(setq ofil (open CSV_name "r")); define "ofile" as open the file
(while (setq rd-pos (read-line ofil));read each line from file
[highlight]
(setq Xcoord (substr rd-pos 1 (vl-string-position (ascii ",") rd-pos)))
(setq Ycoord (substr rd-pos (+ (vl-string-position (ascii ",") rd-pos) 2) (- (strlen rd-pos) (vl-string-position (ascii ",") rd-pos) 1)))
(setq dataL (cons (list (read Xcoord) (read Ycoord)) dataL));add each row to list
[/highlight]
);close while
(close ofil);close file
(setq dataL (reverse dataL));put list in right order
(setq n (length dataL));get the number of points
(setq Pi1 (nth 0 dataL));get initial GPS point
(setq P1 '(0 0));set initial Real point to "0,0"
(setq xi1 (car Pi1));get initial x from first point
(setq yi1 (cadr Pi1));get initial y from first point
(setq cntr 1)
(while (<= cntr n)
(setq P2n (nth cntr dataL));choose second GPS point by cntr
(setq Pxn (car P2n));get x from second point
(setq Pyn (cadr P2n));get y from second point
(setq ang (angle Pi1 P2n));get angle between points
(setq ang (* 180.0 (/ ang pi)));convert angles from radians to degrees
;calculate real distance
(setq ED 6371000);earth diameter in metres
(setq lat1 (* pi (/ xi1 180.0)));lat in rad
(setq lat2 (* pi (/ Pxn 180.0)))
(setq Dlat (* pi(/ (- Pxn xi1) 180.0)))
(setq Dlong (* pi(/ (- Pyn yi1) 180.0))) 
[highlight](setq a (+ (expt (sin (/ 2 Dlat)) 2) (* (* (cos lat1) (cos lat2)) (expt (sin (/ 2 Dlong)) 2))))[/highlight]
(setq c (* 2 (atan (sqrt a) (sqrt (- 1 a)))));calculate earth angle
(setq RD (* ED c));calculate real distance between GPs Points
(setq RD (* RD 1000));convert distance from metres to mm
(setq P2 (polar P1 ang RD));calculate Real P2 by distance and angle
(command ".line" P1 P2 "");draw a line as real distance between GPS points
(setq ent (entlast));select last entoty created
(setq entdata (entget ent));extract info from it
(setq P1 (cdr (assoc 11 entdata)));reset Real P1 as (X2,Y2) from last line
(setq Pi1 P2n);reset GPS Pi1 as next point
(setq xi1 (car Pi1));recalculate x from new point
(setq yi1 (cadr Pi1));recalculate y from new point
;reset counter
(setq cntr (+ 1 cntr)) 
);close while loop
);close lisp

而不是：

copy

...
(setq Pi1 (nth 0 dataL));get initial GPS point
(setq P1 '(0 0));set initial Real point to "0,0"
[color=red](setq xi1 (car Pi1));get initial x from first point
(setq yi1 (cadr Pi1));get initial y from first point[/color]
...

MSasu · 发表于 2022-7-6 07:52:28

问题在于平方根只对正数有意义；至少在一种情况下，论点是否定的。

copy

 
;draw real distances from GPS coordinates
;load points from csv files
(Defun C:gpx ( / Xcoord Ycoord dataL Rpt csv_name ofil rd-pos n Pi1 P1 cntr P2n ang ED lat1 lat2 Dlat Dlong a c RD ent entdata)
;*******************load files***************
(setq dataL '());define an empty list
(setq Rpt '());define an empty list for points
(setq csv_name "");define an empty name for file
(setq csv_name (getfiled "Select Point List .csv file" "" "csv" 2));get the file name 
(setq ofil (open CSV_name "r")); define "ofile" as open the file
   (while (setq rd-pos (read-line ofil));read each line from file
       (setq Xcoord (substr rd-pos 1 (vl-string-position (ascii ",") rd-pos)))
       (setq Ycoord (substr rd-pos (+ (vl-string-position (ascii ",") rd-pos) 2) (- (strlen rd-pos) (vl-string-position (ascii ",") rd-pos) 1)))
       (setq dataL (cons (list (read Xcoord) (read Ycoord)) dataL));add each row to list
   );close while
(close ofil);close file
(setq dataL (reverse dataL));put list in right order
(setq n (length dataL));get the number of points
(setq Pi1 (nth 0 dataL));get initial GPS point
(setq P1 '(0 0));set initial Real point to "0,0"
(setq xi1 (car Pi1));get initial x from first point
(setq yi1 (cadr Pi1));get initial y from first point
(setq Rpt (cons P1 Rpt));add P1 on list as 0,0
(setq cntr 1)
(while (<= cntr n)
   (setq P2n (nth cntr dataL));choose second GPS point by cntr
   (setq Pxn (car P2n));get x from second point
   (setq Pyn (cadr P2n));get y from second point
   (setq ang (angle Pi1 P2n));get angle between points
   (setq ang (* 180.0 (/ ang pi)));convert angles from radians to degrees
;calculate real distance
  (setq ED 6371000);earth diameter in metres
  (setq lat1 (abs(* pi (/ xi1 180.0))));lat in rad
  (setq lat2 (abs(* pi (/ Pxn 180.0))))
  (setq Dlat (abs (* pi(/ (- Pxn xi1) 180.0))));abs added 
  (setq Dlong (abs (* pi(/ (- Pyn yi1) 180.0))));abs added 
  (setq a (+ (expt (sin (/ 2 Dlat)) 2) (* (* (cos lat1) (cos lat2)) (expt (sin (/ 2 Dlong)) 2))))
  (setq c (* 2 (atan (sqrt a) (sqrt (- 1 a)))));calculate earth angle
  (setq RD (* ED c));calculate real distance between GPs Points
  (setq RD (* RD 1000));convert distance from metres to mm
  (setq P2 (polar P1 ang RD));calculate Real P2 by distance and angle
  (setq Rpt (cons P2 Rpt));add next point to the list  
;reset counter
 (setq cntr (+ 1 cntr)) 
);close while loop
(setq Rpt (reverse RPT));put the list in the right order
(command ".line" Rpt "");draw a line as real distance between GPS points
);close lisp

你应该再次验证你的公式-我无法检查它们，因为我不知道它们来自哪里。

（至少第一个）故障制造者数据对是：

Cristi_an_24 · 发表于 2022-7-6 07:56:59

这是哈弗森公式，你可以在这里找到更多细节。

谢谢
克里斯蒂安

Cristi_an_24 · 发表于 2022-7-6 08:00:12

我改变了你家的路。csv-添加了“，”作为分隔符，我在球体上绘制地球和平地面的第一个代码完成了这项工作，但第二个代码给出了错误的结果，因为有许多土壤无法联合，然后与球体相交，形成一个1米厚的外壳。。。

这是我的代码和你的csv-由我修改。。。

您好，M.R。

P、我已经用EvgeniElpanov的三角剖分法来确定区号，但在你们的例子中并没有用，你们必须找到其他方法来确定面积。。。

copy

[s][color=red](while (<= cntr n)[/color][/s]
(while (< cntr n)

copy

[s][color=red](command ".line" Rpt "")[/color][/s]
(command "_.LINE")
(foreach point Rpt
(command "_non" point)
)
(command "")

home road csv。先生txt文件
gpspronsph。lsp
gpsareasph。lsp

MSasu · 发表于 2022-7-6 08:00:31

虽然不是很精确，但我确实得到了ab区域的数据。。。大约75公顷。。。

试试看，这是我的第二个代码：

copy

(ascii "\t")

M.R。
gpsareasph新。lsp

Cristi_an_24 · 发表于 2022-7-6 08:06:19

这很管用，但在这种情况下，我的计算机上的这些GPS坐标不会生成区域对象，因此它无法按预期完成操作。。。但我认为这是最好的解决方案，因为它在计算面积方面是最精确的。。。不过，如果GPS坐标接近赤道，可能无法与曲面形成立体，并且面之间不会相互干扰。。。这仅适用于2012和+（使用SURFSCULPT命令）：

copy