36 #ifndef VIGRA_BLOCKWISE_WATERSHEDS_HXX
37 #define VIGRA_BLOCKWISE_WATERSHEDS_HXX
39 #include "threadpool.hxx"
40 #include "multi_array.hxx"
41 #include "multi_gridgraph.hxx"
42 #include "blockify.hxx"
43 #include "blockwise_labeling.hxx"
44 #include "overlapped_blocks.hxx"
55 namespace blockwise_watersheds_detail
58 template <
class DataArray,
class DirectionsBlocksIterator>
59 void prepareBlockwiseWatersheds(
const Overlaps<DataArray>& overlaps,
60 DirectionsBlocksIterator directions_blocks_begin,
61 BlockwiseLabelOptions
const & options)
63 static const unsigned int N = DataArray::actual_dimension;
64 typedef typename MultiArrayShape<DataArray::actual_dimension>::type Shape;
65 typedef typename DirectionsBlocksIterator::value_type DirectionsBlock;
66 Shape shape = overlaps.shape();
67 vigra_assert(shape == directions_blocks_begin.shape(),
"");
69 MultiCoordinateIterator<DataArray::actual_dimension> itBegin(shape);
70 MultiCoordinateIterator<DataArray::actual_dimension> end = itBegin.getEndIterator();
71 typedef typename MultiCoordinateIterator<DataArray::actual_dimension>::value_type Coordinate;
75 [&](
const int threadId,
const Coordinate iterVal){
77 DirectionsBlock directions_block = directions_blocks_begin[iterVal];
78 OverlappingBlock<DataArray> data_block = overlaps[iterVal];
80 typedef GridGraph<DataArray::actual_dimension, undirected_tag> Graph;
81 typedef typename Graph::NodeIt GraphScanner;
82 typedef typename Graph::OutArcIt NeighborIterator;
84 Graph graph(data_block.block.shape(), options.getNeighborhood());
85 for(GraphScanner node(graph); node != lemon::INVALID; ++node)
87 if(within(*node, data_block.inner_bounds))
89 typedef typename DataArray::value_type Data;
90 Data lowest_neighbor = data_block.block[*node];
92 typedef typename DirectionsBlock::value_type
Direction;
93 Direction lowest_neighbor_direction = std::numeric_limits<unsigned short>::max();
95 for(NeighborIterator arc(graph, *node); arc != lemon::INVALID; ++arc)
97 Shape neighbor_coordinates = graph.target(*arc);
98 Data neighbor_data = data_block.block[neighbor_coordinates];
99 if(neighbor_data < lowest_neighbor)
101 lowest_neighbor = neighbor_data;
102 lowest_neighbor_direction = arc.neighborIndex();
105 directions_block[*node - data_block.inner_bounds.first] = lowest_neighbor_direction;
112 template <
unsigned int N>
113 struct UnionFindWatershedsEquality
117 GridGraph<N, undirected_tag>* graph;
119 template <
class Shape>
120 bool operator()(
unsigned short u,
const unsigned short v,
const Shape& diff)
const
122 static const unsigned short plateau_id = std::numeric_limits<unsigned short>::max();
123 return (u == plateau_id && v == plateau_id) ||
124 (u != plateau_id && graph->neighborOffset(u) == diff) ||
125 (v != plateau_id && graph->neighborOffset(graph->oppositeIndex(v)) == diff);
194 template <
unsigned int N,
class Data,
class S1,
195 class Label,
class S2>
197 MultiArrayView<N, Label, S2> labels,
198 BlockwiseLabelOptions
const & options = BlockwiseLabelOptions())
200 using namespace blockwise_watersheds_detail;
202 typedef typename MultiArrayView<N, Data, S1>::difference_type Shape;
203 Shape shape = data.shape();
204 vigra_precondition(shape == labels.shape(),
"shapes of data and labels do not match");
206 MultiArray<N, unsigned short> directions(shape);
207 Shape block_shape = options.getBlockShapeN<N>();
209 MultiArray<N, MultiArrayView<N, unsigned short> > directions_blocks = blockify(directions, block_shape);
211 Overlaps<MultiArrayView<N, Data, S1> > overlaps(data, block_shape, Shape(1), Shape(1));
212 prepareBlockwiseWatersheds(overlaps, directions_blocks.begin(), options);
213 GridGraph<N, undirected_tag> graph(data.shape(), options.getNeighborhood());
214 UnionFindWatershedsEquality<N> equal = {&graph};
218 template <
unsigned int N,
class Data,
class Label>
220 ChunkedArray<N, Label>& labels,
221 BlockwiseLabelOptions
const & options,
222 ChunkedArray<N, unsigned short>& directions)
224 using namespace blockwise_watersheds_detail;
226 typedef typename ChunkedArray<N, Data>::shape_type Shape;
227 Shape shape = data.shape();
228 vigra_precondition(shape == labels.shape() && shape == directions.shape(),
229 "unionFindWatershedsBlockwise(): shapes of data and labels do not match");
230 Shape chunk_shape = data.chunkShape();
231 vigra_precondition(chunk_shape == labels.chunkShape() && chunk_shape == directions.chunkShape(),
232 "unionFindWatershedsBlockwise(): chunk shapes do not match");
234 Overlaps<ChunkedArray<N, Data> > overlaps(data, data.chunkShape(), Shape(1), Shape(1));
236 prepareBlockwiseWatersheds(overlaps, directions.chunk_begin(Shape(0), shape), options);
238 GridGraph<N, undirected_tag> graph(shape, options.getNeighborhood());
239 UnionFindWatershedsEquality<N> equal = {&graph};
243 template <
unsigned int N,
class Data,
247 ChunkedArray<N, Label>& labels,
248 BlockwiseLabelOptions
const & options = BlockwiseLabelOptions())
250 ChunkedArrayLazy<N, unsigned short> directions(data.shape(), data.chunkShape());
258 #endif // VIGRA_BLOCKWISE_WATERSHEDS_HXX
NeighborCode::Direction Direction
Definition: pixelneighborhood.hxx:321
unsigned int labelMultiArrayBlockwise(...)
Connected components labeling for MultiArrays and ChunkedArrays.
doxygen_overloaded_function(template<...> void separableConvolveBlockwise) template< unsigned int N
Separated convolution on ChunkedArrays.
void parallel_foreach(...)
Apply a functor to all items in a range in parallel.
unsigned int unionFindWatershedsBlockwise(...)
Blockwise union-find watersheds transform for MultiArrays and ChunkedArrays.