Geospatial data representation method based on rhombic triacontahedron equirectangular projection

Abstract

The present application belongs to the technical field of global discrete grid system, and particularly relates to a geographic spatial data representation method based on rhombic triacontahedron equivalent projection. The method first deduces the equivalent projection formula of rhombic triacontahedron and the earth surface, including the forward projection formula and the inverse projection formula, then projects the geographic spatial data onto the rhombic triacontahedron development surface by using the forward projection formula, and constructs a plane map, i.e. the plane representation method of geographic spatial data, and generates a spherical equivalent discrete grid system by using the inverse projection formula, so as to realize the spherical gridding representation of geographic spatial data. The rhombic triacontahedron adopted in the present application has a significantly reduced equivalent grid angle deformation compared with the existing similar schemes, has good geometric properties, has a higher fitting degree to the earth, can effectively improve the expression accuracy of geographic spatial data and reduce data redundancy, and is suitable for statistical analysis and calculation.

Description

Technical Field

[0001] This invention belongs to the field of geographic information technology, specifically relating to a method for representing geospatial data based on the equal-area projection of a rhombic icosahedron. Background Technology

[0002] In geographic information systems, the traditional method of representing geospatial data is based on planar maps, which maps the Earth's surface onto a two-dimensional plane in some way. Based on the changes in geometric properties after projection, it is divided into equal-area projection, conformal projection, etc., to adapt to different application needs.

[0003] With the development of computer technology, it has become possible to describe location-related data in real three-dimensional Earth space. The increasing scale and complexity of geospatial data have also placed higher demands on the methods of its expression, processing, and analysis. As a major implementation method of "Digital Earth," the "Discrete Global Grid System" (DGGS) is a multi-resolution discrete Earth reference model formed by recursively subdividing the entire Earth space. Each unit corresponds to a region on the Earth's surface and is assigned a unique identifier, forming a new data association model with spatial regions as the primary key, providing an optimal solution for the expression, organization, and analysis of geographic information. Among existing DGGS construction methods, DGGS generated based on polyhedral projection methods have superior geometric properties. Due to the difference in geometric properties between a sphere and a plane, it is difficult to directly perform regular grid subdivision on a sphere. By recursively subdividing an ideal polyhedral surface using a specific method to form a multi-resolution regular grid, and then mapping them onto the Earth's surface, a more regular spherical grid can be obtained.

[0004] The construction of DGGS mainly involves five elements: basic polyhedrons, projection methods, polyhedron positioning, subdivision schemes, and coding. Among these, the basic polyhedrons and projection methods are the main factors affecting the deformation of grid system units. Regarding polyhedra, existing methods mainly use Pareto polyhedra, but they have a maximum of only 20 faces, and their fit with the Earth is still not high, resulting in significant distortion of grid units and low accuracy in geospatial data representation.

[0005] In terms of projection methods, the equal volume property of cells is considered a relatively ideal grid characteristic. It can effectively improve the accuracy of geospatial data representation and reduce data redundancy, making it suitable for statistical analysis and calculation. Existing Pareto polyhedral projection (DGGS) mainly uses the Snyder equal-volume polyhedral projection, which suffers from significant angular distortion at the projection boundaries. Furthermore, inverse projection requires iterative solutions, which compromises accuracy and efficiency. However, inverse projection is an essential operation for grid generation and data transformation, and its efficiency directly impacts the efficiency of the grid system. Van Leeuwen and Strebe proposed the "Slice and Dice" equal-volume polyhedral projection, which exhibits smaller angular distortion compared to the Snyder projection. However, they only applied it to Pareto polyhedra and did not provide a formula for inverse projection from a polyhedron to a sphere. Summary of the Invention

[0006] The purpose of this invention is to provide a geospatial data representation method based on the equal-area projection of a rhombic icosahedron, in order to solve the problem of low accuracy of geographic data representation caused by the poor fitting degree when using the Platonic polyhedron for Earth fitting in the prior art.

[0007] To solve the above-mentioned technical problems, the technical solution provided by this invention and the corresponding beneficial effects of the technical solution are as follows:

[0008] The present invention provides a geospatial data representation method based on the rhombic icosahedral equal-area projection, comprising the following steps:

[0009] 1) Derive the formula for the equal-area projection of a rhombic icosahedron onto the Earth's surface using the following method:

[0010] 1.1) Based on the relationship between the rhombic icosahedron used to fit the Earth and the Earth's sphere, determine the position of each vertex of the fitted spherical rhombic icosahedron; connect each vertex of the sphere with a great circle arc on the sphere, so that the area of ​​each rhombic face of the spherical rhombic icosahedron and the arc length of each rhombic face are equal.

[0011] 1.2) Calculate the relevant parameters of the equal-area projection of the rhombic icosahedron; based on the equal-area projection conditions, determine the orthographic projection formula, that is, given the geographic coordinates corresponding to the sphere, calculate the projection coordinates on the surface of the rhombic icosahedron; based on the equal-area projection conditions, determine the inverse projection formula, that is, given the coordinates on the surface of the rhombic icosahedron, calculate the geographic coordinates corresponding to the sphere.

[0012] 2) Using the orthographic projection formula, geospatial data is projected onto the unfolded surface of a rhombic icosahedron to construct a planar map, which is a planar representation method of geospatial data;

[0013] 3) Using the inverse projection formula, a spherical equal-area discrete grid system is generated according to the following method to realize the spherical rasterized representation of geospatial data:

[0014] 3.1) Design a surface subdivision scheme for the rhombic icosahedron, determine the shape of the initial grid unit, and use the initial grid unit as the current layer grid unit; perform the subdivision operation, which means continuing to subdivide downwards under the current layer grid unit; repeat the subdivision operation until the subdivision accuracy requirement is met to obtain the final grid unit.

[0015] 3.2) Using the calculation formula of inverse projection, the final grid cells are mapped to the sphere to generate a spherical grid, thereby completing the construction of the global discrete grid system;

[0016] 3.3) Based on the geographic coordinates and data resolution of the geospatial data, select a grid cell level that is close to the data resolution, use the orthographic projection formula to calculate the grid cell to which it belongs, and assign the data attribute value to the cell to realize the grid rasterization representation of the geospatial data.

[0017] The beneficial effects of the above technical solution are as follows: This method first derives the formula for equal-area projection of a rhombic icosahedron onto the Earth's surface. Based on this, geospatial data can be projected onto the unfolded surface of the rhombic icosahedron to construct a planar map, or a spherical equal-area discrete grid system can be generated to achieve a spherical rasterized representation of geospatial data. The rhombic icosahedron equal-area method used in this invention, compared to similar methods, significantly reduces angular distortion during projection transformation and possesses better geometric properties. The rhombic icosahedron has a high degree of fit to the Earth, enabling the construction of a more uniform grid system, effectively improving the accuracy of geospatial data representation and reducing data redundancy, making it suitable for statistical analysis and calculation.

[0018] Further, the equal-area projection condition mentioned in step 1.2) is the Leeuwen equal-area projection condition; wherein, when determining the calculation formulas for orthographic and inverse projection based on the Leeuwen equal-area projection condition, each rhombus face in the rhombic icosahedron needs to be divided along the short diagonal to obtain 60 isosceles triangle faces, corresponding to 60 spherical triangle faces on the spherical rhombic icosahedron, and then the calculation formulas for orthographic and inverse projection are determined based on the relationship between the spherical triangle faces and the planar triangle faces represented by the Leeuwen equal-area projection condition; point P and point P′ are defined as the corresponding projection points on the spherical triangle face ABC and the planar triangle face A′B′C′, respectively, point B is the acute vertex of the spherical rhombus face, point B′ is the acute vertex of the rhombus face on the rhombic icosahedron, and the great circle arc from point B through point P intersects with the other side of the spherical triangle face ABC. The lines from point B′ on the plane triangle A′B′C′ intersect at point D. The lengths of the great circle arcs from point P to points B and D are x and y, respectively. The lengths of the line segments from point P′ to points B′ and D′ are x′ and y′, respectively.

[0019] The formula for orthographic projection is:

[0020]

[0021] In the formula, (x P′ ,y P′ Let (x) be the coordinates of point P′; D′ Let x, y = 0 be the coordinates of point D′, and x = 0. D′ = q / 2 - qm1, where q is the length of the shorter diagonal in the spherical rhombus, and m1 is the ratio of the areas of spherical triangular faces ABD and BCD. ρ and δ are the angles of vertices B and D in the spherical triangle ABD constructed from the point P to be projected onto the sphere. ρ is calculated using the cross product formula, and δ is calculated using the cosine formula for the spherical triangle ABD. p is the length of the long diagonal in the spherical rhombus.

[0022] The beneficial effects of the above technical solution are as follows: the calculation formula for orthographic projection can be obtained by using the Leeuwen equal area projection condition, thereby realizing the planar map representation based on the rhombus thirty equal area projection and the conversion of existing geospatial data to a grid system.

[0023] Furthermore, the inverse projection formula is:

[0024]

[0025] In the formula, q p Let q be the quaternion coordinates of point P; D Let be the quaternion coordinates of point D, and q B Let q be the quaternion coordinates of point B. A Let A be the quaternion coordinates; x + y be calculated using the 3D coordinates of points B and D; x' and y' be calculated using the coordinates of points within the plane of the triangle, respectively, from point P' to points B' and D', and then the conditions for equal-area projection are applied. The length x of the great circle arc from point P to point B is calculated.

[0026] The beneficial effects of the above technical solution are as follows: the calculation formula for inverse projection can be obtained by using the Leeuwen equal-area projection condition, thereby realizing the generation of spherical grid and the visualization of grid system data. Since no iteration is required, it is more efficient than the existing polyhedral projection method based on iterative solution.

[0027] Further, in step 1.1), the positioning method between the polyhedron and the Earth is as follows: one spherical vertex is located at the North Pole, one spherical vertex is located at the South Pole, the latitudes of the five spherical vertices are 52.622631859° North, 52.622631859° South, 26.565051177° North, 26.565051177° South, 10.812316964° North, and 10.812316964° South; the longitude interval between adjacent vertices at the same latitude is 72°.

[0028] The beneficial effects of the above technical solution are as follows: after projection, the 10 meridians with longitudes of 0°, 36°, 72°, 108°, 154°, and 180° east (or west) all coincide with the diagonals or edge lengths of the rhombus surface, and are symmetrical about the equator. During grid construction, it enables the discretized spherical grid units to be associated with geographical characteristics, which is beneficial for related Earth science research.

[0029] Further, in step 2), the planar map is constructed using the following method: based on the spherical geographic coordinates of the geospatial data, the number of the spherical triangle in which it is located is determined, and the inverse projection formula is used to project it onto the planar triangle. Then, based on the position of the planar triangle on the unfolded surface of the rhombic icosahedron, the planar coordinates of the point are calculated, thereby realizing the construction of the planar map.

[0030] The beneficial effects of the above technical solution are: the projection realizes the equal-area mapping from the sphere to the unfolded plane of the rhombic icosahedron, and has smaller angular deformation compared with existing similar solutions.

[0031] Furthermore, in step 3.1), the initial grid cell shapes include the following three categories:

[0032] The first type uses 30 rhombus faces of a rhombus 30-sided prism as the initial grid unit of the rhombus grid system;

[0033] The second type involves dividing each rhombic face of the rhombic icosahedron along its diagonal to obtain 60 triangular faces, which are then used as the initial grid units of the triangular grid system.

[0034] The third type involves dividing each rhombic face of the icosahedron along its diagonal to obtain 60 triangular faces. These 60 triangular faces are then combined to form 12 pentagonal faces, which serve as the initial grid units for the hexagonal grid system.

[0035] The beneficial effects of the above technical solution are as follows: the present invention can construct a global discrete grid system based on three types of regular units—triangles, quadrilaterals, and hexagons—of equal area, which has good versatility.

[0036] Furthermore, the partitioning operation continues to partition downwards based on the aperture, where the aperture refers to the ratio of the area of ​​the current layer grid cell to the area of ​​the next layer grid cell.

[0037] Furthermore, in the triangular grid system, the apertures of the triangular faces include 4 and 9; in the rhombus grid system, the apertures of the rhombus faces include 4 and 9; and in the hexagonal grid system, the apertures of the hexagonal faces include 3, 4, and 7.

[0038] The beneficial effects of the above technical solutions are as follows: the triangular grid system, the rhomboid grid system and the hexagonal grid system are all equipped with different apertures, and can be subdivided downwards according to the actual accuracy requirements, making them suitable for different accuracy requirements. Attached Figure Description

[0039] Figure 1 This is a schematic diagram of the equal-area forward and reverse projection transformation between a rhombic icosahedron and the Earth's surface;

[0040] Figure 2 This is a schematic diagram of the index markings of the 32 vertices on the unfolded surface of the rhombic icosahedron used in this invention;

[0041] Figure 3(a) is a schematic diagram of the equal-area projection of the spherical triangle surface of a rhombic icosahedron;

[0042] Figure 3(b) is a schematic diagram of the equal-area projection of the plane triangles of the rhombic icosahedron;

[0043] Figure 4 This is a schematic diagram of a planar map on the unfolded surface of a rhombic icosahedron, created using the orthographic projection formula.

[0044] Figure 5 This is a schematic diagram illustrating the local data organization based on this method;

[0045] Figure 6 This is a schematic diagram illustrating the general process of constructing a global discrete grid system using a rhombic icosahedron as the basic polyhedron.

[0046] Figure 7 This is a schematic diagram of a spherical rhombic icosahedron diagram that uses a rhombic icosahedron to fit the Earth, thereby achieving the initial discretization of the Earth.

[0047] Figure 8 This is a schematic diagram of the grid cell distribution on the unfolded surface of the second layer of a rhombic icosahedral hexagonal grid with an aperture of 4.

[0048] Figure 9(a) is a schematic diagram of a triangular spherical equal-area grid after a rhombic 30-sided prism with an aperture of 4 is divided into the third layer of grid units;

[0049] Figure 9(b) is a schematic diagram of a quadrilateral spherical equal-area grid after the rhombic 30-sided octahedron with an aperture of 4 is divided into the third layer of grid units;

[0050] Figure 9(c) is a schematic diagram of a hexagonal spherical equal-area grid after the rhombic icosahedron with an aperture of 4 is divided into the third layer of grid units.

[0051] Figure 10 This is a schematic diagram illustrating the visualization of seawater concentration data in the Arctic region based on a grid system. Detailed Implementation

[0052] This invention constructs an equal-area mapping method between a rhombic icosahedron and the Earth's surface, and based on this, realizes the representation of geospatial data based on a planar equal-area projection map and a spherical discrete grid system. The invention will now be described in detail with reference to the accompanying drawings and embodiments.

[0053] Method Implementation Examples:

[0054] An embodiment of the geospatial data representation method based on the equal-area projection of a rhombic icosahedron according to the present invention, wherein the equal-area projection between the rhombic icosahedron and the Earth's surface is illustrated as follows. Figure 1 As shown. The specific process is as follows:

[0055] Step one: Derive the formula for the equal-area projection transformation between the rhombic icosahedron and the Earth's surface. The specific process includes:

[0056] 1. Determine the positional relationship between the polyhedron used for equal-area mapping of the Earth and the Earth's sphere.

[0057] There are multiple options for positioning parameters. One positioning method for the rhombic icosahedron provided in this embodiment is as follows: Figure 2 The vertices on the unfolded face of the rhombic icosahedron are indexed using numbers 0 to 31. The positioning method used in this embodiment is as follows: place one acute-angled vertex of a face of the rhombic icosahedron at the Earth's North Pole (North Pole). Figure 2 (Vertex 0 in the middle), and the other vertex corresponding to it ( Figure 2 If vertex 31 is placed at the Earth's South Pole, then the geographical coordinates of the 32 vertices are as follows: vertices 1-5, 6-10, and 11-15 have latitudes of 52.622631859°N, 26.565051177°N, and 10.812316964°N, respectively. Adjacent vertices at the same latitude are 72° apart in longitude. Based on symmetry, the geographical coordinates of the other vertices (i.e., their coordinates on the Earth's surface) can be obtained. Connecting the vertices of the sphere with great circle arcs ensures that the area and arc length of each rhomboid face of the spherical rhomboid are equal. Figure 7 The Earth is initially discretized into 30 spherical rhombic faces by connecting the vertices of the sphere with great circular arcs.

[0058] 2. Determine the calculation formulas for orthographic and inverse projection based on the condition of equal area projection.

[0059] Figure 1 This diagram illustrates the orthographic and inverse projections between a rhombic icosahedron and the Earth's sphere. Orthographic projection refers to the process of calculating the projected coordinates on the surface of the polyhedron given the geographic coordinates of the sphere, while inverse projection refers to the process of calculating the geographic coordinates of the sphere given the coordinates of a point on the surface of the polyhedron.

[0060] The equal-area projection of a polyhedron is the equal-area mapping between spherical triangular faces and polyhedron triangular faces. To apply this principle to a rhombic icosahedron, each rhomboid face of the rhombic icosahedron is divided along its short diagonal, resulting in 60 congruent isosceles triangular faces, corresponding to 60 spherical triangular faces on the spherical rhombic icosahedron. Due to symmetry, it is only necessary to study the projection operation between a single spherical triangular face and a planar triangular face.

[0061] Figures 3(a) and 3(b) are schematic diagrams of the spherical and planar triangular faces in equal-area projection, respectively. Points P and P′ are the corresponding projection points on the spherical triangular face ABC and the planar triangular face A′B′C′, respectively. Point B is the acute vertex of the spherical rhombus face, and point B′ is the acute vertex of the rhombus face on the polyhedron. In this embodiment, the Leeuwen equal-area projection segmentation strategy is as follows: the great circle arc that divides the spherical triangular face from the vertex is projected onto the planar triangular face as a straight line passing through the corresponding vertex. The resulting equal-area projection condition is:

[0062] ① The area of ​​a spherical triangle is equal to the area of ​​a planar triangle;

[0063] ② The great circle arc originating from vertex B and passing through projection point P intersects the other side AC at point D; correspondingly, the projected line B′P′ on the corresponding plane intersects A′C′ at point D′. BP and B′P′ divide the spherical triangle and the planar triangle into two regions with corresponding areas of equal ratio, i.e.

[0064] ③ The lengths of the great circle arcs from projection point P to points B and D are x and y, respectively; correspondingly, the lengths of the line segments from P′ to B′ and D′ are x′ and y′, respectively. The positional relationship of the projection points must satisfy...

[0065] Calculate the geometric parameters of the spherical and planar triangular faces based on given conditions for projection calculations. Using the vertex coordinates of the spherical rhombus icosahedron, the arc lengths of the angles and sides of the spherical triangular faces on the unit sphere can be calculated: ∠B = 2π / 5, ∠A = ∠C = π / 3. According to the condition of equal volume projection ①, the surface area of ​​the rhombic icosahedron and the sphere is equal in equal volume projection. The lengths of the edge r, the long diagonal p, and the short diagonal q of the rhombic surface are calculated to be r = 0.68433981476, p = 1.164268433, and q = 0.7195574638, respectively. If a sphere or rhombic icosahedron of other sizes is chosen as the reference, a set of geometric parameters of the sphere and the planar triangular surface can also be calculated. To facilitate planar coordinate calculations, a two-dimensional planar rectangular coordinate system is constructed on the planar triangular surface A′B′C′ shown in Figure 3(b): with A′C′ as the x-axis, the midpoint of A′C′ as the origin, and B′ located in the positive y-axis direction, the coordinates of the three vertices are easily obtained as A′(-q / 2,0), C′(q / 2,0), and B′(0,p / 2).

[0066] Given the geographic coordinates of a point P on a spherical triangle, and its projection onto a planar triangle as P′, the orthographic calculation process for equal-volume projection (orthographic projection) is as follows:

[0067] 1) Determine the position of the projection line l′ containing the projection point P′ on the plane triangle surface based on the equal area condition ②.

[0068] The first proportionality constant can be obtained from the formula for the area of ​​a triangular face on a unit sphere:

[0069]

[0070] In the formula, ρ and δ are the angles of vertices B and D in the spherical triangle ABD constructed from the point P to be projected onto the sphere, respectively; since the latitude and longitude of point P are known, let b, a, and p be the three-dimensional coordinates of points B, A, and P, respectively, and use the cross product formula... Solve for ρ; calculate δ using the cosine formula for the spherical triangular surface ABD:

[0071]

[0072] In the formula, γ is the angle of vertex A of the spherical triangle ABD constructed from the point P to be projected onto the sphere.

[0073] Substituting ρ and δ into equation (1), we can obtain the area ratio m1 of ΔA′B′D′ and ΔB′C′D′, and the coordinates of point D′ (x, y, y). D′ y D′ )for:

[0074]

[0075] 2) Determine the position of the projection point P' on the projection line l' based on the equal area condition ③.

[0076] Given the geographical coordinates of points B and P, calculate the arc length x between the two points; within the spherical triangle ABD, calculate cos(x+y) using the cosine formula:

[0077]

[0078] According to the condition of equal area ③ Since x and cos(x+y) have already been solved, we can obtain Calculate the coordinates (x) of the projection point P′ P′ y P′ The formula is as follows:

[0079]

[0080] The above describes the forward calculation process for the equal-area projection of a rhombic icosahedron onto a sphere.

[0081] Given the planar coordinates of a point P′ on a triangular face of a rhombic icosahedron, the process of calculating the inverse projection of its geographic coordinates onto the spherical triangular face is as follows:

[0082] 1) Find the coordinates of the intersection point D of the great circle arc passing through the projection point on the sphere and the spherical triangle surface.

[0083] For inverse projection, the plane coordinates of point P′ are known, and the scaling factor of the equal area projection condition ② can be easily calculated. Then, ρ+δ can be calculated according to equation (1).

[0084] To solve for ρ and δ, we use the cosine formula for the spherical triangle ABD, i.e., formula (2). Let τ = ρ + δ, and replace δ with τ - ρ on the left side of formula (2). Expanding, we get:

[0085]

[0086] Dividing both sides by cosρ, we get:

[0087]

[0088] The interior angles ρ and δ of the spherical triangle ABD can then be calculated, and the cosine formula can be used to determine the angles.

[0089] because Given that the coordinates of point D are calculated using the spherical linear interpolation (SLERP) formula for quaternions:

[0090]

[0091] In the formula, q D q C q ALet A and C be the quaternion coordinates of the corresponding points on the spherical triangle face. Since the projection coordinates do not require rotation parameters, the three-dimensional coordinates of point D can be calculated using only the three-dimensional coordinates of points A and C.

[0092] 2) Determine the coordinates of the projection point P on the great circle arc BD. This requires using the SLERP interpolation method again. On the great circle arc BD:

[0093]

[0094] In the formula, q B Let B be the quaternion coordinates.

[0095] In the spherical triangular face of Figure 3(a), To find the coordinates of point P, the arc length mentioned above needs to be calculated. Based on the basic projection parameters, the three-dimensional coordinates of points B and D can be easily obtained, and then x+y can be calculated. For x, given the coordinates of points within the plane of the triangle, x′ and y′ can be calculated, and then x can be solved according to the equal-area projection condition ③.

[0096] Substitute x+y, x and y into equation (9) to calculate q. P , which is the coordinate of the projection point P.

[0097] The inverse operations described above are all analytical expressions and do not require iteration.

[0098] Step two: Using the orthographic projection formula, the geospatial data is projected onto the unfolded surface of a rhombic icosahedron to obtain a planar map, such as... Figure 4 The diagram shows a planar representation method for geospatial data based on this approach. The specific process is as follows:

[0099] Based on the geographic coordinates of each point on the sphere, the projection triangle number of that point is calculated. Using the orthographic projection formula, the two-dimensional coordinates of that point on the corresponding plane triangle of the rhombic icosahedron are calculated. Based on the positional relationship of each plane triangle on the unfolded surface, these two-dimensional coordinates are transformed onto the unfolded surface of the rhombic icosahedron. Data from other existing coordinate systems is then converted to the scheme proposed in this invention, achieving data reorganization, such as... Figure 5 As shown.

[0100] Step 3: Using the inverse projection formula, generate a spherical discrete grid system with equal area. The general construction process is as follows: Figure 6 As shown, based on this, a spherical rasterized representation of geospatial data is achieved. The specific process is as follows:

[0101] 1. Design a surface partitioning scheme for the rhombic icosahedron to determine the shape of the initial grid cells, and use the initial grid cells as the current layer grid cells; perform a partitioning operation, which means continuing to partition downwards below the current layer grid cells; repeat the partitioning operation until the required partitioning accuracy is achieved to obtain the final grid cells. The specific process is as follows:

[0102] The 30 faces of the rhombic icosahedron are congruent rhombuses. Its advantage lies in the fact that rhomboid faces are suitable for the division of all regular shapes, including triangles, quadrilaterals and hexagons.

[0103] This embodiment provides the definition of three types of initial grid cells based on the rhombic icosahedral structure, namely the 0th layer grid cells: Figure 2 The 30 rhombic faces form the 0th layer grid unit of the rhombic grid system; each rhombic face is divided into two triangular faces along its diagonal, totaling 60 triangular faces, which form the 0th layer grid unit of the triangular grid system; these 60 triangular faces can form 12 pentagonal faces, which form the 0th layer grid unit of the hexagonal grid system. The pentagonal faces are respectively... Figure 2 The grid centers on 12 vertices with indices 0, 6-10, 21-25, and 31. A notable feature is the hexagonal grid, which at any resolution level consists of hexagonal faces and 12 pentagonal faces, with the area of ​​each pentagonal face being 5 / 6 the area of ​​a hexagonal face.

[0104] After determining the unit shape and initial grid units, the grid is recursively subdivided downwards based on the aperture, which refers to the ratio of the area of ​​the upper and lower layers of units. Commonly used apertures for triangular and quadrilateral faces include 4 and 9, while commonly used apertures for hexagonal faces include 3, 4, and 7.

[0105] Figure 8 The first layer of the hexagonal grid is composed of 30 hexagons and 12 pentagons.

[0106] Based on the distribution structure of the grid cells on the polyhedron, the planar coordinates of the cells are calculated, and then the geographic coordinates of the cells on the sphere can be solved by using the inverse operation of equal area projection. Figure 9(a) , 9(b) 9(c) are schematic diagrams of triangular, quadrilateral, and hexagonal spherical grids after a rhombic 30-sided prism with a aperture of 4 is divided into the third layer of grid units.

[0107] 2. Using the calculation formula of inverse projection, the final grid cells are mapped to the sphere to generate a spherical grid, thereby completing the construction of the global discrete grid system.

[0108] 3. Represent geospatial data using the constructed global discrete grid system. The specific process is as follows: determine the grid cell hierarchy based on the data resolution; generate grid cells within the region; and sample geospatial data into the grid according to location. Figure 10 This diagram illustrates the resampling of Arctic seawater concentration data using an equal-area hexagonal grid system. Traditional latitude and longitude data organization schemes suffer from severe distortion in polar regions. The proposed solution ensures the uniformity of global data sampling, improves data representation accuracy, and simultaneously reduces data storage requirements.

[0109] In summary, the method of this invention has the following characteristics: 1) Based on the rhombic icosahedron, this invention can construct a global discrete grid system with equal area based on three types of regular units: triangles, quadrilaterals, and hexagons, which has good versatility; 2) The angular deformation of the equal area grid constructed by this invention is significantly reduced compared with existing similar schemes, and it has better geometric properties; 3) The projection forward calculation of this invention can realize the conversion of existing geospatial data to a grid system; 4) The inverse projection operation of this invention can realize the generation of spherical grids and the visualization of grid system data, which is more efficient than existing methods.

Claims

1. A method for representing geospatial data based on the equal-area projection of a rhombic icosahedron, characterized in that, Includes the following steps: 1) Derive the formula for the equal-area projection of a rhombic icosahedron onto the Earth's surface using the following method: 1.1) Based on the relationship between the rhombic icosahedron used to fit the Earth and the Earth's sphere, determine the position of each vertex of the fitted spherical rhombic icosahedron; connect each vertex of the sphere with a great circle arc on the sphere, so that the area of ​​each rhombic face of the spherical rhombic icosahedron and the arc length of each rhombic face are equal. 1.2) Calculate the relevant parameters of the equal-area projection of the rhombic icosahedron; based on the equal-area projection conditions, determine the orthographic projection formula, that is, given the geographic coordinates corresponding to the sphere, calculate the projection coordinates on the surface of the rhombic icosahedron; based on the equal-area projection conditions, determine the inverse projection formula, that is, given the coordinates on the surface of the rhombic icosahedron, calculate the geographic coordinates corresponding to the sphere. 2) Using the orthographic projection formula, geospatial data is projected onto the unfolded surface of a rhombic icosahedron to construct a planar map, which is a planar representation method of geospatial data; 3) Using the inverse projection formula, a spherical equal-area discrete grid system is generated according to the following method to realize the spherical rasterized representation of geospatial data: 3.1) Design a surface subdivision scheme for the rhombic icosahedron, determine the shape of the initial grid unit, and use the initial grid unit as the current layer grid unit; perform the subdivision operation, which means continuing to subdivide downwards under the current layer grid unit; repeat the subdivision operation until the subdivision accuracy requirement is met to obtain the final grid unit. 3.2) Using the calculation formula of inverse projection, the final grid cells are mapped to the sphere to generate a spherical grid, thereby completing the construction of the global discrete grid system; 3.3) Based on the geographic coordinates and data resolution of the geospatial data, select a grid cell level that is close to the data resolution, use the orthographic projection formula to calculate the grid cell to which it belongs, and assign the data attribute value to the cell to realize the grid rasterization representation of the geospatial data.

2. The geospatial data representation method based on the rhombic icosahedral equal-area projection according to claim 1, characterized in that, The equal-area projection condition mentioned in step 1.2) is the Leeuwen equal-area projection condition. When determining the calculation formulas for orthographic and inverse projection based on the Leeuwen equal-area projection condition, each rhombus face in the icosahedron must first be divided along its short diagonal to obtain 60 isosceles triangle faces, corresponding to 60 spherical triangle faces on the spherical rhombus icosahedron. Then, the calculation formulas for orthographic and inverse projection are determined based on the relationship between the spherical triangle faces and the planar triangle faces as represented by the Leeuwen equal-area projection condition. Points P and P′ are defined as the corresponding projection points on the spherical triangle face ABC and the planar triangle face A′B′C′, respectively. Point B is the acute vertex of the spherical rhombus face, and point B′ is the acute vertex of the rhombus face on the icosahedron. A great circle arc originating from point B passes through point P and intersects with the other side of the spherical triangle face ABC. The lines from point B′ on the plane triangle A′B′C′ intersect at point D. The lengths of the great circle arcs from point P to points B and D are x and y, respectively. The lengths of the line segments from point P′ to points B′ and D′ are x′ and y′, respectively. The formula for orthographic projection is: In the formula, (x P′ ,y P′ Let (x) be the coordinates of point P′; D′ Let x, y = 0 be the coordinates of point D′, and x = 0. D′ = q / 2 - qm1, where q is the length of the shorter diagonal in the spherical rhombus, and m1 is the ratio of the areas of spherical triangular faces ABD and BCD. ρ and δ are the angles of vertices B and D in the spherical triangle ABD constructed from the point P to be projected onto the sphere. ρ is calculated using the cross product formula, and δ is calculated using the cosine formula for the spherical triangle ABD. p is the length of the long diagonal in the spherical rhombus.

3. The geospatial data representation method based on the equal-area projection of a rhombic icosahedron according to claim 2, characterized in that, The formula for inverse projection is: In the formula, q p Let q be the quaternion coordinates of point P; D Let D be the quaternion coordinates of point D, and q B Let q be the quaternion coordinates of point B. A Let A be the quaternion coordinates; x + y be calculated using the 3D coordinates of points B and D; x' and y' be calculated using the coordinates of points within the plane of the triangle, respectively, from point P' to points B' and D', and then the conditions for equal-area projection are applied. The length x of the great circle arc from point P to point B is calculated.

4. The geospatial data representation method based on the rhombic icosahedral equal-area projection according to claim 1, characterized in that, In step 1.1), the positioning method between the polyhedron and the Earth is as follows: one vertex is located at the North Pole, one vertex is located at the South Pole, the latitudes of the five vertices are 52.622631859° North, 52.622631859° South, 26.565051177° North, 26.565051177° South, 10.812316964° North, and 10.812316964° South; the longitude interval between adjacent vertices at the same latitude is 72°.

5. The geospatial data representation method based on the rhombic icosahedral equal-area projection according to claim 1, characterized in that, In step 2), the following method is used to construct the planar map: Based on the geographic coordinates of the geospatial data, the number of the spherical triangle where it is located is determined. Using the inverse projection formula, it is projected onto the planar triangle. Then, based on the position of the planar triangle on the unfolded surface of the rhombic icosahedron, the planar coordinates of the point are calculated, thereby realizing the construction of the planar map.

6. The geospatial data representation method based on the rhombic icosahedral equal-area projection according to claim 1, characterized in that, In step 3.1), the initial grid cell shapes include the following three categories: The first type uses 30 rhombus faces of a rhombus 30-sided prism as the initial grid unit of the rhombus grid system; The second type involves dividing each rhombic face of the rhombic icosahedron along its diagonal to obtain 60 triangular faces, which are then used as the initial grid units of the triangular grid system. The third type involves dividing each rhombic face of the icosahedron along its diagonal to obtain 60 triangular faces. These 60 triangular faces are then combined to form 12 pentagonal faces, which serve as the initial grid units for the hexagonal grid system.

7. The geospatial data representation method based on the rhombic icosahedral equal-area projection according to claim 6, characterized in that, The partitioning operation continues to partition downwards based on the aperture, where the aperture refers to the ratio of the area of ​​the current layer grid cell to the area of ​​the next layer grid cell.

8. The geospatial data representation method based on the rhombic icosahedral equal-area projection according to claim 7, characterized in that, In a triangular grid system, the apertures of the triangular faces include 4 and 9; in a rhombus grid system, the apertures of the rhombus faces include 4 and 9; and in a hexagonal grid system, the apertures of the hexagonal faces include 3, 4, and 7.

Images