Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Reducing image collection to image returns only one value per band #214

Closed
3 tasks done
DavidDHofmann opened this issue Jan 12, 2022 · 4 comments
Closed
3 tasks done

Comments

@DavidDHofmann
Copy link

  • rgee version: 1.1.2
  • R version: 4.1.2
  • Operating System: x86_64-pc-linux-gnu (64-bit)

At submit an issue, please attached the following information of your rgee session:

  • You have the Python API installed (from terminal):
earthengine -h
  • You can find the credentials file on your system:
library(rgee)
ee_path <- path.expand("~/.config/earthengine/credentials")
file.exists(ee_path)
  • You can run a simple EE command from R:
library(rgee)

# Initialize the Earth Engine module.
ee_Initialize()

# Print metadata for a DEM dataset.
print(ee$Image('USGS/SRTMGL1_003')$getInfo())

Description

Hi there!

I'm trying to download Sentinel-2 data for a specific date. For the date and area of interest, Earth Engine returns an image-collection with two images that overlap in extent. I then try to reduce the collection into a single image using a median-reducer. Ultimately, I clip the reduced image to my area of interest and download the data. Oddly enough, ee_print(reduced) shows that the image only contains one pixel per band. Nevertheless, the image looks fine when displayed on the map. Only when downloaded, the bands contain only one value per band. To me, it seems that the reducer does not apply the median per pixel across bands, but across the entire band. The coordinates of the final image also do not exactly match my specified aoi and I get a warning after download saying:

Warning
In .rasterFromGDAL(x, band = band, objecttype, ...) :
  data seems flipped. Consider using: flip(x, direction='y')

I'm not quite sure what I'm doing wrong. Any suggestions?

What I Did

# Load required packages
library(rgee)
library(terra)

# Init
ee_Initialize()

# Specify an area of interest
aoi <- ee$Geometry$Polygon(
  list(
      c(23.4976149678726, -19.649075296214)
    , c(23.4976149678726, -19.2680589256603)
    , c(23.93394016641, -19.2680589256603)
    , c(23.93394016641, -19.649075296214)
  )
)

# Query sentinel 2 data
collection <- ee$ImageCollection("COPERNICUS/S2_SR")$
  filterDate("2021-04-05", "2021-04-06")$
  filterBounds(aoi)$
  select("B2")

# Check it
ee_print(collection)
────────────────────────────────────────────────────────────── Earth Engine ImageCollection ──
ImageCollection Metadata:
 - Class                      : ee$ImageCollection
 - Number of Images           : 2
 - Number of Properties       : 23
 - Number of Pixels*          : 241142760
 - Approximate size*          : 3.36 GB
Image Metadata (img_index = 0):
 - ID                         : COPERNICUS/S2_SR/20210405T081559_20210405T083938_T34KGD
 - system:time_start          : 2021-04-05 08:45:41
 - system:time_end            : 2021-04-05 08:45:41
 - Number of Bands            : 1
 - Bands names                : B2
 - Number of Properties       : 81
 - Number of Pixels*          : 120571380
 - Approximate size*          : 1.68 GB
Band Metadata (img_band = 'B2'):
 - EPSG (SRID)                : WGS 84 / UTM zone 34S (EPSG:32734)
 - proj4string                : +proj=utm +zone=34 +south +datum=WGS84 +units=m +no_defs
 - Geotransform               : 10 0 699960 0 -10 7900000
 - Nominal scale (meters)     : 10
 - Dimensions                 : 10981 10980
 - Number of Pixels           : 120571380
 - Data type                  : INT
 - Approximate size           : 1.68 GB
 ──────────────────────────────────────────────────────────────────────────────────────────────
 NOTE: (*) Properties calculated considering a constant  geotransform and data type.
# There are two images on this day, let's take the median of the matching bands
reduced <- collection$
  mean()$
  clip(aoi)

# Check again
ee_print(reduced)
# Displayed maps look correct
Map$centerObject(reduced, zoom = 10)
Map$addLayer(reduced$select("B2"), list(max = 1000, min = 0), name = "B2")

image

# Downloaded data shows only one value in the band
ee_as_raster(reduced, dsn = "reduced.tif")
r <- rast("reduced.tif")
show(r)
dimensions  : 1, 1, 1  (nrow, ncol, nlyr)
resolution  : 1, 1  (x, y)
extent      : 23, 24, -20, -19  (xmin, xmax, ymin, ymax)
coord. ref. : lon/lat WGS 84 (EPSG:4326) 
source      : reduced.tif 
name        : B2 
plot(r)

image

@csaybar
Copy link
Collaborator

csaybar commented Jan 14, 2022

Hi @DavidDHofmann you have to define the scale argument:

reduced2 <- ee_as_raster(reduced, dsn = "reduced.tif", scale = 1000)

The warnings appear because GEE retrieved your query using positive y spacing and raster does not automatically support it. I recommend you use ee_as_stars to avoid this problem. evetion/GeoArrays.jl#53

Warning
In .rasterFromGDAL(x, band = band, objecttype, ...) :
  data seems flipped. Consider using: flip(x, direction='y')

@DavidDHofmann
Copy link
Author

Hi @csaybar
Thank you for the swift reply! The code works when providing the scale argument. Does this mean that it is not possible to download the data at its native resolution?

@csaybar
Copy link
Collaborator

csaybar commented Jan 15, 2022

Hi @DavidDHofmann :)

In fact, you can!. The issue is that after applying a reducer collection$mean()$clip(aoi) you will lose the image geotransform properties.

ee_print(collection)

image

ee_print(collection$mean())

image

It is a consequence of GEE lazy evaluation. The client library API does not know if after applying your "selected reducer" (mean) you are conserving the image geometric properties it's something that will be evaluated at the server-side. My advice would be to store the properties before the reducer.

# Load required packages
library(rgee)
library(terra)

# Init
ee_Initialize()

# Specify an area of interest
aoi <- ee$Geometry$Polygon(
  list(
    c(23.4976149678726, -19.649075296214)
    , c(23.4976149678726, -19.2680589256603)
    , c(23.93394016641, -19.2680589256603)
    , c(23.93394016641, -19.649075296214)
  )
)

# Query sentinel 2 data
collection <- ee$ImageCollection("COPERNICUS/S2_SR")$
  filterDate("2021-04-05", "2021-04-06")$
  filterBounds(aoi)$
  select("B2")

collection$first()$projection()$nominalScale()$getInfo()
#10

@DavidDHofmann
Copy link
Author

Ah okey cool, I wasn't aware of this. This also explains why the estimated image size becomes unreasonably small after applying the reudcer. Thank you very much for the clarification. The download works perfectly fine now :)

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants