vba编写kml圆思路

工作原因要用vba编写一个批量生成kml圆的代码，在网上找了很久资料很少，好不容易找到一个还是一直在讲理论，看的晕头转向不说，最后的代码还不是用VB写的，更是看不懂，只能自己想办法了，写出解题过程和代码，供大家参考。
首先在谷歌地球上画一个圆，用txt或excel打开，发现任意一个圆都是由72个不重叠的在圆周上等距的点组成的，因此只需确定这72个点就能画出圆了。画圆时圆心经纬度和半径是已知的。再观察一下这个圆，圆上有四个点是特殊的，点1(J1，W1)，点2(J2，W2)，点3(J3，W3)，点4(J4，W4)，点1和点2的经度与圆心的经度相同，点3和点4的纬度和圆心的纬度相同，如下图：

重点一：同一经度上纬度相差1度，距离相差111195米，设此值为L，因为圆心经纬度(J0，W0)和半径R都是已知的，所以可以算出上下两个特殊点的经纬度，即点1(J0，W0+R/L),点2(J0，W0-R/L)。
重点二：坐标距离公式，在网上可以找到已知两个点经度，求两点的距离公式，如下图：

我们把距离公式中A点看做是圆心，B点看做是圆周上的任意点，那么距离就是已知的R，把公式复制出来，等号左边用R表示如下：
R=6371004*ACOS((SIN(RADIANS(B2))*SIN(RADIANS(D2))+COS(RADIANS(B2))COS(RADIANS(D2))COS(RADIANS(C2-A2))))
重点二，观察上面的距离公式可以发现C2只出现了一次，也就是假设D2为已知时，通过将上面的公式转换，可以等到求C2的一个等式，如下：
C2 = Acos((Cos(R / 6371004) - Sin(Radians(B2)) * Sin(Radians(D2))) / (Cos(Radians(B2)) * Cos(Radians(D2)))) * 57.3 + A2
因此，当知道圆周上任意一点的纬度时，就可以算出他的经度，也就确定了这个点。
重点三，圆周上任意点的纬度怎么算，开始时就说到了，我们要在圆周上取等距的72个点，也就是这72个点与圆心相连形成的夹角为5度，为方便理解以任意一点为例进行解释，如下图（为方便观察5度角画的大了些）：

画出待求点所在纬线，因为任意纬线都和任意经线垂直，所有此线与圆心所在经线垂直，这样就形成了一个直角三角形，通过直角三角形相关公式可以求得L1= R Cos(5)，然后得到L2=R-L1，再通过上面说到的重点一（纬度差与距离的关系），得到待求点的纬度为W0-R Cos(5)/L，再通过重点二说到的距离公式变化出的公式可以得到待求点的经度，综上，圆周上任意一点的经纬度都能求出来了。
在求圆周任意点纬度的时候还需要根据点在圆周不同位置的情况稍微做一些公式调整，不过整体思路是不变的，就不具体说明了。
在编写vb代码的时候cos和sin要注意角度和弧度的转换。
下面是根据上面的思路用VB写的求圆周72点的自定义函数代码，三个参数分别为圆心经度，圆心纬度，半径，返回值为连续的72个经纬度，可以直接放到生成面的kml文件经纬度位置来生成圆。代码写的不好，见笑
在我的资源里有写好的宏文件，可以批量生成点线面圆，可以自定义生成对象的属性如颜色、透明度、图钉样式，在谷歌地球上点击生成的对象可以显示对象的属性卡片，并且这些属性可以自定义，生成的图层最多可以分成三层文件夹，需要的可以下载。

Function 圆周经纬(x, y, r)
Dim j2, w2 '圆周点经纬度
Dim ws, wx '圆最上、最下纬度
Dim jiaodu, hudu
Dim arc
Dim ii '循环变量
Dim wr '纬度间距
Dim yjw, yjw1

wr = 111195 '相同经度上，每相差1纬度时相差的距离，米
arc = 6371229   '米
'r1 = r * Cos(5) '角度的邻边=斜边*cos角度
'r1 = r * Sin(5) '角度的对边=斜边*sin角度

ws = y + r / wr
wx = y - r / wr

'从圆最下偏右6度逆时针，间隔5度
    '1度经纬度
w2 = wx + (r - r * Cos(Application.WorksheetFunction.Radians(1))) / wr
j2 = Application.WorksheetFunction.Acos((Cos(r / 6371004) - Sin(Application.WorksheetFunction.Radians(y)) * Sin(Application.WorksheetFunction.Radians(w2))) / (Cos(Application.WorksheetFunction.Radians(y)) * Cos(Application.WorksheetFunction.Radians(w2)))) * 57.3 + x
yjw = j2 & "," & w2

For ii = 1 To 17
    w2 = wx + (r - r * Cos(Application.WorksheetFunction.Radians(5 * ii + 1))) / wr '邻边,cos和sin内是弧度制
    j2 = Application.WorksheetFunction.Acos((Cos(r / 6371004) - Sin(Application.WorksheetFunction.Radians(y)) * Sin(Application.WorksheetFunction.Radians(w2))) / (Cos(Application.WorksheetFunction.Radians(y)) * Cos(Application.WorksheetFunction.Radians(w2)))) * 57.3 + x
    yjw = yjw & " " & j2 & "," & w2
Next

'从圆最右侧偏上逆时针
    '1度经纬度
w2 = y + (r * Sin(Application.WorksheetFunction.Radians(1))) / wr
j2 = Application.WorksheetFunction.Acos((Cos(r / 6371004) - Sin(Application.WorksheetFunction.Radians(y)) * Sin(Application.WorksheetFunction.Radians(w2))) / (Cos(Application.WorksheetFunction.Radians(y)) * Cos(Application.WorksheetFunction.Radians(w2)))) * 57.3 + x
yjw = yjw & " " & j2 & "," & w2

For ii = 1 To 17
    w2 = y + (r * Sin(Application.WorksheetFunction.Radians(5 * ii + 1))) / wr  '对角边
    j2 = Application.WorksheetFunction.Acos((Cos(r / 6371004) - Sin(Application.WorksheetFunction.Radians(y)) * Sin(Application.WorksheetFunction.Radians(w2))) / (Cos(Application.WorksheetFunction.Radians(y)) * Cos(Application.WorksheetFunction.Radians(w2)))) * 57.3 + x
    yjw = yjw & " " & j2 & "," & w2
Next

'从圆最上侧偏左逆时针
    '1度经纬度
w2 = ws - (r - r * Cos(Application.WorksheetFunction.Radians(1))) / wr
j2 = Application.WorksheetFunction.Acos((Cos(r / 6371004) - Sin(Application.WorksheetFunction.Radians(y)) * Sin(Application.WorksheetFunction.Radians(w2))) / (Cos(Application.WorksheetFunction.Radians(y)) * Cos(Application.WorksheetFunction.Radians(w2)))) * 57.3 + x
j2 = 2 * x - j2
yjw = yjw & " " & j2 & "," & w2

For ii = 1 To 17
    w2 = ws - (r - r * Cos(Application.WorksheetFunction.Radians(5 * ii + 1))) / wr '邻边
    j2 = Application.WorksheetFunction.Acos((Cos(r / 6371004) - Sin(Application.WorksheetFunction.Radians(y)) * Sin(Application.WorksheetFunction.Radians(w2))) / (Cos(Application.WorksheetFunction.Radians(y)) * Cos(Application.WorksheetFunction.Radians(w2)))) * 57.3 + x
    j2 = 2 * x - j2
    yjw = yjw & " " & j2 & "," & w2
Next

'从圆最左侧偏下逆时针
    '1度经纬度
w2 = y - (r * Sin(Application.WorksheetFunction.Radians(1))) / wr
j2 = Application.WorksheetFunction.Acos((Cos(r / 6371004) - Sin(Application.WorksheetFunction.Radians(y)) * Sin(Application.WorksheetFunction.Radians(w2))) / (Cos(Application.WorksheetFunction.Radians(y)) * Cos(Application.WorksheetFunction.Radians(w2)))) * 57.3 + x
j2 = 2 * x - j2
yjw = yjw & " " & j2 & "," & w2

For ii = 1 To 17
    w2 = y - (r * Sin(Application.WorksheetFunction.Radians(5 * ii + 1))) / wr '对角边
    j2 = Application.WorksheetFunction.Acos((Cos(r / 6371004) - Sin(Application.WorksheetFunction.Radians(y)) * Sin(Application.WorksheetFunction.Radians(w2))) / (Cos(Application.WorksheetFunction.Radians(y)) * Cos(Application.WorksheetFunction.Radians(w2)))) * 57.3 + x
    j2 = 2 * x - j2
    yjw = yjw & " " & j2 & "," & w2
Next
圆周经纬 = yjw
End Function