trajectories

[espinielli]

Hi,
I am looking for a decent way to deal with 4D (3D + timestamp) / 7D (4D + speed vector) trajectory data in R since long.
I am super intrigued by the last section in the README where the 4D trajectory is mentioned:

The sphier system explicitly accommodates representations of temporal
dynamics. As the simplest level, a vector trajectory can be represented through
extension of the vertex component to include a t column. The maximum that
the sphier vertex table can contain is the four columns of (x, y, z,
t).

Could you elaborate more on it?

I have an intuition that silicate could help but terminology and R expertise fail me.
To give you a little bit of context my use case is conversion of ADS-B data,
i.e. position reports broadcast by aircraft, to trajectories and then analysis of trajectories
(intersection with airspaces, bundling in flows,...)

Any comments, ideas, clarification is welcom.

[mpadge]

Hi @espinielli, you're definitely in roughly the right place here, and @mdsumner is the man to answer your questions. This sphier repo was a place to work out some ideas on the way to silicate, which should be able to get to what you want, as your intuition rightly suggests. @mdsumner I'll be very keen to hear your answer on this, because i've not yet done any silicate trajectory stuff - what's the status? Workflow? Specific functions / tools?

[mdsumner]

I'd say that silicate can represent such a data set, and round-trip it from whatever format it's in to others - but there's really nothing built-in to work with it.

I'd be really interested to explore specific applications, why don't we get a shared example and start solving some problems?

The flight_tracks data set might give an idea of how I see silicate dealing with tracks, but I haven't pursued it since that was put in as an example data set. See how vertices exist in 4-space, 4D geometry. silicate knows about x_, y_, z_, and t_ but it only preserves those (and any others), the only thing that really does anything with them with anglr::plot3d

 flight_tracks
class       : PATH
type        : Sequential
vertices    : 18942 (4-space)
paths       : 144 (1-space)
coordinates : 18942
crs         : +proj=longlat +datum=WGS84 +no_defs
> plot(flight_tracks)
> SC0(flight_tracks)  ## convert to structural segment form 
class       : SC0
type        : Structural
vertices    : 18942 (4-space)
primitives  : 18798 (2-space)
crs         : +proj=longlat +datum=WGS84 +no_defs

anglr does something with a plot. The next step might be a plot4d() function that would animate this in time. I guess you want something like a discretization of airspace and a way to quantify time spent in 3D air-cells? Please feel free to expand on what you would want.

library(anglr)
 plot3d(flight_tracks); rgl::aspect3d(1, 1, 0.5)

[mdsumner]

oops and I see a print bug,

primitives  : 18798 (2-space)

segment primitives are 1-space because they are just lines, will fix

[mpadge]

@mdsumner Can you please transfer the issue over to the more appropriate repo of silicate, so we can continue the discussion there? I'm excited to see where this goes!

[mdsumner]

(good idea)

[mdsumner]

And my apologies @espinielli moving an issue is a pretty brutal move, I hope you don't mind! Happy to explore discussions

[espinielli]

@mdsumner and @mpadge thank you for your replies and no problem at all to move the issues where it best belong.

My use case is really a practical one: determine the portion of a flight trajectory in an airspace.
With Airspace I have one simple ground-projected polygon (no holes) with one altitude interval.
The example from parse_airspace_prisme in my WIP package {trrrj} gives you an idea:

> library(trrrj)
> library(magrittr)
> es <- system.file("extdata", "airspace_prisme_462.prisme", package = "trrrj")
> readr::read_lines(es) %>%
   parse_airspace_prisme()
Simple feature collection with 589 features and 4 fields
geometry type:  POLYGON
dimension:      XY
bbox:           xmin: 350.6942 ymin: 173.0292 xmax: 364.0546 ymax: 187.7911
epsg (SRID):    4326
proj4string:    +proj=longlat +datum=WGS84 +no_defs
First 10 features:
      unit airspace fl_min fl_max                       geometry
1  LSAZCTA  LSAZESL    125    245 POLYGON ((352.2102 187.0826...
2  LSAZCTA  LSAZESL      0    245 POLYGON ((410.0669 156.5964...
3  LSAZCTA  LSAZESL    105    245 POLYGON ((384.1205 193.2281...
4  LSAZCTA  LSAZESL      0    245 POLYGON ((402.759 165.3511,...
5  LSAZCTA  LSAZESL    125    245 POLYGON ((363.1386 181.2812...
6  LSAZCTA  LSAZWSL    115    245 POLYGON ((304.3561 184.9918...
7  LSAZCTA  LSAZWSL    125    245 POLYGON ((329.2045 178.9664...
8  LSAZCTA  LSAZWSL    125    245 POLYGON ((325.3235 145.3995...
9  LSGGTMA  LSGGSUD      0    195 POLYGON ((247.8045 80.88635...
10 LSAZCTA  LSAZESL    125    245 POLYGON ((362.0952 158.4151...

So I would be interested to find the interceptions of the trajectory with an airspace.
A trajectory then being a bunch of 4D points (in the simple case where I do not have speed vector).

I am thinking of taking segments (consecutive pairs of 4D points) and keep the ones completely internal to the airspace plus finding intersection with airspace using the bisection algorithm for the one segments that have one (and only one) of the 4D points outside.

Any thoughts, again, super welcome.

[espinielli]

I played with your code

library(anglr)
library(rgl)
library(silicate)

plot3d(flight_tracks)
rgl::aspect3d(1, 1, 0.5)
axes3d( edges=c("x--", "y--", "z") )

wit the addition of runways (maybe from osmdata) it would make it further awesome...!

But I am also/more interested in the analysis of the trajectories; once you have the parts belonging
to an airspace you can calculate the flight hours, the number of flights, bundle flows together, determine the level flight portions (which in the terminal area contribute to the fuel burn inefficiency due to not flying continuous descent operations), ...

...

[mdsumner]

I hear you, it's a good example and I hope to spend some time on it - thanks!

[mdsumner]

Hi @espinielli, do you know how to make sense of coordinates like 352.21019662192;187.082627470795?

If you have a reprex with reading an airspace and comparing it to a flight track that would be good. I'll delve into it over time, but those coordinates don't make sense atm.

[mdsumner]

Here's a version of flight tracks that puts on background image fwiw

library(anglr)
library(rgl)
library(silicate)

xy <- sc_coord(flight_tracks)  ## I happen to know this is longlat
im <- ceramic::cc_location(raster::extent(range(xy$x_), range(xy$y_)) + 0.5)
dem <- ceramic::cc_elevation(im)

mesh <- as.mesh3d(dem, image_texture = im)
tracks <- reproj::reproj(flight_tracks, 
                            raster::projection(im))

plot3d(tracks, lit = F)
plot3d(mesh, add = TRUE)
rgl::aspect3d(1, 1, 0.1)
axes3d( edges=c("x--", "y--", "z") )

(ceramic needs a mapbox key but you can replace with any raster-dem + raster-image)

[espinielli]

You found a bug in my code: I am not correctly converting from the internal format (.prisme) to a degree decimal...
Those coordinates are in minutes decimal.
Fixing in {trrrj} right now
Thanks!

[mdsumner]

nice one - I'd very much appreciate a reprex that gets some airspace polygons and some flights ready to go - but another question, how is the altitude component of the airspace/s defined?

Really nothing you can tell me "too simple", I have no idea about these data :)

[espinielli]

In aviation (for airspaces and above a certain altitude) you use flight levels, https://en.wikipedia.org/wiki/Flight_level