[postgis-users] Check If Point Is In New Geometry

Tyler DeWitt tyler.scott.dewitt at gmail.com
Tue Oct 29 08:45:32 PDT 2013


Hi Brent,

I'm putting this back on the list so people can see our continued discussion.

I'm confused by what ST_Dump() does in this example.  Could you explain a little more what it is doing?

I agree with you and Paul that windowing might be my best approach.  I'm looking at how I might accomplish my goals with windowing now.  Do you have much experience windowing? It might be useful to talk about some ideas.

As far as points not extending to the boundaries, I'm not too concerned.  I actually collect 10 points a second and a couple inches here and there is not going to affect me.  I just need to know big picture, "Which zones are a person in?".

Also, do you think it would be useful to turn the segments and line tracks into views instead of tables so that they are updated as someone adds/removes points from the gps_tracks table?

Thanks,
Tyler

On Oct 28, 2013, at 1:22 PM, Brent Wood <pcreso at pcreso.com> wrote:

> Hi Tyler,
> 
> I missed out pasting this bit in the last email:
> 
> # add a record identifier for grouping by later
> psql -d gps -c "alter table segments
>                 add column id  serial  primary key;"
> 
> just before the last SQL.
> 
> This does not do everything you might need, but I figured it might provide some ideas...
> 
> The timestamp of the last point was indeed the same, but it was outside the polygon - to check the spatial & temporal filters work OK.
> 
> I have changed the SQL to generate the segments - it now saves the output from ST_Dump() so each traversal of a polygon is a distinct linestring. This has a few issues in this context as you can see... segments are clipped where they cross each other to creatre only simple linestrings - which is not what you want.
> 
> # get intersections
> psql -d gps -c "create table segments as
>                 SELECT t.user_id,
>                        t.track_date,
>                        p.name,
>                        (ST_Dump(ST_Intersection(t.trackline, p.poly))).geom as segment
>                 FROM trackline t,
>                      polygons p;"
> 
> # add a record identifier for grouping by later
> psql -d gps -c "alter table segments
>                 add column id  serial  primary key;"
> 
> # get time of arrival & no minutes in each polygon
> psql -d gps -c "SELECT s.id,
>                        s.user_id,
>                        s.name,
>                        s.track_date,
>                        min(g.point_time) as arrival_time,
>                        (ST_NPoints(segment)/60.0)::decimal(6,2) as minutes
>                 FROM segments s,
>                      gps_points g
>                 WHERE ST_intersects(g.gps_point,s.segment)
>                 GROUP BY s.id,
>                          s.user_id,
>                          s.name,
>                          s.track_date
>                 ORDER BY arrival_time;"
> 
> 
> I had hoped that using generic Postgis functions in this way might be able to do what you want... but I'm beginning to think they will not, and a function to navigate along the points in the 24hr trackline building up the segment as you go might be more robust. One issue here is that such segments will start & end at the first/last points in each polygon rather than extending to the polygon boundaries. As Paul suggested - windowing might work...
> 
> I have worked with vessel (at sea) gps data in postgis for several years now, and have never trid to do quite what you are :-)
> 
> Brent
> 
> From: Tyler DeWitt <tyler.scott.dewitt at gmail.com>
> To: Brent Wood <pcreso at pcreso.com> 
> Sent: Tuesday, October 29, 2013 3:26 AM
> Subject: Re: [postgis-users] Check If Point Is In New Geometry
> 
> Brent,
> 
> I had to make a slight change to the script:
> 
> # get time of arrival & no minutes in each polygon
> psql -d gps -c "SELECT s.user_id,
>                        s.name,
>                        s.track_date,
>                        min(g.point_time) as arrival_time,
>                        (ST_NPoints(segment)/60.0)::decimal(6,2) as minutes
>                 FROM segments s,
>                      gps_points g
>                 WHERE ST_intersects(g.gps_point,s.segment)
>                 GROUP BY s.segment,
>                          s.user_id,
>                          s.name,
>                          s.track_date
>                 ORDER BY arrival_time;"
> 
> The GROUP BY included a reference to s.id, which doesn't exist.  I think you meant s.segment?
> 
> 
> With that change, the query runs and I get the correct time spent in each zone, but when I tried adding more points to poly1 (after poly2), I am not getting a separate segment, it just adds onto the poly1 segment, so even though the points travel from poly1 to poly2 and back to poly1, I am only seeing 2 segments.
> 
> Here are the points I added:
> insert into gps_points values
>                 (default,
>                  1,
>                  '2013-10-01 00:00:30'::timestamp,
>                  ST_SetSRID(ST_MakePoint(0.01,0.01),4326));
> 
> insert into gps_points values
>                 (default,
>                  1,
>                  '2013-10-01 00:00:31'::timestamp,
>                  ST_SetSRID(ST_MakePoint(0.02,0.02),4326));
> 
> Also, the last gps_point you inserted had the same timestamp as another point, so I changed its time (don't think that has any real affect).
> 
> Thanks,
> Tyler
> 
> On Oct 27, 2013, at 7:46 PM, Brent Wood <pcreso at pcreso.com> wrote:
> 
>> Tyler,
>> 
>> Convert your 24 hours worth of points to a linestring
>> Clip the linestring by your polygons to generate separate lines per polygon per entry tinto the polygon
>> Count the vertices in each linestring to determine how long each period in each polygon was.
>> 
>> The following shell script might give you some ideas... it works fine here, but is a bit simplistic in terms of input data...
>> 
>> Cheers,
>> 
>>    Brent Wood
>> 
>> 
>> #! /bin/bash
>> # gps points to track, clipped by polygons, B Wood, Oct 2013
>> 
>> # create new database to test
>> dropdb gps
>> createdb gps
>> psql -d gps -c "create extension postgis;" 
>> 
>> # create table with some gps points
>> psql -d gps -c "create table gps_points
>>                 (id         serial primary key,
>>                  user_id    integer,
>>                  point_time timestamp,
>>                  gps_point  geometry(point,4326));"
>> 
>> psql -d gps -c "insert into gps_points values
>>                 (default,
>>                  1,
>>                  '2013-10-01 00:00:00'::timestamp,
>>                  ST_SetSRID(ST_MakePoint(0.0,0.0),4326));" 
>> psql -d gps -c "insert into gps_points values
>>                 (default,
>>                  1,
>>                  '2013-10-01 00:00:01'::timestamp,
>>                  ST_SetSRID(ST_MakePoint(0.01,0.01),4326));" 
>> psql -d gps -c "insert into gps_points values
>>                 (default,
>>                  1,
>>                  '2013-10-01 00:00:02'::timestamp,
>>                  ST_SetSRID(ST_MakePoint(0.02,0.02),4326));" 
>> psql -d gps -c "insert into gps_points values
>>                 (default,
>>                  1,
>>                  '2013-10-01 00:00:03'::timestamp,
>>                  ST_SetSRID(ST_MakePoint(0.03,0.03),4326));" 
>> psql -d gps -c "insert into gps_points values
>>                 (default,
>>                  1,
>>                  '2013-10-01 00:00:04'::timestamp,
>>                  ST_SetSRID(ST_MakePoint(0.04,0.04),4326));" 
>> psql -d gps -c "insert into gps_points values
>>                 (default,
>>                  1,
>>                  '2013-10-01 00:00:05'::timestamp,
>>                  ST_SetSRID(ST_MakePoint(0.05,0.05),4326));" 
>> psql -d gps -c "insert into gps_points values
>>                 (default,
>>                  1,
>>                  '2013-10-01 00:00:06'::timestamp,
>>                  ST_SetSRID(ST_MakePoint(0.06,0.06),4326));" 
>> psql -d gps -c "insert into gps_points values
>>                 (default,
>>                  1,
>>                  '2013-10-01 00:00:07'::timestamp,
>>                  ST_SetSRID(ST_MakePoint(0.07,0.07),4326));" 
>> psql -d gps -c "insert into gps_points values
>>                 (default,
>>                  1,
>>                  '2013-10-01 00:00:08'::timestamp,
>>                  ST_SetSRID(ST_MakePoint(0.08,0.08),4326));" 
>> psql -d gps -c "insert into gps_points values
>>                 (default,
>>                  1,
>>                  '2013-10-01 00:00:09'::timestamp,
>>                  ST_SetSRID(ST_MakePoint(0.09,0.09),4326));" 
>> 
>> psql -d gps -c "insert into gps_points values
>>                 (default,
>>                  1,
>>                  '2013-10-01 00:00:10'::timestamp,
>>                  ST_SetSRID(ST_MakePoint(0.10,0.10),4326));" 
>> 
>> psql -d gps -c "insert into gps_points values
>>                 (default,
>>                  1,
>>                  '2013-10-01 00:00:11'::timestamp,
>>                  ST_SetSRID(ST_MakePoint(0.11,0.11),4326));" 
>> 
>> psql -d gps -c "insert into gps_points values
>>                 (default,
>>                  1,
>>                  '2013-10-02 00:00:10'::timestamp,
>>                  ST_SetSRID(ST_MakePoint(0.12,0.12),4326));" 
>> 
>> 
>> # create table & some polygons
>> psql -d gps -c "create table polygons
>>                 ( id     serial primary key,
>>                   name   varchar(10),
>>                   poly   geometry(POLYGON,4326));"
>> 
>> psql -d gps -c "insert into polygons values
>>                 (default,
>>                  'poly1',
>>                  ST_PolygonFromText(
>>   'POLYGON((0.0 0.0, 0.0 0.05, 0.05 0.05, 0.05 0.0, 0.0 0.0))',4326));"
>> 
>> psql -d gps -c "insert into polygons values
>>                 (default,
>>                  'poly2',
>>                  ST_PolygonFromText(
>>   'POLYGON((0.05 0.05, 0.05 0.1, 0.1 0.1, 0.1 0.05, 0.05 0.05))',4326));"
>> 
>> 
>> # query to get trackline
>> psql -d gps -c "create table trackline as 
>>                 SELECT gps.user_id, 
>>                        gps.point_time::date as track_date, 
>>                        ST_MakeLine(gps.gps_point) as trackline
>>                        FROM (SELECT user_id, gps_point, point_time 
>>                              FROM gps_points
>>                              WHERE user_id = 1
>>                                AND point_time::date = '2013-10-01'::date
>>                              ORDER BY point_time) as gps
>>                        GROUP BY gps.user_id,
>>                                 track_date;"
>> 
>> # get intersections
>> psql -d gps -c "create table segments as
>>                 SELECT t.user_id,
>>                        t.track_date,
>>                        p.name,
>>                        ST_Intersection(t.trackline, p.poly) as segment
>>                 FROM trackline t,
>>                      polygons p;"
>> 
>> # get time of arrival & no minutes in each polygon
>> psql -d gps -c "SELECT s.user_id,
>>                        s.name,
>>                        s.track_date,
>>                        min(g.point_time) as arrival_time,
>>                        (ST_NPoints(segment)/60.0)::decimal(6,2) as minutes
>>                 FROM segments s,
>>                      gps_points g
>>                 WHERE ST_intersects(g.gps_point,s.segment)
>>                 GROUP BY s.id,
>>                          s.user_id,
>>                          s.name,
>>                          s.track_date
>>                 ORDER BY arrival_time;"
>> 
>> 
>> 
>> 
>> 
>> From: Tyler DeWitt <tyler.scott.dewitt at gmail.com>
>> To: postgis-users at lists.osgeo.org 
>> Sent: Monday, October 28, 2013 7:52 AM
>> Subject: [postgis-users] Check If Point Is In New Geometry
>> 
>> I have a table of geometries (with a geom column of type MultiPolygon Z).
>> 
>> I collect real time location data (1 point a second) and store that in a tracked_points table (with a geom column of type PointZ, and user_id column of type int).
>> 
>> I'd like to know when a user (tracked by their point) enters a new/different geometry.
>> 
>> I can use ST_Contains(geometries.geom, tracked_points.geom) to figure out which geometry a user is currently in, but I'd like to find an efficient way to check a 24 hour period and say "User 1 went from geometry A to geometry B to Geometry C and back to Geometry B".  I could find how many points a user has in each geometry, and therefore how long they spent in each geometry, but I don't know about the path the user took.
>> 
>> 
>> Thanks,
>> Tyler
>> 
>> _______________________________________________
>> postgis-users mailing list
>> postgis-users at lists.osgeo.org
>> http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users
>> 
>> 
> 
> 
> 

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/postgis-users/attachments/20131029/004dc117/attachment.html>


More information about the postgis-users mailing list