search

google / Xee

An Xarray extension for Google Earth Engine
Apache License 2.0
240 stars 28 forks source link

Error when passing in ee.Projection() as parameter for open_dataset() #158

Open mhutchinson-woolpert opened 4 months ago

mhutchinson-woolpert commented 4 months ago

Good afternoon,

I've tried to reproduce one of the examples from the website, but rather using a different dataset (Sentinel-2). The goal was to open the dataset using the same projection as the satellite image. So far these have been tried:

# ic = ee.ImageCollection(ee.Image('COPERNICUS/S2_SR_HARMONIZED/20230728T221531_20230728T221534_T05WPU'))
ic = ee.ImageCollection('COPERNICUS/S2_SR_HARMONIZED').filter(ee.Filter.eq('PRODUCT_ID', 'S2A_MSIL2A_20230728T221531_N0509_R115_T05WPU_20230729T014233') )

rect = ee.Geometry.Rectangle(-148.46380797863654, 70.30003598701754, -148.39544465908497, 70.32308553183186)
# rect_projected = rect.transform(ic.first().select(0).projection(), 1)

xr_ds = xarray.open_dataset(ic, 
                            engine='ee',
                            projection=ic.first().select('B4').projection(),
                            geometry=rect)

Version 0.0.12 of Xee is being used. Note that other data has been loaded successfully using Xee (different parameters, different datasets etc.), so it's not an issue related to installation or dependencies.

The error can be seen below:

image

naschmitz commented 4 months ago

I took a quick look. Xee doesn't do a great job exporting data in this image's native projection. Earth Engine batch export tools are the way to go for something like that.

I think you've hit a bug. We calculate the specified region's bounds in EPSG:32605 projection and it ends up being a complex multi-polygon (which Xee isn't equipped to deal with). You could consider using a different projection for the output (EPSG:4326 for example).

warframeAlpha commented 1 month ago

I got the same issue, so I tried to void passing a geometry.

jasper_bbox = ee.Geometry.Rectangle(-118.19, 52.00, -117.92, 53.00)
lc9 = ee.ImageCollection('LANDSAT/LC09/C02/T1_L2') \
.filterBounds(jasper_bbox)\
.filterMetadata('CLOUD_COVER', 'less_than', 1)
crs = lc9.first().select('SR_B2').projection().getInfo()['crs']
lc9 = lc9.map(lambda img: img.reproject(crs, scale = 30))

lc9_clip = lc9.map(lambda img: img.clip(jasper_bbox))

ds = xr.open_dataset(
    lc9,
    engine='ee',
    crs = crs,
    scale = 30
)

But no matter whether I pass lc9 or lc9_clip, I got the same dimension (time: 6 X: 12299 Y: 310967). Do you know any workaround?

BTW, I also tried to convert the projection of the geometry to UTM and found it returns a polygon with 19 vertices. Got the same Value error in return.