7

Compute intersection size with HANA spatial data

 1 year ago
source link: https://blogs.sap.com/2023/04/07/compute-intersection-size-with-hana-spatial-data/
Go to the source link to view the article. You can view the picture content, updated content and better typesetting reading experience. If the link is broken, please click the button below to view the snapshot at that time.
neoserver,ios ssh client
April 7, 2023 4 minute read

Compute intersection size with HANA spatial data

Imagine you were given a map divided into blocks, with a lake marked on it. And now you’re supposed to decide to which block this lake belongs.

Map%20with%20blocks%20and%20lake

Map with blocks and lake

For purposes of this post, let’s assume we should select the block with the biggest overlap with the lake itself. In this case, both the blocks and the lake are represented as ST_GEOMETRY types within our database. And we are accessing our database using the AMDPs.

The HANA database allows native handling of spatial data (see here for HELP page).

ST_INTERSECTS

To start us off, let’s select just the blocks that have some kind of intersection with our lake. For this end, we can use the ST_Intersects( ) method. It is available for all ST_Geometry types.

(Help page)

If an intersection exists, this will return 1. Otherwise, 0.

If (BLOCKLOCATION.ST_Intersection(:lake) = 1) 
Then 
  -- this Blocklocation does have an intersection with our lake
End if;

So, we can get all the blocks that have some intersection with our lake. Now we need to sort them by their size.

ST_INTERSECTION

We can get an ST_Geometry representation of the intersection using the ST_Intersection( ) method.(Help page)

It returns an empty geometry if there is no intersection for the given two geometries. Otherwise, it return a geometry that describes the intersection.

Intersection = Blocklocation.st_intersection(:lake);

ST_AREA

Since we do have the geometry representation of the intersection, it should be easy to compute it’s area, right? Just use ST_Area( ) method on the intersection?

Intersection_area = intersection.ST_Area( )

Well, not so quickly. The ST_Area method can be used with geometries of types ST_Polygon and ST_MultiPolygon.

What are the possible results of the ST_Intersection method? Well, any ST_Geometry. Which can include ST_Point, ST_Line and ST_GeometryCollection. And these will cause issues. So let’s look at these in detail.

ST_DIMENSION

ST_Point and ST_Line

It can happen that the whole intersection is actually just a point/line where the geometries touch. But they do not actually overlap. So, how can we filter these out? For example, by using the ST_Dimension( ) method.

(Help page)

If (intersection_area.st_dimension( ) = 0 or intersection_area.st_dimension( ) = 1)
Then
  -- these are points and lines, not interesting for us now
End if;

ST_GEOMETRYTYPE

ST_GeometryCollection

Another possibility is, that there can be multiple intersections. Like in this picture.

st_map2-1.png

In this case, the value would be of type ST_GeometryCollection. It’s dimension would then be equal to the highest dimension of its parts. So, if our collection contains a line and a polygon, the dimension would be 2.

How can we recognize that we’re dealing with a collection? We can use the ST_GeometryType( ) method.

If (intersection_area.st_geometryType( ) = 'ST_GeometryCollection')
Then
  -- now we know that our intersection area is a collection
End if;

So, how do we calculate the area of a collection? One option is to iterate over all the geometries within the collection and add all their respective areas.

ST_NUMGEOMETRIES

For iteration, it would be nice to know the number of geometries in the collection. Fortunately, there is a method for this, the ST_Numgeometries method.

Help Page

And once we know the overall number, we can iterate over the collection and retrieve the individual geometries using the ST_GeometryN method.

Help page

for index_coll in 1..intsec.st_numgeometries(  )
do
  -- now we can iterate over the collection and access the Nth geometry each iteration like this
  Geom = intsec.st_geometryN( :index_coll );
End for;

Putting it all together

So, how do we put this all together?

Let’s return to our 4 blocks and a lake from the beginning and say we have:

  • Lake — the geometry representing the lake
  • Blocks — a table that has the geometries representing our blocks (geometry), plus their numbers (blocknumber)
-- we will use theselater  for iterating
DECLARE index_block, index_coll int;

-- here we get all the blocks with some intersection
Tmp_blocks = SELECT blocknumber,
  b.geometry.st_intersecion(:lake) AS intersection
  FROM blocks b
  WHERE b.geometry.st_intersects(:lake) = 1;

-- we get the intersections with actual area and for polygons, we calculate the area straight
Tmp_areas = SELECT blocknumber,
  Intersection,
  Intersection.st_dimension( ) AS dimension,
  Intersection.st_geometrytype( ) AS geomtype,
  CASE intsec.st_geometrytype( )
     WHEN 'ST_Polygon' THEN intersection.st_area(  )
     WHEN 'ST_MultiPolygon' THEN inertsection.st_area(  )
     ELSE 0
  END AS intsec_area
  FROM :tmp_blocks
  WHERE intersection.st_dimension( ) = 2;

-- Now let's iterate over the collections and sum the areas of their parts
FOR index_block IN 1..record_count(tmp_areas)
  DO
  -- check that we're dealing with a collection
  IF ( :tmp_areas.geomtype[ :index_block ] = 'ST_GeometryCollection'  )
    THEN
    DECLARE intsec st_geometry = :tmp_areas.intersection[:index_block];
    -- iterate over the collection
    FOR index_coll IN 1..intsec.st_numgeometries(  )
      DO
      -- check if this geometry has an area
      IF (:intsec.st_geometryN( :index_coll ).st_dimension(  ) = 2)
        THEN
        -- add it's area to the sum
        tmp_areas.intsec_area[:index_block] = :tmp_areas.intsec_area[:index_block] + 
        :intsec.st_geometryN(:index_coll).st_area(  );
      END IF;
    END FOR;
END FOR;

And what do we have now in our tmp_areas table? A list of blocks that have an intersection with our lake (that is not just a border touch) and for each we have a intsec_area value, in which is the sum of the overlaps of this block with our lake. So now we just select the one with biggest value and we’re done!

Conclusion

In this post, I have provided a very brief overview on a possibility how to calculate the areas of intersections of geometries using HANA native spatial data handling capabilities. I went over the basics of usage of the relevant methods, like ST_Intersects, ST_Intersection, ST_Area, ST_GeometryType, ST_Dimension, ST_Numgeometries and ST_GeometryN.


About Joyk


Aggregate valuable and interesting links.
Joyk means Joy of geeK