查看单个帖子
旧 2009-04-26, 03:11 PM   #1
yang686526
高级会员
 
注册日期: 06-11
帖子: 14579
精华: 1
现金: 224494 标准币
资产: 234494 标准币
yang686526 向着好的方向发展
默认 [求助]高斯正算朋友写的,总算不对,帮忙看下

[求助]高斯正算朋友写的,总算不对,帮忙看下
www.dimcax.com
[求助]高斯正算朋友写的,总算不对,帮忙看下
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;sub function;;;;;;;;;;;;;;;;;;;
;;;;;;;;;expain:
;;;;;;;; 本函数用于根据大地坐标使用高斯投影正算法计算国家平面坐标
;;;;;;;; 型参:
;;;;;;;; l ,大地经度坐标值
;;;;;;;; b ,大地纬度坐标值
;;;;;;;; fn_flags, 三度带分幅 还是 六度带分幅标志
;;;;;;;; 函数返回值:
;;;;;;;; 一个包含平面坐标 x,y 的表!
;;;;******************************************************************
(defun tkzs_gszs (lb fn_flags
/ cos_dbl_b a0
a4 a6 a3
a5 p xl
n x y
x_y l b
l_d b_d
)

;;;;;;; the first step: translate 度为弧度;;;;;;;;;;;;;;;;;;;;

(setq l (* (car lb)
(/ pi 180)
)
)
(setq b (* (cadr lb)
(/ pi 180)
)
)
(setq l_d (car lb))
(setq b_d (cadr lb))
(setq cos_dbl_b (* (cos b)
(cos b)
)
)
(setq a0 (- 32140.404 (* cos_dbl_b
(- 135.3302 (* cos_dbl_b
(- 0.7092 (* 0.0040 cos_dbl_b))
)
)
)
)
)

;;;;;;;;

(setq a4 (- (* cos_dbl_b
(+ 0.25 (* 0.00252 cos_dbl_b))
)
0.04166
)
);;end setq
;;;;;;;;;
(setq a6 (* cos_dbl_b
(- (* 0.166 cos_dbl_b)
0.084
)
)
);;end setq
;;;;;;;;;
(setq a3 (- (* cos_dbl_b
(+ 0.3333333 (* 0.001123 cos_dbl_b))
)
0.1666667
)
);;;;end setq
;;;;;;;;;
(setq a5 (- 0.0083
(* cos_dbl_b
(- 0.1667
(* cos_dbl_b
(+ 0.1968
(* 0.0040 cos_dbl_b)
)
)
)
)
)
)
;;;end setq
;;;;;;;;;;
(if (= fn_flags "0")
(progn ;;;;;;;;;假如为 6 度带
(setq fixed_l (fix l_d))
(setq l0 (- (* (1+ (fix (/ fixed_l 6)))
6
)
3
)
);;end setq
)
(progn
(setq fixed_l (fix l_d))
(setq l0 (* (1+ (fix (/ fixed_l 3)))
3.0
)
);;end setq
(setq l0 (* (atoi(rtos (/ l_d 3.0) 2 0))
3
)
);;end setq
)
);;end if
(setq p (* (/ 180 pi) 3600.0))

(setq xl (/ (* (- l_d l0) 3600.0) p))
(setq n (* 0.004 cos_dbl_b cos_dbl_b cos_dbl_b cos_dbl_b))
(setq n (+ n (* 108.973 cos_dbl_b cos_dbl_b)))
(setq n (+ 6399698.902 n))
(setq n (- n (* 21562.267 cos_dbl_b) (* 0.612 cos_dbl_b cos_dbl_b cos_dbl_b)))
;;;;;;;;;
(setq x (- (/ (* 6367558.4969 (* b_d 3600.0)) p)
(* (sin b)
(cos b)
(- a0 (* xl xl n (+ 0.5
(* (+ a4 (* xl xl a6)) xl xl)
)
)
)
)
)
);;end setq
(setq y (+ 500000
(* xl
n
(cos b)
(+ 1
(* xl xl (+ a3 (* a5 xl xl)))
)
)
)
)
;;;;;;;将函数返回值列表;;;;;;;;;;;;;
;(setq x_y (list y (- x 3000000)))
(setq x_y (list y x))
);;end def
它山之石,可以攻玉。
学以致用,用而知之。
d
yang686526离线中   回复时引用此帖
GDT自动化论坛(仅游客可见)