running_max_vertical2d.cpp 4.13 KB
Newer Older
Edwin Carlinet's avatar
Edwin Carlinet committed
1
2
3
4
5
6
7
8
#include <mln/morpho/experimental/private/dilation_vertical_block2d.hpp>

#include <memory>

namespace mln::morpho::experimental::details
{


Edwin Carlinet's avatar
Edwin Carlinet committed
9
10
  void vertical_running_max_algo_base_t::running_max_block2d(std::byte* __restrict f, std::byte* __restrict g,
                                                             std::byte* __restrict h, std::ptrdiff_t f_byte_stride,
Edwin Carlinet's avatar
Edwin Carlinet committed
11
                                                             std::ptrdiff_t g_byte_stride, std::ptrdiff_t h_byte_stride,
Edwin Carlinet's avatar
Edwin Carlinet committed
12
                                                             mln::experimental::box2d roi, int k, bool use_extension)
Edwin Carlinet's avatar
Edwin Carlinet committed
13
  {
Edwin Carlinet's avatar
Edwin Carlinet committed
14
15
16
17
    int x0     = roi.x();
    int y0     = roi.y();
    int width  = roi.width();
    int height = roi.height();
Edwin Carlinet's avatar
Edwin Carlinet committed
18

Edwin Carlinet's avatar
Edwin Carlinet committed
19
    assert(width <= this->get_block_width());
Edwin Carlinet's avatar
Edwin Carlinet committed
20
21
22
23
24
25
26
27
28
29
30

    const int alpha = 2 * k + 1;

    {
      int chunk_start = use_extension ? -k : 0;
      int rem         = use_extension ? (height + 2 * k) : height;

      for (; rem > 0; chunk_start += alpha, rem -= alpha)
      {
        int chunk_size = std::min(rem, alpha);

Edwin Carlinet's avatar
Edwin Carlinet committed
31
32
33
34
35
36
37
        // Copy the into the tile$
        if (m_tile_loader)
        {
          mln::experimental::box2d region(x0, y0 + chunk_start, width, chunk_size);
          m_tile_loader->load_tile(f + chunk_start * f_byte_stride, f_byte_stride, region);
        }

Edwin Carlinet's avatar
Edwin Carlinet committed
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
        // Forward pass
        // Compute g[x] = Max f(y), y ∈ [α * ⌊x / α⌋ : x]
        this->partial_sum_block2d(f + chunk_start * f_byte_stride, //
                                  g + chunk_start * g_byte_stride, //
                                  width, chunk_size, f_byte_stride, g_byte_stride);

        // Backward pass
        // Compute h[x] = Max f(y) y ∈ [x : α * (⌊x/α⌋+1))
        this->partial_sum_block2d(f + (chunk_start + chunk_size - 1) * f_byte_stride, //
                                  h + (chunk_start + chunk_size - 1) * h_byte_stride, //
                                  width, chunk_size, -f_byte_stride, -h_byte_stride);
      }
    }


    // Compute local maximum out[x] = Max f(y) y ∈ [x-k:x+k]
    // out[x] = Max   (Max f[x-k:b),  Max f[b:x+k]) with b = α.⌈(x-k)/α⌉ = α.⌊(x+k)/α⌋
    //        = Max( h[x-k], g[x+k] )
    {
Edwin Carlinet's avatar
Edwin Carlinet committed
57
58
      const int kBlockHeight = 16;
      for (int y = 0; y < height; y += kBlockHeight)
Edwin Carlinet's avatar
Edwin Carlinet committed
59
      {
Edwin Carlinet's avatar
Edwin Carlinet committed
60
61
62
63
64
65
66
67
68
69
70
71
72
73
        int hroi = std::min(kBlockHeight, height - y);
        for (int i = 0; i < hroi; ++i)
        {
          this->apply_sup(h + (y + i - k) * h_byte_stride, //
                          g + (y + i + k) * g_byte_stride, //
                          f + (y + i) * f_byte_stride, width);
        }

        // Write the tile
        if (m_tile_writer)
        {
          mln::experimental::box2d region(x0, y0 + y, width, hroi);
          m_tile_writer->write_tile(f + y * f_byte_stride, f_byte_stride, region);
        }
Edwin Carlinet's avatar
Edwin Carlinet committed
74
75
76
77
78
      }
    }
  }


Edwin Carlinet's avatar
Edwin Carlinet committed
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
  void vertical_running_max_algo_base_t::execute(mln::experimental::box2d roi, int k, bool use_extension, bool vertical)
  {
    int            kBlockWidth    = this->get_block_width();
    auto           sz             = this->get_sample_size();
    std::ptrdiff_t kBlockByteSize = sz * kBlockWidth;

    const int x0    = (vertical) ? roi.x() : roi.y();
    const int y0    = (vertical) ? roi.y() : roi.x();
    const int width = (vertical) ? roi.width() : roi.height();
    const int height = (vertical) ? roi.height() : roi.width();


    std::byte* f = (std::byte*)std::malloc(kBlockByteSize * (height + 2 * k));
    std::byte* g = (std::byte*)std::malloc(kBlockByteSize * (height + 2 * k));
    std::byte* h = (std::byte*)std::malloc(kBlockByteSize * (height + 2 * k));

    std::byte* f_shifted = f + k * kBlockByteSize;
    std::byte* g_shifted = g + k * kBlockByteSize;
    std::byte* h_shifted = h + k * kBlockByteSize;

    for (int x = 0; x < width; x += kBlockWidth)
    {
      int w = std::min(kBlockWidth, width - x);


      mln::experimental::box2d region(x0 + x, y0, w, height);
      this->running_max_block2d(f_shifted, g_shifted, h_shifted, //
                                kBlockByteSize, kBlockByteSize, kBlockByteSize, region, k, use_extension);
    }

    std::free(f);
    std::free(g);
    std::free(h);
  }

Edwin Carlinet's avatar
Edwin Carlinet committed
114
115

} // namespace mln::morpho::experimental::details