locked
How to create a LineString or Polygon from Points? RRS feed

  • Question

  •  

    I have several geography Points.  What's the best way to create a LineString or Polygon out of them?  It seems like there should be a MultiPoint.ToLineString() or MultiPoint.ToPolygon methods.  I've seen that other SQL spatial databases have a function like LineString.AddPoint that would also do this.

     

    For now I'm using Parse and the Lat and Long functions for Point but this seems to be pretty ugly and inefficient.

     

    Anyone have a better way?

     

    Geometry has some shortcuts like using STConvexHull but this only works for LineString if you have only two points and it only works for Polygon if all your points are on the perimeter.

     

    Thanks,

    Matt

    Wednesday, May 14, 2008 6:29 AM

Answers

  • Hi Folks,

    This might be of interest to folks on this thread.

    Cheers,
    -Isaac
    Saturday, May 31, 2008 1:53 AM

All replies

  • Its not clear from what is here what dimensions the polygon should have, if any. For example, do you want a LineString that passes through all the points? Does that curve close on itself? Or is a rectangle that encompasses all of the points good enough? If the latter, STEnvelope should fill the bill nicely, albeit not its not a very 'tight' fit except at the extremes:

     

    Code Snippet

    use scratch

    go

    drop table dbo.mpenner;

    go

    create table dbo.mpenner(shape geometry);

    go

    insert into dbo.mpenner values (

    geometry::STGeomFromText('MULTIPOINT ((8.023739416 7.89803918),(0.562674104 3.959565114),(-4.658036269 0.233642215),(-0.340279814 5.888079549),(-0.501298487 5.738512107),(-9.009625164 1.622385096),(0.91769339 7.706463706),(6.187174516 3.688827171),(-2.905915841 -1.962295621),(2.846264994 -5.236967242))',0))

    go

    select shape.STEnvelope().ToString(),shape.STEnvelope().STBoundary().ToString()

    from dbo.mpenner;

     

     

    This doesn't sound like what you want though, so can you be more precise as to what it is you are really looking for?

    Wednesday, May 14, 2008 12:40 PM
    Answerer
  • The problem is that there is not a single unique linestring or polygon which can be created from a bunch of points - it matters what order the points are defined in.
    It is obvious that a linestring from (0 0) to (5 10) to (25 30) is different from the linestring that goes from (0 0) to (25 30) to (5 10). However, it is also true that a polygon whose exterior ring is defined as (0 0), (0 5), (5 5), (5 0), (0 0) is different from one defined as (0 0), (5 0), (5 5), (0 5), (0 0). This significance of ring ordering is especially important for polygons using the geography datatype you mention, in order to resolve ambiguity of which side of the ring represents the area contained by the polygon....

    So, to answer your question, exactly how do you have your geography points stored? Are they stored as distinct POINT items in a geography column (in which case, is there another column which could be used to order those points in a logical way?) Or are they defined in a single MULTIPOINT element? In which case - are they already listed in the same order as would be used to create the linestring or polygon you are looking for?

    STConvexHull(), STEnvelope() etc. can be used to define broad shapes which describe the overall extent of the geometry, i.e. the spread of multiple poiint objects, but I don't think they will do exactly what you're looking for. However, if you can explain a little more about your data structure, there might be solutions involving manipulating the binary representation of the point instances that will be more efficient than Parse(), Lat(), and Long()

    Wednesday, May 14, 2008 4:21 PM
    Answerer
  •  

    I guess I should be more specific.  I would like a method to create a line from a series of ordered points or a polygon from the same.

     

    From your reponses I see that MultiPoint.ToLineString() or MultiPoint.ToPolygon methods don't make sense.

     

    However, a function similar to LineStringFromWKT like LineStringFromPoints('Point(lat long), Point(lat long), ...') would be helpful.

     

    Something where I could provide a set of my own ordered points, either through a text string, or better yet, from a SQL statement that would be great.

     

    These points are already in a SQL table so retrieving them in an ordered fashion is not a problem.

    Wednesday, May 14, 2008 4:47 PM
  • Okay, so, let say that the order of points doesn't really matter for a LineString, or, more correctly, that you are able to order them appropriately from your source.

     

    Considering a setup like this:

     

    Code Snippet

    use scratch

    go

    drop table mpenner.points;

    drop schema mpenner;

    drop aggregate dbo.BuildLineStringSrid4326

    drop assembly BuildLineString

    go

    create schema mpenner

    go

    create table mpenner.points(ID int identity(1,1) primary key,descr varchar(50),lat float,lon float);

    go

    set nocount on;

    insert into mpenner.points values

    ('Joe Foss Airfield - Sioux Falls',43.5846,-96.7421),

    ('Chicago O''Hare International Airport',41.9788,-87.9089),

    ('Boston Logan Airport',42.3614,-71.0135),

    ('Los Angles International Airport',33.9441,-118.408),

    ('Denver International Airport',39.8581,-104.694),

    ('Epply Airfield - Omaha',41.303,-95.9014),

    ('Hartsfiled Airport - Atlanta',33.6402,-84.4267),

    ('Orlando International Airport',28.431,-81.3091),

    ('New Orleans International Airport',29.9934,-90.258),

    ('Dallas-Ft.Worth International Airport',32.8991,-97.0483);

    go

     

     

    We can take advantage of the improvements in SQL Server 2008 User-Defined SQL/CLR aggregation feature to build your line string. Note that this is somewhat of a naive implementation -- I really would prefer to working with WKB than WKT for efficiency reasons and code could be adopted to that end. What we can do is write a UDA that takes in your points (the ability to accept more than one parameter in a UDA is new in 2008) and returns an Geometry instance. The code for that looks like this:

     

    Code Snippet
    using System;
    using System.Data;
    using System.Data.SqlTypes;
    using Microsoft.SqlServer.Server;
    using Microsoft.SqlServer.Types;
    using System.Text;
    namespace DM.Examples {
      [Serializable]
      [SqlUserDefinedAggregate(Format.UserDefined, IsInvariantToNulls = true,MaxByteSize=-1)]
      public struct LineStringBuilder : IBinarySerialize {
        StringBuilder sb;
        public void Init() {
          sb = new StringBuilder();
        }
        [SqlMethod(OnNullCall = false)]
        public void Accumulate(SqlDouble Lat, SqlDouble Lon) {
          sb.Append(Lat);
          sb.Append(' ');
          sb.Append(Lon);
          sb.Append(',');
        }
        public void Merge(LineStringBuilder Group) {
          this.sb.Append(Group.sb);
        }
        public SqlGeometry Terminate() {
          SqlChars c = new SqlChars("LINESTRING(" + sb.ToString().TrimEnd(',') + ')');
          return SqlGeometry.STLineFromText(c, 4326);
        }
        public void Read(System.IO.BinaryReader r) {
          sb = new StringBuilder(r.ReadString());
        }
        public void Write(System.IO.BinaryWriter w) {
          w.Write(sb.ToString());
        }
      }
    };

     

     

    Next we need to catalog with compiled asssembly into SQL Server. The deploy feature in VS90 (and with the VS90 SP1 beta) has issues with this, so I'll do it manually here:

     

    Code Snippet

    create assembly BuildLineString

    from '[elided]\BuildLineStrring\bin\Debug\BuildLineStrring.dll'

    go

    create aggregate dbo.BuildLineStringSrid4326(@lat float,@lon float)

    returns geometry

    external name BuildLineString.[DM.Examples.LineStringBuilder];

    go

     

     

    Then we just need to call it. I'll assume the default ordering of the points is appropriate here.

     

    Code Snippet

    select dbo.BuildLineStringSrid4326(lat,lon).ToString() from mpenner.points;

     

    Hopefully that'll give you a nudge in the right direction.

     

    Wednesday, May 14, 2008 6:24 PM
    Answerer
  • Alternatively, you could use a cursor to loop through the recordset, extracting the binary values of each point and then building up a new WKB definition as follows:

    First, create some test data:
    Code Snippet

    CREATE TABLE #BigListOfPoints (
      geometry geometry
    )

    INSERT INTO #BigListOfPoints VALUES
    (geometry::STPointFromText('POINT (51.48642 -0.12367)',4326)),
    (geometry::STPointFromText('POINT (51.48631 -0.12062)',4326)),
    (geometry::STPointFromText('POINT (51.48653 -0.11917)',4326)),
    (geometry::STPointFromText('POINT (51.48667 -0.11822)',4326)),
    (geometry::STPointFromText('POINT (51.48694 -0.11653)',4326)),
    (geometry::STPointFromText('POINT (51.48703 -0.11508)',4326)),
    (geometry::STPointFromText('POINT (51.48722 -0.11284)',4326)),
    (geometry::STPointFromText('POINT (51.48874 -0.11125)',4326))


    Now, use a cursor to build a binary variable @Coordinatestream containing a stream of each of the co-ordinate values from their WKB representation. Since a point contains exactly two co-ordinates, we know that this involves taking 16 bytes of data, starting at byte 6.
    Code Snippet

    DECLARE @Point geometry
    DECLARE @CoordinateStream varbinary(max)
    DECLARE @i int
    SET @i = 0

    DECLARE GeomCursor CURSOR FOR SELECT geometry FROM #BigListOfPoints WITH (NOLOCK);
    OPEN GeomCursor;
    FETCH NEXT FROM GeomCursor INTO @Point;
    WHILE @@FETCH_STATUS = 0
        BEGIN
            IF (@i=0) BEGIN
                SET @CoordinateStream = SUBSTRING(@Point.STAsBinary(),6,16)
            END
            ELSE BEGIN
                SET @CoordinateStream = @CoordinateStream + SUBSTRING(@Point.STAsBinary(),6,16)
            END;
            SET @i = @i+1
            FETCH NEXT FROM GeomCursor INTO @Point;
        END;
    CLOSE GeomCursor;
    DEALLOCATE GeomCursor;


    To build a WKB linestring from this 'coordinatestream', you can do the following:
    Code Snippet

    SELECT geometry::STGeomFromWKB(
     0x01 -- Byte order is big-endian (since this is what STAsBinary() outputs)
     + 0x02000000 -- This is a 4-byte value showing that we are building a line
     + CONVERT(varbinary,REVERSE(CONVERT(varbinary,@i))) -- This is the number of points in the line which we get from our counter. We have to REVERSE it to get big-endian
     + @CoordinateStream -- This is the stream of point co-ordinates
    ,0).STAsText()

    Result is:
    LINESTRING (51.48642 -0.12367, 51.48631 -0.12062, 51.48653 -0.11917, 51.48667 -0.11822, 51.48694 -0.11653, 51.48703 -0.11508, 51.48722 -0.11284, 51.48874 -0.11125)

    Or, to build a polygon:
    Code Snippet

    SELECT geometry::STGeomFromWKB(
     0x01 -- Again, sepcify big-endian byte order
     + 0x03000000 -- This is a 4-byte value showing that we are building a polygon
     + 0x01000000 -- The number of rings in the polygon. Let's keep it simple and assume we only have an exterior ring
     + CONVERT(varbinary,REVERSE(CONVERT(varbinary,@i+1))) -- There is one more point in the polygon definition than in the original set of points, since we repeat the start/end point
     + @CoordinateStream + SUBSTRING(@CoordinateStream,1,16)-- We take the co-ordinate stream, and repeat the first point co-ordinates to close the ring
    ,0).STAsText()

    Result is:
    POLYGON ((51.48642 -0.12367, 51.48631 -0.12062, 51.48653 -0.11917, 51.48667 -0.11822, 51.48694 -0.11653, 51.48703 -0.11508, 51.48722 -0.11284, 51.48874 -0.11125, 51.48642 -0.12367))

    In order to change the order in which the points are defined, you simply have to change the SELECT statement in the cursor.

    I should also mention that for more information about how WKB is formed, you should see Simon Sabin's excellent blog post at: http://www.sqlskills.com/blogs/simon/2007/12/31/SQLServer2008SpatialUsingItInNetCode.aspx

    t.



    Wednesday, May 14, 2008 8:09 PM
    Answerer
  • Hi Folks,

    This might be of interest to folks on this thread.

    Cheers,
    -Isaac
    Saturday, May 31, 2008 1:53 AM
  • Hey Isaac,

     

    Thanks!  This is great!  This would certainly solve this problem in a great abstract and reusable way.  Would also solve many common problems such as outputting spatial objects as a GeoRSS feed for Virtual Earth.

     

    Thanks!

    Matt

    Saturday, May 31, 2008 7:15 AM
  • Excellent! I look forward to the next RC so we can start trying this out - hopefully there will be some really interesting classes which people can start developing and sharing to solve the sort of problems you mention...

    Will this be demoed at teched? And if so, is there any chance you will make the presentation available afterwards for those of us who can't make it across the atlantic?

    thanks.
    Saturday, May 31, 2008 8:23 PM
    Answerer
  • Hi,

     

    Yes, I will be demoing this a bit at TechEd.  I don't have a great track record for publishing presentations, but heckle me and I'll try to do so.  Smile

     

    And RC0 shouldn't be very far off, so you shall be playing with it soon.

     

    Cheers,

    -Isaac

     

    Saturday, May 31, 2008 11:53 PM
  •  

    I'm officially starting the heckling. 

     

    I won't be able to make it to TechEd this year so I'm definitely interested in any form of the presentation I can get my hands on.

     

    Looking forward to RC0!

     

    Thanks for the great communication Isaac.

    Matt

    Sunday, June 1, 2008 1:47 AM
  • ahem, heckle, heckle...
    Wednesday, September 28, 2011 10:33 AM
  • Hi all!

    First of all thank you for the solution you proposed it totally fits what I have to do. I know this is an old subject, but I am facing some problems with it...

    I have a table filled with coordinate like that :

    beam_id;sat_name;gain;lon;sat
    1;906z41downlink;-3.00;125.091;-20.655
    1;906z41downlink;-3.00;125.299;-20.632
    1;906z41downlink;-3.00;126.144;-20.593
    1;906z41downlink;-3.00;127.337;-20.494
    1;906z41downlink;-3.00;134.841;-19.253
    (... d autres points avec le même beam_id ...)
    2;906z41downlink;-1.00;141.187;24.476
    2;906z41downlink;-1.00;141.298;24.487
    2;906z41downlink;-1.00;141.443;24.501
    2;906z41downlink;-1.00;142.351;24.596
    2;906z41downlink;-1.00;144.550;25.075
    2;906z41downlink;-1.00;144.529;25.332
    (... other points with the same beam_id ...)
    3;906z41downlink;-4.00;105.037;-20.751
    3;906z41downlink;-4.00;104.821;-20.626
    3;906z41downlink;-4.00;104.821;-20.626
    (...again and again...)

    In that way I have to create as much polygon as the beam_id value is at the end of the TABLE.

    I have modified a little bit your code tanoshimi. When the beam_id changes it means I have to create a new polygon. But I am not sure I am doing it right because I have a FormatException saying that the "WKB is not valid".

    What did I do wrong?

    Thank you very much for you help!! :)

    -- Table filled with coordinate
    CREATE TABLE sat_beam (
    	beam_id int NOT NULL,
    	sat_ref char(50) NOT NULL,
    	gain real,	
    	longitude real,
    	latitude real
    );
     
    -- Add coordinates to sat_beam
    ALTER TABLE sat_beam ADD coordinate geography;
     
    UPDATE sat_beam
       SET coordinate = geography::STPointFromText('POINT(' + CAST(longitude AS VARCHAR(20)) + ' ' + 
                        CAST(latitude AS VARCHAR(20)) + ')', 4326);
    GO
    ----------------------------
    --- Create polygon perimeters with sat_beal Point data
    DECLARE @Point geography
    DECLARE @CoordinateStream varbinary(max)
    DECLARE @i int
    DECLARE @valb int
    DECLARE @Bid int
    DECLARE @RefS char(50)
    DECLARE @Gain real
    DECLARE @PrevRefS char(50)
    DECLARE @PrevGain real
    DECLARE @Cpt int
    SET @i = 0
    SET @valb = 2 -- to compare with beam_id
    SET @Cpt = 0
    DECLARE GeomCursor CURSOR FOR SELECT beam_id,sat_ref,gain,coordinate FROM sat_beam WITH (NOLOCK);
    OPEN GeomCursor;
    FETCH NEXT FROM GeomCursor INTO @Bid,@RefS,@Gain,@Point;
    WHILE @@FETCH_STATUS = 0
        BEGIN
    		-- First time the cursor look for the table
    		IF (@i=0) BEGIN
                SET @CoordinateStream = SUBSTRING(@Point.STAsBinary(),6,16)
            END
    		-- Change of beam ID value
    		IF @valb = @Bid BEGIN
    			INSERT INTO ListOfBeam ( nb_points ,sat_ref, gain, stream ) VALUES ( @i-@Cpt,@PrevRefS,@PrevGain,@CoordinateStream )
    			
    			-- Init and update for new beam
    			SET @CoordinateStream = NULL
    			SET @CoordinateStream = SUBSTRING(@Point.STAsBinary(),6,16)
    			SET @valb = @valb + 1
    			SET @Cpt = @i
    		END
    		ELSE BEGIN
    			SET @CoordinateStream = @CoordinateStream + SUBSTRING(@Point.STAsBinary(),6,16)
    		END;
    		SET @i = @i+1
    		SET @PrevRefS = @RefS
    		SET @PrevGain = @Gain
    		FETCH NEXT FROM GeomCursor INTO @Bid,@RefS,@Gain,@Point;
        END;
    CLOSE GeomCursor;
    DEALLOCATE GeomCursor;
    GO
    -- To build a WKB polygon from this 'coordinatestream'
    ALTER TABLE ListOfBeam ADD perimeter geography;
    GO
    -- Insert the created beam perimeter to ListOfBeam
    SELECT geography::STPolyFromWKB(
     0x01
     + 0x03000000
     + 0x01000000
     + CONVERT(varbinary,REVERSE(CONVERT(varbinary,[nb_points]+1)))
     + stream + SUBSTRING(stream,1,16)
     ,4326).STAsText()


    • Edited by Slevihn Tuesday, April 16, 2013 4:41 PM
    Tuesday, April 16, 2013 2:35 PM
  • I achieved to solve my problem myself!! :)

    It appears that the definition of my input files of the beams (to create polygons) are variously defined. The first point is not always repeated at the end and that creates wrong POLYGON definition.

    So the solution above works!! I have just modified a few things to optimize my work on the created polygons.

    Thanks again for what you proposed many years from now; it is still very usefull!! ;)

    Louis

    Wednesday, April 17, 2013 6:05 PM