The STs* of PostGIS
* saints?* spatial types?
preso/
balticgit-2022/
where's that bus at?
(supposedly)(data by the Estonian Road Administration)
length: ~ 3.5 km
Based on current time calculate fraction of time between trip start and end and then
st_lineinterpolatepoint
the same amount
length: ~ 3.5 km, total time: 2 minutes
(but with intermediate stops this bus will never be on time)
A to B: 1.8 km, B to A: 1.7 km
A to B: 30 secs, B to A: 90 secs
(add impedance to traveltime at stops)
Adding acceleration/deceleration to transit location calculation.
(see here for details)
(see here for details)
accelerate / decelarate: 8 secs, stop 3 + 3 secs
Red: interpolate whole trip. Blue: interpolate using stops. Green: interpolate using stops and impeded time
example dashboard and
API
the projection tango
madagascatorI think I've earned my geonerd points today.
— ⚛⚛⚛Iván Sánchez⚛⚛⚛ (@RealIvanSanchez) April 29, 2022
(re: https://t.co/84PGLpHLDZ ) pic.twitter.com/MZGwI76zSG
select
*
from (
select
gid, (st_dumppoints(geom)).*
from
ne_10m_admin_0_countries ) f
;
st_dumppoints for all vertices
select
st_project(
st_setsrid(st_point(-30,0), 4326)::geography,
st_distance(st_setsrid(st_point(0,90),4326)::geography, geom::geography),
st_azimuth(st_setsrid(st_point(0,90),4326)::geography, geom::geography) - radians(125)
)::geometry(point, 4326) as geom,
row_number() over()::int as oid, gid, path
from (
select
gid, (st_dumppoints(geom)).*
from
ne_10m_admin_0_countries ) f
where
st_y(geom) >-88
;
st_project
a new point (
geography
) with the old distance and azimuth (minus some whatever constant)
Based on NaturalEarth
country boundaries, and proccessed in PostGIS using
this SQL
select
*
from (
select gid, (st_dumppoints(geom)).*
from ne_10m_admin_0_countries
) f
;
st_dumppoints for all vertices
select
st_buffer(
st_setsrid(st_point(0,0),4326)::geography,
st_distance(st_setsrid(st_point(0,90),4326)::geography, geom::geography) / pi()
)::geometry(polygon, 4326)
from (
select gid, (st_dumppoints(geom)).*
from ne_10m_admin_0_countries
) f
;
st_project
will not be of use here, let's pretend we're on a flat surface so
st_buffer instead
select
st_rotate(
st_exteriorring(
st_buffer(
st_setsrid(st_point(0,0),4326)::geography,
st_distance(st_setsrid(st_point(0,90),4326)::geography, geom::geography) / pi()
)::geometry(polygon, 4326)
),
radians(90)
)
from (
select gid, (st_dumppoints(geom)).*
from ne_10m_admin_0_countries
) f
;
Take the resulting st_exteriorring
and st_rotate it 90 degrees
counter-clockwise
select
st_scale(
st_rotate(
st_exteriorring(
st_buffer(
st_setsrid(st_point(0,0),4326)::geography,
st_distance(st_setsrid(st_point(0,90),4326)::geography, geom::geography) / pi()
)::geometry(polygon, 4326)
),
radians(90)
), 1.70, 1.70
)
from (
select gid, (st_dumppoints(geom)).*
from ne_10m_admin_0_countries
) f
where st_y(geom) >-85 and st_y(geom) < 85
;
select
st_lineinterpolatepoint(
st_scale(
st_rotate(
st_exteriorring(
st_buffer(
st_setsrid(st_point(0,0),4326)::geography,
st_distance(st_setsrid(st_point(0,90),4326)::geography, geom::geography) / pi()
)::geometry(polygon, 4326)
),
radians(90)
), 1.70, 1.70
),
st_azimuth(st_setsrid(st_point(0, 90),4326)::geography, geom::geography) / (2*pi())
) as geom, gid, path
from (
select gid, (st_dumppoints(geom)).*
from ne_10m_admin_0_countries
) f
where st_y(geom) >-85 and st_y(geom) < 85
;
And finally st_lineinterpolatepoint
from the start of the exterior line to the fraction defined by
st_azimuth between the
North Pole and the original vertice.
Based on NaturalEarth
country boundaries, and proccessed in PostGIS using
this SQL
#30DayMapChallenge
(on GitHub and Twitter)(data courtesy of Estonian Environment Agency)
Unknown Pleasures by Joy Division
(and air temperature forecast)
image: beyondvinyl.co.uk
select
st_force3d(
st_point(
x,
y + temperature * 1000
)
)
from (
select temperature, valid_from, valid_to, x, y from phenomen
) f
;
create a point (st_point)
shifting its y coordinate to the north by coefficent of temperature, and st_force3d.
select
st_setsrid(
st_translate(
st_force3d(
st_point(
x,
y + temperature * 1000
)
), 0, 0, temperature
), 3301
)
from (
select temperature, valid_from, valid_to, x, y from phenomen
) f
;
select
st_makeline(
array_agg(
st_setsrid(
st_translate(
st_force3d(
st_point(
x,
y + temperature * 1000
)
), 0, 0, temperature
), 3301
) order by x
)
)
from (
select temperature, valid_from, valid_to, x, y from phenomen
) f
group by y, valid_from, valid_to;
st_makeline of all points on the same y
select
st_chaikinsmoothing(
st_makeline(
array_agg(
st_setsrid(
st_translate(
st_force3d(
st_point(
x,
y + temperature * 1000
)
), 0, 0, temperature
), 3301
) order by x
)
), 5
) as geom, y, valid_from, valid_to
from (
select temperature, valid_from, valid_to, x, y from phenomen
) f
group by y, valid_from, valid_to;
apply st_chaikinsmoothing on the lines
(sure, why not but for that we'll need to ...)
with segments as (
select
oid, st_dumppoints(geom) as pt
from
chaikinlines_temperature
)
select
*
from
segments
;
st_dumppoints
for all points of st_chaikinsmoothing-ed
lines.
with segments as (
select
st_makeline(
lag((pt).geom, 1, null) over(
partition by y, valid_from order by valid_from, y, (pt).path
),
(pt).geom
) as geom, valid_from, valid_to
from (
select
oid, st_dumppoints(geom) as pt
from
chaikinlines_temperature
) as dumps
)
select
*
from
segments
where
geom is not null;
st_makeline
of two consequtive points on the same y
with segments as (
select
st_makeline(
lag((pt).geom, 1, null) over(
partition by y, valid_from order by valid_from, y, (pt).path
),
(pt).geom
) as geom, valid_from, valid_to
from (
select
oid, st_dumppoints(geom) as pt
from
chaikinlines_temperature
) as dumps
)
select
geom, st_z(st_lineinterpolatepoint(geom, 0.5)) as temperature, valid_from, valid_to
from
segments
where
geom is not null;
take the midpoint of the segment with
st_lineinterpolatepoint
and its (interpolated)
st_z as temperature value
Polyps
(and wind speed + direction forecast)
- eccentricity for wind speed
- turn angle for wind direction
select
*
from (
select
st_transform(
st_project(
grid::geography, coalesce(wind_speed_mps,0)*500, radians(coalesce(wind_direction_deg, 0))-pi()
)::geometry(point, 4326),
3301) as geom, wind_speed_mps, wind_direction_deg, valid_from, valid_to
from
phenomen
) f;
shift grid point off its location by a factor of speed in the wind direction using
st_project and then
st_transform to
epsg:3301
select
st_ellipse(
st_x(geom), st_y(geom),
2500, 2500 + (coalesce(wind_speed_mps,0)*1000), 0, 16
)
from (
select
st_transform(
st_project(
grid::geography, coalesce(wind_speed_mps,0)*500, radians(coalesce(wind_direction_deg, 0))-pi()
)::geometry(point, 4326),
3301) as geom, wind_speed_mps, wind_direction_deg, valid_from, valid_to
from
phenomen
) f;
create a st_ellipse
in the previously shifted position with the longer axis as a factor of wind speed
select
st_rotate(
st_setsrid(
st_ellipse(
st_x(geom), st_y(geom),
2500, 2500 + (coalesce(wind_speed_mps,0)*1000), 0, 16
),
3301
),
-1*(radians(wind_direction_deg) + pi()),
geom
) as geom, wind_speed_mps, wind_direction_deg, valid_from, valid_to
from (
select
st_transform(
st_project(
grid::geography, coalesce(wind_speed_mps,0)*500, radians(coalesce(wind_direction_deg, 0))-pi()
)::geometry(point, 4326),
3301) as geom, wind_speed_mps, wind_direction_deg, valid_from, valid_to
from
phenomen
) f;
st_rotate
the ellipse in the desired direction anchoring on the shifted position.
(me going full-on estonian)
as well useful in real life
proceed to the directors cut
tonis.kardi@gmail.com
https://tkardi.ee