Properly extending GIL to support multispectral imagery

[mpettit2253]

I am working to create extensions to GIL to support multispectral and hyperspectral imagery, and I might be misunderstanding the GIL extension pattern.

My first project is WorldView 2 and 3. These imagers product 1, 3, 4, 8 and 12 channel imagery at a variety of wavelengths and processing levels. It is important for me to distinguish the specific channels and bit-depths (ie, 8-channel visual and near-infrared data is different from the 8-channel short-wave infrared). It looks like I need to define a nearly identical homogeneous_color_base with each possible number of channels (8, 12 and lots, lots more when I get to the imagers from MODIS and others). While tedious, this can be done in my own source file and injected into the boost::gil::detail namespace. However, it also seems like I need to modify the actual definition of planar_pixel_iterator to have a constructor that takes 8 or 12, etc. different channels. When I go in and make changes inside the GIL codebase, I can get it to work nicely, but this seems to really break the encapsulation of the library. I failed to fork GIL itself, and so I'm struggling to go in and make those changes again after I pulled the latest GIL.

Here is my custom color space

namespace worldview
{
    struct pan_color_t {};           // 450  - 800  nanometers

    struct coastal_t {};             //  400 -  450  nanometers
    struct blue_t {};                //  450 -  510  nanometers
    struct green_t {};               //  510 -  580  nanometers
    struct yellow_t {};              //  585 -  625  nanometers
    struct red_t {};                 //  630 -  690  nanometers
    struct red_edge_t {};            //  705 -  745  nanometers
    struct near_infrared_t {};       //  770 -  895  nanometers
    struct near_infrared_2_t {};     //  860 - 1040 nanometers

    struct swir_1_t {};              // 1195 - 1225 nanometers
    struct swir_2_t {};              // 1550 - 1590 nanometers
    struct swir_3_t {};              // 1640 - 1680 nanometers
    struct swir_4_t {};              // 1710 - 1750 nanometers
    struct swir_5_t {};              // 2145 - 2185 nanometers
    struct swir_6_t {};              // 2185 - 2225 nanometers
    struct swir_7_t {};              // 2235 - 2285 nanometers
    struct swir_8_t {};              // 2295 - 2365 nanometers

    struct desert_clouds_t {};       //  405 -  420 nanometers
    struct aerosol_1_t {};           //  459 -  509 nanometers
    struct cavis_green_t {};         //  525 -  585 nanometers
    struct aerosol_2_t {};           //  620 -  670 nanometers
    struct water_1_t {};             //  845 -  885 nanometers
    struct water_2_t {};             //  897 -  927 nanometers
    struct water_3_t {};             //  930 -  965 nanometers
    struct ndvi_swir_t {};           // 1220 - 1252 nanometers
    struct cirrus_t {};              // 1350 - 1410 nanometers
    struct snow_t {};                // 1620 - 1680 nanometers
    struct aerosol_3_t {};           // 2105 - 2245 nanometers
    struct aerosol_3_p_t {};         // 2105 - 2245 nanometers


    using pan_t = boost::mp11::mp_list<pan_color_t>;

    using msi_all_t = boost::mp11::mp_list<
        coastal_t, 
        blue_t, 
        green_t, 
        yellow_t,
        red_t, 
        red_edge_t, 
        near_infrared_t, 
        near_infrared_2_t>;

    using msi_1_t = boost::mp11::mp_list<
        blue_t, 
        green_t, 
        red_t, 
        near_infrared_t>;

    using msi_2_t = boost::mp11::mp_list<
        coastal_t, yellow_t, red_edge_t, near_infrared_2_t>;

    using swir_t = boost::mp11::mp_list<
        swir_1_t, swir_2_t, swir_3_t, swir_4_t, swir_5_t, swir_6_t, swir_7_t, swir_8_t>;

    using cavis_t = boost::mp11::mp_list<
        desert_clouds_t, aerosol_1_t, cavis_green_t, aerosol_2_t,
        water_1_t, water_2_t, water_3_t, ndvi_swir_t,
        cirrus_t, snow_t, aerosol_3_t, aerosol_3_p_t>;

    using pan_layout_t = boost::gil::layout<pan_t>;

    using msi_all_layout_t = boost::gil::layout<msi_all_t>;

    using msi_1_layout_t = boost::gil::layout<msi_1_t>;

	using msi_2_layout_t = boost::gil::layout<msi_2_t>;

	using swir_layout_t = boost::gil::layout<swir_t>;

	using cavis_layout_t = boost::gil::layout<cavis_t>;

    using pan_16_pixel_t = boost::gil::pixel<uint16_t, pan_layout_t>;

    using msi_all_16_pixel_t = boost::gil::pixel<uint16_t, msi_all_layout_t>;

    using msi_1_16_pixel_t = boost::gil::pixel<uint16_t, msi_1_layout_t>;

	using msi_2_16_pixel_t = boost::gil::pixel<uint16_t, msi_2_layout_t>;

	using swir_16_pixel_t = boost::gil::pixel<uint16_t, swir_layout_t>;

	using cavis_16_pixel_t = boost::gil::pixel<uint16_t, cavis_layout_t>;

	using pan_16_planar_image_t = boost::gil::image<pan_16_pixel_t, true, std::allocator<unsigned char>>;

	using msi_all_16_planar_image_t = boost::gil::image<msi_all_16_pixel_t, true, std::allocator<unsigned char>>;

    using msi_1_16_planar_image_t = boost::gil::image<msi_1_16_pixel_t, true, std::allocator<unsigned char>>;

	using msi_2_16_planar_image_t = boost::gil::image<msi_2_16_pixel_t, true, std::allocator<unsigned char>>;

	using swir_16_planar_image_t = boost::gil::image<swir_16_pixel_t, true, std::allocator<unsigned char>>;

	using cavis_16_planar_image_t = boost::gil::image<cavis_16_pixel_t, true, std::allocator<unsigned char>>;

    // Add all of the correct layouts, pixels, images and view types.

    template <typename ChannelPtr>
    inline
        auto planar_msi_all_view(std::size_t width, std::size_t height,
            ChannelPtr c, ChannelPtr b, ChannelPtr g, ChannelPtr y,
            ChannelPtr r, ChannelPtr re, ChannelPtr n, ChannelPtr n2,
            std::ptrdiff_t rowsize_in_bytes)
        -> typename boost::gil::type_from_x_iterator<boost::gil::planar_pixel_iterator<ChannelPtr, msi_all_t> >::view_t
    {
        using pixel_iterator_t = boost::gil::planar_pixel_iterator<ChannelPtr, msi_all_t>;
        using view_t = typename boost::gil::type_from_x_iterator<pixel_iterator_t>::view_t;
        using locator_t = typename view_t::locator;

        locator_t loc(pixel_iterator_t(c, b, g, y, r, re, n, n2), rowsize_in_bytes);
        return view_t(width, height, loc);
    }

    template <typename ChannelPtr>
    inline
        auto planar_msi_1_view(std::size_t width, std::size_t height,
            ChannelPtr b, ChannelPtr g, ChannelPtr r, ChannelPtr n,
            std::ptrdiff_t rowsize_in_bytes)
        -> typename boost::gil::type_from_x_iterator<boost::gil::planar_pixel_iterator<ChannelPtr, msi_1_t> >::view_t
    {
        using pixel_iterator_t = boost::gil::planar_pixel_iterator<ChannelPtr, msi_1_t>;
        using view_t = typename boost::gil::type_from_x_iterator<pixel_iterator_t>::view_t;
        using locator_t = typename view_t::locator;

        locator_t loc(pixel_iterator_t(b, g, r, n), rowsize_in_bytes);
        return view_t(width, height, loc);
    }

    template <typename ChannelPtr>
    inline
        auto planar_msi_2_view(std::size_t width, std::size_t height,
            ChannelPtr c, ChannelPtr y, ChannelPtr r, ChannelPtr n,
            std::ptrdiff_t rowsize_in_bytes)
        -> typename boost::gil::type_from_x_iterator<boost::gil::planar_pixel_iterator<ChannelPtr, msi_2_t> >::view_t
    {
        using pixel_iterator_t = boost::gil::planar_pixel_iterator<ChannelPtr, msi_2_t>;
        using view_t = typename boost::gil::type_from_x_iterator<pixel_iterator_t>::view_t;
        using locator_t = typename view_t::locator;

        locator_t loc(pixel_iterator_t(c, y, r, n), rowsize_in_bytes);
        return view_t(width, height, loc);
    }

    template <typename ChannelPtr>
    inline
        auto planar_swir_view(std::size_t width, std::size_t height,
            ChannelPtr s1, ChannelPtr s2, ChannelPtr s3, ChannelPtr s4,
            ChannelPtr s5, ChannelPtr s6, ChannelPtr s7, ChannelPtr s8,
            std::ptrdiff_t rowsize_in_bytes)
        -> typename boost::gil::type_from_x_iterator<boost::gil::planar_pixel_iterator<ChannelPtr, swir_t> >::view_t
    {
        using pixel_iterator_t = boost::gil::planar_pixel_iterator<ChannelPtr, swir_t>;
        using view_t = typename boost::gil::type_from_x_iterator<pixel_iterator_t>::view_t;
        using locator_t = typename view_t::locator;

        locator_t loc(pixel_iterator_t(s1, s2, s3, s4, s5, s6, s7, s8), rowsize_in_bytes);
        return view_t(width, height, loc);
    }

    template <typename ChannelPtr>
    inline
        auto planar_cavis_view(std::size_t width, std::size_t height,
            ChannelPtr dc, ChannelPtr a1, ChannelPtr g, ChannelPtr a2,
            ChannelPtr w1, ChannelPtr w2, ChannelPtr w3, ChannelPtr ndvi_swir,
			ChannelPtr cirrus, ChannelPtr snow, ChannelPtr a3, ChannelPtr a3p,
            std::ptrdiff_t rowsize_in_bytes)
        -> typename boost::gil::type_from_x_iterator<boost::gil::planar_pixel_iterator<ChannelPtr, cavis_t> >::view_t
    {
        using pixel_iterator_t = boost::gil::planar_pixel_iterator<ChannelPtr, cavis_t>;
        using view_t = typename boost::gil::type_from_x_iterator<pixel_iterator_t>::view_t;
        using locator_t = typename view_t::locator;

        locator_t loc(pixel_iterator_t(dc, a1, g, a2, w1, w2, w3, ndvi_swir, cirrus, snow, a3, a3p), rowsize_in_bytes);
        return view_t(width, height, loc);
    }



} // namespace worldview

Any guidance on this is appreciated.

[mloskot]

Briefly, yet simply, don't fear from 'injecting' your extension into boost::gil namespace. It's not an uncommon practice and I wouldn't die hard to avoid it. For example, see https://github.com/bogado/gil_threaded for extension example hosted outside GIL codebase.

You can also consider contributing your extension to GIL codebase.

[mpettit2253]

Thanks. Claude the AI gave me a suggestion to write my own color base, planar references and iterators, etc. using variadic techniques. I'll follow your suggestion for my simple case of 8 and 12 channel types, and I'll experiment with the AI suggestion below to see if I can support not only the variable number of channels but the variable uncompressed and compressed pixel encoding and tiling schemes that I have to support. I also have to support partially resident image views. The images are in the 1-20 GB range, and sometimes the images themselves are split into multiple files - but it is still the same imaging session with the same camera-to-world transformation parameters.

namespace boost { namespace gil {

// Generic homogeneous color base for N channels
template<typename ChannelValue, std::size_t N>
struct dynamic_homogeneous_color_base {
std::array<ChannelValue, N> channels;

// Generic constructor
template<typename... Channels>
dynamic_homogeneous_color_base(Channels... init_channels) 
    : channels{init_channels...} {
    static_assert(sizeof...(Channels) == N, "Must initialize all channels");
}

};

// Generic planar pixel reference
template<typename ChannelValue, std::size_t N>
class dynamic_planar_pixel_reference {
std::array<ChannelValue*, N> channel_ptrs;

public:
template<typename... Ptrs>
dynamic_planar_pixel_reference(Ptrs... ptrs)
: channel_ptrs{ptrs...} {
static_assert(sizeof...(Ptrs) == N, "Must provide all channel pointers");
}
};

// Generic planar pixel iterator
template<typename ChannelValue, std::size_t N>
class dynamic_planar_pixel_iterator {
std::array<ChannelValue*, N> channel_ptrs;
// Add iterator implementation
};

}} // namespace boost::gil

[sdebionne]

Remember that the code was first released in the pre-c++11 area and while we try to modernize it, there is a lot of room for improvement. Like the one the AI suggested: make the template parameters variadic for a number of class, especially the one that handle planar representation. So a PR in this direction will be very much appreciated.

About the declaration of the typedefs for you plannar types specialization, you can take inspiration from the macros in typedefs.hpp.

I am a bit worried about the large number of types that you will need for your application and the potential code bloat associated and compile time performance. Is there any advantage of having the number of channels encoded in the type? Will you be using the dynamic_extension too?

About "partially resident image views", do you have in mind loading partial images in memory and process them by block? I know that some I/O drivers support that (e.g. tiif/geotiff), so the framework supports reading a sub-part of an image.

[mloskot]

By the way, a long ago I prototyped GIL IO extension for geospatial imagery library GDAL.
Here is my branch:
https://github.com/mloskot/boost-gil-workshop/tree/gil-io2-gdal/

I'm not sure what image file formats are needed by @mpettit2253 but if GDAL can read them, then I'd try piggy-backing on GDAL for I/O operations.

p.s. Although it's been a while since I touched that extension, I still think it would be a useful addition to GIL.