swisstopo / topo-satromo

Erdbeobachtungs-Satellitendaten fürs Trockenheitsmonitoring (SATROMO)
BSD 3-Clause "New" or "Revised" License
14 stars 3 forks source link

Remove 0 / no-data value for a Google Drive export for areas outside of of the export geometry's bounding box #56

Closed davidoesch closed 3 months ago

davidoesch commented 8 months ago

Is your feature request related to a problem? Please describe. When exporting the image, I use the bounding box of the geometry.

The no-data values within the bounding box are correct 9999) but I cannot figure out how to set the no-data value of the pixels outside of the bounding box; these seem to default to 0.

image

Describe the solution you'd like

davidoesch commented 5 months ago

it is not fixed yet

davidoesch commented 3 months ago

needs to be fixed in the satromo_publisher step at https://github.com/swisstopo/topo-satromo/blob/30e4e884f27bb841fcaff1a7c0c29c5ffbeb7360/satromo_publish.py#L671 as seprate util function?, since in the GEE assets it is already masked correctly The approach is: to get from metafile the asset ID the get the masked /non masked area

`` asset_id = 'projects/satromo-int/assets/COL_S2_SR_HARMONIZED_SWISS/S2-L2A_mosaic_2024-07-05T104021_bands-10m'

Load the asset.

asset = ee.Image(asset_id)

Select the first band of the image and convert it to an integral image.

first_band = asset.select(0).toInt()

Create a binary mask.

binary_mask = first_band.mask()

Invert the mask to get the non-masked areas.

non_masked = binary_mask.Not()

Convert the non-masked areas to vectors.

vectors = non_masked.reduceToVectors(geometryType='polygon', maxPixels=1e9, bestEffort=True)

Convert the vectors to GeoJSON.

geojson = vectors.getInfo()

Extract the features from GeoJSON

features = geojson['features']

Extract geometries and properties

geometries = [shape(feature['geometry']) for feature in features] properties = [feature['properties'] for feature in features]

Create a GeoDataFrame

gdf = gpd.GeoDataFrame(properties, geometry=geometries)

Define output shapefile path

output_shapefile = 'output.shp'

Define the CRS

crs = {'init': 'epsg:2056'}

Export to Shapefile with specified CRS

gdf.to_file(output_shapefile, driver='ESRI Shapefile', crs=crs)

print(f"Shapefile exported successfully to {output_shapefile}")

``

davidoesch commented 3 months ago

We solve thsi differently: We precompute for each orbit a buffer https://code.earthengine.google.com/279232b086bcfb5d456fe26eda6b09c6 the apply this buffer when creating the mosaic in satromo_publisher: Advantage: we save EECU hours, and we alawys have the orginal data masked in GEE assets. Disadvantage: with orbit shifts we will miss some tiles

davidoesch commented 3 months ago

fixed on dev