white water border in pixel classification

[krd57]

When running SDS_shoreline.extract_shorelines the images have a white water/water border over the entire image. I believe this is due to the way the 0 values at the shore are being masked out in SDS_preprocess.preprocess_single. I have figured out a fix that switches the order of the logic to remove the 0s at the border.

Here is the detection image with the border:

vs with the border masked out:

[2320sharon]

Hi @krd57
I've encountered this same bug with coastsat for the S2 imagery. Can you share your fix to this issue?

Thanks for reporting this bug!

[krd57]

@2320sharon

If you change lines 264-267 in SDS_preprocess.py to:

im_zeros = np.zeros(im_nodata.shape).astype(bool)
im_zeros = np.logical_or(np.isin(im_ms[:,:,1],0), im_zeros) # Green
im_zeros = np.logical_or(np.isin(im_ms[:,:,3],0), im_zeros) # NIR
im_zeros = np.logical_or(np.isin(im_ms[:,:,4],0), im_zeros) # SWIR

it should fix it. I have not encountered any issues with this change in subsequent processes but have not tested everything. Let me know if that works for you!

[2320sharon]

Thank you so much @krd57
I'll give this a solution a shot and see if it works. Hopefully it fixes that issue with the S2 imagery and doesn't cause additional artifact with normal images.

[2320sharon]

Hey @krd57
I tried your solution on a few different locations were there was that issue with the ring around the S2 imagery and your solution is working nicely. You might want to submit a pull request with your fix.

Thank you for your solution

[kvos]

hi @krd57, thanks for looking into this, it's definitely a good idea to fix that.
To provide some context, the band on the edge is due to the different resolution of the B,G,R,NIR at 10m and SWIR1 at 20m in Sentinel-2. The bands have different resolution AND footprint. When I reproject the 20 m band onto the 10m pixel grid (with bilinear interpolation in GDAL), it creates that band of 0 SWIR1 values on the edge. Then the classifier does some funky stuff there as the gradients are high so it thinks that it is whitewater...
So if I understood your solution, you are adding the pixels that have 0 intensity in the SWIR1 to the nodata mask, is that right? While it may work here, I don't know how general this solution is, is a 0 SWIR1 intensity really unlikely? Another solution could be to pad the edges, something like:

padding = 5 # pad 5 pixels on each edge, should be more than enough
im_edges = np.zeros(cloud_mask.shape).astype(bool)
im_edges[pad:-pad, pad:-pad] = True
im_nodata = np.logical_or(im_nodata,im_edges)

You can try this as well see if it works.

[2320sharon]

def find_edge_padding(im_band):
    # Assuming non-data values are zeros. Adjust the condition if needed.
    is_data = im_band != 0

    # Function to find padding for one edge
    def find_edge_data(is_data_along_edge):
        for idx, has_data in enumerate(is_data_along_edge):
            if has_data:
                return idx
        return len(is_data_along_edge)  # Return full length if no data found

    # Calculate padding for each side
    top_padding = find_edge_data(np.any(is_data, axis=1))
    bottom_padding = find_edge_data(np.any(is_data, axis=1)[::-1])
    left_padding = find_edge_data(np.any(is_data, axis=0))
    right_padding = find_edge_data(np.any(is_data, axis=0)[::-1])

    return top_padding, bottom_padding, left_padding, right_padding

def pad_edges(im_swir, im_nodata) -> np.ndarray:
    top_pad, bottom_pad, left_pad, right_pad = find_edge_padding(im_swir)

    # Apply this padding to your masks or other arrays as needed
    im_nodata[:top_pad, :] = True
    im_nodata[-bottom_pad:, :] = True
    im_nodata[:, :left_pad] = True
    im_nodata[:, -right_pad:] = True
    return im_nodata

Hey @kvos
Thank you for explaining how the wrapping process creates these artifacts. Its helpful to understand where and why these artifacts are created. Your solution with the padding is on the right track, but it needs account for the fact that those edge artifacts aren't always the same size. I created these functions using your idea and these functions find the edges of the image that are only 0s and adds those edges to the no data mask. This should make the solution work for images of different sizes.

[2320sharon]

Here is the same image but with the edge solution I implemented. (please ignore the reference shoreline buffer plotted on the image that's related to some additional code I wrote)

[kvos]

looks good thanks @2320sharon , nicely written code. Feel free to make a pull request, just add it in SDS_preprocess.preprocess_single() inside an if satnamme == 'S2'. Does this also happen in Landsat? it could since all the bands are getting pan-sharpened/down-sampled to the 15m pixel grid using that same gdal.warp function.

[2320sharon]

I'll submit a PR with this fix later today after I test this solution a bit more.
I'll make sure to add it to SDS_preprocess.preprocess_single() inside an if satnamme == 'S2'.
Strangely enough I've never encountered an issue with this with any of the LandSat imagery. Only the S2 imagery.

[2320sharon]

Update

I found out why the code wasn't working correctly for all S2 imagery. If the S2 imagery didn't have any 0's at the bottom or right side of the image, the find_edge_padding(im_swir) would return a 0 for the bottom_pad and right_pad. When the the 0 was passed to a line like im_nodata[-bottom_pad:, :] = True it would result in im_nodata[0:, :] = True which would set all the rows and columns to True resulting in the entire no_data masking being set to True, which was incorrect.

The solution is to check for if the bottom or right side padding is greater than 0, then the code works great.

New Code

def pad_edges(im_swir, im_nodata) -> np.ndarray:
    top_pad, bottom_pad, left_pad, right_pad = find_edge_padding(im_swir)
    # Apply this padding to your masks or other arrays as needed

    # if bottom pad is 0 the entire image gets set to True
    if bottom_pad >0:
        im_nodata[-bottom_pad:, :] = True
    # if right pad is 0 the entire image gets set to True
    if right_pad >0:
        im_nodata[:, -right_pad:] = True
    
    im_nodata[:, :left_pad] = True
    im_nodata[:top_pad, :] = True
    return im_nodata

[kvos]

great Sharon, feel free to PR directly to master. Thanks.