Within the snappy framework, I would like to apply a mask to a sentinel-1 image using a wkt string.
I know this can be done with a Land-Sea-Mask, but the geometry needs to be added to the product first.
Adding the geometry can be done directly using the Import-Vector, but it only allows a shp. from file (please correct me if i’m wrong).
To use a wkt string I have tried the example script in snappy_geo_roi.py
(see below for example to try-need to change the IO file path of course, and the file name is original so people want they can search for it here: https://scihub.copernicus.eu/dhus/#/home)
However, the vector it creates is empty (see attached pic) . So when I try and run the Land-Sea-Mask I get an error:
“RuntimeError: org.esa.snap.core.gpf.OperatorException: expression: Undefined symbol ‘shape’. due to Undefined symbol ‘shape’.”
I know this error also appears if the image is out of bounds, but im confident it is within the extent (I drew it manually in the toolbox):
As a temporary workaround, I dont mind if i have to write the wkt to a shapefile, and then import it via ‘Import-Vector’ operation (export to shapefile is possible from the he Sentinel GUI).
If there is another way to import the wkt, please let me know! thanks in advance!
Test script:
> import sys
> import os
> from snappy import (ProgressMonitor, VectorDataNode,
> WKTReader, ProductIO, PlainFeatureFactory,
> SimpleFeatureBuilder, DefaultGeographicCRS,
> ListFeatureCollection, FeatureUtils, WKTReader, HashMap, GPF)
> test_folder = r'C:\Users\myName\Documents\sentinel_1'
> inputFileName = 'S1B_IW_GRDH_1SDV_20161205T061407_20161205T061432_003257_0058D8_C74B.zip'
> wellKnownText = 'POLYGON ((-0.3788970856103323 51.80800809372078, -0.376076940574716 51.80800809372078,
> -0.376076940574716 51.80640687188327, -0.3788970856103323 51.80640687188327,
> -0.3788970856103323 51.80800809372078)) '
> geom_name = 'shape'
> full_path_input = os.path.join(test_folder,inputFileName)
> geometry = WKTReader().read(wellKnownText)
> product = ProductIO.readProduct(full_path_input)
> wktFeatureType = PlainFeatureFactory.createDefaultFeatureType(DefaultGeographicCRS.WGS84)
> featureBuilder = SimpleFeatureBuilder(wktFeatureType)
> wktFeature = featureBuilder.buildFeature(geom_name)
> wktFeature.setDefaultGeometry(geometry)
> newCollection = ListFeatureCollection(wktFeatureType)
> newCollection.add(wktFeature)
> #productFeatures = FeatureUtils.clipFeatureCollectionToProductBounds(newCollection, product, None, ProgressMonitor.NULL)
> node = VectorDataNode(geom_name, wktFeatureType)
> product.getVectorDataGroup().add(node)
> #vdGroup = product.getVectorDataGroup()
> #for i in range(vdGroup.getNodeCount()):
> # print('Vector data = ', vdGroup.get(i).getName())
> #
> #maskGroup = product.getMaskGroup()
> #for i in range(maskGroup.getNodeCount()):
> # print('Mask = ', maskGroup.get(i).getName())
> target_name = 'test_empty_geom.dim'
> """Can manually check the outout file in sentinel toolbox for missing 'wkt coords'"""
> target_name = 'test_empty_geom.dim'
> full_path_vector = os.path.join(test_folder,target_name)
> ProductIO.writeProduct(product, full_path_vector, 'BEAM-DIMAP')
> print('Done with test to add product vector')