CN113204607A - Vector polygon rasterization method for balancing area, topology and shape features - Google Patents

Abstract

The invention relates to a vector diagram rasterization method for balancing area, topology and shape characteristics, which repairs the topology change of a grid region according to the real topology of a vector polygon; the grid pixels of the polygonal narrow channel, the convex part and the concave part are reserved, so that the loss of shape characteristics can be reduced; and comparing the area difference between the adjacent grid areas and the original polygon, and reconcile the rasterization areas. The experimental result shows that compared with the conventional algorithm, the algorithm can effectively retain the area, shape and topological characteristics of the polygon, and the total rasterization error is obviously reduced.

Description

Vector polygon rasterization method for balancing area, topology and shape features

Technical Field

The invention relates to a vector polygon rasterization method for balancing area, topology and shape characteristics, and belongs to the technical field of methods for converting vector data into raster data.

Technical Field

Vectors and grids are the basic format of spatial data in Geographic Information Systems (GIS) (Maguire et al, 1991). The fundamental differences between these two data structures (Mattikali, 1995) lead to a lack of versatility in many data processing and analysis tools, and thus often face vector and grid interconversion problems. Due to the advantages of simple structure, simple and convenient data acquisition and analysis method and the like, the raster data is generally more suitable for scientific research than vector data (Peuquet, 1984; Goodchild, 2011). In recent years, the application scene of raster data is greatly expanded by a remote sensing technology. In particular, GIS-based applications such as spatial analysis and spatial modeling (Demers, 2001; Kuhnert et al, 2005), such as resource and environment regulation (Hong et al, 2016), land utilization/land cover change simulation (Mustafa et al, 2017), have led to deeper research on vector-grid transformation. The internal pixel values are generally derived from the polygon during rasterization of the polygon, and the key problem is the determination and assignment of boundary pixels. Algorithms such as boundary algebra (Feito,1995), boundary tracing (Wu,1998), etc. have been developed to find boundary pixels. Since each grid pixel can only be assigned one attribute value, rules must be made to trade off multiple attribute values when border pixels overlap multiple polygons (Wehde et al, 1980). The conventional selection and selection rules of the pixels comprise a center point method and a maximum area method, wherein the center point method is based on a polygon containing the center point of the pixel, and the maximum area method is based on the polygon with the maximum overlapping area with the pixel to determine the value of the pixel.

Discrete pixels of raster data cannot express smooth boundaries of polygons, resulting in vector-raster conversion being a lossy process (Van Der Knaap, 1992). Rasterization errors are mainly manifested as changes in the perimeter, area, shape, location, topology, etc. of polygons (Carver & Brunsdon, 1994; Galton, 2003; Shortridge, 2004). Although reducing the pixel size can effectively reduce various rasterization errors, it can result in a significant increase in data volume and processing time (Bettinger et al, 1996; Bai et al, 2011), so in practical applications it is not as small as possible for grid pixels. Analyzing and evaluating rasterization errors is an important task in the GIS field, where area errors are dominant (Wade et al, 2003). Bregt et al (1991) propose a double-transform-based area error analysis method, which finds a linear positive correlation with a boundary index (boundary length of a polygon per unit area). Kam et al (2000) found that area errors are different due to different rules of various GIS software for selecting grid pixels, which causes problems for spatial analysis depending on the area and number of pixels. Zhou et al (2006) minimize area errors using domain compensation principles based on the spatial relationship of polygons to their boundary pixels. Liao et al (2012) have established a functional relationship between the area error and the arc length density, the polygon density, and the grid pixel size by analyzing the influence factors of the rasterization error. Wei et al (2019) perform a multi-scale rasterization test based on a grid purity index, and an area error analysis result shows that the method is obviously superior to a maximum area method.

The vector-grid conversion also significantly affects the shape and topological features of the polygons. Congalton (1988) studied the spatial complexity and error patterns of 3 landscapes of arable land, grassland and woodland and found that polygons with a large perimeter-area ratio were easily lost or severely deformed in the grid representation. Bettinger et al (1996) noted that when the grid pixel size is too coarse for polygons, the shape and topology of the latter can change significantly, e.g., a polygon is broken due to the loss of a narrow channel, and two adjacent polygons are separated. Congalton (1997) explored the area and shape changes of multi-scale rasterized polygons, indicating that smaller grid pixels are required to retain the shape of the polygon than merely maintaining the polygon. Galton (2003) analyzes the influence of rasterization on the characteristics of polygon area, shape, topology, position and direction in detail, wherein the shape loss is discussed according to symmetry and convexity, and the topology difference is judged according to the region connectivity. Zhou et al (2017) compared different types of topology changes and their causes of formation and then proposed a topology preserving rasterization algorithm.

Conventional rasterization algorithms based on simple center point methods or maximum area methods are still widely applied (Longley et al, 2005), and a large number of rasterization error studies show that the spatial features of polygons are seriously lost in the conversion process. A few scholars have been working on preserving the area or topology of polygons, ignoring the fact that this method can cause greater damage to other features. The area, the shape and the topology are important spatial features of a polygon, and the invention aims to develop an algorithm for balancing and retaining the three features in vector-grid conversion so as to achieve the aim of comprehensively reducing rasterization errors. The accuracy and performance of the algorithm under various grid pixel sizes are evaluated and compared with the situation of repairing a single error. The area change of each type of land in the conventional rasterization result and the repairing result thereof is analyzed. And finally, discussing the influence of different pixel sizes and vector data on the algorithm, some application scenes of the algorithm and the future work.

Disclosure of Invention

The technical problem to be solved by the invention is as follows: the vector polygon rasterization method overcomes the defects of the prior art and provides a vector polygon rasterization method for balancing area, topology and shape characteristics.

In order to solve the technical problems, the technical scheme provided by the invention is as follows: a vector polygon rasterization method for balancing area, topology and shape features converts polygons P1 and P2 into region sets R ═ { R ═ R₁, R₂,...,R_nAnd region R2, performing the following steps to fix the rasterization error:

step 1, obtaining a polygonal boundary pixel as a candidate pixel to be classified;

step 2, if

Or n>1, performing shape restoration, preferentially ensuring that grid areas of the polygon P1 exist and are communicated to obtain an area R1, and specifically performing the following steps:

obtaining boundary pixel lost by polygon P1 as candidate set P_o；

If it is

I.e. no grid area is present, from the candidate set P_oTaking out a pixel to be allocated to P1;

if n is>1, namely the grid area is not connected, searching the pixels at the disconnected part of the polygon P1 and taking the pixels from P_oTaken out and distributed to P1;

identifying missing convex portions of polygon P1, from P_oThe pixel of the part is taken out and allocated to P1;

thus, a region R1 corresponding to the polygon P1 is obtained;

step 3, if the region R1 is separated from the region R2 and the polygon P1 is adjacent to the polygon P2, or the region R1 is adjacent to the region R2 and the polygon P1 is separated from the polygon P2, repairing the topological relation of the region, specifically including the following steps:

if the region R1 is separated from the region R2 and the polygon P1 is adjacent to the polygon P2, the lost pixels at the polygon common boundary are obtained as the candidate set P_o(ii) a The pixels occupied by the polygon P1 and the polygon P2 are respectively used as starting points and stopping points, and the shortest path idea is adopted in P_oSearching for pixels and assigning to polygon P1 or polygon P2 overlapping therewith;

if region R1 is adjacent to region R2 and polygon P1 is separated from polygon P2, find the pixel pairs that result in region adjacency as candidate set P_c＝{(p₁₁,p₁₂),(p₂₁,p₂₂),...,(p_n1,p_n2) H, mixing P_cOne of each pair of pixels in the polygon is assigned to a polygon other than the polygon P1, the polygon P2;

step 4, if the long and narrow convex part of the polygon is lost, the boundary pixel of the part is allocated to the polygon again;

and 5, combining the area of the polygon with the size of the grid pixel, and coordinating the areas of the areas R1 and R2 by using the rest boundary pixels, so that the difference between the area of the polygon and the area of the grid region is not more than half of the area of the pixel.

The invention repairs the topology change of the grid area according to the real topology of the vector polygon. And the grid pixels of the polygonal narrow channel, the convex part and the concave part are reserved, so that the loss of shape characteristics can be reduced. And comparing the area difference between the adjacent grid areas and the original polygon, and reconcile the rasterization areas. Experimental results show that compared with the conventional method, the method disclosed by the invention can effectively retain the area, shape and topological characteristics of the polygon, and the total rasterization error is obviously reduced.

Drawings

The invention will be further described with reference to the accompanying drawings.

Fig. 1 is a plot of experimental data (s ═ 10m), where (a) is a vector data plot and (b) is a grid data plot obtained using a balancing strategy.

Fig. 2 is a comparison graph (s ═ 10m) of the rasterization results, where (a) is a vector polygon; (b) is an initial rasterization result; (c) - (e) is a graph of the results of the three conventional strategies; (f) and (5) processing result graph of the balance strategy.

Detailed Description

The technical route and the operation steps of the present invention will be more clearly understood from the following detailed description of the present invention with reference to the accompanying drawings.

The invention relates to a vector polygon rasterization method for balancing area, topology and shape characteristics, which converts polygons P1 and P2 into region sets R ═ { R ═ R respectively₁,R₂,...,R_nAnd region R2, performing the following steps to fix the rasterization error:

step 1, obtaining a polygonal boundary pixel as a candidate pixel to be classified;

step 2, if

Or n>1, performing shape restoration, preferentially ensuring that grid areas of the polygon P1 exist and are communicated to obtain an area R1, and specifically performing the following steps:

obtaining boundary pixel lost by polygon P1 as candidate set P_o；

If it is

I.e. no grid area is present, from the candidate set P_oTaking out a pixel to be allocated to P1;

if n is>1, namely the grid area is not connected, searching the pixels at the disconnected part of the polygon P1 and taking the pixels from P_oTaken out and distributed to P1;

identifying missing convex portions of polygon P1, from P_oThe pixel of the part is taken out and allocated to P1;

thus, a region R1 corresponding to the polygon P1 is obtained;

step 3, if the region R1 is separated from the region R2 and the polygon P1 is adjacent to the polygon P2, or the region R1 is adjacent to the region R2 and the polygon P1 is separated from the polygon P2, repairing the topological relation of the region, specifically including the following steps:

if the region R1 is separated from the region R2 and the polygon P1 is adjacent to the polygon P2, the lost pixels at the polygon common boundary are obtained as the candidate set P_o(ii) a The pixels occupied by the polygon P1 and the polygon P2 are respectively used as starting points and stopping points, and the shortest path idea is adopted in P_oSearching for pixels and assigning to polygon P1 or polygon P2 overlapping therewith;

if region R1 is adjacent to region R2 and polygon P1 is separated from polygon P2, find the pixel pairs that result in region adjacency as candidate set P_c＝{(p₁₁,p₁₂),(p₂₁,p₂₂),...,(p_n1,p_n2) H, mixing P_cOne of each pair of pixels in the polygon is assigned to a polygon other than the polygon P1, the polygon P2;

step 4, if the long and narrow convex part of the polygon is lost, the boundary pixel of the part is allocated to the polygon again;

and 5, combining the area of the polygon with the size of the grid pixel, and coordinating the areas of the areas R1 and R2 by using the rest boundary pixels, so that the difference between the area of the polygon and the area of the grid region is not more than half of the area of the pixel.

The specific execution process of the step 5 is as follows:

a1, calculating the area difference between the polygon P and the region R obtained after rasterization, and determining the number n of pixels needing to be reclassified according to the area difference, wherein adjacent pixels with the same attribute value form a region R;

in the formula, S_pIs the area of the polygon P, S_RIs the area of the region R, and s is the area of the grid pixel;

a2, obtaining boundary image elements of the polygon P, and taking the image elements belonging to the region R as a candidate set P_aOther pixels as candidate set P_o；

a3, if n>0, from the candidate set P_aGet anda pixel is distributed to other polygons; if n is<0, from the candidate set P_oTaking out a pixel to be allocated to the polygon P;

a4, repeating the step a3 until the average area error of the region R and its adjacent regions is no longer reduced or the candidate set P_oIs empty.

After rasterization is finished, rasterization error evaluation is carried out, an area error AE, a shape error SE and a topology error TE are respectively calculated by the following formula,

wherein A is_v、A_rIs the area before and after polygon rasterization, max (A)_v,A_r) Is expressed as taking A_vAnd A_rThe greater of the two; c_v、C_rCompactness before and after rasterization of polygons respectively, | C if polygons disappear or are broken into a plurality of small regions after conversion_r-C_v|＝1.0；

N_rFor the number of pairs or polygons of the actually occurring topological change, N_vPolygon pairs or polygon numbers that are topological changes that may occur;

the compactness calculation formula of the polygon is as follows:

in the formula, A and P are the area and the perimeter of the polygon respectively, and the value range of C is (0, 1).

The total number of polygon pairs or polygons with possible topology changes is obtained according to the following formula:

where i and j represent two polygons, d (i, j) represents the shortest distance between the polygon i and the boundary of the polygon j, and s is the area of the grid pixel.

If the error of the initial rasterization result is AE₀、SE₀And TE₀Error in repair result is AE₁、SE₁And TE₁The overall repair capacity OFC of the algorithm was evaluated using the following formula:

the test data is the land utilization data of the village street in Changzhou city in China, and is derived from the investigation result of land utilization change in Changzhou city in 2018. The vector data is in the format of ESRI Shapefile, and the geographical location is shown in FIG. 1. There are 14,000 polygons in the layer, for a total area of 7,462 hectares. Including 6 land utilization types: cultivated land, garden land, water area, construction land, transportation land and other land. Table 1 lists the number of polygons, the average area, and the average compactness for different types of land.

TABLE 1 number of polygons, average area and average compactness for different types of land

The algorithm is implemented using Python programming language. The Geospatial Data Abstraction Library (GDAL) provides a series of spatial classes and spatial manipulation functions for reading, writing, converting and processing vector and raster Data. And calling an existing vector-grid conversion algorithm to obtain an initial rasterization result. Each polygon in the vector data has a unique id value, and the value is given to a grid pixel in the rasterization process. And calculating the optimal number of reclassified pixels required for minimizing the area error of the polygon and the area of the initial grid region of the polygon. To efficiently query the adjoining polygons, a quadtree spatial index (manolooulos et al, 2018) is designed and applied against the experimental data. Based on the basic idea of the Bresemhan algorithm (Gaol,2013), all the pixels that a polygon boundary passes through are calculated. And acquiring the adjacency state of the grid region by searching the boundary pixel and the 4-neighbor thereof for identifying the generated topology change. In addition, when the grid pixels are classified again, new topological changes are prevented from being introduced.

The algorithm runs on a Windows 10 platform and is provided with an Intel (R) core (TM) i7-6700 processor (with the clock frequency of 3.4GHz), an 8GB memory and a 1TB hard disk.

Control experiments were performed at various rasterization scales. An initial rasterization result is obtained based on the center point rule. And repairing the corresponding type of rasterization error by utilizing three conventional strategies (hereinafter, referred to as an area strategy, a shape strategy and a topology strategy) and a balance strategy respectively. These repair results were compared and analyzed in detail. Meanwhile, the accuracy and performance of the algorithm are evaluated, and the area change of different types of land after the error is repaired.

Taking a 10m grid pixel as an example, the vector polygon of fig. 1(a) is rasterized. The three conventional strategies and the balance strategy obtain four repaired rasterization results. FIG. 1(b) is the result of the balancing strategy.

As shown in fig. 2, the same position is selected in each rasterization result for comparing the effects of different strategies on the spatial features of the polygon and its neighbors. Fig. 2(a) is a vector layer, and the shape of the polygon P2 is approximately triangular. FIG. 2(b) is the initial rasterization result with polygon P1 transformed into two unconnected regions with an area error of 2.72%; the region R2 is approximately rectangular in shape. In fig. 2(c), the area variation is reduced, and the area error of the polygon P1 is reduced to 0.02%, but the regions R1 and R4 are adjacent, and R5 and R6 are separated. In fig. 2(d), region R2 retains one lobe of polygon P2, but it is contiguous with R3; the area errors of P1 and P2 increase. The topological variation is repaired in fig. 2(e), but the area error of the polygon P1 is still 2.61%, and the inclined narrow passage becomes vertical. FIG. 2(f) shows the result of the balancing strategy, where the area error of the polygon P1 is reduced to 0.02%; the narrow passage and lobe losses of P1 and P2 are reduced; the topology of the region is consistent with that of the polygon. The comparison shows that the balance strategy has better comprehensive effect on the reserved area, the shape and the topology.

The rasterization test was carried out using six pixel sizes of 5, 10, 15, 20, 25, 30 m. Table 2 records the area error AE, shape error SE, and topology error TE of the rasterization results obtained with the conventional strategy and the equilibrium strategy. An index OFC reflecting the repair capability of each strategy for the rasterization error is calculated.

The results show that the rasterization error is sensitive to the pixel size. As the picture element increases, AE, SE, and TE all increase significantly. The conventional strategy has a good repairing effect on the single type rasterization error. However, the area strategy causes the increase of shape error and topology error, and is especially obvious when the pixel size is small. The shape strategy or the topology strategy helps to reduce shape errors and topology errors, but introduces new area errors. The balancing strategy can reduce these three rasterization errors simultaneously. The OFC index shows that as the pixel increases, the repair capacity of the area strategy and the shape strategy increases, and the repair capacity of the topology strategy and the balance strategy gradually decreases. The main reason is that too coarse grid pixels severely lose the spatial features of the polygons. Many topological changes are difficult to repair and the interplay between different errors is strong. In summary, although the rasterization error is always present, the results of the balancing strategy are clearly superior to the other strategies. When the pixel size is less than 10m, more than half of the error can be repaired.

TABLE 2 rasterization errors for conventional and equilibrium strategies at different raster pixel sizes

Vector-to-grid conversion destroys spatial features such as the area, shape, and topology of polygons, thereby introducing various types of rasterization errors. The invention designs three conventional strategies for independently repairing the area error, the shape error and the topological error, and develops a balance strategy capable of simultaneously repairing the three errors on the basis. Under different grid pixel sizes, the balancing strategy can more effectively repair various types of grid errors. Conventional strategies repair one type of error while increasing the other type of error. When the pixel size is less than 10m, more than 50% of the rasterization errors in the result of the balancing strategy are repaired, which is obviously better than that of the conventional strategy, but the former requires a longer execution time.

In addition to the above embodiments, the present invention may have other embodiments. All technical solutions formed by adopting equivalent substitutions or equivalent transformations fall within the protection scope of the claims of the present invention.

Claims (5)

1. A vector polygon rasterization method for balancing area, topology and shape features converts polygons P1 and P2 into region sets R ═ R respectively₁,R₂,...,R_nAnd region R2, performing the following steps to fix the rasterization error:

step 1, obtaining a polygonal boundary pixel as a candidate pixel to be classified;

step 2, if

Or n>1, performing shape restoration, preferentially ensuring that grid areas of the polygon P1 exist and are communicated to obtain an area R1, and specifically performing the following steps:

obtaining boundary pixel lost by polygon P1 as candidate set P_o；

If it is

I.e. no grid area is present, from the candidate set P_oTaking out a pixel to be allocated to P1;

if n is>1, namely the grid area is not connected, searching the pixels at the disconnected part of the polygon P1 and taking the pixels from P_oTaken out and distributed to P1;

identifying missing convex portions of polygon P1, from P_oThe pixel of the part is taken out and allocated to P1;

thus, a region R1 corresponding to the polygon P1 is obtained;

step 3, if the region R1 is separated from the region R2 and the polygon P1 is adjacent to the polygon P2, or the region R1 is adjacent to the region R2 and the polygon P1 is separated from the polygon P2, repairing the topological relation of the region, specifically including the following steps:

if the region R1 is separated from the region R2 and the polygon P1 is adjacent to the polygon P2, the lost pixels at the polygon common boundary are obtained as the candidate set P_o(ii) a The pixels occupied by the polygon P1 and the polygon P2 are respectively used as starting points and stopping points, and the shortest path idea is adopted in P_oSearching for pixels and assigning to polygon P1 or polygon P2 overlapping therewith;

if region R1 is adjacent to region R2 and polygon P1 is separated from polygon P2, find the pixel pairs that result in region adjacency as candidate set P_c＝{(p₁₁,p₁₂),(p₂₁,p₂₂),...,(p_n1,p_n2) H, mixing P_cOne of each pair of pixels in the polygon is assigned to a polygon other than the polygon P1, the polygon P2;

step 4, if the long and narrow convex part of the polygon is lost, the boundary pixel of the part is allocated to the polygon again;

and 5, combining the area of the polygon with the size of the grid pixel, and coordinating the areas of the areas R1 and R2 by using the rest boundary pixels, so that the difference between the area of the polygon and the area of the grid region is not more than half of the area of the pixel.

2. The vector polygon rasterization method of claim 1 for balancing area, topology and shape features, wherein: the specific execution process of the step 5 is as follows:

a1, calculating the area difference between the polygon P and the region R obtained after rasterization, and determining the number n of pixels needing to be reclassified according to the area difference, wherein adjacent pixels with the same attribute value form a region R;

in the formula, S_pIs the area of the polygon P, S_RIs the area of the region R, and s is the area of the grid pixel;

a2, obtaining boundary image elements of the polygon P, and taking the image elements belonging to the region R as a candidate set P_aOther pixels as candidate set P_o；

a3, if n>0, from the candidate set P_aTaking out a pixel to be distributed to other polygons; if n is<0, from the candidate set P_oTaking out a pixel to be allocated to the polygon P;

a4, repeating the step a3 until the average area error of the region R and its adjacent regions is no longer reduced or the candidate set P_oIs empty.

3. The vector polygon rasterization method for balancing area, topology and shape features as recited in claim 1, wherein: after rasterization is finished, rasterization error evaluation is carried out, an area error AE, a shape error SE and a topology error TE are respectively calculated by the following formula,

wherein A is_v、A_rBefore or after rasterization of polygonsArea, max (A)_v,A_r) Is expressed as taking A_vAnd A_rThe greater of the two; c_v、C_rCompactness before and after rasterization of polygons respectively, | C if polygons disappear or are broken into a plurality of small regions after conversion_r-C_v|＝1.0；N_rFor the number of pairs or polygons of the actually occurring topological change, N_vPolygon pairs or polygon numbers that are topological changes that may occur;

if the error of the initial rasterization result is AE₀、SE₀And TE₀Error in repair result is AE₁、SE₁And TE₁The overall repair capacity OFC of the algorithm was evaluated using the following formula:

4. a vector polygon rasterization method to balance area, topology and shape features according to claim 3 and characterized in that: the compactness calculation formula of the polygon is as follows:

in the formula, A and P are the area and the perimeter of the polygon respectively, and the value range of C is (0, 1).

5. A vector polygon rasterization method to balance area, topology and shape features according to claim 3 and characterized in that: the total number of polygon pairs or polygons with possible topology changes is obtained according to the following formula:

where i and j represent two polygons, d (i, j) represents the shortest distance between the polygon i and the boundary of the polygon j, and s is the area of the grid pixel.

Images