2 * CompactExponentialModel.icpp
4 * Created on: Aug 19, 2019
8 namespace ModelFitting {
10 template<typename ImageType>
11 CompactExponentialModel<ImageType>::CompactExponentialModel(
13 std::shared_ptr<BasicParameter> i0, std::shared_ptr<BasicParameter> k,
14 std::shared_ptr<BasicParameter> x_scale, std::shared_ptr<BasicParameter> y_scale,
15 std::shared_ptr<BasicParameter> rotation, double width, double height,
16 std::shared_ptr<BasicParameter> x, std::shared_ptr<BasicParameter> y,
17 std::tuple<double, double, double, double> transform)
18 : CompactModelBase<ImageType>(x_scale, y_scale, rotation, width, height, x, y, transform),
19 m_sharp_radius_squared(float(sharp_radius * sharp_radius)),
23 template<typename ImageType>
24 double CompactExponentialModel<ImageType>::getValue(double x, double y) const {
25 ExponentialModelEvaluator model_eval;
26 model_eval.transform = getCombinedTransform(1.0);
27 model_eval.i0 = m_i0->getValue();
28 model_eval.k = m_k->getValue();
29 model_eval.max_r_sqr = 1e30;
31 auto area_correction = (1.0 / fabs(m_jacobian[0] * m_jacobian[3] - m_jacobian[1] * m_jacobian[2]));
32 return model_eval.evaluateModel(x, y) * area_correction;
35 template<typename ImageType>
36 ImageType CompactExponentialModel<ImageType>::getRasterizedImage(double pixel_scale, std::size_t size_x, std::size_t size_y) const {
37 //std::cout << "]] " << getX() << " " << getY() << "\n";
38 using Traits = ImageTraits<ImageType>;
40 if (size_x % 2 == 0 || size_y % 2 == 0) {
41 throw Elements::Exception() << "Rasterized image dimensions must be odd numbers "
42 << "but got (" << size_x << ',' << size_y << ")";
45 ImageType image = Traits::factory(size_x, size_y);
47 auto combined_tranform = getCombinedTransform(pixel_scale);
49 ExponentialModelEvaluator model_eval;
50 model_eval.transform = getCombinedTransform(pixel_scale);
51 model_eval.i0 = m_i0->getValue();
52 model_eval.k = m_k->getValue();
53 model_eval.max_r_sqr = getMaxRadiusSqr(size_x, size_y, combined_tranform);
55 float area_correction = (1.0 / fabs(m_jacobian[0] * m_jacobian[3] - m_jacobian[1] * m_jacobian[2])) * pixel_scale * pixel_scale;
57 for (int x=0; x<(int)size_x; ++x) {
58 int dx = x - size_x / 2;
59 for (int y=0; y<(int)size_y; ++y) {
60 int dy = y - size_y / 2;
61 if (dx*dx + dy*dy < m_sharp_radius_squared) {
62 Traits::at(image, x, y) = adaptiveSamplePixel(model_eval, dx, dy, 8, 1.01) * area_correction;
64 Traits::at(image, x, y) = model_eval.evaluateModel(dx, dy) * area_correction;