Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
73 changes: 40 additions & 33 deletions src/spurt/workflows/emcf/_bulk_offset.py
Original file line number Diff line number Diff line change
Expand Up @@ -52,41 +52,48 @@ def get_bulk_offsets(
olap_counts.append(grp["c1"].shape[0])
offsets.append(grp["stats"][:, 5])

# Get number of bands to merge
nbands = len(offsets[0])
ntiles = tiledata.ntiles
counts = np.zeros(ntiles)

# Add counts of valid pixels in each tile
arr = stack.read_temporal_coherence(np.s_[:, :]) > stack.temp_coh_threshold
for ii, tile in enumerate(tiledata.tiles):
counts[ii] = np.sum(arr[tile.space])

# If integer solver requested
logger.info(f"Solving for bulk offsets with method: {mrg_settings.method}")
if mrg_settings.bulk_method == "integer":
bulk_offset: np.ndarray = np.zeros((nbands, ntiles), dtype=np.int32)
obj: np.ndarray = np.zeros(nbands, dtype=np.int32)

# Solve band by band
for ii in range(nbands):
logger.info(f"Computing bulk offsets for band {ii + 1}")
off: np.ndarray = np.array([x[ii] for x in offsets])
bulk_offset[ii, :], obj[ii] = _solve_int_offsets(
overlaps, labels, off, counts
)
# Use L2 method
elif mrg_settings.bulk_method == "L2":
# These can be floating point
off = np.zeros((len(overlaps), nbands))
for ii in range(len(overlaps)):
off[ii, :] = offsets[ii]

bulk_offset, obj = _solve_l2_min(overlaps, off, ntiles, counts)
if offsets:
# Get number of bands to merge
nbands = len(offsets[0])
ntiles = tiledata.ntiles
counts = np.zeros(ntiles)

# Add counts of valid pixels in each tile
arr = stack.read_temporal_coherence(np.s_[:, :]) > stack.temp_coh_threshold
for ii, tile in enumerate(tiledata.tiles):
counts[ii] = np.sum(arr[tile.space])

# If integer solver requested
logger.info(f"Solving for bulk offsets with method: {mrg_settings.method}")
if mrg_settings.bulk_method == "integer":
bulk_offset: np.ndarray = np.zeros((nbands, ntiles), dtype=np.int32)
obj: np.ndarray = np.zeros(nbands, dtype=np.int32)

# Solve band by band
for ii in range(nbands):
logger.info(f"Computing bulk offsets for band {ii + 1}")
off: np.ndarray = np.array([x[ii] for x in offsets])
bulk_offset[ii, :], obj[ii] = _solve_int_offsets(
overlaps, labels, off, counts
)
# Use L2 method
elif mrg_settings.bulk_method == "L2":
# These can be floating point
off = np.zeros((len(overlaps), nbands))
for ii in range(len(overlaps)):
off[ii, :] = offsets[ii]

bulk_offset, obj = _solve_l2_min(overlaps, off, ntiles, counts)

else:
errmsg = f"Unsupported bulk offset method {mrg_settings.bulk_method}"
raise RuntimeError(errmsg)

else:
errmsg = f"Unsupported bulk offset method {mrg_settings.bulk_method}"
raise RuntimeError(errmsg)
# Write empty file as availability of this file marks completion of
# this stage
bulk_offset = np.empty(0)
obj = np.empty(0)

# Write HDF5 file with bulk offsets
with h5py.File(str(gen_settings.offsets_filename), "w") as fid:
Expand Down
13 changes: 10 additions & 3 deletions src/spurt/workflows/emcf/_merge.py
Original file line number Diff line number Diff line change
Expand Up @@ -95,11 +95,18 @@ def merge_tiles(
# Initialize each tile for the right band with bulk offset
for jj, tile in tiles.items():
tile.reset_band_index(idx)
tile.add_correction(np.float32(-2 * np.pi * bulk_offsets[idx, jj]))
if bulk_offsets.size != 0:
tile.add_correction(np.float32(-2 * np.pi * bulk_offsets[idx, jj]))
else:
# This is to set tile metadata
tile.add_correction(np.float32(0.0))
tile.increment_correction()

# Adjust the tiles
_adjust_tiles(tiles, overlap_map, gen_settings, max_degree, debug_stats=False)
if bulk_offsets.size != 0:
# Adjust the tiles
_adjust_tiles(
tiles, overlap_map, gen_settings, max_degree, debug_stats=False
)

# Write file to output band-by-band
for ii in idx:
Expand Down