操作矢量数据

1. 多边形的分解和组装

当部分直线和面片被组合到元素用NaN分隔的很大的矢量中以后很难识别。可以用polysplit函数把这些多边形或直线矢量分解成基本图形单元,该函数将列矢量作为输入变量。[大谦MATLAB,dqmatlab点com]

下面的例子用polysplit函数提取和连接多边形或直线段。

输入两个列矢量形式的以NaN为界的数组。

code.matlab
>> lat=[45.6 -23.47 78 NaN 43.9 -67.14 90 -89]';
>> long=[13 -97.45 165 NaN 0 -114.2 -18 0]';

用polysplit函数创建两个单元数组latc和lonc。

code.matlab
>> [latc,lonc]=polysplit(lat,long)
latc=
    [3x1 double]    [4x1 double]
lonc=
    [3x1 double]    [4x1 double]

查看单元数组的内容。

code.matlab
>> [latc{1} lonc{1}]
ans =
                      45.6                        13
                    -23.47                    -97.45
                        78                       165
>> [latc{2} lonc{2}]
ans =
                      43.9                         0
                    -67.14                    -114.2
                        90                      -18
                       -89                        0

注意,每个单元数组元素包含原直线的一段。

进行反操作,使用polyjoin函数。

code.matlab
>> [lat2,lon2]=polyjoin(latc,lonc);

连接后的线段与最初的lat和lon数组等价。

code.matlab
>> [lat long] == [lat2 lon2]
ans =
     1     1
     1     1
     1     1
     0     0
     1     1
     1     1
     1     1
     1     1

可以测试两套数据是否完全相等,键入下面的命令行。

code.matlab
>> isequalwithequalnans(lat,lat2) & isequalwithequalnans(long,lon2)
ans =
     1

2. 将直线段连接成多边形

将具有公共端点的直线段进行连接是与直线段有关的一个通用操作。polymerge命令比较保存在纬度数据矢量和经度数据矢量中的线段端点数据,以便识别精确匹配或距离在一定范围内的端点,然后将匹配的直线段进行连接。重复此过程,直到没有公用点可以找到为止。有两个必需的变量是纬度矢量和经度矢量。

下面的例子通过将直线段连接成多边形演示了这个过程。

创建表示坐标值的列矢量。

code.matlab
>> lat=[3 2 NaN 1 2 NaN 5 6 NaN 3 4]';
>> lon=[13 12 NaN 11 12 NaN 15 16 NaN 13 14]';

将精确匹配的线段连接起来。

code.matlab
>> [latm,lonm]=polymerge(lat,lon)
ans =
     5    15
     6    16
   NaN   NaN
     1    11
     2    12
     2    12
     3    13
     3    13
     4    14

原来的4条线段变成了两条。

polymerge函数还有可选的第3个变量,它设置允许不精确匹配的距离容限。第4个变量使用户可以指定函数输出矢量还是单元数组。

3. 地理插值

可以用interpm函数在点之间进行线性填充。考虑一系列纬度和经度点,对它们进行线性插值处理,使得在每个方向上的距离都不超过1°。

code.matlab
>> lats=[1 2 4 5]; longs=[1 3 4 5]; maxdiff=1;
>> [newlats,newlongs]=interpm(lats,longs,maxdiff)
newlats =
    1.0000
    1.5000
    2.0000
    3.0000
    4.0000
    5.0000
newlongs =
    1.0000
    2.0000
    3.0000
    3.5000
    4.0000
    5.0000

在原来的lats数据中,2和4之间相差2,所以,newlats和newlongs两个新变量中插入了(3, 3.5)这个点。类似地,原来的longs数据中,1和3之间也相差2,所以newlats和newlongs两个变量中插入了(1.5, 2)这个点。现在,newlats和newlongs中再没有相邻的点的距离比maxdiff还大了。

用intrplat函数进行插值,可以选择用线性、样条、三次、恒向线或大圆弧等方法进行插值。利用下面给定的数据进行计算,求与经度7.3°对应的纬度,用线性、大圆和恒向线3种方法进行计算。

code.matlab
>> longs=[1 3 4 9 13]; lats=[57 68 60 65 56]; newlong=7.3;
>> newlat=intrplat(longs,lats,newlong,'linear')
newlat =
    63.3000
>> newlat=intrplat(longs,lats,newlong,'gc')
newlat =
    63.5029
>> newlat=intrplat(longs,lats,newlong,'rh')
newlat =
    63.3937

intrplon函数的功能与intrplat的相似,只不过它求的是经度。

4. 矢量相交

Mapping工具箱提供了一系列函数对矢量数据进行求交计算。这些函数还可以计算任意矢量数据的交。

下面的例子计算圆心在(0°, 0°),半径为1250海里的小圆与圆心在(5ºN, 30ºE),半径为2500公里的小圆的交点。

code.matlab
>> [lat,long]=scxsc(0,0,nm2deg(1250),5,30,km2deg(2500))
lat =
    17.7487 -12.9839
long =
    11.0624 16.4170

其他功能相似的命令包括rhxrh、gcxgc和gcxsc等,其中,rhxrh函数计算恒向线的交点,gcxgc函数计算大圆的交点,gcxsc函数计算大圆和小圆的交点。

5. 多边形的面积

可以用函数areaint计算多边形格式矢量数据的地理面积。

下面用usalo数据集计算美国大陆的面积。

code.matlab
>> load usalo
>> earthradius=almanac('earth','radius');
>> area=areaint(uslat,uslon,earthradius)
area =
1.0e+06 *
        7.9256
        0.0035
        0.0004