| 查看: 1362 | 回复: 19 | |||
| 【奖励】 本帖被评价5次，作者3s_studio增加金币 4.5 个 | |||
| 当前主题已经存档。 | |||
| 当前只显示满足指定条件的回帖，点击这里查看本话题的所有回帖 | |||
[资源]
【分享】由经纬度坐标求距离的VB程序（附源代码）
|
|||
|
经纬度坐标的地图中求距离的VB6.0程序，距离单位为M，利用大地正反算白塞尔公式编制而成。 下载： 程序及源代码下载 [ Last edited by 3s_studio on 2008-12-13 at 09:38 ] |
回复此楼» 猜你喜欢
长江大学石油工程学院有调剂名额
已经有0人回复
-
首次在隐伏矿体上方发现金、硫酸铅、硝酸铅、三氧化钨等含金属纳米微粒
已经有8人回复
-
地质学论文润色/翻译怎么收费?
已经有154人回复
-
求助：有没有大神有流体包裹体计算软件
已经有0人回复
-
地区已中，感谢虫友们！
已经有109人回复
|
更详细源码请看http://bbs.i-location.org/thread-294-1-1.html /// /// 贝赛尔大地问题正解 /// /// 已知点纬度 /// 已知点经度 /// 大地线长 /// 大地方位角 /// 待求点纬度 /// 待求点经度 /// 大地反方位角 /// 参考椭球长半轴 /// 参考椭球扁率倒数 public static void Bessel_PntSA_Pnt(double B1, double L1, double S, double A1, out double B2, out double L2, out double A2, double a, double f) { double ee = (2 * f - 1) / f / f; //第一偏心率的平方 double ee2 = ee / (1 - ee); //第二偏心率的平方 double u1 = Math.Atan(Math.Sqrt(1 - ee) * Math.Tan(B1 * Math.PI / 180)); double m = Math.Cos(u1) * Math.Sin(A1); m = Math.Atan(m / Math.Sqrt(1 - m * m)); double M_ = Math.Atan(Math.Tan(u1) / Math.Cos(A1)); if (M_ < 0) M_ += Math.PI; double KK = ee2 * Math.Pow(Math.Cos(m), 2); double alpha = Math.Sqrt(1 + ee2) * (1 - KK / 4 + 7 * KK * KK / 64 - 15 * Math.Pow(KK, 3) / 256) / a; double beta = KK / 4 - KK * KK / 8 + 37 * Math.Pow(KK, 3) / 512; double gama = KK * KK * (1 - KK) / 128; double sigma, temp; sigma = alpha * S; temp = 0; while (Math.Abs(temp - sigma) > 0.0000000001) { temp = sigma; sigma = alpha * S + beta * Math.Sin(sigma) * Math.Cos(2 * M_ + sigma) + gama * Math.Sin(2 * sigma) * Math.Cos(4 * M_ + 2 * sigma); } A2 = Math.Atan(Math.Tan(m) / Math.Cos(M_ + sigma)); if (Math.Cos(M_) > 0) A2 += Math.PI; if (m < 0 && Math.Cos(M_) < 0) A2 += 2 * Math.PI; double u2 = Math.Atan(-Math.Cos(A2) * Math.Tan(M_ + sigma)); double lamda1; lamda1 = Math.Atan(Math.Sin(u2) * Math.Tan(A2)); if (Math.Cos(M_) < 0) lamda1 += Math.PI; if (m < 0 && Math.Cos(M_) > 0) lamda1 += 2 * Math.PI; double lamda2; lamda2 = Math.Atan(Math.Sin(u1) * Math.Tan(A1)); if (m > 0) { if (Math.Cos(M_ + sigma) < 0) lamda2 += Math.PI; else if (Math.Sin(M_ + sigma) < 0) lamda2 += 2 * Math.PI; } else { if (Math.Cos(2 * Math.PI - M_ - sigma) < 0) lamda2 += Math.PI; else if (Math.Sin(M_ + sigma) < 0) lamda2 += 2 * Math.PI; } B2 = Math.Atan(Math.Sqrt(1 + ee2) * Math.Tan(u2)) * 180 / Math.PI; KK = ee * Math.Pow(Math.Cos(m), 2); alpha = ee / 2 + ee * ee / 8 + Math.Pow(ee, 3) / 16 - ee * (1 + ee) * KK / 16 + 3 * ee * KK * KK / 128; gama = ee * KK * KK / 256; L2 = lamda1 + lamda2 - Math.Sin(m) * (alpha * sigma + beta * Math.Sin(sigma) * Math.Cos(2 * M_ + sigma) + gama * Math.Sin(2 * sigma) * Math.Cos(4 * M_ + 2 * sigma)); L2 = L1 + L2 * 180 / Math.PI; } /// /// 由两点大地坐标求解大地方位角 /// /// 起点纬度 /// 起点经度 /// 末点纬度 /// 末点经度 /// 参考椭球长半轴 /// 参考椭球扁率倒数 /// public static double Bessel_BL_A(double B1, double L1, double B2, double L2, double a, double f) { if (L1 == L2) return 0; double ee = (2 * f - 1) / f / f; //第一偏心率的平方 double u1 = Math.Atan(Math.Sqrt(1 - ee) * Math.Tan(B1 * Math.PI / 180)); double u2 = Math.Atan(Math.Sqrt(1 - ee) * Math.Tan(B2 * Math.PI / 180)); double dL = (L2 - L1) * Math.PI / 180; double sigma = Math.Sin(u1) * Math.Sin(u2) + Math.Cos(u1) * Math.Cos(u2) * Math.Cos(dL); sigma = Math.Atan(Math.Sqrt(1 - sigma * sigma) / sigma); if (sigma <= 0) sigma += Math.PI; double m = Math.Cos(u1) * Math.Cos(u2) * Math.Sin(dL) / Math.Sin(sigma); m = Math.Atan(m / Math.Sqrt(1 - m * m)); double KK = ee * Math.Pow(Math.Cos(m), 2); double alpha = ee / 2 + ee * ee / 8 + Math.Pow(ee, 3) / 16 - ee * (1 + ee) * KK / 16 + 3 * ee * KK * KK / 128; double lamda = dL + alpha * sigma * Math.Sin(m); sigma += Math.Sin(m) * (alpha * sigma * Math.Sin(m)); m = Math.Cos(u1) * Math.Cos(u2) * Math.Sin(lamda) / Math.Sin(sigma); m = Math.Atan(m / Math.Sqrt(1 - m * m)); double A = Math.Atan(Math.Sin(lamda) / (Math.Cos(u1) * Math.Tan(u2) - Math.Sin(u1) * Math.Cos(lamda))); if (A <= 0) A += Math.PI; if (m <= 0) A += Math.PI; double M_ = Math.Atan(Math.Sin(u1) * Math.Tan(A) / Math.Sin(m)); if (M_ <= 0) M_ += Math.PI; KK = ee * Math.Pow(Math.Cos(m), 2); alpha = ee / 2 + ee * ee / 8 + Math.Pow(ee, 3) / 16 - ee * (1 + ee) * KK / 16 + 3 * ee * KK * KK / 128; double beta = ee * (1 + ee) * KK / 16 - ee * KK * KK / 32; lamda = dL + Math.Sin(m) * (alpha * sigma + beta * Math.Sin(sigma) * Math.Cos(2 * M_ + sigma)); A = Math.Atan(Math.Sin(lamda) / (Math.Cos(u1) * Math.Tan(u2) - Math.Sin(u1) * Math.Cos(lamda))); if (A <= 0) A += Math.PI; if (m <= 0) A += Math.PI; return A * 180 / Math.PI; } /// /// 由大地坐标计算大地线长 /// /// 起点纬度 /// 起点经度 /// 末点纬度 /// 末点经度 /// 参考椭球长半轴 /// 参考椭球扁率倒数 /// public static double Bessel_BL_S(double B1, double L1, double B2, double L2, double a, double f) { if (L1 == L2) return Math.Abs(MeridianLength(B1, a, f) - MeridianLength(B2, a, f)); double ee = (2 * f - 1) / f / f; //第一偏心率的平方 double ee2 = ee / (1 - ee); //第二偏心率的平方 double u1 = Math.Atan(Math.Sqrt(1 - ee) * Math.Tan(B1 * Math.PI / 180)); double u2 = Math.Atan(Math.Sqrt(1 - ee) * Math.Tan(B2 * Math.PI / 180)); double dL = (L2 - L1) * Math.PI / 180; double sigma = Math.Sin(u1) * Math.Sin(u2) + Math.Cos(u1) * Math.Cos(u2) * Math.Cos(dL); sigma = Math.Atan(Math.Sqrt(1 - sigma * sigma) / sigma); if (sigma <= 0) sigma += Math.PI; double m = Math.Cos(u1) * Math.Cos(u2) * Math.Sin(dL) / Math.Sin(sigma); m = Math.Atan(m / Math.Sqrt(1 - m * m)); double KK = ee * Math.Pow(Math.Cos(m), 2); double alpha = ee / 2 + ee * ee / 8 + Math.Pow(ee, 3) / 16 - ee * (1 + ee) * KK / 16 + 3 * ee * KK * KK / 128; double lamda = dL + alpha * sigma * Math.Sin(m); sigma += Math.Sin(m) * (alpha * sigma * Math.Sin(m)); m = Math.Cos(u1) * Math.Cos(u2) * Math.Sin(lamda) / Math.Sin(sigma); m = Math.Atan(m / Math.Sqrt(1 - m * m)); double A = Math.Atan(Math.Sin(lamda) / (Math.Cos(u1) * Math.Tan(u2) - Math.Sin(u1) * Math.Cos(lamda))); if (A <= 0) A += Math.PI; if (m <= 0) A += Math.PI; double M_ = Math.Atan(Math.Sin(u1) * Math.Tan(A) / Math.Sin(m)); if (M_ <= 0) M_ += Math.PI; KK = ee * Math.Pow(Math.Cos(m), 2); alpha = ee / 2 + ee * ee / 8 + Math.Pow(ee, 3) / 16 - ee * (1 + ee) * KK / 16 + 3 * ee * KK * KK / 128; double beta = ee * (1 + ee) * KK / 16 - ee * KK * KK / 32; lamda = dL + Math.Sin(m) * (alpha * sigma + beta * Math.Sin(sigma) * Math.Cos(2 * M_ + sigma)); sigma = Math.Sin(u1) * Math.Sin(u2) + Math.Cos(u1) * Math.Cos(u2) * Math.Cos(lamda); sigma = Math.Atan(Math.Sqrt(1 - sigma * sigma) / sigma); if (sigma < 0) sigma += Math.PI; KK = ee2 * Math.Pow(Math.Cos(m), 2); alpha = Math.Sqrt(1 + ee2) / a * (1 - KK / 4 + 7 * KK * KK / 64 - 15 * Math.Pow(KK, 3) / 256); beta = KK / 4 - KK * KK / 8 + 37 * Math.Pow(KK, 3) / 512; double gama = KK * KK * (1 - KK) / 128; double S = (sigma - beta * Math.Sin(sigma) * Math.Cos(2 * M_ + sigma) - gama * Math.Sin(2 * gama) * Math.Cos(4 * M_ + 2 * sigma)) / alpha; return S; } [ Last edited by moonlight3988 on 2008-7-29 at 21:06 ] |
12楼2008-06-26 12:49:33
2楼2008-06-24 09:53:42
3楼2008-06-24 09:59:28
4楼2008-06-24 10:50:34