What is this PostGIS query doing to show great circle connections?

Question

I'm trying to draw a flight map in CartoDB. My starting point is a single lat long point. We'll pretend it is this:

ST_GeomFromText('POINT(-71.060316 48.432044)', 4326)

And I want draw lines between that point and the points in my CartoDB table. I had a hard time even getting started, but someone pointed me to this map:

http://bl.ocks.org/andrewxhill/raw/8406976/

Which is fancier than what I need (I don't need to redraw everytime you click someplace) but includes the basic query:

SELECT ST_Transform(ST_Segmentize(ST_MakeLine(
   ST_Transform(the_geom, 953027),
   ST_Transform(
          (SELECT ST_PointOnSurface(the_geom) 
          FROM lower_48_zips 
          WHERE cartodb_id = "+data.cartodb_id+"), 953027)), 100000), 4326) 
   the_geom FROM us_ski_areas WHERE mindist < 1000

Someone also pointed me towards an example which uses a very similar query:

UPDATE experimental.airroutes
SET the_geom =
(SELECT ST_Transform(ST_Segmentize(ST_MakeLine(
       ST_Transform(a.the_geom, 953027),
       ST_Transform(b.the_geom, 953027)
     ), 100000 ), 4326 )
 FROM experimental.airports a, experimental.airports b
 WHERE a.id = airroutes.source_id 
   AND b.id = airroutes.dest_id
);

But I'm having trouble pulling either of these apart. I'm decently conversant in MySQL but totally new to PostGIS and I would love help taking apart what is happening here and replacing source_id with a specific point.

I thought I could do something like

UPDATE my_table
SET the_geom =
(SELECT ST_Transform(ST_Segmentize(ST_MakeLine(
       ST_Transform(my_table.the_geom, 953027),
       ST_Transform(ST_GeomFromText('POINT(-71.060316 48.432044)', 4326), 953027)
     ), 100000 ), 4326 )
 FROM my_table
);

But I know I need a "where" clause in there and I'm not even sure what it should be.

Answer 1

Short answer first, you are super close. Try this instead,

UPDATE my_table
SET the_geom =
   ST_Transform(
     ST_Segmentize(
       ST_MakeLine(
         ST_Transform(my_table.the_geom, 953027),
         ST_Transform(CDB_LatLng(48.432044, -71.060316), 953027)
       ), 
       100000
     ), 
    4326 
   )

This would update your table from being points, to being a line from the original point to the point at -71, 48. I used the CDB_LatLng helper function. You didn't need to As_Text function, because that just turns geometries into human readable text.

Now the longer is that each of these functions has a documentation page, you can find those here,

ST_Transform: http://postgis.net/docs/ST_Transform.html

ST_Segmentize: http://postgis.net/docs/ST_Segmentize.html

ST_MakeLine: http://postgis.net/docs/ST_MakeLine.html

You are basically reprojecting the reference map to make it so that the straightest line between points follows more closely a curved globe (actually in this case I think it is conic, but minor detail). Next, you segmentize so that you have waypoints all along that line that follow the shortest path across the curved map. Then, your reproject back to WGS84. If you just grabbed the start and end points in the curved world, when you reproject back it would still just be a straight line. By segmentizing you get the shortest line plus waypoints along that line that use the curved map. When you reproject back to WGS84 then, it will still use the same way points which will now appear curved.

Hope that makes some sense.

Answer 2

The examples given, and the previous answer all depend on taking the line into a planar projection and then segmentizing. As of PostGIS 2.1 (which I think CartoDB is on) you can do a general-purpose job of it using ST_Segmentize(geography), so like this:

SELECT 
  ST_Segmentize(
  ST_MakeLine(
    ST_SetSRID(ST_MakePoint(-71.060316, 48.432044),4326), 
    pt
  )::geography
  , 10000
  )::geometry
FROM mytable;

In words,

• make a central point in SRS 4326
• make a line between that point and each point in the table
• convert that line into the geography type
• segmentize in geography space (along the great circle route)
• convert that detailed geography back into geometry