About Earth Timelapse

Page Summary

outlined_flag

• Earth Timelapse provides a zoomable video view of how our planet has changed since 1984, showcasing different types of planetary change.

• The Timelapse videos page features preview thumbnails with buttons to view or download videos in various 2D and 3D formats, including 4K resolution and options with or without labels.

• Attribution is required for using Earth Timelapse, citing "Google Earth Timelapse (Google, Landsat, Copernicus)".

• Earth Timelapse is licensed under a Creative Commons Attribution 4.0 International License.

The Earth Timelapse dataset provides a four-decade visual record of planetary change spanning from 1984 to 2022. Synthesized by Google from petabytes of spaceborne observations acquired by NASA/USGS Landsat missions (Landsat 4, 5, 7, 8, and 9) and ESA's Copernicus Sentinel-2 constellation, the dataset is offered as annual mosaics for visual interpretation.

The annual mosaics are globally complete, gap-filled, and water-masked visual basemaps that power the public Google Earth Timelapse interactive viewer. Pre-1999 observation voids are smoothly reconstructed using pixel-wise temporal linear regression interpolation across adjacent valid observation years. Open ocean waters are replaced with styled NOAA ETOPO1 shaded-relief bathymetry and global water masking (using Hansen GFC and MOD44W datasets). These data are suitable for: visual interpretation, educational storytelling, base mapping, or custom video exports.

Spatial Resolution Eras

The Earth Timelapse collection spans two distinct spatial resolution eras:

• 1984–2014 (30 m / pixel): Synthesized from USGS/NASA Landsat 4, 5, 7, and 8 observations. Images are gridded in Web Mercator (EPSG:3857) at a 30-meter resolution.
• 2015–Present (19.11 m / pixel): Synthesized by fusing ESA Copernicus Sentinel-2A/2B MultiSpectral Instrument (MSI) observations (10 m / 20 m native) with Landsat 8 and 9. Images are gridded in Web Mercator (EPSG:3857) at ~19.11-meter resolution.

The annual mosaics cover all global landmasses and coastal zones from approximately -82.6°S to 83.69°N latitude. Synthesizing over four decades of Earth observation data requires harmonizing disparate satellite constellations across operational lifespans, sensor modalities, orbital geometries, and atmospheric perturbations.

The Earth Timelapse pipeline ingests millions of individual scenes from the USGS/NASA Landsat program (Landsat 4, 5, 7, 8, and 9) and the European Space Agency (ESA) Copernicus Sentinel-2 constellation. To transform raw and Top-of-Atmosphere (TOA) observations into consistent annual basemaps, each scene is processed with:

• Phenological and seasonal scene filtering to optimize peak vegetation growth and minimize transient seasonal effects.
• Missing channel restoration and sensor anomaly mitigation (e.g., Landsat 8 TIRS saturation, Landsat 7 SLC-off artifacts).
• Multi-sensor cloud, haze, and shadow quality scoring combining Landsat simpleCloudScore, heuristic HSV space cloud recovery, and Sentinel-2 Cloud Score+.
• Panchromatic normalization and structural pseudo-pansharpening to increase the apparent resolution.
• MODIS BRDF-calibrated radiometric normalization using multi-scale low-pass spatial filtering to eliminate scene boundaries while preserving local surface textures.
• Water-conscious edge sharpening to avoid coastal shoreline bleaching.
• Temporal 16-bit median compositing to reject transient atmospheric noise, shadows, and ephemeral artifacts.
• Temporal linear regression interpolation (for Global Final assets) to reconstruct missing observations in pre-1999 mosaics.
• Global ocean bathymetry modeling and water masking combining NOAA ETOPO1, Hansen Global Forest Change datamask, and MODIS water masks.
• Global color balancing, polar ice white-point preservation, multi-scale local contrast enhancement (LCE) for visual enhancement.

Source Data & Sensor Specifications

The Earth Timelapse collection leverages the following spaceborne Earth observation missions:

Satellite and Sensor Specifications

Table 1. Core optical sensor systems integrated into the Earth Timelapse annual mosaic production pipeline.

Auxiliary and Reference Datasets

1. MODIS Global BRDF / Surface Reflectance Baselines:
   • Pre-computed multi-year surface reflectance composites derived from Terra and Aqua MODIS (MODIS/006/MCD43A4 and MOD09GA/MYD09GA) at 500-meter resolution serve as the planetary radiometric reference target for broad-scale illumination and haze harmonization.
2. NOAA ETOPO1 Global Relief:
   • 1 Arc-Minute global topography and bathymetry utilized to generate realistic shaded-relief ocean floor terrain and depth contouring in the global water-masked product.
3. Hansen Global Forest Change (UMD/hansen/global_forest_change_2015):
   • The datamask layer (distinguishing land, permanent water bodies, and coastal zones) serves as the primary land and water separation boundary.
4. MODIS Global Land/Water Mask (MODIS/MOD44W/MOD44W_005_2000_02_24):
   • Used in conjunction with Hansen GFC to isolate terrestrial surfaces and prevent coastal water-mask bleeding.
5. Sentinel-2 Cloud Score+ (GOOGLE/CLOUD_SCORE_PLUS/V1/S2_HARMONIZED):
   • Pixel-level quality assessment providing cs (clear sky confidence) and cs_cdf (cumulative distribution clear sky probability) metrics for Sentinel-2 MSI data.
6. Global Cloud Climatology Statistics:
   • Empirically derived 25th-percentile and minimum cloud score baselines that adapt cloud-masking thresholds in persistently cloudy tropical zones.

Spatial and Temporal Filtering

The global annual mosaics balance the need to acquire enough valid observations in cloud-prone tropical regions while preventing winter snow contamination and extreme solar zenith angles in high-latitude regions.

Latitudinal and Phenological Windowing

• Northern Hemisphere Polar / High Latitude Zone (>60°N to 83.69°N):
  • Direct filtering by Day of Year (DOY): Observations are restricted to DOY 150 through 270 (late May through September). This targets the peak vegetative growing season (maximum greenness) and minimizes seasonal snow, ice cover, and long terrain shadows caused by low solar elevation angles (solar elevation ≤ 0°).
• Temperate and Southern Hemisphere Zone (57°S to 60°N):
  • Utilizes the full calendar year (January 1 to December 31) to maximize scene availability.

Processing Methodology: Masked Annual Mosaics

The masked annual mosaic pipeline generates un-interpolated global composites.

Missing Channel Recovery & Sensor Preparation

Before atmospheric correction and quality scoring, raw landsat data are standardized:

• Thermal Band Reconstruction: During operational periods where Landsat 8 TIRS thermal bands were degraded or uncalibrated, placeholder zero-variance thermal channels are synthesized to satisfy the internal interface requirements of automated Landsat cloud scoring algorithms.
• Synthetic Panchromatic Generation for Landsat 4/5: Landsat 4 and 5 TM sensors lack an optical panchromatic channel (Band 8). A synthetic pseudo-panchromatic band is constructed per scene:

Here \(I_{\text{raw}}(\cdot)\) is an image function which provides raw pixels and might take one or more bands as parameters. For example, \(I_{\text{raw}}(Pan)\) returns the raw panchromatic band and \(I_{\text{raw}}(RGB)\) returns the raw RGB bands. mean() is the mean function.

Calibration to Top-of-Atmosphere Reflectance (TOA)

Landsat raw scenes are processed to TOA using ee.Algorithms.Landsat.TOA.

ee.Algorithms.Landsat.TOA is an image function which converts Landsat raw data to Landsat TOA data (\(I_{\text{TOA}}\)).

Multi-Sensor Quality and Cloud Masking

Orbital Swath Boundary Tapering

A distance-decay soft alpha mask helps minimize hard seams along Landsat orbital tracks:

where \(\text{Mask}_{\text{MODIS}}\) is the mask derived from MODIS surface reflectance annual composites (MOD09GA/MYD09GA), \(\text{Mask}_{\text{raw}}\) is the mask of a Landsat raw scene, min() is the minimization function, and \(\text{Gaussian}\) is a Gaussian spatial convolution with specified standard deviation.

Landsat 7 SLC-off Water Masking

Following the failure of the Scan Line Corrector (SLC) on Landsat 7 in May 2003, Landsat 7 scenes contain linear data gaps, which can lead to severe artifacts over open water. Landsat 7 pixels over open water are masked out using a blurred MODIS water mask, prioritizing Landsat 4/5/8/9 or Sentinel-2 depending on the year of mosaic production.

Landsat Cloud Scoring & Anomaly Recovery

Landsat scenes are given a per-pixel \(\text{Cloud}\) score (in [0, 100]) using ee.Algorithms.Landsat.simpleCloudScore. Two additional refinement stages are applied:

1. Undetected Cloud Recovery in HSV Space: High-reflectance, low-saturation cloud anomalies that receive an erroneous score of Cloud = 0 are detected by transforming \(I_{\text{TOA}}(RGB)\) into Hue-Saturation-Value (HSV) space and updating the \(\text{Cloud}\) mask:

1. Adaptive Climatological Thresholding: The threshold for discarding cloudy pixels (i.e. \(\text{Cloud} \gt \text{Threshold}_{\text{cloud}}\)) is dynamically calculated from pre-computed spatial cloud statistics:

where \(\text{Cloud}_{\text{p25}}\) is the per pixel 25th percentile \(\text{Cloud}\) score clamped to be between 25 and 65. The threshold is raised over verified bright land surfaces (desert sand dunes, salt flats):

Sentinel-2 Cloud Scoring

Sentinel-2 Harmonized scenes are paired with GOOGLE/CLOUD_SCORE_PLUS/V1/S2_HARMONIZED:

• Pixel Filtering: Retained when clear-sky confidence cs ≥ 0.60 and cumulative clear-sky probability cs_cdf ≥ 0.70.

Panchromatic Normalization and Structural Fusion

To leverage the 15-meter panchromatic resolution of Landsat 7/8/9 without spectral distortion across sensor transitions, panchromatic normalization is applied to the TOA data:

where \(\mu_{\text{3x3}}\) is the 3x3 pixel neighborhood mean function.

Radiometric Normalization against MODIS BRDF Baselines

Non-linear Contrast Scaling: Top-of-Atmosphere reflectance is scaled to 8-bit display range ([0,255]) with gamma adjustment:

where \(\text{Visualize}\) is the ee.Image.visualize() function. Hereafter, 8-bit RGB color space is implied unless noted otherwise.

Overcorrection Rejection Mask: Protects real rapid land-cover change from over-smoothing:

where \(I_{\text{BRDF}}\) produces 8-bit RGB data from MODIS BRDF-Adjusted Reflectance (MCD43A4) composites.

Low-Pass Spatial Field Extraction:

Haze Subtraction and Reflectance Harmonization:

Terrestrial Land Spatial Sharpening

Applies Laplacian edge convolution to a combination of \(I_{\text{scaled}}\) and \(I_{\text{normalized}}\) exclusively on terrestrial land surfaces.

Temporal Median Compositing

For each calendar year \(Y\), the normalized, cloud-masked image collection is reduced into a single multi-band composite using a 16-bit precision median reducer:

The resulting image represents a raw annual mosaic, featuring red, green, blue, and cloud bands with un-interpolated valid pixel masks. Note that time (\(t\)) is a parameter to the image function.

Post-Processing: Global Gap-Filled & Water-Masked Assets

To produce the seamless, globally complete basemap utilized in the interactive Earth Timelapse viewer and the annual mosaics, the raw annual composites are post-processed.

Temporal Linear Regression Gap-Filling (Pre-1999 Mosaics)

Prior to the launch of Landsat 7 in 1999, global satellite acquisition was constrained by historical downlink infrastructure, onboard tape recorder limitations, and persistent cloud cover. Consequently, pre-1999 annual mosaics contain significant spatial voids, particularly across Central Africa, Southeast Asia, Siberia, and the Amazon.

To eliminate distracting gray voids while preserving temporal transitions, the pipeline implements a pixel-wise temporal linear regression interpolation:

Temporal Anchor Search

For each target year Y, the algorithm searches the multi-decadal collection to find the most recent, valid, non-masked pixel observation before (\(t_{\text{before}} \le Y\)) and after (\(t_{\text{after}} \ge Y\)).

Linear Regression Fit

An ordinary least squares regression model is evaluated per pixel across the bounding temporal observations with independent variables [1, t] and dependent variables [Red, Green, Blue]:

Simulated Value Evaluation

The linear trajectory is evaluated at target year Y:

Observation Layering

The annual pixels from year Y are mosaicked with the interpolated pixels:

This ensures that wherever real satellite observations exist for that year, they are preserved , while historical data gaps are interpolated across time.

High-Latitude Glacial Infill

Landsat orbital inclination limits observations near the extreme polar caps (>82.6°N). For interior regions of Northern Greenland and Arctic ice shelves lacking optical coverage, a calibrated, normalized multi-year baseline is blended to maintain clean, seamless polar basemaps.

Global Ocean Bathymetry & Water Masking

In raw annual mosaics, ocean surfaces contain sun-glint, cloud shadows, and transient wave artifacts. The final asset replaces open ocean waters with a global, shaded-relief bathymetric basemap:

1. Bathymetric Shading: Topographic elevation from the NOAA ETOPO1 model is styled using an ocean depth color ramp:
   • Deep Abyssal Plains (-5000 m): Deep Navy (#000927)
   • Continental Slope (-1000 m): Slate Blue (#000E3A)
   • Continental Shelf (-100 m): Azure Navy (#000E3B)
   • Coastline (0 m): Royal Cobalt (#001146)
2. Hillshading & High-Pass Gaussian Enhancement: Combined with analytical hillshading (ee.Terrain.hillshade) and 5000 m Gaussian unsharp masking to emphasize oceanic trenches, mid-ocean ridges, and seamounts.
3. Land/Water Boundary Harmonization:
   • Combines the Hansen Global Forest Change datamask (distinguishing land from ocean) with the MODIS water mask (MOD44W).
   • Paints inland lakes and enclosed seas (e.g., Caspian Sea, Great Lakes, Lake Baikal, Aral Sea) to preserve natural water color dynamics.
   • Applies a 3-stage morphological circular kernel reduction and 3000 m Gaussian blur to smoothly taper the boundary between shallow coastal waters and offshore bathymetry without hard clipping.

Global Color Balancing & Local Contrast Enhancement (LCE)

1. HSV Value Boost & Gamma Balancing: Channel gammas ($\gamma_R = 0.98, \gamma_G = 1.00, \gamma_B = 1.04$), followed by a 10% value channel boost in HSV color space.
2. Polar Ice White-Point Correction: Flags high-reflectance polar surfaces (Greenland, Antarctica, alpine ice sheets) and balances them to pure neutral white ([255, 255, 255]).
3. Multi-Scale LCE: Enhance contrast and terrain sharpness for different feature types using a combination of Laplacian convolution, mean filtering, and Gaussian unsharp masking.

Table 2. Detailed technical comparison between Raw and Global Final Earth Timelapse collections.

Limitations and Analytical Considerations

1. Dual Resolution Eras (30 m vs. 19.11 m): Users performing multi-decadal time-series analysis must account for the resolution transition in 2015. Composite images from 1984–2014 are gridded at 30.0 meters per pixel, whereas images from 2015 onwards are gridded at 19.11 meters per pixel due to the integration of Sentinel-2 MSI data.
2. Interpolated Pixels (Pre-1999): In the Global Final collection, missing pixels in pre-1999 mosaics are temporally interpolated. In regions undergoing sudden land-use transitions (e.g., rapid deforestation or reservoir construction) during a multi-year gap, the interpolated pixels will depict a gradual linear transition rather than an abrupt discrete event.
3. Altered Band Ratios: Radiometric normalization optimizes visual consistency across scenes against MODIS BRDF targets. While relative spatial patterns are preserved, derived spectral indexes (e.g., NDVI, EVI) differ from Level-2 Surface Reflectance values.
4. Phenological Mixing: High-latitude northern regions (>60°N) represent mid-summer conditions (DOY 150–270), whereas temperate and tropical zones represent annual composite medians.

Attribution

Google Earth Timelapse (Google, Landsat, Copernicus)