文本处理，求助，怎么计算笛卡尔积

kwunkau

有一个文体a：经度 纬度 坐标点
112.475842        23.060001        A
112.486325        23.064556        B
112.466711        23.088511        C
112.454637        23.085901        D
112.479413        23.100844        E
112.454296        23.095857        F
112.417969        23.088477        G
112.500282        23.069834        H
112.500498        23.104442        I
112.446222        23.083907        J
112.436018        23.062628        K
112.43079        23.078846        L
112.421788        23.082244        M
112.511062        23.078412        N



现在A点要与其他点计算经纬度的距离，同理，也计算其他的点，经纬度距离公式可以随便找一个的，贴上excel的公式：
=6371004*ACOS(1-(POWER((SIN((90-A2)*PI()/180)*COS(B2*PI()/180)-SIN((90-$F$2)*PI()/180)*COS($E$2*PI()/180)),2)+POWER((SIN((90-A2)*PI()/

希望得出的结果如下：

经度        纬度        坐标点A        坐标点B        距离(M)
112.475842        23.060001        A        A        0
112.486325        23.064556        B        A        1181.636768
112.466711        23.088511        C        A        1580.854665
112.454637        23.085901        D        A        2602.063314
112.479413        23.100844        E        A        1781.157223
112.454296        23.095857        F        A        2839.188799
112.417969        23.088477        G        A        6547.774843
112.500282        23.069834        H        A        2749.59611
112.500498        23.104442        I        A        3330.030006
112.446222        23.083907        J        A        3446.620678
112.436018        23.062628        K        A        4429.635023
112.43079        23.078846        L        A        5073.083759
112.421788        23.082244        M        A        6084.284147
112.511062        23.078412        N        A        3993.837677

.
.
.
112.421788        23.082244        M        N            ...
112.511062        23.078412        N      N         ...


各位有什么好的办法吗？或者简单的输出格式，不一定按照上面的，只有计算除自身的坐标点与其他的距离，这样构成了笛卡尔乘积.

~~~新人求助，在这谢谢大家了

jason680

回复 1# kwunkau

$ awk 'BEGIN{pi=atan2(0,-1);OFS="\t"}function abs(n){return n>0?n:-n}function acos(x){return atan2(sqrt(abs(1-x*x)),x)}function r(n){return n/180*pi}function rad(n){return n/180*pi}function dist(p1,p2){lon1=lon[p1];lon2=lon[p2];lat1=lat[p1];lat2=lat[p2];d=6371004*acos(sin(r(lat1))*sin(r(lat2))+cos(r(lat1))*cos(r(lat2))*cos(r(lon2-lon1)));if(p1==p2)d=0;print lon1,lat1,p1,p2,d}{a[++c]=$3;lon[$3]=$1;lat[$3]=$2}END{for(n=1;n<=c;++n)for(m=1;m<=c;++m){dist(a[m],a[n])}}' a
112.475842    23.060001    A    A    0
112.486325    23.064556    B    A    1186.08
112.466711    23.088511    C    A    3304.92
...
112.43079    23.078846    L    N    8211.63
112.421788    23.082244    M    N    9142.18
112.511062    23.078412    N    N    0

kwunkau

回复 2# jason680

谢谢，学习下，太精髓了