23 #ifndef LIBSEGMENTATION_SRC_PROFILESAMPLER_INL_HH 24 #define LIBSEGMENTATION_SRC_PROFILESAMPLER_INL_HH 40 template<
class T,
int Dim>
42 : _numSamples(numSamples), _samplingDistUm(samplingDistUm)
45 template<
class T,
int Dim>
47 blitz::Array<T,Dim>
const &data,
48 blitz::TinyVector<double,Dim>
const &elSize,
49 blitz::TinyVector<double,Dim>
const ¢erUm,
50 blitz::TinyVector<double,Dim>
const &positionUm)
const 53 blitz::TinyVector<double,Dim> samplingStepUm(positionUm - centerUm);
57 samplingStepUm *= _samplingDistUm /
58 std::sqrt(
blitz::dot(samplingStepUm, samplingStepUm));
63 blitz::TinyVector<double,Dim> samplingOffsetUm(
64 0.5 * samplingStepUm * _numSamples);
65 blitz::TinyVector<double,Dim> samplingStartUm(
66 positionUm - samplingOffsetUm);
69 blitz::Array<T, 1> edgeSamples(_numSamples);
70 for (
int i = 0; i < _numSamples; ++i)
72 blitz::TinyVector<double,Dim> samplePosUm(
73 samplingStartUm + samplingStepUm * static_cast<double>(i));
74 blitz::TinyVector<double,Dim> samplePosPx(samplePosUm / elSize);
75 edgeSamples(i) = ip.
get(data, samplePosPx);
79 blitz::Array<T,1> dEdge(_numSamples - 1);
80 for (
int i = 0; i < _numSamples - 1; ++i)
81 dEdge(i) = edgeSamples(i + 1) - edgeSamples(i);
84 T normalizer = blitz::max(blitz::abs(dEdge));
87 if(normalizer > 0) dEdge /= normalizer;
92 template<
class T,
int Dim>
94 blitz::Array<T,Dim>
const &data,
95 blitz::TinyVector<double,Dim>
const &elSize,
96 blitz::TinyVector<double,Dim>
const ¢erUm,
97 std::vector< blitz::TinyVector<double,Dim> >
const &candidates,
98 blitz::Array<T,2> &profiles,
101 profiles.resize(candidates.size(), _numSamples - 1);
103 ptrdiff_t progressCnt = 0;
104 ptrdiff_t progressUpdate = 0;
105 if (progress != NULL)
106 progressUpdate = candidates.size() /
110 #pragma omp parallel for 112 for (ptrdiff_t i = 0; i < static_cast<ptrdiff_t>(candidates.size()); ++i)
114 if (progress != NULL)
117 if (progressCnt % progressUpdate == 0)
119 static_cast<int>(progressCnt / progressUpdate +
127 sample(data, elSize, centerUm, candidates[i]);
131 template<
class T,
int Dim>
133 blitz::Array<T,Dim>
const &data,
134 blitz::TinyVector<double,Dim>
const &elSize,
135 blitz::TinyVector<double,Dim>
const ¢erUm,
136 blitz::Array<T,Dim>
const &candidates, T threshold,
137 blitz::Array<T,2> &profiles,
138 std::vector< blitz::TinyVector<double,Dim> > &usedCandidates,
139 blitz::TinyVector<double,Dim>
const &lbUm,
140 blitz::TinyVector<double,Dim>
const &ubUm,
143 usedCandidates.clear();
145 blitz::TinyVector<ptrdiff_t,Dim> lb(lbUm / elSize);
146 blitz::TinyVector<ptrdiff_t,Dim> ub(ubUm / elSize);
147 for (
int d = 0; d < Dim; ++d)
149 if (lb(d) < 0) lb(d) = 0;
150 if (ub(d) < 0) ub(d) = 0;
151 if (lb(d) >= data.extent(d)) lb(d) = data.extent(d) - 1;
152 if (ub(d) >= data.extent(d)) ub(d) = data.extent(d) - 1;
154 blitz::RectDomain<Dim> roi(lb, ub);
156 typename blitz::Array<T, Dim>::const_iterator it;
157 for (it = candidates(roi).begin(); it != candidates(roi).end(); ++it)
159 if (*it < threshold)
continue;
160 blitz::TinyVector<double,Dim> positionUm(it.position() + lb);
161 positionUm *= elSize;
162 usedCandidates.push_back(positionUm);
164 sample(data, elSize, centerUm, usedCandidates, profiles, progress);
168 #endif //LIBSEGMENTATION_SRC_PROFILESAMPLER_INL_HH virtual bool isAborted() const =0
DataT get(blitz::Array< DataT, Dim > const &data, blitz::TinyVector< double, Dim > const &pos) const
Get the Array value at the given subpixel position.
ProfileSampler(int numSamples, T samplingDistUm)
The LinearInterpolator class provides sub-pixel access to blitz++ Arrays using linear interpolation...
virtual bool updateProgress(int progress)=0
int taskProgressMax() const
int taskProgressMin() const
Classes and functions for sub-pixel Array access with different interpolation strategies.
bool all(blitz::TinyMatrix< bool, NRows, NColumns > const &matrix)
all() reduction for boolean blitz::TinyMatrix.
double dot(blitz::Array< DataT, 1 > const &vec1, blitz::Array< DataT, 1 > const &vec2)
Calculation of the inner product (dot product) of the given vectors.
blitz::Array< T, 1 > sample(blitz::Array< T, Dim > const &data, blitz::TinyVector< double, Dim > const &elSize, blitz::TinyVector< double, Dim > const ¢erUm, blitz::TinyVector< double, Dim > const &positionUm) const