[ VIGRA Homepage | Function Index | Class Index | Namespaces | File List | Main Page ]

blockwise_watersheds.hxx |

1 /************************************************************************/
2 /* */
3 /* Copyright 2013-2014 by Martin Bidlingmaier and Ullrich Koethe */
4 /* */
5 /* This file is part of the VIGRA computer vision library. */
6 /* The VIGRA Website is */
7 /* http://hci.iwr.uni-heidelberg.de/vigra/ */
8 /* Please direct questions, bug reports, and contributions to */
9 /* ullrich.koethe@iwr.uni-heidelberg.de or */
10 /* vigra@informatik.uni-hamburg.de */
11 /* */
12 /* Permission is hereby granted, free of charge, to any person */
13 /* obtaining a copy of this software and associated documentation */
14 /* files (the "Software"), to deal in the Software without */
15 /* restriction, including without limitation the rights to use, */
16 /* copy, modify, merge, publish, distribute, sublicense, and/or */
17 /* sell copies of the Software, and to permit persons to whom the */
18 /* Software is furnished to do so, subject to the following */
19 /* conditions: */
20 /* */
21 /* The above copyright notice and this permission notice shall be */
22 /* included in all copies or substantial portions of the */
23 /* Software. */
24 /* */
25 /* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND */
26 /* EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES */
27 /* OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND */
28 /* NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT */
29 /* HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, */
30 /* WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING */
31 /* FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR */
32 /* OTHER DEALINGS IN THE SOFTWARE. */
33 /* */
34 /************************************************************************/
35
36 #ifndef VIGRA_BLOCKWISE_WATERSHEDS_HXX
37 #define VIGRA_BLOCKWISE_WATERSHEDS_HXX
38
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"
45
46 #include <limits>
47
48 namespace vigra
49 {
50
51 /** \addtogroup SeededRegionGrowing
52 */
53 //@{
54
55 namespace blockwise_watersheds_detail
56 {
57
58 template <class DataArray, class DirectionsBlocksIterator>
59 void prepareBlockwiseWatersheds(const Overlaps<DataArray>& overlaps,
60 DirectionsBlocksIterator directions_blocks_begin,
61 BlockwiseLabelOptions const & options)
62 {
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(), "");
68
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;
72
73 parallel_foreach(options.getNumThreads(),
74 itBegin,end,
75 [&](const int threadId, const Coordinate iterVal){
76
77 DirectionsBlock directions_block = directions_blocks_begin[iterVal];
78 OverlappingBlock<DataArray> data_block = overlaps[iterVal];
79
80 typedef GridGraph<DataArray::actual_dimension, undirected_tag> Graph;
81 typedef typename Graph::NodeIt GraphScanner;
82 typedef typename Graph::OutArcIt NeighborIterator;
83
84 Graph graph(data_block.block.shape(), options.getNeighborhood());
85 for(GraphScanner node(graph); node != lemon::INVALID; ++node)
86 {
87 if(within(*node, data_block.inner_bounds))
88 {
89 typedef typename DataArray::value_type Data;
90 Data lowest_neighbor = data_block.block[*node];
91
92 typedef typename DirectionsBlock::value_type Direction;
93 Direction lowest_neighbor_direction = std::numeric_limits<unsigned short>::max();
94
95 for(NeighborIterator arc(graph, *node); arc != lemon::INVALID; ++arc)
96 {
97 Shape neighbor_coordinates = graph.target(*arc);
98 Data neighbor_data = data_block.block[neighbor_coordinates];
99 if(neighbor_data < lowest_neighbor)
100 {
101 lowest_neighbor = neighbor_data;
102 lowest_neighbor_direction = arc.neighborIndex();
103 }
104 }
105 directions_block[*node - data_block.inner_bounds.first] = lowest_neighbor_direction;
106 }
107 }
108 }
109 );
110 }
111
112 template <unsigned int N>
113 struct UnionFindWatershedsEquality
114 {
115 // FIXME: this graph object shouldn't be necessary, most functions (and state) of graph are not used
116 // this probably needs some refactoring in GridGraph
117 GridGraph<N, undirected_tag>* graph;
118
119 template <class Shape>
120 bool operator()(unsigned short u, const unsigned short v, const Shape& diff) const
121 {
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);
126 }
127
128 struct WithDiffTag
129 {};
130 };
131
132 } // namespace blockwise_watersheds_detail
133
134 /*************************************************************/
135 /* */
136 /* unionFindWatershedsBlockwise */
137 /* */
138 /*************************************************************/
139
140 /** \brief Blockwise union-find watersheds transform for MultiArrays and ChunkedArrays.
141
142 <b> Declaration:</b>
143
144 \code
145 namespace vigra { namespace blockwise {
146 template <unsigned int N, class Data, class S1,
147 class Label, class S2>
148 Label
149 unionFindWatershedsBlockwise(MultiArrayView<N, Data, S1> data,
150 MultiArrayView<N, Label, S2> labels,
151 BlockwiseLabelOptions const & options);
152
153 template <unsigned int N, class Data, class Label>
154 Label
155 unionFindWatershedsBlockwise(const ChunkedArray<N, Data>& data,
156 ChunkedArray<N, Label>& labels,
157 BlockwiseLabelOptions const & options = BlockwiseLabelOptions());
158
159 // provide temporary directions storage
160 template <unsigned int N, class Data, class Label>
161 Label
162 unionFindWatershedsBlockwise(const ChunkedArray<N, Data>& data,
163 ChunkedArray<N, Label>& labels,
164 BlockwiseLabelOptions const & options,
165 ChunkedArray<N, unsigned short>& temporary_storage);
166 }}
167 \endcode
168
169 The resulting labeling is equivalent to a labeling by \ref watershedsUnionFind, that is,
170 the components are the same but may have different ids.
171 If \a temporary_storage is provided, this array is used for intermediate result storage.
172 Otherwise, a newly created \ref vigra::ChunkedArrayLazy is used.
173
174 Return: the number of labels assigned (=largest label, because labels start at one)
175
176 <b> Usage: </b>
177
178 <b>\#include </b> <vigra/blockwise_watersheds.hxx><br>
179 Namespace: vigra
180
181 \code
182 Shape3 shape = Shape3(10);
183 Shape3 chunk_shape = Shape3(4);
184 ChunkedArrayLazy<3, int> data(shape, chunk_shape);
185 // fill data ...
186
187 ChunkedArrayLazy<3, size_t> labels(shape, chunk_shape);
188
189 unionFindWatershedsBlockwise(data, labels, IndirectNeighborhood);
190 \endcode
191 */
192 doxygen_overloaded_function(template <...> unsigned int unionFindWatershedsBlockwise)
193
194 template <unsigned int N, class Data, class S1,
195 class Label, class S2>
196 Label unionFindWatershedsBlockwise(MultiArrayView<N, Data, S1> data,
197 MultiArrayView<N, Label, S2> labels,
198 BlockwiseLabelOptions const & options = BlockwiseLabelOptions())
199 {
200 using namespace blockwise_watersheds_detail;
201
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");
205
206 MultiArray<N, unsigned short> directions(shape);
207 Shape block_shape = options.getBlockShapeN<N>();
208
209 MultiArray<N, MultiArrayView<N, unsigned short> > directions_blocks = blockify(directions, block_shape);
210
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};
215 return labelMultiArrayBlockwise(directions, labels, options, equal);
216 }
217
218 template <unsigned int N, class Data, class Label>
219 Label unionFindWatershedsBlockwise(const ChunkedArray<N, Data>& data,
220 ChunkedArray<N, Label>& labels,
221 BlockwiseLabelOptions const & options,
222 ChunkedArray<N, unsigned short>& directions)
223 {
224 using namespace blockwise_watersheds_detail;
225
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");
233
234 Overlaps<ChunkedArray<N, Data> > overlaps(data, data.chunkShape(), Shape(1), Shape(1));
235
236 prepareBlockwiseWatersheds(overlaps, directions.chunk_begin(Shape(0), shape), options);
237
238 GridGraph<N, undirected_tag> graph(shape, options.getNeighborhood());
239 UnionFindWatershedsEquality<N> equal = {&graph};
240 return labelMultiArrayBlockwise(directions, labels, options, equal);
241 }
242
243 template <unsigned int N, class Data,
244 class Label>
245 inline Label
246 unionFindWatershedsBlockwise(const ChunkedArray<N, Data>& data,
247 ChunkedArray<N, Label>& labels,
248 BlockwiseLabelOptions const & options = BlockwiseLabelOptions())
249 {
250 ChunkedArrayLazy<N, unsigned short> directions(data.shape(), data.chunkShape());
251 return unionFindWatershedsBlockwise(data, labels, options, directions);
252 }
253
254 //@}
255
256 } // namespace vigra
257
258 #endif // VIGRA_BLOCKWISE_WATERSHEDS_HXX
vigra::MultiArrayView::difference_type
MultiArrayShape< actual_dimension >::type difference_type
Definition: multi_array.hxx:687
vigra
Definition: accessor.hxx:43
vigra::FourNeighborhood::Direction
NeighborCode::Direction Direction
Definition: pixelneighborhood.hxx:321
vigra::doxygen_overloaded_function
doxygen_overloaded_function(template<... > void separableConvolveBlockwise) template< unsigned int N
Separated convolution on ChunkedArrays.
vigra::labelMultiArrayBlockwise
unsigned int labelMultiArrayBlockwise(...)
Connected components labeling for MultiArrays and ChunkedArrays.
vigra::MultiArrayShape::type
TinyVector< MultiArrayIndex, N > type
Definition: multi_shape.hxx:250
vigra::parallel_foreach
void parallel_foreach(...)
Apply a functor to all items in a range in parallel.
vigra::unionFindWatershedsBlockwise
unsigned int unionFindWatershedsBlockwise(...)
Blockwise union-find watersheds transform for MultiArrays and ChunkedArrays.

© Ullrich Köthe (ullrich.koethe@iwr.uni-heidelberg.de) Heidelberg Collaboratory for Image Processing, University of Heidelberg, Germany | html generated using doxygen and Python vigra 1.11.0