Skip to content

Instantly share code, notes, and snippets.

@pecard
Last active February 20, 2021 11:14
Show Gist options
  • Save pecard/3b8424677ee43ca69785e3f52a1e916b to your computer and use it in GitHub Desktop.
Save pecard/3b8424677ee43ca69785e3f52a1e916b to your computer and use it in GitHub Desktop.
Wind trajectory with ECMWF ERA-5 hourly wind data

GEE Code script

/**
 * @license
 * Copyright 2020 Google LLC.
 * SPDX-License-Identifier: Apache-2.0
 * 
 * Generates wind trajectory maps and animations from ECMWF ERA-5 hourly wind
 * data. A set of random points are drawn within an area of interest; a path is
 * generated for each point by iteratively moving it according to wind speed and
 * direction for a defined number of steps. Each step (segment) represents a
 * given duration in minutes. Wind conditions vary across space (~11 km pixels)
 * and time (1-hour cadance).
 
 Original Code provided by Justin Braaten [here](https://code.earthengine.google.com/d327f1117f6feafcf8dfcf961510ef39) and
 
[here](https://code.earthengine.google.com/917e5be419e988dfdd29fae77b960445)
 
 */


// #############################################################################
// ### INPUTS ###
// #############################################################################

// region to include in animation AOI
var AOI = ee.Geometry.BBox(-9.3617, 38.5565, -8.6256, 39.2405);   // region to sample initialization points from
// region to sample initialization points from ANIMATION_REGION
var ANIMATION_REGION = ee.Geometry.BBox(-9.3198, 38.5196, -8.6633, 9.1486);  // region to include in animation

var TARGET_TIME = ee.Date('2020-10-01T06:00:00');  // initialization time
var SEG_DURATION = 2;  // minutes per wind path segment
var N_SEGMENTS = 50;  // number wind path segments 
var N_POINTS = 150;  // number of sample points
var CRS = 'EPSG:3763';  // animation projection
var MAX_DIM = 600;  // max animation size; reduce if too many pixel error occurs
var FPS = 16;  // animation frames per second

// #############################################################################
// Filter the wind collection by time; add date range as an image property.
var endTime = TARGET_TIME.advance(N_SEGMENTS*SEG_DURATION, 'minutes');
var wind = ee.ImageCollection('ECMWF/ERA5_LAND/HOURLY')
  .filter(ee.Filter.date(TARGET_TIME, endTime))
  .map(function(img) {
    var start = img.get('system:time_start');
    var end = img.get('system:time_end');
    return img.set('dateRange', ee.DateRange(start, end));
  })
  .select(['u_component_of_wind_10m', 'v_component_of_wind_10m'],
          ['east_wind', 'north_wind']);

// Get wind dataset scale for use in sampling.
var scale = wind.first().projection().nominalScale();

// Sample random points from the wind dataset within the AOI.
var pts = wind.first().sample({
  region: AOI,
  scale: scale,
  numPixels: N_POINTS,
  geometries: true
});

// Get a series of segment labels to define iteration.
var segmentSeq = ee.List.sequence(1, N_SEGMENTS);
var segmentInfo = segmentSeq.map(function(segment) {
    segment = ee.Number(segment).subtract(1);
    var currentTime = TARGET_TIME.advance(ee.Number(SEG_DURATION).multiply(segment), 'minutes');
    var thisWind = wind.filter(ee.Filter.dateRangeContains({
      leftField: 'dateRange',
      rightValue: currentTime
    })).first();
    return ee.List([segment, thisWind]);
  });

// Converts meters per minutes to degrees per minutes.
function metersToDegrees(meters, minutes) {
  return ee.Number(meters)
    .multiply(60 * minutes)  // m/s > m/n minutes
    .divide(1000)  // m/n minutes > km/n minutes
    .divide(111);  // km/n minutes > degree/n minutes
}

// For each random point, generate a path from segments.
var pathsList = pts.toList(pts.size()).map(function(pt) {
  // Initialize a path segment list.
  pt = ee.Feature(pt);
  var initSegment = ee.Feature(
    ee.Geometry.LineString([pt.geometry(), pt.geometry()]), {segment: 0});
  var segmentList = ee.List([initSegment]);

  // Define a function to move a point and append path segments to path list.
  function moveIt(segment, segList) {
    segment = ee.List(segment);
    var segmentLabel = ee.Number(segment.get(0));
    var segmentWind = ee.Image(segment.get(1)).unmask(0);  // set ocean as 0
    
    // Get the current point location.
    segList = ee.List(segList);
    var currentLoc = ee.Feature(segList.get(-1)).geometry();
    var locT1 = ee.Geometry.Point(ee.List(currentLoc.coordinates()).get(1));
    
    // Get wind attributes for current location.
    var samp = segmentWind.sample({
      region: locT1,
      scale: scale,
      numPixels: 1
    });
    
    // Get x and y delta.
    var pt = samp.first();
    var east = metersToDegrees(pt.getNumber('east_wind'), SEG_DURATION);
    var north = metersToDegrees(pt.getNumber('north_wind'), SEG_DURATION);
    
    // Determine the new location.
    var locT2x = ee.Number(ee.List(locT1.coordinates()).get(0)).add(east);
    var locT2y = ee.Number(ee.List(locT1.coordinates()).get(1)).add(north);
    
    // Build a line segment.
    var locT2 = ee.Geometry.Point([locT2x, locT2y]);
    var line = ee.Geometry.LineString([locT1, locT2]);
    var segFeature = ee.Feature(line, {segment: segmentLabel});
    
    // Add the segment to this point's segment list.
    return segList.add(segFeature);
  }
  
  // Iterate over the segment sequence to generate segments for point path.
  var segments = ee.List.sequence(1, N_SEGMENTS);
  var path = segmentInfo.iterate(moveIt, segmentList);
  return path;
});

// Convert segment list to FeatureCollection; filter out initial segment.
var pathsFc = ee.FeatureCollection(pathsList.flatten())
  .filter(ee.Filter.gt('segment', 0));



// #############################################################################
// ### DISPLAY TRAJECTORIES PATHS TO THE MAP ###
// #############################################################################

// Define a function to paint vector to raster and visualize.
function paintFc(fc, paintParams, visParams) {
  paintParams.featureCollection = fc;
  return ee.Image().byte()
    .paint(paintParams)
    .visualize(visParams);
}

// Define visualization parameters and display wind paths.
var visParamsMap = {
  palette: ['3742FA', 'FFFFFF'],
  min: 1,
  max: N_SEGMENTS
};

var paintParamsMap = {
  featureCollection: null,
  color: 'segment',
  width: 1.5
};

// Define hillshade visualization layer.
var dem = ee.Image('MERIT/DEM/v1_0_3');
var hillshade = ee.Terrain.hillshade(dem.multiply(15))
  .visualize({min: 0, max: 800, opacity: 0.5});
  
// Add administrative boundaries to the hillshade layer.
var adminUnitsFc = ee.FeatureCollection('FAO/GAUL_SIMPLIFIED_500m/2015/level1');
var adminUnitsImg = ee.Image().byte()
  .paint(adminUnitsFc, 1, 2)
  .visualize({palette: ['D3D3D3']});
hillshade = adminUnitsImg.blend(hillshade);

// Paint paths to an image; color segments according to sequence label.
var pathsImg = paintFc(pathsFc, paintParamsMap, visParamsMap);

Map.centerObject(AOI, 6);
Map.addLayer(hillshade);
Map.addLayer(pathsImg, null, 'Wind Paths', true, 0.7);

// Define a dark base map.
Map.setOptions('Dark Map', {'Dark Map': darkMap()});
function darkMap() {
  return [{"elementType":"geometry","stylers":[{"color":"#212121"}]},{"elementType":"labels.icon","stylers":[{"visibility":"off"}]},{"elementType":"labels.text.fill","stylers":[{"color":"#757575"}]},{"elementType":"labels.text.stroke","stylers":[{"color":"#212121"}]},{"featureType":"administrative","elementType":"geometry","stylers":[{"color":"#757575"}]},{"featureType":"administrative.country","elementType":"labels.text.fill","stylers":[{"color":"#9e9e9e"}]},{"featureType":"administrative.land_parcel","stylers":[{"visibility":"off"}]},{"featureType":"administrative.locality","elementType":"labels.text.fill","stylers":[{"color":"#bdbdbd"}]},{"featureType":"poi","elementType":"labels.text.fill","stylers":[{"color":"#757575"}]},{"featureType":"poi.park","elementType":"geometry","stylers":[{"color":"#181818"}]},{"featureType":"poi.park","elementType":"labels.text.fill","stylers":[{"color":"#616161"}]},{"featureType":"poi.park","elementType":"labels.text.stroke","stylers":[{"color":"#1b1b1b"}]},{"featureType":"road","elementType":"geometry.fill","stylers":[{"color":"#2c2c2c"}]},{"featureType":"road","elementType":"labels.text.fill","stylers":[{"color":"#8a8a8a"}]},{"featureType":"road.arterial","elementType":"geometry","stylers":[{"color":"#373737"}]},{"featureType":"road.highway","elementType":"geometry","stylers":[{"color":"#3c3c3c"}]},{"featureType":"road.highway.controlled_access","elementType":"geometry","stylers":[{"color":"#4e4e4e"}]},{"featureType":"road.local","elementType":"labels.text.fill","stylers":[{"color":"#616161"}]},{"featureType":"transit","elementType":"labels.text.fill","stylers":[{"color":"#757575"}]},{"featureType":"water","elementType":"geometry","stylers":[{"color":"#000000"}]},{"featureType":"water","elementType":"labels.text.fill","stylers":[{"color":"#3d3d3d"}]}];
}



// #############################################################################
// ### ANIMATE PARTICLE PATHS ###
// #############################################################################

// Define parameters for painting features to images.
var paintParamsGif = {
  featureCollection: null,
  color: 'segment',
  width: 1.5
};

// Define parameters for visualizing paths.
var visParamsGif = {
  palette: ['#DFFF00'],
  min: 1,
  max: N_SEGMENTS
};

// Make an ImageCollection of segments with a fade filter overlay.
var fadeCol = ee.ImageCollection.fromImages(segmentSeq.map(function(segment) {
  // Visualize this segment.
  var thisSegFc = pathsFc.filter(ee.Filter.eq('segment', segment));
  var thisSegImg = paintFc(thisSegFc, paintParamsGif, visParamsGif);
  
  // Visualize segments <= to this segment; used to define fade filter mask.
  var theseSegsFc = pathsFc.filter(ee.Filter.lte('segment', segment));
  var theseSegsImg = paintFc(theseSegsFc, paintParamsGif, visParamsGif);

  // Make a fade filter; use background, which is hillshade in this case.
  var fadeFilter = hillshade.updateMask(theseSegsImg.select(0).mask())
    .visualize({opacity: 0.1});
  
  // Overlay the fade filter on the current segment.
  return thisSegImg.blend(fadeFilter).set('segment', segment);
}));

// Make an image collection of progressively composited faded segments.
var vidCol = ee.ImageCollection.fromImages(segmentSeq.map(function(segment) {
  // Mosaic faded segments up to and including the current segment.
  var segmentsMosaic = fadeCol.filter(ee.Filter.lte('segment', segment))
    .sort('year')
    .mosaic();
  
  // Paint the current segment to an image and visualize.
  var thisSegFc = pathsFc.filter(ee.Filter.eq('segment', segment));
  var thisSegImg = paintFc(thisSegFc, paintParamsGif, visParamsGif);
  
  // Blend all the components.
  return hillshade.blend(segmentsMosaic).blend(thisSegImg).set('segment', segment);
}));

// Define animation parameters and display the GIF to the console.
var videoArgs = {
  dimensions: MAX_DIM, // Max dimension (pixels), min dimension is proportionally scaled.
  region: ANIMATION_REGION,
  framesPerSecond: 16,
  crs: CRS
};
print(ui.Thumbnail(vidCol, videoArgs));
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment