Applied-GeoSolutions / gips

Geospatial Image Processing System
GNU General Public License v3.0
17 stars 5 forks source link

did I/ should I invert the HLS mask #502

Open bhbraswell opened 5 years ago

bhbraswell commented 5 years ago

I am wondering if the hls mask is inverted from the others (e.g. landsat). Such that when I apply the mask, I get something like the attached that leaves the clouds, and eliminates the good stuff. Alternatively is there a switch I can use to flip the mask so that the 1s are the good stuff?

Or possibly I breathed too many fumes on my morning commute. I think it's likely I'm missing something.

Any input on this appreciated. Thanks as always -- Rob

clouds

ra-tolson commented 5 years ago

The "1" value in the mask means "this pixel has clouds":

https://github.com/Applied-GeoSolutions/gips/blob/15d157702a21750094c140b92fe59fe0e4318691/gips/data/hls/hls.py#L221

I'm pretty sure that's by design, though I can't remember the specifics.

Edit: hey github permalinks to code render pretty good

bhbraswell commented 5 years ago

Hey thanks Tom. I see that and it makes sense in a way, but I'm just wondering for the purposes of using gips_mask if I need to manually flip 1 and 0, which I can totally do, or make a separate cloud product, whatever you guys think is the best way to go (if a solution doesn't already exist).

ircwaves commented 5 years ago

So, the reason for it being 1==cloud and 0==nodata is for mosaicking in the gdal_merge fashion. In the case where two scenes overlap there are redundant pixels. The cloudmask arranged this way results in a pixel being labeled cloud if it is identified as cloud in either image. Which is the default (cloudy-data-averse) GIPS stance. I believe this is consistent across drivers...but I'm checking.

bhbraswell commented 5 years ago

Thanks Ian, so I might be doing something double bass-ackwards as my grandad used to say. But, when I invoked gips_mask -p ndvi --pmask cloudmask on hls export, it turned all the non-cloud pixels into missing value, and left the clouds.

On top of that, when I tweaked the cloud mask product in the driver such that the 0s and 1s are reversed, then run gips_mask using that mask, I get a good result in that I see land instead of clouds.

ircwaves commented 5 years ago

Yeah. To apply a gips *cloudmask.tif (from any driver), you'll need to invert the mask. This flip was added ~8 months ago to improve optical data fusion.

I thought there used to be a gips_mask option to --invert, but don't see it now. It would need to be wired in something like this (UNTESTED):

diff --git a/gips/scripts/mask.py b/gips/scripts/mask.py
index 27a94b0..027a3cc 100644
--- a/gips/scripts/mask.py
+++ b/gips/scripts/mask.py
@@ -39,6 +39,7 @@ def main():
     group = parser.add_argument_group('masking options')
     group.add_argument('--filemask', help='Mask all files with this static mask', default=None)
     group.add_argument('--pmask', help='Mask files with this corresponding product', nargs='*', default=[])
+    group.add_argument('--invert', help='Invert the masks from corresponding products', nargs='*', default=[])
     h = 'Write mask to original image instead of creating new image'
     group.add_argument('--original', help=h, default=False, action='store_true')
     h = 'Overwrite existing files when creating new'
@@ -76,7 +77,11 @@ def main():
                         img.AddMask(mask_file[0])
                         meta = basename(args.filemask) + ' '
                     for mask in available_masks:
-                        img.AddMask(inv[date].open(mask)[0])
+                        mask_img = inv[date].open(mask)[0]
+                        if mask in args.invert:
+                            mask_img = mask_img.BXOR(1)
+                            meta += 'inverted-'
+                        img.AddMask(mask_img)
                         meta = meta + basename(inv[date][mask]) + ' '
                     if meta != '':
                         if args.original:
ra-tolson commented 5 years ago

If you want to submit a PR implementing that @bhbraswell, let me know if you'd like some help figuring out how to write a system test for it.