Commit 34d7f4a0 authored by Gaetano Raffaele's avatar Gaetano Raffaele
Browse files

ENH: new light version of lsgrm without aggregation.

No related merge requests found
Showing with 35 additions and 8 deletions
+35 -8
...@@ -76,9 +76,13 @@ def grm_process_tile(input_image, params : LSGRMParams, tile_width, tile_height, ...@@ -76,9 +76,13 @@ def grm_process_tile(input_image, params : LSGRMParams, tile_width, tile_height,
op = regionprops(seg['array'][:,:,0].astype(np.int32)) op = regionprops(seg['array'][:,:,0].astype(np.int32))
to_del = [] to_del = []
for o in op: for o in op:
if not ((tie_lines[0] <= o.bbox[0] < tie_lines[2]) if not (tie_lines[0] <= o.bbox[0] < tie_lines[2]
and tie_lines[1] <= o.bbox[1] < tie_lines[3]): and tie_lines[1] <= o.bbox[1] < tie_lines[3]):
to_del.append(o.label) # This one to check if a potentially deleted object trespass the overlap area,
# in which case it is kept.
if not (2*params.margin <= o.bbox[2] < tile_height and
2*params.margin <= o.bbox[3] < tile_width):
to_del.append(o.label)
seg['array'][np.isin(seg['array'], to_del)] = 0 seg['array'][np.isin(seg['array'], to_del)] = 0
seg['array'], nlab = label(seg['array'], background=0, connectivity=1, return_num=True) seg['array'], nlab = label(seg['array'], background=0, connectivity=1, return_num=True)
seg['array'] = np.ascontiguousarray(seg['array']).astype(float) seg['array'] = np.ascontiguousarray(seg['array']).astype(float)
...@@ -158,23 +162,29 @@ def get_ls_seg_parameter(input_image, roi=None, margin=0, n_proc=None, memory=No ...@@ -158,23 +162,29 @@ def get_ls_seg_parameter(input_image, roi=None, margin=0, n_proc=None, memory=No
return in_img, n_tiles_x, n_tiles_y, nominal_tw, nominal_th return in_img, n_tiles_x, n_tiles_y, nominal_tw, nominal_th
def lsgrm_light(input_image, params : LSGRMParams, out_seg, n_proc=None, memory=None, roi=None, force_parallel=False): def lsgrm_light(input_image, params : LSGRMParams, out_seg, n_proc=None, memory=None,
roi=None, force_parallel=False, remove_tiles_if_useless=True):
# Check output file type # Check output file type
ext = os.path.splitext(out_seg)[-1] ext = os.path.splitext(out_seg)[-1]
if ext in ['.tif']: if ext == '.vrt':
mode = 'vrt'
elif ext in ['.tif']:
mode = 'raster' mode = 'raster'
elif ext in ['.shp', '.gpkg', 'gml']: elif ext in ['.shp', '.gml', '.gpkg']:
mode = 'vector' mode = 'vector'
elif ext in ['.vrt']:
mode = 'vrt'
else: else:
raise ValueError('Output type {} not recognized/supported.'.format(ext)) raise ValueError('Output type {} not supported.'.format(ext))
in_img, n_tiles_x, n_tiles_y, nominal_tw, nominal_th = get_ls_seg_parameter(input_image, roi, params.margin, n_proc, in_img, n_tiles_x, n_tiles_y, nominal_tw, nominal_th = get_ls_seg_parameter(input_image, roi, params.margin, n_proc,
memory, force_parallel) memory, force_parallel)
print('[INFO] Using a layout of {} x {} tiles.'.format(n_tiles_x,n_tiles_y))
print('[INFO] Nominal tile size w/margin: {} x {} pixels'.format(nominal_tw+2*params.margin,
nominal_th+2*params.margin))
if (n_tiles_x == 1 and n_tiles_y == 1): if (n_tiles_x == 1 and n_tiles_y == 1):
# Fallback to classical GRM # Fallback to classical GRM
print('[INFO] Fallback to one-tile GRM.')
grm = otb.Registry.CreateApplication('GenericRegionMerging') grm = otb.Registry.CreateApplication('GenericRegionMerging')
grm.SetParameterInputImage('in', in_img.GetParameterOutputImage('out')) grm.SetParameterInputImage('in', in_img.GetParameterOutputImage('out'))
grm.SetParameterFloat('threshold', params.threshold) grm.SetParameterFloat('threshold', params.threshold)
...@@ -198,9 +208,26 @@ def lsgrm_light(input_image, params : LSGRMParams, out_seg, n_proc=None, memory= ...@@ -198,9 +208,26 @@ def lsgrm_light(input_image, params : LSGRMParams, out_seg, n_proc=None, memory=
tile.nodata = 0 tile.nodata = 0
cumul_label += l cumul_label += l
print('[INFO] Output mode {}.'.format(mode))
if mode == 'vrt': if mode == 'vrt':
vrtopt = gdal.BuildVRTOptions(separate=False, srcNodata=0, VRTNodata=0) vrtopt = gdal.BuildVRTOptions(separate=False, srcNodata=0, VRTNodata=0)
gdal.BuildVRT(out_seg, [x[0] for x in out_tiles], options=vrtopt) gdal.BuildVRT(out_seg, [x[0] for x in out_tiles], options=vrtopt)
elif mode in ['raster', 'vector']:
mos = otb.Registry.CreateApplication('Mosaic')
mos.SetParameterStringList('il', [x[0] for x in out_tiles])
mos.SetParameterInt('nodata', 0)
if mode == 'raster':
mos.SetParameterOutputImagePixelType('out', otb.ImagePixelType_uint32)
mos.SetParameterString('out', out_seg)
mos.ExecuteAndWriteOutput()
elif mode == 'vector':
mos.Execute()
vec = otb.Registry.CreateApplication('SimpleVectorization')
vec.SetParameterInputImage('in', mos.GetParameterOutputImage('out'))
vec.SetParameterString('out', out_seg)
vec.ExecuteAndWriteOutput()
if remove_tiles_if_useless:
[os.remove(x[0]) for x in out_tiles]
return out_seg return out_seg
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment