1 /*M///////////////////////////////////////////////////////////////////////////////////////
2 //
3 // IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
4 //
5 // By downloading, copying, installing or using the software you agree to this license.
6 // If you do not agree to this license, do not download, install,
7 // copy or use the software.
8 //
9 //
10 // Intel License Agreement
11 // For Open Source Computer Vision Library
12 //
13 // Copyright (C) 2000, Intel Corporation, all rights reserved.
14 // Third party copyrights are property of their respective owners.
15 //
16 // Redistribution and use in source and binary forms, with or without modification,
17 // are permitted provided that the following conditions are met:
18 //
19 // * Redistribution's of source code must retain the above copyright notice,
20 // this list of conditions and the following disclaimer.
21 //
22 // * Redistribution's in binary form must reproduce the above copyright notice,
23 // this list of conditions and the following disclaimer in the documentation
24 // and/or other materials provided with the distribution.
25 //
26 // * The name of Intel Corporation may not be used to endorse or promote products
27 // derived from this software without specific prior written permission.
28 //
29 // This software is provided by the copyright holders and contributors "as is" and
30 // any express or implied warranties, including, but not limited to, the implied
31 // warranties of merchantability and fitness for a particular purpose are disclaimed.
32 // In no event shall the Intel Corporation or contributors be liable for any direct,
33 // indirect, incidental, special, exemplary, or consequential damages
34 // (including, but not limited to, procurement of substitute goods or services;
35 // loss of use, data, or profits; or business interruption) however caused
36 // and on any theory of liability, whether in contract, strict liability,
37 // or tort (including negligence or otherwise) arising in any way out of
38 // the use of this software, even if advised of the possibility of such damage.
39 //
40 //M*/
41
42 #ifndef _CV_GCGRAPH_H_
43 #define _CV_GCGRAPH_H_
44
45 template <class TWeight> class GCGraph
46 {
47 public:
48 GCGraph();
49 GCGraph( unsigned int vtxCount, unsigned int edgeCount );
50 ~GCGraph();
51 void create( unsigned int vtxCount, unsigned int edgeCount );
52 int addVtx();
53 void addEdges( int i, int j, TWeight w, TWeight revw );
54 void addTermWeights( int i, TWeight sourceW, TWeight sinkW );
55 TWeight maxFlow();
56 bool inSourceSegment( int i );
57 private:
58 class Vtx
59 {
60 public:
61 Vtx *next; // initialized and used in maxFlow() only
62 int parent;
63 int first;
64 int ts;
65 int dist;
66 TWeight weight;
67 uchar t;
68 };
69 class Edge
70 {
71 public:
72 int dst;
73 int next;
74 TWeight weight;
75 };
76
77 std::vector<Vtx> vtcs;
78 std::vector<Edge> edges;
79 TWeight flow;
80 };
81
82 template <class TWeight>
GCGraph()83 GCGraph<TWeight>::GCGraph()
84 {
85 flow = 0;
86 }
87 template <class TWeight>
GCGraph(unsigned int vtxCount,unsigned int edgeCount)88 GCGraph<TWeight>::GCGraph( unsigned int vtxCount, unsigned int edgeCount )
89 {
90 create( vtxCount, edgeCount );
91 }
92 template <class TWeight>
~GCGraph()93 GCGraph<TWeight>::~GCGraph()
94 {
95 }
96 template <class TWeight>
create(unsigned int vtxCount,unsigned int edgeCount)97 void GCGraph<TWeight>::create( unsigned int vtxCount, unsigned int edgeCount )
98 {
99 vtcs.reserve( vtxCount );
100 edges.reserve( edgeCount + 2 );
101 flow = 0;
102 }
103
104 template <class TWeight>
addVtx()105 int GCGraph<TWeight>::addVtx()
106 {
107 Vtx v;
108 memset( &v, 0, sizeof(Vtx));
109 vtcs.push_back(v);
110 return (int)vtcs.size() - 1;
111 }
112
113 template <class TWeight>
addEdges(int i,int j,TWeight w,TWeight revw)114 void GCGraph<TWeight>::addEdges( int i, int j, TWeight w, TWeight revw )
115 {
116 CV_Assert( i>=0 && i<(int)vtcs.size() );
117 CV_Assert( j>=0 && j<(int)vtcs.size() );
118 CV_Assert( w>=0 && revw>=0 );
119 CV_Assert( i != j );
120
121 if( !edges.size() )
122 edges.resize( 2 );
123
124 Edge fromI, toI;
125 fromI.dst = j;
126 fromI.next = vtcs[i].first;
127 fromI.weight = w;
128 vtcs[i].first = (int)edges.size();
129 edges.push_back( fromI );
130
131 toI.dst = i;
132 toI.next = vtcs[j].first;
133 toI.weight = revw;
134 vtcs[j].first = (int)edges.size();
135 edges.push_back( toI );
136 }
137
138 template <class TWeight>
addTermWeights(int i,TWeight sourceW,TWeight sinkW)139 void GCGraph<TWeight>::addTermWeights( int i, TWeight sourceW, TWeight sinkW )
140 {
141 CV_Assert( i>=0 && i<(int)vtcs.size() );
142
143 TWeight dw = vtcs[i].weight;
144 if( dw > 0 )
145 sourceW += dw;
146 else
147 sinkW -= dw;
148 flow += (sourceW < sinkW) ? sourceW : sinkW;
149 vtcs[i].weight = sourceW - sinkW;
150 }
151
152 template <class TWeight>
maxFlow()153 TWeight GCGraph<TWeight>::maxFlow()
154 {
155 const int TERMINAL = -1, ORPHAN = -2;
156 Vtx stub, *nilNode = &stub, *first = nilNode, *last = nilNode;
157 int curr_ts = 0;
158 stub.next = nilNode;
159 Vtx *vtxPtr = &vtcs[0];
160 Edge *edgePtr = &edges[0];
161
162 std::vector<Vtx*> orphans;
163
164 // initialize the active queue and the graph vertices
165 for( int i = 0; i < (int)vtcs.size(); i++ )
166 {
167 Vtx* v = vtxPtr + i;
168 v->ts = 0;
169 if( v->weight != 0 )
170 {
171 last = last->next = v;
172 v->dist = 1;
173 v->parent = TERMINAL;
174 v->t = v->weight < 0;
175 }
176 else
177 v->parent = 0;
178 }
179 first = first->next;
180 last->next = nilNode;
181 nilNode->next = 0;
182
183 // run the search-path -> augment-graph -> restore-trees loop
184 for(;;)
185 {
186 Vtx* v, *u;
187 int e0 = -1, ei = 0, ej = 0;
188 TWeight minWeight, weight;
189 uchar vt;
190
191 // grow S & T search trees, find an edge connecting them
192 while( first != nilNode )
193 {
194 v = first;
195 if( v->parent )
196 {
197 vt = v->t;
198 for( ei = v->first; ei != 0; ei = edgePtr[ei].next )
199 {
200 if( edgePtr[ei^vt].weight == 0 )
201 continue;
202 u = vtxPtr+edgePtr[ei].dst;
203 if( !u->parent )
204 {
205 u->t = vt;
206 u->parent = ei ^ 1;
207 u->ts = v->ts;
208 u->dist = v->dist + 1;
209 if( !u->next )
210 {
211 u->next = nilNode;
212 last = last->next = u;
213 }
214 continue;
215 }
216
217 if( u->t != vt )
218 {
219 e0 = ei ^ vt;
220 break;
221 }
222
223 if( u->dist > v->dist+1 && u->ts <= v->ts )
224 {
225 // reassign the parent
226 u->parent = ei ^ 1;
227 u->ts = v->ts;
228 u->dist = v->dist + 1;
229 }
230 }
231 if( e0 > 0 )
232 break;
233 }
234 // exclude the vertex from the active list
235 first = first->next;
236 v->next = 0;
237 }
238
239 if( e0 <= 0 )
240 break;
241
242 // find the minimum edge weight along the path
243 minWeight = edgePtr[e0].weight;
244 CV_Assert( minWeight > 0 );
245 // k = 1: source tree, k = 0: destination tree
246 for( int k = 1; k >= 0; k-- )
247 {
248 for( v = vtxPtr+edgePtr[e0^k].dst;; v = vtxPtr+edgePtr[ei].dst )
249 {
250 if( (ei = v->parent) < 0 )
251 break;
252 weight = edgePtr[ei^k].weight;
253 minWeight = MIN(minWeight, weight);
254 CV_Assert( minWeight > 0 );
255 }
256 weight = fabs(v->weight);
257 minWeight = MIN(minWeight, weight);
258 CV_Assert( minWeight > 0 );
259 }
260
261 // modify weights of the edges along the path and collect orphans
262 edgePtr[e0].weight -= minWeight;
263 edgePtr[e0^1].weight += minWeight;
264 flow += minWeight;
265
266 // k = 1: source tree, k = 0: destination tree
267 for( int k = 1; k >= 0; k-- )
268 {
269 for( v = vtxPtr+edgePtr[e0^k].dst;; v = vtxPtr+edgePtr[ei].dst )
270 {
271 if( (ei = v->parent) < 0 )
272 break;
273 edgePtr[ei^(k^1)].weight += minWeight;
274 if( (edgePtr[ei^k].weight -= minWeight) == 0 )
275 {
276 orphans.push_back(v);
277 v->parent = ORPHAN;
278 }
279 }
280
281 v->weight = v->weight + minWeight*(1-k*2);
282 if( v->weight == 0 )
283 {
284 orphans.push_back(v);
285 v->parent = ORPHAN;
286 }
287 }
288
289 // restore the search trees by finding new parents for the orphans
290 curr_ts++;
291 while( !orphans.empty() )
292 {
293 Vtx* v2 = orphans.back();
294 orphans.pop_back();
295
296 int d, minDist = INT_MAX;
297 e0 = 0;
298 vt = v2->t;
299
300 for( ei = v2->first; ei != 0; ei = edgePtr[ei].next )
301 {
302 if( edgePtr[ei^(vt^1)].weight == 0 )
303 continue;
304 u = vtxPtr+edgePtr[ei].dst;
305 if( u->t != vt || u->parent == 0 )
306 continue;
307 // compute the distance to the tree root
308 for( d = 0;; )
309 {
310 if( u->ts == curr_ts )
311 {
312 d += u->dist;
313 break;
314 }
315 ej = u->parent;
316 d++;
317 if( ej < 0 )
318 {
319 if( ej == ORPHAN )
320 d = INT_MAX-1;
321 else
322 {
323 u->ts = curr_ts;
324 u->dist = 1;
325 }
326 break;
327 }
328 u = vtxPtr+edgePtr[ej].dst;
329 }
330
331 // update the distance
332 if( ++d < INT_MAX )
333 {
334 if( d < minDist )
335 {
336 minDist = d;
337 e0 = ei;
338 }
339 for( u = vtxPtr+edgePtr[ei].dst; u->ts != curr_ts; u = vtxPtr+edgePtr[u->parent].dst )
340 {
341 u->ts = curr_ts;
342 u->dist = --d;
343 }
344 }
345 }
346
347 if( (v2->parent = e0) > 0 )
348 {
349 v2->ts = curr_ts;
350 v2->dist = minDist;
351 continue;
352 }
353
354 /* no parent is found */
355 v2->ts = 0;
356 for( ei = v2->first; ei != 0; ei = edgePtr[ei].next )
357 {
358 u = vtxPtr+edgePtr[ei].dst;
359 ej = u->parent;
360 if( u->t != vt || !ej )
361 continue;
362 if( edgePtr[ei^(vt^1)].weight && !u->next )
363 {
364 u->next = nilNode;
365 last = last->next = u;
366 }
367 if( ej > 0 && vtxPtr+edgePtr[ej].dst == v2 )
368 {
369 orphans.push_back(u);
370 u->parent = ORPHAN;
371 }
372 }
373 }
374 }
375 return flow;
376 }
377
378 template <class TWeight>
inSourceSegment(int i)379 bool GCGraph<TWeight>::inSourceSegment( int i )
380 {
381 CV_Assert( i>=0 && i<(int)vtcs.size() );
382 return vtcs[i].t == 0;
383 }
384
385 #endif
386