https://traffic-viz.github.io/scenarios/windfield.html

traffic Logo
[                    ] 
  * Installation
  * Quickstart
  * Core structure
  * Sources of data
  * Algorithms
  * Exporting and storing data
  * Advanced usage
  * Example Gallery
  * Scenarios and use cases
      + Is this a plane?
      + A Mode S based wind field
      + Calibration flights for VOR and ILS systems
      + Geomagnetic declination in Mode S data
      + Sun elevation along a flight
      + Impact of COVID-19 on worldwide aviation
      + Flight density heatmaps
  * Publications

traffic

  *  >>
  * Scenarios and use cases >>
  * A Mode S based wind field
  * View page source

---------------------------------------------------------------------

A Mode S based wind fieldP

ADS-B data broadcasted by aircraft contains information about
groundspeed and true track angle of its trajectory. When proper
requests are sent by a Secondary surveillance radar, aircraft also
send more information with true airspeed and heading angle in
specific BDS5,0 and BDS6,0 messages (see Junzi Sun's website)

Groundspeeds and true track angles are derived from the GNSS
positions whereas true airspeed is computed with traditional onboard
instruments like Pitot tubes.

Here, aircraft behave as a distributed network of moving sensors and
researchers 1 2 have been recommending methods to derive wind fields
from Mode S data. The following method is a very basic approach to
compute wind.

Please note I am not a meteorology specialist, so feel free to
improve this page where it should.


The example above is wind averaged between 25degW and 55degE and between
32degN and 65degN, from FL200 and above on February 23th 2019, between
14:00 and 16:30 UTC.

from traffic.core import Traffic

t = Traffic.from_file("<fill here>")

t_extended = (
    traffic
    # download and decode EHS messages (DF 20/21)
    .query_ehs()
    # resample/interpolate one sample per second
    .resample("1s")
    # median filters
    .filter(altitude=23, track=53, heading=53, groundspeed=53, TAS=53)
    # wind triangle computation
    .compute_wind()
    # median filter
    .filter(wind_u=53, wind_v=53)
    # resample one sample per minute
    .resample("1T")
    # do not use multiprocessing to avoid denial of service
    .eval(desc="preprocessing")
)

# t_extended.to_pickle("wind_backup.pkl")

The result of this computation is a set of trajectories: each
aircraft yields one point per minute with a 4D-position (timestamp,
latitude, longitude, altitude) and a wind vector decomposed along a
zonal speed (wind_u) and a meridional speed (wind_v).

We then use ipyleaflet to display the wind field. The Velocity widget
requires two 2D matrices of zonal and meridional components of the
wind. The following method rounds lat/lon coordinates to the closest
integer and average wind in each resulting cell.

import xarray as xr

def compute_grid(traffic: Traffic) -> xr.Dataset:

    avg = (
        traffic
        # remove NaN values, just in case
        .query("wind_u == wind_u")
        # prepare coordinates for the 4d-grid, also remove NaN in wind
        .assign(
            # round coordinates to the closest .33 latitude/longitude
            lat_=lambda df: (3 * df.latitude.round(0)) / 3,
            lon_=lambda df: (3 * df.longitude.round(0)) / 3,
            # This basic version averages on all altitudes/timeranges
            # but it is easy to use the following fields to display
            # wind fields in particular time ranges and altitude levels.
            alt_=lambda df: (3e-3 * df.altitude).round(0) / 3e-3,
            hour=lambda df: df.timestamp.dt.round("h"),
        )
        # compute the average wind
        .data[["wind_u", "wind_v", "lat_", "lon_"]]
        .groupby(["lat_", "lon_"])
        .mean()
    )

    # Unstack then fill the holes where possible (2D interpolation)
    u = avg[["wind_u"]].unstack().interpolate().values
    v = avg[["wind_v"]].unstack().interpolate().values

    return xr.Dataset(
        data_vars={
            "u_wind": xr.DataArray(u, coords=avg.index.levels),
            "v_wind": xr.DataArray(v, coords=avg.index.levels),
        }
    )

The following is a basic rendering delegated to ipyleaflet library.

from ipyleaflet import Map, Velocity, basemaps

# t_extended = Traffic.from_file("wind_backup.pkl")

map_ = Map(
    center=(52, 15),
    zoom=4,
    interpolation="nearest",
    basemap=basemaps.CartoDB.DarkMatter,
)

wind = Velocity(
    data=compute_grid(t_extended),
    zonal_speed="u_wind",
    meridional_speed="v_wind",
    latitude_dimension="lat_",
    longitude_dimension="lon_",
    velocity_scale=0.002,
    max_velocity=150,
)

map_.add_layer(wind)

map_

1
    Hurter, C., R. Alligier, D. Gianazza, S. Puechmorel, G.
    Andrienko, and N. Andrienko.
    << Wind Parameters Extraction from Aircraft Trajectories >>.
    Computers, Environment and Urban Systems 47 (2014): 28-43.
    https://doi.org/10.1016/j.compenvurbsys.2014.01.005.
2
    Sun, Junzi, Huy Vu, Joost Ellerbroek, and Jacco Hoekstra.
    << Ground-Based Wind Field Construction from Mode-S and ADS-B Data
    with a Novel Gas Particle Model >>, 2017, 9.

Next Previous
---------------------------------------------------------------------

(c) Copyright 2018, Xavier Olive.

Built with Sphinx using a theme provided by Read the Docs.