3 Value at a point

Author

Andres Patrignani and Joaquin Peraza

Published

January 26, 2024

Perhaps, one of the simplest and most basic operations in Google Earth Engine is to retrieve information for a single point. Typical examples include elevation, mean annual temperature, and soil properties. In other words, this type of operations is mostly intended at retrieving information that does not change with time.

# Import modules
import ee
from datetime import datetime
from pprint import pprint

# Authenticate
ee.Authenticate()

# Initialize the library.
ee.Initialize()

Tip

While we need to run ee.Initialize() everytime we start a notebook or restart the Python kernel, we only need to authenticate once, so a handy tip is to mute the ee.Initialize() line with a # once you are done. Note that after several hours or days of inactivity you will need to authenticate again.

Example 1: Elevation

In the following example we will retrieve the elevation for a specific point on Earth.

Product: USGS/SRTMGL1_003

# Define geographic coordinates
lat = 39.186512 # This is y
lon = -96.576844 # This is x

# Convert coordinates into a Point geometry following the x,y notation
point = ee.Geometry.Point([lon, lat])

# Explore point geometry
point.getInfo()

{'type': 'Point', 'coordinates': [-96.576844, 39.186512]}

Note

In this tutorial we will use the .getInfo() method to print the underlying data and metadata of multiple objects, so that you can see the output and become more familiar with GEE data structures and the geoJSON format. This information is also valubale for debugging code errors. However, you don’t typically want to clutter your notebooks with these long outputs, so we will skip some of these steps in future tutorials.

# Load digital elevation model (DEM) from Shuttle Radar Topography Mission (SRTM) 
dem = ee.Image('USGS/SRTMGL1_003')

# Obtain some information about the DEM layer (output is long!)
pprint(dem.getInfo())

# Retrieve the elevation value: This step will get us close to the answer, but we are not there yet
pprint(dem.sample(point, 1).getInfo())

{'columns': {'elevation': 'Short'},
 'features': [{'geometry': None,
               'id': '0',
               'properties': {'elevation': 317},
               'type': 'Feature'}],
 'properties': {'band_order': ['elevation']},
 'type': 'FeatureCollection'}

# Retrieve the elevation value: This step will get us even closer to the answer, but we are not there yet
dem.sample(point, 1).first().getInfo()

{'type': 'Feature',
 'geometry': None,
 'id': '0',
 'properties': {'elevation': 317}}

# Retrieve the elevation value: This step will get us the correct value
elev = dem.sample(point, 1).first().getNumber('elevation').getInfo()
print(f'Elevation: {elev} m')

Elevation: 317 m

Example 2: Weather variables

Obtain the long-term mean annual air temperature and precipitation for a specific location

Product: For more information about bands and units visit: WorldClim BIO

# Define geographic coordinates
lat = 39.186512 # This is y
lon = -96.576844 # This is x

# Convert coordinates into a Point geometry following the x,y notation
point = ee.Geometry.Point(lon, lat)

# Load WorldClim BIO dataset
dataset = ee.Image('WORLDCLIM/V1/BIO')

# Access metadata of the product (output is long!)
dataset.getInfo()

# Get long-term mean annual air temperature
T_mean = dataset.select('bio01').sample(point,1).first().getNumber('bio01').multiply(0.1).getInfo()
print(f'Mean annual air temperature: {round(T_mean,1)} Celsius')

Mean annual air temperature: 12.2 Celsius

# Get long-term annual precipitation
P_annual = dataset.select('bio12').sample(point,1).first().get('bio12').getInfo()
print(f'Mean annual precipitation: {round(P_annual,1)} mm')

Mean annual precipitation: 857 mm

Example 3: Reference ET

In this example we will retrieve daily values of reference evapotranspiration for a point.

Product: For more information visit the description in Google Earth Engine of GRIDMET

# Define point
point = ee.Geometry.Point([-96.576844, 39.186512])

# Obtain GRIDMET dataset for a specific period. End date is excluded from the call
start_date = ee.Date('2022-07-01')
end_date = ee.Date('2022-07-02')
day_of_interest = ee.Date('2022-03-15')
dataset = ee.ImageCollection('IDAHO_EPSCOR/GRIDMET').filterDate(start_date, end_date)

# Information about image collection
pprint(dataset.getInfo())

# Use the get methog to retrieve a specific property
dataset.get('description').getInfo()

# Information about first image within the collection (output is long!)
dataset.first().getInfo()

# Information about feature collection
dataset.first().sample(point,1).getInfo() # Sample at 1 meter resolution

{'type': 'FeatureCollection',
 'columns': {'bi': 'Float',
  'erc': 'Float',
  'eto': 'Float',
  'etr': 'Float',
  'fm100': 'Float',
  'fm1000': 'Float',
  'pr': 'Float',
  'rmax': 'Float',
  'rmin': 'Float',
  'sph': 'Float',
  'srad': 'Float',
  'th': 'Float',
  'tmmn': 'Float',
  'tmmx': 'Float',
  'vpd': 'Float',
  'vs': 'Float'},
 'properties': {'band_order': ['pr',
   'rmax',
   'rmin',
   'sph',
   'srad',
   'th',
   'tmmn',
   'tmmx',
   'vs',
   'erc',
   'eto',
   'bi',
   'fm100',
   'fm1000',
   'etr',
   'vpd']},
 'features': [{'type': 'Feature',
   'geometry': None,
   'id': '0',
   'properties': {'bi': 0,
    'erc': 20,
    'eto': 4.5,
    'etr': 5.300000190734863,
    'fm100': 14.899999618530273,
    'fm1000': 17.799999237060547,
    'pr': 43.70000076293945,
    'rmax': 90.9000015258789,
    'rmin': 60.20000076293945,
    'sph': 0.014340000227093697,
    'srad': 249.39999389648438,
    'th': 114,
    'tmmn': 293.70001220703125,
    'tmmx': 300.1000061035156,
    'vpd': 0.7300000190734863,
    'vs': 2.799999952316284}}]}

# Information about feature
eto = dataset.first().sample(point,1).first().getNumber('eto').getInfo()
print(f'The grass reference ET is {eto} mm')

The grass reference ET is 4.400000095367432 mm

Alternative solution: Access the image directly rather than the collection.

ee.Image('IDAHO_EPSCOR/GRIDMET/20220701').sample(point,1).first().getNumber('eto').getInfo()

Get dataset timestamps

# Obtain the START time of the dataset to check that it matches our request.
# Response is in milliseconds since 1-Jan-1970
dataset.select('eto').first().getNumber('system:time_start').getInfo()

1656655200000

# Obtain the END time of the dataset to check that it matches our request.
# Response is in milliseconds since 1-Jan-1970
dataset.select('eto').first().getNumber('system:time_end').getInfo()

1656741600000

# Find the datetime of the serial date numbers
# Input in "fromtimestamp()" has to be in seconds, so we divide by 1000
print('Start time:', datetime.fromtimestamp(1656655200000/1000).strftime('%Y-%m-%d %H:%M:%S.%f'))
print('End time:', datetime.fromtimestamp(1656741600000/1000).strftime('%Y-%m-%d %H:%M:%S.%f'))

Start time: 2022-07-01 01:00:00.000000
End time: 2022-07-02 01:00:00.000000

Example 4: Retrieve soil properties for a given location

Product: SoilGrids

Source: https://samapriya.github.io/awesome-gee-community-datasets/projects/isric/

# Load the SoilGrids dataset from GEE
soil_grids = ee.Image("projects/soilgrids-isric/sand_mean")

# Define a point geometry
lat = 37.839154
lon = -99.101594
point = ee.Geometry.Point(lon,lat)
sand = soil_grids.sample(point,250).first().getNumber('sand_0-5cm_mean').multiply(0.1).getInfo()

print(f'The percentage of sand at ({lat},{lon}) is: {round(sand)} %')

The percentage of sand at (37.839154,-99.101594) is: 53 %