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

tensorutilities.hxx
1/************************************************************************/
2/* */
3/* Copyright 2002-2004 by 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_TENSORUTILITIES_HXX
37#define VIGRA_TENSORUTILITIES_HXX
38
39#include <cmath>
40#include "utilities.hxx"
41#include "mathutil.hxx"
42#include "multi_shape.hxx"
43
44namespace vigra {
45
46/** \addtogroup TensorImaging Tensor Image Processing
47*/
48//@{
49
50/********************************************************/
51/* */
52/* vectorToTensor */
53/* */
54/********************************************************/
55
56/** \brief Calculate the tensor (outer) product of a 2D vector with itself.
57
58 This function is useful to transform vector images into a tensor representation
59 that can be used as input to tensor based processing and analysis functions
60 (e.g. tensor smoothing). The input pixel type must be vectors of length 2, whereas
61 the output must contain vectors of length 3 which will represent the tensor components
62 in the order t11, t12 (== t21 due to symmetry), t22.
63
64 <b>Note:</b> In order to account for the left-handedness of the image coordinate system,
65 the second tensor component (t12) can be negated by setting <tt>negateComponent2 = false</tt>.
66 Angles will then be interpreted counter-clockwise rather than clockwise. By default,
67 this behavior is switched off.
68
69 <b> Declarations:</b>
70
71 pass 2D array views:
72 \code
73 namespace vigra {
74 template <class T1, class S1,
75 class T2, class S2>
76 void
77 vectorToTensor(MultiArrayView<2, T1, S1> const & src,
78 MultiArrayView<2, T2, S2> dest,
79 bool negateComponent2 = false);
80 }
81 \endcode
82
83 \deprecatedAPI{vectorToTensor}
84 pass \ref ImageIterators and \ref DataAccessors :
85 \code
86 namespace vigra {
87 template <class SrcIterator, class SrcAccessor,
88 class DestIterator, class DestAccessor>
89 void vectorToTensor(SrcIterator sul, SrcIterator slr, SrcAccessor src,
90 DestIterator dul, DestAccessor dest,
91 bool negateComponent2 = false);
92 }
93 \endcode
94 use argument objects in conjunction with \ref ArgumentObjectFactories :
95 \code
96 namespace vigra {
97 template <class SrcIterator, class SrcAccessor,
98 class DestIterator, class DestAccessor>
99 void vectorToTensor(triple<SrcIterator, SrcIterator, SrcAccessor> s,
100 pair<DestIterator, DestAccessor> d,
101 bool negateComponent2 = false);
102 }
103 \endcode
104 \deprecatedEnd
105
106 <b> Usage:</b>
107
108 <b>\#include</b> <vigra/tensorutilities.hxx><br/>
109 Namespace: vigra
110
111 \code
112 MultiArray<2, float> img(w,h);
113 MultiArray<2, TinyVector<float, 2> > gradient(w,h);
114 MultiArray<2, TinyVector<float, 3> > tensor(w,h);
115 ...
116
117 gaussianGradient(img, gradient, 2.0);
118 vectorToTensor(gradient, tensor);
119 \endcode
120
121 \deprecatedUsage{vectorToTensor}
122 \code
123 FImage img(w,h);
124 FVector2Image gradient(w,h);
125 FVector3Image tensor(w,h);
126
127 gaussianGradient(srcImageRange(img), destImage(gradient), 2.0);
128 vectorToTensor(srcImageRange(gradient), destImage(tensor));
129 \endcode
130 \deprecatedEnd
131*/
132doxygen_overloaded_function(template <...> void vectorToTensor)
133
134template <class SrcIterator, class SrcAccessor,
135 class DestIterator, class DestAccessor>
138 bool negateComponent2)
139{
140 vigra_precondition(src.size(sul) == 2,
141 "vectorToTensor(): input image must have 2 bands.");
142 vigra_precondition(dest.size(dul) == 3,
143 "vectorToTensor(): output image must have 3 bands.");
144
145 int w = slr.x - sul.x;
146 int h = slr.y - sul.y;
147
148 for(int y=0; y<h; ++y, ++sul.y, ++dul.y)
149 {
150 typename SrcIterator::row_iterator s = sul.rowIterator();
151 typename SrcIterator::row_iterator send = s + w;
152 typename DestIterator::row_iterator d = dul.rowIterator();
154 {
155 for(; s < send; ++s, ++d)
156 {
157 dest.setComponent(sq(src.getComponent(s, 0)), d, 0);
158 dest.setComponent(-src.getComponent(s, 0)*src.getComponent(s, 1), d, 1);
159 // ^ negative sign to turn left-handed into right-handed coordinates
160 dest.setComponent(sq(src.getComponent(s, 1)), d, 2);
161 }
162 }
163 else
164 {
165 for(; s < send; ++s, ++d)
166 {
167 dest.setComponent(sq(src.getComponent(s, 0)), d, 0);
168 dest.setComponent(src.getComponent(s, 0)*src.getComponent(s, 1), d, 1);
169 dest.setComponent(sq(src.getComponent(s, 1)), d, 2);
170 }
171 }
172 }
173}
174
175template <class SrcIterator, class SrcAccessor,
176 class DestIterator, class DestAccessor>
177inline
178void vectorToTensor(SrcIterator sul, SrcIterator slr, SrcAccessor src,
179 DestIterator dul, DestAccessor dest)
180{
181 vectorToTensor(sul, slr, src, dul, dest, false);
182}
183
184template <class SrcIterator, class SrcAccessor,
185 class DestIterator, class DestAccessor>
186inline void
187vectorToTensor(triple<SrcIterator, SrcIterator, SrcAccessor> s,
188 pair<DestIterator, DestAccessor> d,
189 bool negateComponent2)
190{
191 vectorToTensor(s.first, s.second, s.third, d.first, d.second, negateComponent2);
192}
193
194template <class SrcIterator, class SrcAccessor,
195 class DestIterator, class DestAccessor>
196inline void
197vectorToTensor(triple<SrcIterator, SrcIterator, SrcAccessor> s,
198 pair<DestIterator, DestAccessor> d)
199{
200 vectorToTensor(s.first, s.second, s.third, d.first, d.second, false);
201}
202
203template <class T1, class S1,
204 class T2, class S2>
205inline void
206vectorToTensor(MultiArrayView<2, T1, S1> const & src,
207 MultiArrayView<2, T2, S2> dest,
208 bool negateComponent2 = false)
209{
210 vigra_precondition(src.shape() == dest.shape(),
211 "vectorToTensor(): shape mismatch between input and output.");
212 vectorToTensor(srcImageRange(src), destImage(dest), negateComponent2);
213}
214
215/********************************************************/
216/* */
217/* tensorEigenRepresentation */
218/* */
219/********************************************************/
220
221/** \brief Calculate eigen representation of a symmetric 2x2 tensor.
222
223 This function turns a 3-band image representing the tensor components
224 t11, t12 (== t21 due to symmetry), t22 into the a 3-band image holding the eigen
225 representation e1, e2, and angle, where e1 > e2. When the tensor is
226 defined in a left-handed coordinate system (the default on images), the angle will
227 then be given in clockwise orientation, starting at the x-axis. Otherwise, it
228 will be given in counter-clockwise orientation.
229
230 <b> Declarations:</b>
231
232 pass 2D array views:
233 \code
234 namespace vigra {
235 template <class T1, class S1,
236 class T2, class S2>
237 void
238 tensorEigenRepresentation(MultiArrayView<2, T1, S1> const & src,
239 MultiArrayView<2, T2, S2> dest);
240 }
241 \endcode
242
243 \deprecatedAPI{tensorEigenRepresentation}
244 pass \ref ImageIterators and \ref DataAccessors :
245 \code
246 namespace vigra {
247 template <class SrcIterator, class SrcAccessor,
248 class DestIterator, class DestAccessor>
249 void tensorEigenRepresentation(SrcIterator sul, SrcIterator slr, SrcAccessor src,
250 DestIterator dul, DestAccessor dest);
251 }
252 \endcode
253 use argument objects in conjunction with \ref ArgumentObjectFactories :
254 \code
255 namespace vigra {
256 template <class SrcIterator, class SrcAccessor,
257 class DestIterator, class DestAccessor>
258 void tensorEigenRepresentation(triple<SrcIterator, SrcIterator, SrcAccessor> s,
259 pair<DestIterator, DestAccessor> d);
260 }
261 \endcode
262 \deprecatedEnd
263
264 <b> Usage:</b>
265
266 <b>\#include</b> <vigra/tensorutilities.hxx><br/>
267 Namespace: vigra
268
269 \code
270 MultiArray<2, TinyVector<float, 3> > tensor(w,h),
271 eigen(w,h);
272
273 tensorEigenRepresentation(tensor, eigen);
274 \endcode
275
276 \deprecatedUsage{tensorEigenRepresentation}
277 \code
278 FVector3Image tensor(w,h);
279 FVector3Image eigen(w,h);
280
281 tensorEigenRepresentation(srcImageRange(tensor), destImage(eigen));
282 \endcode
283 \deprecatedEnd
284*/
285doxygen_overloaded_function(template <...> void tensorEigenRepresentation)
286
287template <class SrcIterator, class SrcAccessor,
288 class DestIterator, class DestAccessor>
291{
292 vigra_precondition(src.size(sul) == 3,
293 "tensorEigenRepresentation(): input image must have 3 bands.");
294 vigra_precondition(dest.size(dul) == 3,
295 "tensorEigenRepresentation(): output image must have 3 bands.");
296
297 int w = slr.x - sul.x;
298 int h = slr.y - sul.y;
299
300 for(int y=0; y<h; ++y, ++sul.y, ++dul.y)
301 {
302 typename SrcIterator::row_iterator s = sul.rowIterator();
303 typename SrcIterator::row_iterator send = s + w;
304 typename DestIterator::row_iterator d = dul.rowIterator();
305 for(; s < send; ++s, ++d)
306 {
307 typedef typename
308 NumericTraits<typename SrcAccessor::component_type>::RealPromote TmpType;
309 TmpType d1 = src.getComponent(s,0) + src.getComponent(s,2);
310 TmpType d2 = src.getComponent(s,0) - src.getComponent(s,2);
311 TmpType d3 = TmpType(2.0) * src.getComponent(s,1);
313
314 dest.setComponent(0.5 * (d1 + d4), d, 0); // large EV
315 dest.setComponent(0.5 * (d1 - d4), d, 1); // small EV
316 if(d2==0.0 && d3==0.0)
317 {
318 dest.setComponent(0, d, 2); // orientation
319 }
320 else
321 {
322 dest.setComponent(0.5 * VIGRA_CSTD::atan2(d3, d2), d, 2); // orientation
323 }
324 }
325 }
326}
327
328template <class SrcIterator, class SrcAccessor,
329 class DestIterator, class DestAccessor>
330inline void
331tensorEigenRepresentation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
332 pair<DestIterator, DestAccessor> dest)
333{
334 tensorEigenRepresentation(src.first, src.second, src.third, dest.first, dest.second);
335}
336
337template <class T1, class S1,
338 class T2, class S2>
339inline void
340tensorEigenRepresentation(MultiArrayView<2, T1, S1> const & src,
341 MultiArrayView<2, T2, S2> dest)
342{
343 vigra_precondition(src.shape() == dest.shape(),
344 "tensorEigenRepresentation(): shape mismatch between input and output.");
345 tensorEigenRepresentation(srcImageRange(src), destImage(dest));
346}
347
348/********************************************************/
349/* */
350/* tensorTrace */
351/* */
352/********************************************************/
353
354/** \brief Calculate the trace of a 2x2 tensor.
355
356 This function turns a 3-band image representing the tensor components
357 t11, t12 (== t21 due to symmetry), t22 into the a 1-band image holding the
358 tensor trace t11 + t22.
359
360 <b> Declarations:</b>
361
362 pass 2D array views:
363 \code
364 namespace vigra {
365 template <class T1, class S1,
366 class T2, class S2>
367 void
368 tensorTrace(MultiArrayView<2, T1, S1> const & src,
369 MultiArrayView<2, T2, S2> dest);
370 }
371 \endcode
372
373 \deprecatedAPI{tensorTrace}
374 pass \ref ImageIterators and \ref DataAccessors :
375 \code
376 namespace vigra {
377 template <class SrcIterator, class SrcAccessor,
378 class DestIterator, class DestAccessor>
379 void tensorTrace(SrcIterator sul, SrcIterator slr, SrcAccessor src,
380 DestIterator dul, DestAccessor dest);
381 }
382 \endcode
383 use argument objects in conjunction with \ref ArgumentObjectFactories :
384 \code
385 namespace vigra {
386 template <class SrcIterator, class SrcAccessor,
387 class DestIterator, class DestAccessor>
388 void tensorTrace(triple<SrcIterator, SrcIterator, SrcAccessor> s,
389 pair<DestIterator, DestAccessor> d);
390 }
391 \endcode
392 \deprecatedEnd
393
394 <b> Usage:</b>
395
396 <b>\#include</b> <vigra/tensorutilities.hxx><br/>
397 Namespace: vigra
398
399 \code
400 MultiArray<2, TinyVector<float, 3> > tensor(w,h);
401 MultiArray<2, float> trace(w,h);
402
403 tensorTrace(tensor, trace);
404 \endcode
405
406 \deprecatedUsage{tensorTrace}
407 \code
408 FVector3Image tensor(w,h);
409 FImage trace(w,h);
410
411 tensorTrace(srcImageRange(tensor), destImage(trace));
412 \endcode
413 \deprecatedEnd
414*/
415doxygen_overloaded_function(template <...> void tensorTrace)
416
417template <class SrcIterator, class SrcAccessor,
418 class DestIterator, class DestAccessor>
421{
422 vigra_precondition(src.size(sul) == 3,
423 "tensorTrace(): input image must have 3 bands.");
424
425 int w = slr.x - sul.x;
426 int h = slr.y - sul.y;
427
428 for(int y=0; y<h; ++y, ++sul.y, ++dul.y)
429 {
430 typename SrcIterator::row_iterator s = sul.rowIterator();
431 typename SrcIterator::row_iterator send = s + w;
432 typename DestIterator::row_iterator d = dul.rowIterator();
433 for(; s < send; ++s, ++d)
434 {
435 dest.set(src.getComponent(s,0) + src.getComponent(s,2), d);
436 }
437 }
438}
439
440template <class SrcIterator, class SrcAccessor,
441 class DestIterator, class DestAccessor>
442inline void
443tensorTrace(triple<SrcIterator, SrcIterator, SrcAccessor> src,
444 pair<DestIterator, DestAccessor> dest)
445{
446 tensorTrace(src.first, src.second, src.third, dest.first, dest.second);
447}
448
449template <class T1, class S1,
450 class T2, class S2>
451inline void
452tensorTrace(MultiArrayView<2, T1, S1> const & src,
453 MultiArrayView<2, T2, S2> dest)
454{
455 vigra_precondition(src.shape() == dest.shape(),
456 "tensorTrace(): shape mismatch between input and output.");
457 tensorTrace(srcImageRange(src), destImage(dest));
458}
459
460/********************************************************/
461/* */
462/* tensorToEdgeCorner */
463/* */
464/********************************************************/
465
466/** \brief Decompose a symmetric 2x2 tensor into its edge and corner parts.
467
468 This function turns a 3-band image representing the tensor components
469 t11, t12 (== t21 due to symmetry), t22 into the a 2-band image holding
470 the tensor's edgeness (difference of the tensor's
471 eigenvalues) and orientation, and a 1-band image representing its corner part
472 (equal to the twice the small eigen value). The original tensor must be
473 positive definite and defined in a right-handed coordinate system (e.g.
474 the tensor resulting from \ref boundaryTensor()).
475
476 <b> Declarations:</b>
477
478 pass 2D array views:
479 \code
480 namespace vigra {
481 template <class T1, class S1,
482 class T21, class S21,
483 class T22, class S22>
484 void
485 tensorToEdgeCorner(MultiArrayView<2, T1, S1> const & src,
486 MultiArrayView<2, T21, S21> edge,
487 MultiArrayView<2, T22, S22> corner);
488 }
489 \endcode
490
491 \deprecatedAPI{tensorToEdgeCorner}
492 pass \ref ImageIterators and \ref DataAccessors :
493 \code
494 namespace vigra {
495 template <class SrcIterator, class SrcAccessor,
496 class DestIterator1, class DestAccessor1,
497 class DestIterator2, class DestAccessor2>
498 void tensorToEdgeCorner(SrcIterator sul, SrcIterator slr, SrcAccessor src,
499 DestIterator1 edgeul, DestAccessor1 edge,
500 DestIterator2 cornerul, DestAccessor2 corner);
501 }
502 \endcode
503 use argument objects in conjunction with \ref ArgumentObjectFactories :
504 \code
505 namespace vigra {
506 template <class SrcIterator, class SrcAccessor,
507 class DestIterator1, class DestAccessor1,
508 class DestIterator2, class DestAccessor2>
509 void tensorToEdgeCorner(triple<SrcIterator, SrcIterator, SrcAccessor> s,
510 pair<DestIterator1, DestAccessor1> edge,
511 pair<DestIterator2, DestAccessor2> corner);
512 }
513 \endcode
514 \deprecatedEnd
515
516 <b> Usage:</b>
517
518 <b>\#include</b> <vigra/tensorutilities.hxx><br/>
519 Namespace: vigra
520
521 \code
522 MultiArray<2, TinyVector<float, 3> > tensor(w,h);
523 MultiArray<2, TinyVector<float, 2> > edgePart(w,h);
524 MultiArray<2, float> cornerPart(w,h);
525 ...
526
527 tensorTrace(tensor, edgePart, cornerPart);
528 \endcode
529
530 \deprecatedUsage{tensorToEdgeCorner}
531 \code
532 FVector3Image tensor(w,h);
533 FVector2Image edgePart(w,h);
534 FImage cornerPart(w,h);
535
536 tensorTrace(srcImageRange(tensor), destImage(edgePart), destImage(cornerPart));
537 \endcode
538 \deprecatedEnd
539*/
540doxygen_overloaded_function(template <...> void tensorToEdgeCorner)
541
542template <class SrcIterator, class SrcAccessor,
543 class DestIterator1, class DestAccessor1,
544 class DestIterator2, class DestAccessor2>
548{
549 vigra_precondition(src.size(sul) == 3,
550 "tensorToEdgeCorner(): input image must have 3 bands.");
551 vigra_precondition(edge.size(edgeul) == 2,
552 "tensorToEdgeCorner(): edge image must have 2 bands.");
553
554 int w = slr.x - sul.x;
555 int h = slr.y - sul.y;
556
557 for(int y=0; y<h; ++y, ++sul.y, ++edgeul.y, ++cornerul.y)
558 {
559 typename SrcIterator::row_iterator s = sul.rowIterator();
560 typename SrcIterator::row_iterator send = s + w;
561 typename DestIterator1::row_iterator e = edgeul.rowIterator();
562 typename DestIterator2::row_iterator c = cornerul.rowIterator();
563 for(; s < send; ++s, ++e, ++c)
564 {
565 typedef typename
566 NumericTraits<typename SrcAccessor::component_type>::RealPromote TmpType;
567 TmpType d1 = src.getComponent(s,0) + src.getComponent(s,2);
568 TmpType d2 = src.getComponent(s,0) - src.getComponent(s,2);
569 TmpType d3 = 2.0 * src.getComponent(s,1);
571
572 edge.setComponent(d4, e, 0); // edgeness = difference of EVs
573 if(d2 == 0.0 && d3 == 0.0)
574 {
575 edge.setComponent(0.0, e, 1); // orientation
576 }
577 else
578 {
579 edge.setComponent(0.5 * VIGRA_CSTD::atan2(d3, d2), e, 1); // orientation
580 }
581 corner.set(d1 - d4, c); // cornerness = 2 * small EV
582 }
583 }
584}
585
586template <class SrcIterator, class SrcAccessor,
587 class DestIterator1, class DestAccessor1,
588 class DestIterator2, class DestAccessor2>
589inline void
590tensorToEdgeCorner(triple<SrcIterator, SrcIterator, SrcAccessor> src,
591 pair<DestIterator1, DestAccessor1> edge,
592 pair<DestIterator2, DestAccessor2> corner)
593{
594 tensorToEdgeCorner(src.first, src.second, src.third,
595 edge.first, edge.second, corner.first, corner.second);
596}
597
598template <class T1, class S1,
599 class T21, class S21,
600 class T22, class S22>
601inline void
602tensorToEdgeCorner(MultiArrayView<2, T1, S1> const & src,
603 MultiArrayView<2, T21, S21> edge,
604 MultiArrayView<2, T22, S22> corner)
605{
606 vigra_precondition(src.shape() == edge.shape(),
607 "tensorToEdgeCorner(): shape mismatch between input and output.");
608 tensorToEdgeCorner(srcImageRange(src),
609 destImage(edge), destImage(corner));
610}
611
612//@}
613
614} // namespace vigra
615
616#endif /* VIGRA_TENSORUTILITIES_HXX */
Class for a single RGB value.
Definition rgbvalue.hxx:128
size_type size() const
Definition tinyvector.hxx:913
void vectorToTensor(...)
Calculate the tensor (outer) product of a 2D vector with itself.
void tensorTrace(...)
Calculate the trace of a 2x2 tensor.
void tensorToEdgeCorner(...)
Decompose a symmetric 2x2 tensor into its edge and corner parts.
void tensorEigenRepresentation(...)
Calculate eigen representation of a symmetric 2x2 tensor.
FixedPoint16< IntBits, OverflowHandling > hypot(FixedPoint16< IntBits, OverflowHandling > v1, FixedPoint16< IntBits, OverflowHandling > v2)
Length of hypotenuse.
Definition fixedpoint.hxx:1640

© 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.12.1