# Image Rotation and Transpose Functions Implementation in C++

This is a follow-up question for Gaussian Fisheye Image Generator Implementation in C++ and An Updated Multi-dimensional Image Data Structure with Variadic Template Functions in C++. I am trying to implement image rotation function in this post. The example output:

Image Input Output

The experimental implementation

• rotate, rotate_detail and rotate_degree template functions implementation (in file image_operations.h)

rotate_detail template function performs image rotation operation between 0° to 90°.

namespace TinyDIP
{
//  rotate_detail template function implementation
//  rotate_detail template function performs image rotation between 0° to 90°
template<arithmetic ElementT, std::floating_point FloatingType = double>
constexpr static auto rotate_detail(const Image<ElementT>& input, FloatingType radians)
{
if (input.getDimensionality()!=2)
{
throw std::runtime_error("Unsupported dimension!");
}
//  if 0° rotation case
{
return input;
}
//  if 90° rotation case
if(radians == std::numbers::pi_v<long double> / 2.0)
{
Image<ElementT> output(input.getHeight(), input.getWidth());
for (std::size_t y = 0; y < input.getHeight(); ++y)
{
for (std::size_t x = 0; x < input.getWidth(); ++x)
{
output.at(input.getHeight() - y - 1, x) =
input.at(x, y);
}
}
return output;
}

FloatingType half_width = static_cast<FloatingType>(input.getWidth()) / 2.0;
FloatingType half_height = static_cast<FloatingType>(input.getHeight()) / 2.0;
FloatingType new_width = 2 *
std::hypot(half_width, half_height) *
FloatingType new_height = 2 *
std::hypot(half_width, half_height) *

Image<ElementT> output(input.getWidth(), input.getHeight());
for (std::size_t y = 0; y < input.getHeight(); ++y)
{
for (std::size_t x = 0; x < input.getWidth(); ++x)
{

FloatingType distance_x = x - half_width;
FloatingType distance_y = y - half_height;
FloatingType distance = std::hypot(distance_x, distance_y);
FloatingType angle = std::atan2(distance_y, distance_x) + radians;

FloatingType width_ratio = new_width / static_cast<FloatingType>(input.getWidth());
FloatingType height_ratio = new_height / static_cast<FloatingType>(input.getHeight());
FloatingType distance_weight = (input.getWidth() > input.getHeight())?(1 / height_ratio):(1 / width_ratio);
FloatingType new_distance = distance * distance_weight;
FloatingType new_distance_x = new_distance * std::cos(angle);
FloatingType new_distance_y = new_distance * std::sin(angle);
output.at(
static_cast<std::size_t>(new_distance_x + half_width),
static_cast<std::size_t>(new_distance_y + half_height)) =
input.at(x, y);
}
}
return output;
}

//  rotate template function implementation
template<arithmetic ElementT, std::floating_point FloatingType = double>
constexpr static auto rotate(const Image<ElementT>& input, FloatingType radians)
{
if (input.getDimensionality()!=2)
{
throw std::runtime_error("Unsupported dimension!");
}
auto output = input;
while(radians >= 2 * std::numbers::pi_v<long double>)
{
}
while(radians > std::numbers::pi_v<long double> / 2.0)
{
output = rotate_detail(output, std::numbers::pi_v<long double> / 2.0);
}
return output;
}

//  rotate template function implementation
template<typename ElementT, class FloatingType = double>
requires ((std::same_as<ElementT, RGB>) || (std::same_as<ElementT, HSV>))
constexpr static auto rotate(const Image<ElementT>& input, FloatingType radians)
{
if (input.getDimensionality()!=2)
{
throw std::runtime_error("Unsupported dimension!");
}
return apply_each(input, [&](auto&& planes) { return rotate(planes, radians); });
}

//  rotate template function implementation
template<arithmetic ElementT, std::integral T = int>
constexpr static auto rotate(const Image<ElementT>& input, T radians)
{
if (input.getDimensionality()!=2)
{
throw std::runtime_error("Unsupported dimension!");
}
}

//  rotate_degree template function implementation
template<typename ElementT, class T = double>
constexpr static auto rotate_degree(const Image<ElementT>& input, T degrees)
{
if (input.getDimensionality()!=2)
{
throw std::runtime_error("Unsupported dimension!");
}
return rotate(input, static_cast<double>(degrees) * std::numbers::pi_v<long double> / 180.0);
}
}

• transpose template function implementation (in file image_operations.h)

namespace TinyDIP
{
//  transpose template function implementation
template<typename ElementT>
constexpr static auto transpose(const Image<ElementT>& input)
{
if (input.getDimensionality()!=2)
{
throw std::runtime_error("Unsupported dimension!");
}
Image<ElementT> output(input.getHeight(), input.getWidth());
for (std::size_t y = 0; y < input.getHeight(); ++y)
{
for (std::size_t x = 0; x < input.getWidth(); ++x)
{
output.at(y, x) =
input.at(x, y);
}
}
return output;
}
}

• Image class implementation (in file image.h)

namespace TinyDIP
{
template <typename ElementT>
class Image
{
public:
Image() = default;

template<std::same_as<std::size_t>... Sizes>
Image(Sizes... sizes): size{sizes...}, image_data((1 * ... * sizes))
{}

template<std::same_as<int>... Sizes>
Image(Sizes... sizes)
{
size.reserve(sizeof...(sizes));
(size.push_back(sizes), ...);
image_data.resize(
std::reduce(
std::ranges::cbegin(size),
std::ranges::cend(size),
std::size_t{1},
std::multiplies<>()
)
);
}

template<std::ranges::input_range Range,
std::same_as<std::size_t>... Sizes>
Image(const Range& input, Sizes... sizes):
size{sizes...}, image_data(begin(input), end(input))
{
if (image_data.size() != (1 * ... * sizes)) {
throw std::runtime_error("Image data input and the given size are mismatched!");
}
}

Image(std::vector<ElementT>&& input, std::size_t newWidth, std::size_t newHeight)
{
size.reserve(2);
size.emplace_back(newWidth);
size.emplace_back(newHeight);
if (input.size() != newWidth * newHeight)
{
throw std::runtime_error("Image data input and the given size are mismatched!");
}
image_data = std::move(input);              //  Reference: https://stackoverflow.com/a/51706522/6667035
}

Image(const std::vector<std::vector<ElementT>>& input)
{
size.reserve(2);
size.emplace_back(input[0].size());
size.emplace_back(input.size());
for (auto& rows : input)
{
image_data.insert(image_data.end(), std::ranges::begin(input), std::ranges::end(input));    //  flatten
}
return;
}

//  at template function implementation
template<typename... Args>
constexpr ElementT& at(const Args... indexInput)
{
return const_cast<ElementT&>(static_cast<const Image &>(*this).at(indexInput...));
}

//  at template function implementation
//  Reference: https://codereview.stackexchange.com/a/288736/231235
template<typename... Args>
constexpr ElementT const& at(const Args... indexInput) const
{
checkBoundary(indexInput...);
constexpr std::size_t n = sizeof...(Args);
if(n != size.size())
{
throw std::runtime_error("Dimensionality mismatched!");
}
std::size_t i = 0;
std::size_t stride = 1;
std::size_t position = 0;

auto update_position = [&](auto index) {
position += index * stride;
stride *= size[i++];
};
(update_position(indexInput), ...);

return image_data[position];
}

constexpr std::size_t count() const noexcept
{
return std::reduce(std::ranges::cbegin(size), std::ranges::cend(size), 1, std::multiplies());
}

constexpr std::size_t getDimensionality() const noexcept
{
return size.size();
}

constexpr std::size_t getWidth() const noexcept
{
return size[0];
}

constexpr std::size_t getHeight() const noexcept
{
return size[1];
}

constexpr auto getSize() noexcept
{
return size;
}

std::vector<ElementT> const& getImageData() const noexcept { return image_data; }      //  expose the internal data

void print(std::string separator = "\t", std::ostream& os = std::cout) const
{
if(size.size() == 1)
{
for(std::size_t x = 0; x < size[0]; ++x)
{
//  Ref: https://isocpp.org/wiki/faq/input-output#print-char-or-ptr-as-number
os << +at(x) << separator;
}
os << "\n";
}
else if(size.size() == 2)
{
for (std::size_t y = 0; y < size[1]; ++y)
{
for (std::size_t x = 0; x < size[0]; ++x)
{
//  Ref: https://isocpp.org/wiki/faq/input-output#print-char-or-ptr-as-number
os << +at(x, y) << separator;
}
os << "\n";
}
os << "\n";
}
else if (size.size() == 3)
{
for(std::size_t z = 0; z < size[2]; ++z)
{
for (std::size_t y = 0; y < size[1]; ++y)
{
for (std::size_t x = 0; x < size[0]; ++x)
{
//  Ref: https://isocpp.org/wiki/faq/input-output#print-char-or-ptr-as-number
os << +at(x, y, z) << separator;
}
os << "\n";
}
os << "\n";
}
os << "\n";
}
}

//  Enable this function if ElementT = RGB
void print(std::string separator = "\t", std::ostream& os = std::cout) const
requires(std::same_as<ElementT, RGB>)
{
for (std::size_t y = 0; y < size[1]; ++y)
{
for (std::size_t x = 0; x < size[0]; ++x)
{
os << "( ";
for (std::size_t channel_index = 0; channel_index < 3; ++channel_index)
{
//  Ref: https://isocpp.org/wiki/faq/input-output#print-char-or-ptr-as-number
os << +at(x, y).channels[channel_index] << separator;
}
os << ")" << separator;
}
os << "\n";
}
os << "\n";
return;
}

Image<ElementT>& setAllValue(const ElementT input)
{
std::fill(std::ranges::begin(image_data), std::ranges::end(image_data), input);
return *this;
}

friend std::ostream& operator<<(std::ostream& os, const Image<ElementT>& rhs)
{
const std::string separator = "\t";
rhs.print(separator, os);
return os;
}

Image<ElementT>& operator+=(const Image<ElementT>& rhs)
{
check_size_same(rhs, *this);
std::transform(std::ranges::cbegin(image_data), std::ranges::cend(image_data), std::ranges::cbegin(rhs.image_data),
std::ranges::begin(image_data), std::plus<>{});
return *this;
}

Image<ElementT>& operator-=(const Image<ElementT>& rhs)
{
check_size_same(rhs, *this);
std::transform(std::ranges::cbegin(image_data), std::ranges::cend(image_data), std::ranges::cbegin(rhs.image_data),
std::ranges::begin(image_data), std::minus<>{});
return *this;
}

Image<ElementT>& operator*=(const Image<ElementT>& rhs)
{
check_size_same(rhs, *this);
std::transform(std::ranges::cbegin(image_data), std::ranges::cend(image_data), std::ranges::cbegin(rhs.image_data),
std::ranges::begin(image_data), std::multiplies<>{});
return *this;
}

Image<ElementT>& operator/=(const Image<ElementT>& rhs)
{
check_size_same(rhs, *this);
std::transform(std::ranges::cbegin(image_data), std::ranges::cend(image_data), std::ranges::cbegin(rhs.image_data),
std::ranges::begin(image_data), std::divides<>{});
return *this;
}

friend bool operator==(Image<ElementT> const&, Image<ElementT> const&) = default;

friend bool operator!=(Image<ElementT> const&, Image<ElementT> const&) = default;

friend Image<ElementT> operator+(Image<ElementT> input1, const Image<ElementT>& input2)
{
return input1 += input2;
}

friend Image<ElementT> operator-(Image<ElementT> input1, const Image<ElementT>& input2)
{
return input1 -= input2;
}

friend Image<ElementT> operator*(Image<ElementT> input1, ElementT input2)
{
return multiplies(input1, input2);
}

friend Image<ElementT> operator*(ElementT input1, Image<ElementT> input2)
{
return multiplies(input2, input1);
}

#ifdef USE_BOOST_SERIALIZATION

void Save(std::string filename)
{
const std::string filename_with_extension = filename + ".dat";
//    Reference: https://stackoverflow.com/questions/523872/how-do-you-serialize-an-object-in-c
std::ofstream ofs(filename_with_extension, std::ios::binary);
boost::archive::binary_oarchive ArchiveOut(ofs);
//    write class instance to archive
ArchiveOut << *this;
//    archive and stream closed when destructors are called
ofs.close();
}

#endif
private:
std::vector<std::size_t> size;
std::vector<ElementT> image_data;

template<typename... Args>
void checkBoundary(const Args... indexInput) const
{
constexpr std::size_t n = sizeof...(Args);
if(n != size.size())
{
throw std::runtime_error("Dimensionality mismatched!");
}
std::size_t parameter_pack_index = 0;
auto function = [&](auto index) {
if (index >= size[parameter_pack_index])
throw std::out_of_range("Given index out of range!");
parameter_pack_index = parameter_pack_index + 1;
};

(function(indexInput), ...);
}

#ifdef USE_BOOST_SERIALIZATION
friend class boost::serialization::access;
template<class Archive>
void serialize(Archive& ar, const unsigned int version)
{
ar& size;
ar& image_data;
}
#endif

};

template<typename T, typename ElementT>
concept is_Image = std::is_same_v<T, Image<ElementT>>;
}

#endif


The usage of rotate function:

std::string file_path = "InputImages/1";
bmp1 = rotate(bmp1, 1.0);
TinyDIP::bmp_write("test", bmp1);


The usage transpose function:

template<typename ElementT>
void transposeTest(const size_t size = 5)
{
std::cout << "Test with 2D image:\n";
auto image2d = TinyDIP::Image<ElementT>(size, size);
image2d.setAllValue(1);
image2d.at(0, 1) = 3;
image2d.print();
TinyDIP::transpose(image2d).print();

return;
}


TinyDIP on GitHub

All suggestions are welcome.

The summary information:

# Avoid trigonometric functions where possible

Trigonometric functions like atan2(), sin() and cos() are very slow. You do three of those operations per pixel, which is going to be the major performance issue with your code.

There is a way to avoid doing any per-pixel trigonometry though, and that is by decomposing the rotation into three shear operations, like for example shown in this StackOverflow question.

# Loop over the output pixels

You loop over the input pixels, and remap their coordinates to the output. However, you will then have a risk that you miss some output pixels, as some might be skipped due to floating point rounding. Maybe you get lucky most of the time, but to avoid the problem, loop over the output pixels, and remap their coordinates to the input to find out which color to copy.

Even better, if you use the shearing technique then you can avoid this issue completely.

# Handle negative rotation angles

rotate() seems to only check for radians being too large in the positive direction, but doesn't check for potentially large negative values.

I also recommend you use std::fmod() or something similar to truncate the angle here, instead of using a while-loop.

# Avoid rotating up to four times

If the angle is between $$\1\frac{1}{2} \pi\$$ and $$\2 \pi\$$, you are performing four rotations. You should be able to rotate everything in one go (or if you use shear operations, using not more than 3 shears).

If your algorithm really can't handle rotations larger than 90 degrees, then at least implement a fast algorithm to do 90, 180 and 270 degree rotations, as these are trivial to implement.

• "decomposing the rotation into three shear operations, like for example shown in this StackOverflow question." - Matt Parker did a video on this very recently: youtube.com/watch?v=1LCEiVDHJmc. No affiliation, just a demonstration I thought was neat. Commented Mar 26 at 19:11
• Another advantage of the three sheer technique is that you interpolate in 1D only, which is much simpler than interpolating in 2D when sampling a pixel at an arbitrary location in the image. Commented Mar 27 at 1:23
• @G. Sliepen Thank you for letting me know the shearing technique. It looks good to deal with image rotation with it. :) Commented Mar 27 at 5:48

G. Sliepen already pointed out the three sheer technique, the need to iterate over the output pixels rather than the input, and the expense of doing large rotations 90 degrees at the time. There are a few things to add here:

1. The three sheer technique is not only cheaper, it also makes adding interpolation very easy, as all you do is shift lines of pixels over by a certain distance.

2. You do a lot of computation repeatedly inside the loop, you should move these outside the loop. For example, rotated coordinates are computed by multiplying the input coordinates by a rotation matrix, which is fixed for the whole operation. So inside the loop all you have is four multiplications and two additions to rotate each coordinate. Though of course you should iterate over the output image, and so you use the inverse rotation matrix (or negate the angle) to find the input coordinate corresponding to each output pixel.

3. Indexing using .at(x,y) is quite expensive. To access an arbitrary element it is great, but in a loop that iterates over each pixel, it is much cheaper to add one to the pointer for one pixel to get to the next pixel. You also don’t want to continuously test for indices being inside the bounds, you already know they are.

4. requires ((std::same_as<ElementT, RGB>) || (std::same_as<ElementT, HSV>)) means that you need to rewrite this function (and all the other ones in the library that do a similar thing) if you ever decide to add a color space. I would suggest making this more generic for all color images. Maybe have a base class for both RGB and HSV, that you can use in the requirement?

PS: looking up the definition of apply_each in your repo (it’s not included in the question), I see:

template<class T = RGB, class F, class... Args>
requires (std::same_as<T, RGB>)
constexpr static auto apply_each(Image<T> input, F operation, Args&&... args)
{…}


That requirement makes T not a template parameter. Prefer:

template<class F, class... Args>
constexpr static auto apply_each(Image<RGB> input, F operation, Args&&... args)
{…}


Just because it’s so much simpler. Also:

template<class T = RGB>
requires (std::same_as<T, RGB>)
constexpr static auto getRplane(Image<T> input)
{…}


This is not a template at all, just do:

    constexpr static auto getRplane(Image<RGB> input)
{…}

• > (Image<RGB> input...) makes error: 'input' has incomplete type error. Commented Mar 27 at 4:46
• I think it's called "partial template specialization" and it is not allowed. Commented Mar 27 at 4:53
• @JimmyHu Do you mean in apply_each? You'd have to update each of the templates in the same way Commented Mar 27 at 6:34
• Yes, I already updated each of apply_each templates in the same way but the error is still occurred. Commented Mar 27 at 6:59
• @JimmyHu The templated image processing libraries I know do something like Image<vec3<float>> for color images, without specifying the color space. In your case I was thinking of struct Color {}; struct RGB : Color { uint8 value[3]; }; and so on. Then you can still do Image<RGB>, but in your templated functions you can restrict to T being derived from Color. There might be other solutions. But you do want to avoid listing all the color spaces your library knows for every function that applies generically to color images. Commented Mar 28 at 13:35