36 #ifndef VIGRA_SYMMETRY_HXX
37 #define VIGRA_SYMMETRY_HXX
39 #include "utilities.hxx"
40 #include "numerictraits.hxx"
41 #include "stdimage.hxx"
42 #include "convolution.hxx"
143 template <
class SrcIterator,
class SrcAccessor,
144 class DestIterator,
class DestAccessor>
147 DestIterator dul, DestAccessor ad,
150 vigra_precondition(scale > 0.0,
151 "radialSymmetryTransform(): Scale must be > 0");
153 int w = slr.x - sul.x;
154 int h = slr.y - sul.y;
156 if(w <= 0 || h <= 0)
return;
159 NumericTraits<typename SrcAccessor::value_type>::RealPromote TmpType;
161 typedef BasicImage<TmpType> TmpImage;
162 typedef typename TmpImage::Iterator TmpIterator;
166 IImage orientationCounter(w,h);
167 TmpImage magnitudeAccumulator(w,h);
170 destImage(gx), destImage(gy),
173 orientationCounter.init(0);
174 magnitudeAccumulator.init(NumericTraits<TmpType>::zero());
176 TmpIterator gxi = gx.upperLeft();
177 TmpIterator gyi = gy.upperLeft();
179 for(y=0; y<h; ++y, ++gxi.y, ++gyi.y)
181 typename TmpIterator::row_iterator gxr = gxi.rowIterator();
182 typename TmpIterator::row_iterator gyr = gyi.rowIterator();
184 for(
int x = 0; x<w; ++x, ++gxr, ++gyr)
189 if(magnitude < NumericTraits<TmpType>::epsilon()*10.0)
192 int dx = NumericTraits<int>::fromRealPromote(scale *
VIGRA_CSTD::cos(angle));
193 int dy = NumericTraits<int>::fromRealPromote(scale *
VIGRA_CSTD::sin(angle));
198 if(xx >= 0 && xx < w && yy >= 0 && yy < h)
200 orientationCounter(xx, yy) += 1;
201 magnitudeAccumulator(xx, yy) += detail::RequiresExplicitCast<TmpType>::cast(magnitude);
207 if(xx >= 0 && xx < w && yy >= 0 && yy < h)
209 orientationCounter(xx, yy) -= 1;
210 magnitudeAccumulator(xx, yy) -= detail::RequiresExplicitCast<TmpType>::cast(magnitude);
215 int maxOrientation = 0;
216 TmpType maxMagnitude = NumericTraits<TmpType>::zero();
220 for(
int x = 0; x<w; ++x)
224 if(o > maxOrientation)
236 for(
int x = 0; x<w; ++x)
238 double o = (double)orientationCounter(x, y) / maxOrientation;
239 magnitudeAccumulator(x, y) = detail::RequiresExplicitCast<TmpType>::cast(o * o * magnitudeAccumulator(x, y) / maxMagnitude);
243 gaussianSmoothing(srcImageRange(magnitudeAccumulator), destIter(dul, ad), 0.25*scale);
246 template <
class SrcIterator,
class SrcAccessor,
247 class DestIterator,
class DestAccessor>
250 triple<SrcIterator, SrcIterator, SrcAccessor> src,
251 pair<DestIterator, DestAccessor> dest,
255 dest.first, dest.second,
FixedPoint16< 2, OverflowHandling > atan2(FixedPoint16< IntBits, OverflowHandling > y, FixedPoint16< IntBits, OverflowHandling > x)
Arctangent. Accuracy better than 1/3 degree (9 significant bits).
Definition: fixedpoint.hxx:1640
void gaussianGradient(...)
Calculate the gradient vector by means of a 1st derivatives of Gaussian filter.
linalg::TemporaryMatrix< T > sin(MultiArrayView< 2, T, C > const &v)
void radialSymmetryTransform(...)
Find centers of radial symmetry in an image.
void gaussianSmoothing(...)
Perform isotropic Gaussian convolution.
FFTWComplex< R >::NormType abs(const FFTWComplex< R > &a)
absolute value (= magnitude)
Definition: fftw3.hxx:1002
linalg::TemporaryMatrix< T > cos(MultiArrayView< 2, T, C > const &v)
SquareRootTraits< FixedPoint< IntBits, FracBits > >::SquareRootResult sqrt(FixedPoint< IntBits, FracBits > v)
square root.
Definition: fixedpoint.hxx:616
BasicImage< Int32 > IImage
Definition: stdimage.hxx:116