WMS - Getting data clipped to a polygon

Apps / Services
1

Hello,

I’m getting started accessing the OGC API through Python to download Sentinel data.
I would like to know if it is possible to make a request for data clipped to a polygon instead of a rectangular bounding box.

It’s possible to do in the EO browser, but I haven’t found a way to do it directly from the API.

I can of course download get the image cut to the bounding box, and then clip it, but being able to directly download only the data I need would save me some computer ressources, and possibly some processing units.

Thank you in advance for your help.

1 Like
2

Hello,

The way I set a polygon in my Python OGC requests is by using the custom_url_params argument in my request. You can find more information here.

You will still need to set a bounding box, but you can add your polygon as an optional parameter. Here is a code snippet in Python:

# Set bounding box
bbox_loc = BBox(bbox=[7.3994, 43.7190, 7.4531, 43.7564], crs=CRS.WGS84)

# Polygon WKT (pasted polygon WKT, you can import a file using shapely for example)
geometryWKT = "POLYGON((7.41877555847168 43.72459139811755,7.410106658935547 43.7294293330051,7.414312362670898 43.732158252099964,7.413797378540039 43.734700996616425,7.423152923583985 43.74102638042922,7.432508468627931 43.735693257872555,7.41877555847168 43.72459139811755))"

# WMS request 
request = WmsRequest(data_source=DataSource.SENTINEL2_L1C,
layer='1_TRUE_COLOR',
data_folder=data_folder,
bbox=bbox_loc,
time='latest', 
width=512,
config=config,
custom_url_params={CustomUrlParam.GEOMETRY: geometryWKT}
 )

When making a request for data clipped to a polygon, it’s a good idea to activate the transparency for the areas outside the polygon in the Evalscript that you are calling. Here I am requesting the layer 1_TRUE_COLOR from my dashboard, which has the following evalscript:

//VERSION=3
function setup() {
  return {
    input: ["B02", "B03", "B04", "dataMask"],
    output: { bands: 4 }
  };
}

function evaluatePixel(sample) {
  return [2.5 * sample.B04, 2.5 * sample.B03, 2.5 * sample.B02, sample.dataMask];
}
4 Likes
3

Thank you very much for taking the time to write such a comprehensive answer, this is exactly what I was looking for, and you explained it very clearly. I hope I can get enough experience to one day give this kind of contribution to the forum.

3 Likes
4

Is there an option to use this also inside SentinelHubRequest without using WMS? for example, adding this parameter to the request?

5

You don’t need this option in SentinelHubRequest, because you can directly specify a Geometry (see options).

The Geometry can be initialised with any of the following geometry representations:

  • shapely.geometry.Polygon or shapely.geometry.MultiPolygon
  • A GeoJSON dictionary with (multi)polygon coordinates
  • A WKT string with (multi)polygon coordinates

e.g.

geometry = Geometry(geometry={"type":"Polygon","coordinates":[[[12.448527,41.92361],[12.480804,41.927186],[12.529563,41.919651],[12.478744,41.922972],[12.448527,41.92361]]]}, crs=CRS.WGS84)
1 Like
6

In SentinelHubRequest, how is possible send more than one polygon? I try create a Multipolygon but it seems not to be the way. What is the right way?

7

Hi,

Sentinelhub supports multipolygons. In the documentation you can see information about the use of a Geometry in SentinelHubRequest:

It currently supports polygons and multipolygons. It can be initialized with any of the following geometry representations:

  • shapely.geometry.Polygon or shapely.geometry.MultiPolygon
  • A GeoJSON dictionary with (multi)polygon coordinates
  • A WKT string with (multi)polygon coordinates

How did you specify the parameter in your request?

9

That’s a great idea, Maxim. So, how can I do this transparency to NDVI? I’m trying doing as follow but always return a single band

//VERSION=3
function setup() {
  return {
    input: [{
      bands:["B04", "B08", "dataMask"],
    }],
    output: {
      id: "default",
      bands: 2,
    }
  }
}
function evaluatePixel(sample, scene) {
    
    var ndvi = (sample.B08 - sample.B04) / (sample.B08 + sample.B04)
     if (sample.dataMask == 1)  {
    return [ndvi, sample.dataMask ]
  } else {
    return [-1, -1, -1]}
    if (ndvi<-0) return [1.65,0,0.38] 
    else if (ndvi<0.1) return [0.87,0.85,0.61]
    else if (ndvi<0.2) return [0.8,0.78,0.51]
    else if (ndvi<0.3) return [0.74,0.72,0.42]
    else if (ndvi<0.4) return [0.69,0.76,0.38]
    else if (ndvi<0.5) return [0.64,0.8,0.35]
    else if (ndvi<0.6) return [0.57,0.75,0.32]
    else if (ndvi<0.7) return [0.5,0.7,0.28]
    else if (ndvi<0.8) return [0.44,0.64,0.25]
    else if (ndvi<0.9) return [0,1.04,0.55] 
    else return [sample.dataMask]
}
10

There are a couple of inconsistencies in the evalscript you posted, but these can be easily fixed.

  1. The lines starting from if (ndvi<-0) statement are never called because your first if/else statement catches all situations.
  2. Your return an array of 2 values in your first if statement, and 3 values in your else, but have set your output bands to 2 in the setup function (see documentation).

It depends what you want to return: the NDVI values or a visualisation of NDVI. If you want to return a visualisation of NDVI with only data showing (where dataMask = 1), you would want to set your output bands to 4 and return 4 values in your evaluatePixel function, where the first three correspond to RGB values, and the last one is the dataMask band (because 0 = noData, and will be interpreted as transparent).

1 Like