Merge pull request #149 from betonr/master

fix generate cube sentinel with 0.8.0 version -  get nodata values by params

Author: betonr

cube_builder_aws/cube_builder_aws/controller.py

@@ -468,7 +468,8 @@ def get_cube_status(self, cube_name):
 
 
     def start_harmonization_process(self, params):
-        _ = prepare_harm(self, params['scenes'], params['bucket_dst'], params['bucket_angles'])
+        _ = prepare_harm(self, params['scenes'], params['bucket_dst'], 
+                         params['bucket_angles'], params['satellite'])
 
         return dict(
             message='Harmonization processing started with succesfully'
@@ -534,10 +535,17 @@ def start_process(self, params):
 
         bands_list = []
         bands_ids_list = {}
+        quality_nodata = 0
+        nodata = 0
         for band in bands:
             if band.name.upper() not in [i['common_name'].upper() for i in indexes]:
                 bands_list.append(band.name)
                 bands_ids_list[band.id] = band.name
+                if band.name == quality_band:
+                    quality_nodata = band.nodata
+                else:
+                    nodata = band.nodata
+
             elif band._metadata and band._metadata.get('expression') and band._metadata['expression'].get('value'):
                 meta = deepcopy(band._metadata)
                 meta['data_type'] = band.data_type
@@ -564,7 +572,7 @@ def start_process(self, params):
         formatted_version = format_version(cube_infos.version)
         not_started = prepare_merge(self, cube_infos.name, cube_infos_irregular.name, collections, satellite,
             bands_list, bands_ids_list, bands_ql_list, float(bands[0].resolution_x), 
-            float(bands[0].resolution_y), int(bands[0].nodata), crs, quality_band, functions, formatted_version, 
+            float(bands[0].resolution_y), crs, nodata, quality_nodata, quality_band, functions, formatted_version, 
             params.get('force', False), mask, bands_expressions=bands_expressions, 
             indexes_only_regular_cube=params.get('indexes_only_regular_cube'), 
             landsat_harmonization=landsat_harmonization)

cube_builder_aws/cube_builder_aws/maestro.py

@@ -275,14 +275,15 @@ def harmonization(self, activity):
 # MERGE
 ###############################
 def prepare_merge(self, datacube, irregular_datacube, datasets, satellite, bands, bands_ids, 
-                  quicklook, resx, resy, nodata, crs, quality_band, functions, version,
+                  quicklook, resx, resy, crs, nodata, quality_nodata, quality_band, functions, version,
                   force=False, mask=None, bands_expressions=dict(), indexes_only_regular_cube=False,
                   landsat_harmonization=None):
     services = self.services
 
     # Build the basics of the merge activity
     activity = {}
     activity['action'] = 'merge'
+    activity['bucket_name'] = services.bucket_name
     activity['datacube'] = datacube
     activity['irregular_datacube'] = irregular_datacube
     activity['version'] = version
@@ -294,9 +295,9 @@ def prepare_merge(self, datacube, irregular_datacube, datasets, satellite, bands
     activity['quicklook'] = quicklook
     activity['resx'] = resx
     activity['resy'] = resy
-    activity['nodata'] = nodata
     activity['srs'] = crs
-    activity['bucket_name'] = services.bucket_name
+    activity['nodata'] = int(nodata)
+    activity['quality_nodata'] = int(quality_nodata)
     activity['quality_band'] = quality_band
     activity['functions'] = functions
     activity['force'] = force
@@ -549,7 +550,7 @@ def merge_warped(self, activity):
         ymax = float(activity['ymax'])
         dist_x = float(activity['dist_x'])
         dist_y = float(activity['dist_y'])
-        nodata = int(activity['nodata'])
+        nodata = int(activity['quality_nodata']) if is_quality_band else int(activity['nodata'])
 
         shape = activity.get('shape', None)
         if shape:
@@ -570,14 +571,13 @@ def merge_warped(self, activity):
 
         is_sentinel_landsat_quality_fmask = ('LANDSAT' in satellite or satellite == 'SENTINEL-2') and \
                                             (is_quality_band and activity_mask['nodata'] != 0)
-        source_nodata = 0
+        source_nodata = activity.get('source_nodata', 0)
 
         # Quality band is resampled by nearest, other are bilinear
         if is_quality_band:
             resampling = Resampling.nearest
 
-            nodata = activity_mask['nodata']
-            source_nodata = nodata
+            source_nodata = source_nodata if activity.get('source_nodata') else activity_mask['nodata']
 
             raster = numpy.zeros((numlin, numcol,), dtype=numpy.uint16)
             raster_merge = numpy.full((numlin, numcol,), dtype=numpy.uint16, fill_value=source_nodata)
@@ -616,96 +616,98 @@ def merge_warped(self, activity):
             new_url = url
             if landsat_harmonization:
                 if not is_quality_band:
-                    bucket_src = HARMONIZATION['landsat']['bucket_src']
-                    new_url = new_url.replace(bucket_src, landsat_harmonization['bucket_dst'])
+                    key_path = new_url.replace(HARMONIZATION['landsat']['bucket_src'], '')
+                    if services.s3_file_exists(bucket_name=landsat_harmonization['bucket_dst'], key=key_path, request_payer=True):
+                        bucket_src = HARMONIZATION['landsat']['bucket_src'].replace('s3://', '')
+                        new_url = new_url.replace(bucket_src, landsat_harmonization['bucket_dst'])
 
-                from rasterio.session import AWSSession
+            from rasterio.session import AWSSession
 
-                aws_session = AWSSession(services.session, requester_pays=True)
-                with rasterio.Env(aws_session, AWS_SESSION_TOKEN=""):
+            aws_session = AWSSession(services.session, requester_pays=True)
+            with rasterio.Env(aws_session, AWS_SESSION_TOKEN=""):
 
-                    with rasterio.open(new_url) as src:
+                with rasterio.open(new_url) as src:
 
-                        kwargs = src.meta.copy()
+                    kwargs = src.meta.copy()
+                    kwargs.update({
+                        'width': numcol,
+                        'height': numlin
+                    })
+                    if not shape:
                         kwargs.update({
-                            'width': numcol,
-                            'height': numlin
+                            'crs': activity['srs'],
+                            'transform': transform
                         })
-                        if not shape:
-                            kwargs.update({
-                                'crs': activity['srs'],
-                                'transform': transform
-                            })
-
-                        if src.profile['nodata'] is not None:
-                            source_nodata = src.profile['nodata']
-                        
-                        elif 'source_nodata' in activity:
-                            source_nodata = activity['source_nodata']
 
-                        elif 'LANDSAT' in satellite and not is_quality_band:
-                            source_nodata = nodata if src.profile['dtype'] == 'int16' else 0
+                    if src.profile['nodata'] is not None:
+                        source_nodata = src.profile['nodata']
+                    
+                    elif 'source_nodata' in activity:
+                        source_nodata = activity['source_nodata']
+
+                    elif 'LANDSAT' in satellite and not is_quality_band:
+                        source_nodata = nodata if src.profile['dtype'] == 'int16' else 0
+
+                    elif 'CBERS' in satellite and not is_quality_band:
+                        source_nodata = nodata
 
-                        elif 'CBERS' in satellite and not is_quality_band:
-                            source_nodata = nodata
+                    kwargs.update({
+                        'nodata': source_nodata
+                    })
 
+                    if is_quality_band and activity_mask.get('bits'):
                         kwargs.update({
-                            'nodata': source_nodata
+                            'blockxsize': 2048,
+                            'blockysize': 2048,
+                            'tiled': True
                         })
 
-                        if is_quality_band and activity_mask.get('bits'):
-                            kwargs.update({
-                                'blockxsize': 2048,
-                                'blockysize': 2048,
-                                'tiled': True
-                            })
-
-                        with MemoryFile() as memfile:
-                            with memfile.open(**kwargs) as dst:
-                                if shape:
-                                    raster = src.read(1)
-                                else:
-                                    reproject(
-                                        source=rasterio.band(src, 1),
-                                        destination=raster,
-                                        src_transform=src.transform,
-                                        src_crs=src.crs,
-                                        dst_transform=transform,
-                                        dst_crs=activity['srs'],
-                                        src_nodata=source_nodata,
-                                        dst_nodata=nodata,
-                                        resampling=resampling)
-
-                                if not is_quality_band or is_sentinel_landsat_quality_fmask:
-                                    valid_data_scene = raster[raster != nodata]
-                                    raster_merge[raster != nodata] = valid_data_scene.reshape(numpy.size(valid_data_scene))
-
-                                    valid_data_scene = None
-                                else:
-                                    factor = raster * raster_mask
-                                    raster_merge = raster_merge + factor
+                    with MemoryFile() as memfile:
+                        with memfile.open(**kwargs) as dst:
+                            if shape:
+                                raster = src.read(1)
+                            else:
+                                reproject(
+                                    source=rasterio.band(src, 1),
+                                    destination=raster,
+                                    src_transform=src.transform,
+                                    src_crs=src.crs,
+                                    dst_transform=transform,
+                                    dst_crs=activity['srs'],
+                                    src_nodata=source_nodata,
+                                    dst_nodata=nodata,
+                                    resampling=resampling)
+
+                            if not is_quality_band or is_sentinel_landsat_quality_fmask:
+                                valid_data_scene = raster[raster != nodata]
+                                raster_merge[raster != nodata] = valid_data_scene.reshape(numpy.size(valid_data_scene))
+
+                                valid_data_scene = None
+                            else:
+                                factor = raster * raster_mask
+                                raster_merge = raster_merge + factor
 
-                                    raster_mask[raster != nodata] = 0
+                                raster_mask[raster != nodata] = 0
 
-                                if build_provenance:
-                                    raster_masked = numpy.ma.masked_where(raster == nodata, raster)
+                            if build_provenance:
+                                raster_masked = numpy.ma.masked_where(raster == nodata, raster)
 
-                                    where_valid = numpy.invert(raster_masked.mask)
-                                    raster_provenance[where_valid] = datasets.index(platforms[url])
+                                where_valid = numpy.invert(raster_masked.mask)
+                                raster_provenance[where_valid] = datasets.index(platforms[url])
 
-                                    where_valid = None
-                                    raster_masked = None
+                                where_valid = None
+                                raster_masked = None
 
-                                if template is None:
-                                    template = dst.profile
+                            if template is None:
+                                template = dst.profile
 
-                                    template['driver'] = 'GTiff'
+                                template['driver'] = 'GTiff'
 
-                                    raster_blocks = list(dst.block_windows())
+                                raster_blocks = list(dst.block_windows())
 
-                                    if not is_quality_band:
-                                        template.update({'dtype': 'int16'})
-                                        template['nodata'] = nodata
+                                if not is_quality_band:
+                                    template.update({'dtype': 'int16'})
+                                    template['nodata'] = nodata
 
 
         raster = None
@@ -734,7 +736,7 @@ def merge_warped(self, activity):
             else:
                 raster_merge = raster_merge.astype(numpy.uint8)
 
-            nodata = activity_mask['nodata']
+            nodata = activity['quality_nodata']
 
             # Save merged image on S3
             create_cog_in_s3(services, template, key, raster_merge, bucket_name, nodata=nodata)
@@ -779,7 +781,7 @@ def next_blend(services, mergeactivity):
     blendactivity = {}
     blendactivity['action'] = 'blend'
     for key in ['datasets', 'satellite', 'bands', 'quicklook', 'srs', 'functions', 'bands_ids',
-                'tileid', 'start', 'end', 'dirname', 'nodata', 'bucket_name', 'quality_band',
+                'tileid', 'start', 'end', 'dirname', 'nodata', 'bucket_name', 'quality_band', 'quality_nodata',
                 'internal_bands', 'force', 'version', 'datacube', 'irregular_datacube', 'mask',
                 'bands_expressions', 'indexes_only_regular_cube', 'empty_file', 'landsat_harmonization']:
         blendactivity[key] = mergeactivity.get(key, '')
@@ -950,9 +952,7 @@ def blend(self, activity):
     band = activity['band']
     numscenes = len(activity['scenes'])
 
-    nodata = int(activity.get('nodata', -9999))
-    if band == activity['quality_band']:
-        nodata = activity_mask['nodata']
+    nodata = int(activity.get('nodata', -9999)) if band != activity['quality_band'] else activity['quality_nodata']
 
     # Check if band ARDfiles are in activity
     for date_ref in activity['scenes']:
@@ -1010,6 +1010,8 @@ def blend(self, activity):
                 resolution = 10
                 mask_tuples.append((100. * efficacy / resolution, key))
 
+            mask_tuples = sorted(mask_tuples, reverse=True)
+
         provenance_merge_map = dict()
         build_provenance = activity.get('internal_band') == PROVENANCE_NAME
         build_clear_observation = activity.get('internal_band') == CLEAR_OBSERVATION_NAME

cube_builder_aws/cube_builder_aws/services.py

@@ -516,7 +516,9 @@ def search_STAC(self, activity):
             _ = stac['instance'].catalog
 
             filter_opts['collections'] = [stac['collection']]
-            filter_opts['query'] = dict(collections=[stac['collection']])
+            # TODO: STAC https://landsatlook.usgs.gov/  does not conform to specification
+            if '.usgs.gov' in stac['url']:
+                filter_opts['query'] = dict(collections=[stac['collection']])
             items = stac['instance'].search(filter=filter_opts)
 
             res = self._parse_stac_result(items, stac['collection'], bands, harm_bands_map)
@@ -561,11 +563,14 @@ def create_bucket(self, name, requester_pay=True):
             return False
         return True
 
-    def s3_file_exists(self, bucket_name=None, key=''):
+    def s3_file_exists(self, bucket_name=None, key='', request_payer=False):
         try:
             if not bucket_name:
                 bucket_name = self.bucket_name
-            return self.S3client.head_object(Bucket=bucket_name, Key=key)
+            if request_payer:
+                return self.S3client.head_object(Bucket=bucket_name, Key=key, RequestPayer='requester')
+            else:
+                return self.S3client.head_object(Bucket=bucket_name, Key=key)
         except ClientError:
             return False