3D Inverse Discrete Cosine Transformation Implementation in C++

This is a follow-up question for 3D Discrete Cosine Transformation Implementation in C++. In this post, I am trying to follow user17732522's answer to update dct3_detail and dct3 template functions, and propose idct3_detail and idct3 template functions. The formula of 3D Inverse Discrete Cosine Transformation is as follows.

The 3D inverse discrete cosine transformation $$x(n_{1}, n_{2}, n_{3})$$ of size $$N_{1} \times N_{2} \times N_{3}$$ is

$$$$\begin{split} {x(n_{1}, n_{2}, n_{3})} = \sum_{{k_1 = 0}}^{N_1 - 1} \sum_{{k_2 = 0}}^{N_2 - 1} \sum_{{k_3 = 0}}^{N_3 - 1} \epsilon_{k_{1}} \epsilon_{k_{2}} \epsilon_{k_{3}} X(k_{1}, k_{2}, k_{3}) \\ \times \cos({\frac {\pi}{2N_{1}} (2n_{1} + 1)k_{1}}) \\ \times \cos({\frac {\pi}{2N_{2}} (2n_{2} + 1)k_{2}}) \\ \times \cos({\frac {\pi}{2N_{3}} (2n_{3} + 1)k_{3}}) \end{split} \label{3DIDCTMainFormula}$$$$

where

$$$$\begin{split} n_{1} = 0, 1, \dots, N_{1} - 1 \\ n_{2} = 0, 1, \dots, N_{2} - 1 \\ n_{3} = 0, 1, \dots, N_{3} - 1 \\ \epsilon_{k_{i}} = \begin{cases} \frac{1}{\sqrt{2}} & \text{for k_{i} = 0} \\ 1 & \text{otherwise} \end{cases} i = 1, 2, 3 \end{split} \label{3DIDCTMainFormulaDetail}$$$$

The experimental implementation

• non-static idct3_detail template function implementation:

template<std::floating_point ElementT = double, std::floating_point OutputT = ElementT>
Image<OutputT> idct3_detail(const std::vector<Image<ElementT>>& input, const std::size_t plane_index)
{
auto N1 = static_cast<OutputT>(input[0].getWidth());
auto N2 = static_cast<OutputT>(input[0].getHeight());
auto N3 = input.size();
auto output = Image<OutputT>(input[plane_index].getWidth(), input[plane_index].getHeight());
for (std::size_t y = 0; y < output.getHeight(); ++y)
{
for (std::size_t x = 0; x < output.getWidth(); ++x)
{
OutputT sum{};
for (std::size_t inner_z = 0; inner_z < N3; ++inner_z)
{
auto plane = input[inner_z];
for (std::size_t inner_y = 0; inner_y < plane.getHeight(); ++inner_y)
{
for (std::size_t inner_x = 0; inner_x < plane.getWidth(); ++inner_x)
{
auto l1 = (std::numbers::pi_v<OutputT> / (2 * N1) * (2 * x + 1) * static_cast<OutputT>(inner_x));
auto l2 = (std::numbers::pi_v<OutputT> / (2 * N2) * (2 * y + 1) * static_cast<OutputT>(inner_y));
auto l3 = (std::numbers::pi_v<OutputT> / (2 * static_cast<OutputT>(N3)) * (2 * static_cast<OutputT>(plane_index) + 1) * static_cast<OutputT>(inner_z));
OutputT alpha1 = (inner_x == 0) ? (std::numbers::sqrt2_v<OutputT> / 2) : (OutputT{1.0});
OutputT alpha2 = (inner_y == 0) ? (std::numbers::sqrt2_v<OutputT> / 2) : (OutputT{1.0});
OutputT alpha3 = (inner_z == 0) ? (std::numbers::sqrt2_v<OutputT> / 2) : (OutputT{1.0});
sum += alpha1 * alpha2 * alpha3 * static_cast<OutputT>(plane.at(inner_x, inner_y)) *
std::cos(l1) * std::cos(l2) * std::cos(l3);
}
}
}
output.at(x, y) = sum;
}
}
return output;
}

• non-static idct3 template function implementation:

template<std::floating_point ElementT = double, std::floating_point OutputT = ElementT>
std::vector<Image<OutputT>> idct3(const std::vector<Image<ElementT>>& input)
{
std::vector<Image<OutputT>> output;
output.resize(input.size());
for (std::size_t i = 0; i < input.size(); ++i)
{
output[i] = idct3_detail<ElementT, OutputT>(input, i);
}
return output;
}

• The updated version dct3_detail template function implementation:

• static keyword removed

• const std::size_t type plane_index parameter

• Update input parameter with type const std::vector<Image<ElementT>>&.

• Use OutputT template parameter in output type.

• Use std::numbers::pi_v<OutputT> instead of using std::numbers::pi

• About the suggestion of swapping ElementT and OutputT template parameters, I want to make OutputT follows the input ElementT if it isn't specified, i.e. std::floating_point OutputT = ElementT. This makes OutputT needs to be placed after ElementT. As far as I know, template<std::floating_point OutputT = ElementT, std::floating_point ElementT = double> this way doesn't work because ElementT hasn't been declared, it can't be assigned as the default value to OutputT.

It seems that you don't really access input in dct3_detail except as input[plane_index]. The only exception is input[0].getWidth(), but it appears that the result is required to be independent of the index anyway, right? If so, don't pass input and plane_index to dct3_detail. Just pass the element you currently want to use as a reference:

constexpr Image<OutputT> dct3_detail(const Image<ElementT>& input_element)
{
// input[plane_index] -> input_element
// input[0] -> input_element
}

//...

for (auto& input_element : input)
{
output.push_back(dct3_detail(input_element));
}


As far as I know, for calculating each output element in 3D DCT cube, the full input 3D DCT cube would be used. From the mathematics aspect, as the formula presented:

$$$$\begin{split} {X(k_{1}, k_{2}, k_{3})} = {\frac {8}{N_{1} N_{2} N_{3}}} \epsilon_{k_{1}} \epsilon_{k_{2}} \epsilon_{k_{3}} \sum_{{n_1 = 0}}^{N_1 - 1} \sum_{{n_2 = 0}}^{N_2 - 1} \sum_{{n_3 = 0}}^{N_3 - 1} x(n_{1}, n_{2}, n_{3}) \\ \times \cos({\frac {\pi}{2N_{1}} (2n_{1} + 1)k_{1}}) \\ \times \cos({\frac {\pi}{2N_{2}} (2n_{2} + 1)k_{2}}) \\ \times \cos({\frac {\pi}{2N_{3}} (2n_{3} + 1)k_{3}}) \end{split} \label{eq:3DDCTMainFormula}$$$$

The term $${X(k_{1}, k_{2}, k_{3})}$$ represents each element in the transformed output from 3D DCT calculation. From the aspect of programming, there is input[inner_z] used in the inner loops. The whole transformed output can be constructed by planes. dct3_detail template function plays the role of calculating each output plane. By the way, if the dimension of input DCT cube is large, dct3 may take long execution time. Instead of continuous computation, each plane can be calculated (saved) separately with dct3_detail function. If there is any misunderstanding, please let me know.

template<std::floating_point ElementT = double, std::floating_point OutputT = ElementT>
Image<OutputT> dct3_detail(const std::vector<Image<ElementT>>& input, const std::size_t plane_index)
{
auto N1 = static_cast<OutputT>(input[0].getWidth());
auto N2 = static_cast<OutputT>(input[0].getHeight());
auto N3 = input.size();
auto alpha1 = (plane_index == 0) ? (std::numbers::sqrt2_v<OutputT> / 2) : (OutputT{1.0});
auto output = Image<OutputT>(input[plane_index].getWidth(), input[plane_index].getHeight());
for (std::size_t y = 0; y < output.getHeight(); ++y)
{
OutputT alpha2 = (y == 0) ? (std::numbers::sqrt2_v<OutputT> / 2) : (OutputT{1.0});
for (std::size_t x = 0; x < output.getWidth(); ++x)
{
OutputT sum{};
OutputT alpha3 = (x == 0) ? (std::numbers::sqrt2_v<OutputT> / 2) : (OutputT{1.0});
for (std::size_t inner_z = 0; inner_z < N3; ++inner_z)
{
auto plane = input[inner_z];
for (std::size_t inner_y = 0; inner_y < plane.getHeight(); ++inner_y)
{
for (std::size_t inner_x = 0; inner_x < plane.getWidth(); ++inner_x)
{
auto l1 = (std::numbers::pi_v<OutputT> / (2 * N1) * (2 * static_cast<OutputT>(inner_x) + 1) * x);
auto l2 = (std::numbers::pi_v<OutputT> / (2 * N2) * (2 * static_cast<OutputT>(inner_y) + 1) * y);
auto l3 = (std::numbers::pi_v<OutputT> / (2 * static_cast<OutputT>(N3)) * (2 * static_cast<OutputT>(inner_z) + 1) * static_cast<OutputT>(plane_index));
sum += static_cast<OutputT>(plane.at(inner_x, inner_y)) *
std::cos(l1) * std::cos(l2) * std::cos(l3);
}
}
}
output.at(x, y) = 8 * alpha1 * alpha2 * alpha3 * sum / (N1 * N2 * N3);
}
}
return output;
}

• The updated version dct3 template function implementation:

• static keyword removed

• Update input parameter with type const std::vector<Image<ElementT>>&.

• Use OutputT template parameter in output type.

template<std::floating_point ElementT = double, std::floating_point OutputT = ElementT>
std::vector<Image<OutputT>> dct3(const std::vector<Image<ElementT>>& input)
{
std::vector<Image<OutputT>> output;
output.resize(input.size());
for (std::size_t i = 0; i < input.size(); ++i)
{
output[i] = dct3_detail<ElementT, OutputT>(input, i);
}
return output;
}


Full Testing Code

The full tests for dct3 and idct3 template functions:

#include <algorithm>
#include <array>
#include <cassert>
#include <chrono>
#include <cmath>
#include <concepts>
#include <exception>
#include <execution>
#include <fstream>
#include <functional>
#include <iostream>
#include <iterator>
#include <numbers>
#include <numeric>
#include <ranges>
#include <string>
#include <type_traits>
#include <utility>
#include <vector>

using BYTE = unsigned char;

struct RGB
{
BYTE channels[3];
};

using GrayScale = BYTE;

namespace TinyDIP
{
#define is_size_same(x, y) {assert(x.getWidth() == y.getWidth()); assert(x.getHeight() == y.getHeight());}

//  Reference: https://stackoverflow.com/a/58067611/6667035
template <typename T>
concept arithmetic = std::is_arithmetic_v<T>;

//  recursive_depth function implementation
template<typename T>
constexpr std::size_t recursive_depth()
{
return 0;
}

template<std::ranges::input_range Range>
constexpr std::size_t recursive_depth()
{
return recursive_depth<std::ranges::range_value_t<Range>>() + 1;
}

//  recursive_invoke_result_t implementation
template<std::size_t, typename, typename>
struct recursive_invoke_result { };

template<typename T, typename F>
struct recursive_invoke_result<0, F, T> { using type = std::invoke_result_t<F, T>; };

template<std::size_t unwrap_level, typename F, template<typename...> typename Container, typename... Ts>
requires (std::ranges::input_range<Container<Ts...>> &&
requires { typename recursive_invoke_result<unwrap_level - 1, F, std::ranges::range_value_t<Container<Ts...>>>::type; })
struct recursive_invoke_result<unwrap_level, F, Container<Ts...>>
{
using type = Container<typename recursive_invoke_result<unwrap_level - 1, F, std::ranges::range_value_t<Container<Ts...>>>::type>;
};

template<std::size_t unwrap_level, typename F, typename T>
using recursive_invoke_result_t = typename recursive_invoke_result<unwrap_level, F, T>::type;

template<std::size_t, typename, typename, typename...>

template<typename F, class...Ts1, template<class...>class Container1, typename... Ts>
{
using type = Container1<std::invoke_result_t<F,
std::ranges::range_value_t<Container1<Ts1...>>,
std::ranges::range_value_t<Ts>...>>;
};

template<std::size_t unwrap_level, typename F, class...Ts1, template<class...>class Container1, typename... Ts>
requires (  std::ranges::input_range<Container1<Ts1...>> &&
unwrap_level - 1,
F,
std::ranges::range_value_t<Container1<Ts1...>>,
std::ranges::range_value_t<Ts>...>::type; })                //  The rest arguments are ranges
{
using type = Container1<
unwrap_level - 1,
F,
std::ranges::range_value_t<Container1<Ts1...>>,
std::ranges::range_value_t<Ts>...
>::type>;
};

template<std::size_t unwrap_level, typename F, typename T1, typename... Ts>

template<typename OutputIt, typename NAryOperation, typename InputIt, typename... InputIts>
OutputIt transform(OutputIt d_first, NAryOperation op, InputIt first, InputIt last, InputIts... rest) {
while (first != last) {
*d_first++ = op(*first++, (*rest++)...);
}
return d_first;
}

//  recursive_transform for the multiple parameters cases (the version with unwrap_level)
template<std::size_t unwrap_level = 1, class F, class Arg1, class... Args>
constexpr auto recursive_transform(const F& f, const Arg1& arg1, const Args&... args)
{
if constexpr (unwrap_level > 0)
{
static_assert(unwrap_level <= recursive_depth<Arg1>(),
"unwrap level higher than recursion depth of input");
transform(
std::inserter(output, std::ranges::end(output)),
[&f](auto&& element1, auto&&... elements) { return recursive_transform<unwrap_level - 1>(f, element1, elements...); },
std::ranges::cbegin(arg1),
std::ranges::cend(arg1),
std::ranges::cbegin(args)...
);
return output;
}
else
{
return f(arg1, args...);
}
}

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

Image(const std::size_t width, const std::size_t height):
width(width),
height(height),
image_data(width * height) { }

Image(const std::size_t width, const std::size_t height, const ElementT initVal):
width(width),
height(height),
image_data(width * height, initVal) {}

Image(std::vector<ElementT> input, std::size_t newWidth, std::size_t newHeight):
width(newWidth),
height(newHeight)
{
if (input.size() != newWidth * newHeight)
{
throw std::runtime_error("Image data input and the given size are mismatched!");
}
image_data = std::move(input);
}

constexpr ElementT& at(const unsigned int x, const unsigned int y)
{
checkBoundary(x, y);
return image_data[y * width + x];
}

constexpr ElementT const& at(const unsigned int x, const unsigned int y) const
{
checkBoundary(x, y);
return image_data[y * width + x];
}

constexpr std::size_t getWidth() const
{
return width;
}

constexpr std::size_t getHeight() const noexcept
{
return height;
}

constexpr auto getSize() noexcept
{
return std::make_tuple(width, height);
}

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

void print(std::string separator = "\t", std::ostream& os = std::cout) const
{
for (std::size_t y = 0; y < height; ++y)
{
for (std::size_t x = 0; x < width; ++x)
{
//  Ref: https://isocpp.org/wiki/faq/input-output#print-char-or-ptr-as-number
os << +at(x, y) << separator;
}
os << "\n";
}
os << "\n";
return;
}

//  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 < height; ++y)
{
for (std::size_t x = 0; x < width; ++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;
}

friend std::ostream& operator<<(std::ostream& os, const Image<ElementT>& rhs)
{
const std::string separator = "\t";
for (std::size_t y = 0; y < rhs.height; ++y)
{
for (std::size_t x = 0; x < rhs.width; ++x)
{
//  Ref: https://isocpp.org/wiki/faq/input-output#print-char-or-ptr-as-number
os << +rhs.at(x, y) << separator;
}
os << "\n";
}
os << "\n";
return os;
}

Image<ElementT>& operator+=(const Image<ElementT>& rhs)
{
assert(rhs.width == this->width);
assert(rhs.height == this->height);
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)
{
assert(rhs.width == this->width);
assert(rhs.height == this->height);
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)
{
assert(rhs.width == this->width);
assert(rhs.height == this->height);
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)
{
assert(rhs.width == this->width);
assert(rhs.height == this->height);
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;
}

Image<ElementT>& operator=(Image<ElementT> const& input) = default;  //  Copy Assign

Image<ElementT>& operator=(Image<ElementT>&& other) = default;       //  Move Assign

Image(const Image<ElementT> &input) = default;                       //  Copy Constructor

Image(Image<ElementT> &&input) = default;                            //  Move Constructor

private:
std::size_t width;
std::size_t height;
std::vector<ElementT> image_data;

void checkBoundary(const size_t x, const size_t y) const
{
assert(x < width);
assert(y < height);
}

};

template<std::floating_point ElementT = double, std::floating_point OutputT = ElementT>
Image<OutputT> dct3_detail(const std::vector<Image<ElementT>>& input, const std::size_t plane_index)
{
auto N1 = static_cast<OutputT>(input[0].getWidth());
auto N2 = static_cast<OutputT>(input[0].getHeight());
auto N3 = input.size();
auto alpha1 = (plane_index == 0) ? (std::numbers::sqrt2_v<OutputT> / 2) : (OutputT{1.0});
auto output = Image<OutputT>(input[plane_index].getWidth(), input[plane_index].getHeight());
for (std::size_t y = 0; y < output.getHeight(); ++y)
{
OutputT alpha2 = (y == 0) ? (std::numbers::sqrt2_v<OutputT> / 2) : (OutputT{1.0});
for (std::size_t x = 0; x < output.getWidth(); ++x)
{
OutputT sum{};
OutputT alpha3 = (x == 0) ? (std::numbers::sqrt2_v<OutputT> / 2) : (OutputT{1.0});
for (std::size_t inner_z = 0; inner_z < N3; ++inner_z)
{
auto plane = input[inner_z];
for (std::size_t inner_y = 0; inner_y < plane.getHeight(); ++inner_y)
{
for (std::size_t inner_x = 0; inner_x < plane.getWidth(); ++inner_x)
{
auto l1 = (std::numbers::pi_v<OutputT> / (2 * N1) * (2 * static_cast<OutputT>(inner_x) + 1) * x);
auto l2 = (std::numbers::pi_v<OutputT> / (2 * N2) * (2 * static_cast<OutputT>(inner_y) + 1) * y);
auto l3 = (std::numbers::pi_v<OutputT> / (2 * static_cast<OutputT>(N3)) * (2 * static_cast<OutputT>(inner_z) + 1) * static_cast<OutputT>(plane_index));
sum += static_cast<OutputT>(plane.at(inner_x, inner_y)) *
std::cos(l1) * std::cos(l2) * std::cos(l3);
}
}
}
output.at(x, y) = 8 * alpha1 * alpha2 * alpha3 * sum / (N1 * N2 * N3);
}
}
return output;
}

template<std::floating_point ElementT = double, std::floating_point OutputT = ElementT>
std::vector<Image<OutputT>> dct3(const std::vector<Image<ElementT>>& input)
{
std::vector<Image<OutputT>> output;
output.resize(input.size());
for (std::size_t i = 0; i < input.size(); ++i)
{
output[i] = dct3_detail<ElementT, OutputT>(input, i);
}
return output;
}

template<std::floating_point ElementT = double, std::floating_point OutputT = ElementT>
Image<OutputT> idct3_detail(const std::vector<Image<ElementT>>& input, const std::size_t plane_index)
{
auto N1 = static_cast<OutputT>(input[0].getWidth());
auto N2 = static_cast<OutputT>(input[0].getHeight());
auto N3 = input.size();
auto output = Image<OutputT>(input[plane_index].getWidth(), input[plane_index].getHeight());
for (std::size_t y = 0; y < output.getHeight(); ++y)
{
for (std::size_t x = 0; x < output.getWidth(); ++x)
{
OutputT sum{};
for (std::size_t inner_z = 0; inner_z < N3; ++inner_z)
{
auto plane = input[inner_z];
for (std::size_t inner_y = 0; inner_y < plane.getHeight(); ++inner_y)
{
for (std::size_t inner_x = 0; inner_x < plane.getWidth(); ++inner_x)
{
auto l1 = (std::numbers::pi_v<OutputT> / (2 * N1) * (2 * x + 1) * static_cast<OutputT>(inner_x));
auto l2 = (std::numbers::pi_v<OutputT> / (2 * N2) * (2 * y + 1) * static_cast<OutputT>(inner_y));
auto l3 = (std::numbers::pi_v<OutputT> / (2 * static_cast<OutputT>(N3)) * (2 * static_cast<OutputT>(plane_index) + 1) * static_cast<OutputT>(inner_z));
OutputT alpha1 = (inner_x == 0) ? (std::numbers::sqrt2_v<OutputT> / 2) : (OutputT{1.0});
OutputT alpha2 = (inner_y == 0) ? (std::numbers::sqrt2_v<OutputT> / 2) : (OutputT{1.0});
OutputT alpha3 = (inner_z == 0) ? (std::numbers::sqrt2_v<OutputT> / 2) : (OutputT{1.0});
sum += alpha1 * alpha2 * alpha3 * static_cast<OutputT>(plane.at(inner_x, inner_y)) *
std::cos(l1) * std::cos(l2) * std::cos(l3);
}
}
}
output.at(x, y) = sum;
}
}
return output;
}

template<std::floating_point ElementT = double, std::floating_point OutputT = ElementT>
std::vector<Image<OutputT>> idct3(const std::vector<Image<ElementT>>& input)
{
std::vector<Image<OutputT>> output;
output.resize(input.size());
for (std::size_t i = 0; i < input.size(); ++i)
{
output[i] = idct3_detail<ElementT, OutputT>(input, i);
}
return output;
}
}

template<typename ElementT>
void print3(std::vector<TinyDIP::Image<ElementT>> input)
{
for (std::size_t i = 0; i < input.size(); i++)
{
input[i].print();
std::cout << "*******************\n";
}
}

void idct3Test()
{
std::size_t N1 = 10, N2 = 10, N3 = 10;
std::vector<TinyDIP::Image<double>> test_input;
for (std::size_t z = 0; z < N3; z++)
{
test_input.push_back(TinyDIP::Image<double>(N1, N2));
}
for (std::size_t z = 1; z <= N3; z++)
{
for (std::size_t y = 1; y <= N2; y++)
{
for (std::size_t x = 1; x <= N1; x++)
{
test_input[z - 1].at(y - 1, x - 1) = x * 100 + y * 10 + z;
}
}
}
print3(test_input);

auto dct3_output = TinyDIP::dct3(test_input);
print3(dct3_output);

auto idct3_output = TinyDIP::idct3(dct3_output);
print3(TinyDIP::recursive_transform(
[](auto&& input1, auto&& input2)
{
return input1 - input2;
}, idct3_output, test_input));
}

int main()
{
auto start = std::chrono::system_clock::now();
idct3Test();
auto end = std::chrono::system_clock::now();
std::chrono::duration<double> elapsed_seconds = end - start;
std::time_t end_time = std::chrono::system_clock::to_time_t(end);
std::cout << "Computation finished at " << std::ctime(&end_time) << "elapsed time: " << elapsed_seconds.count() << '\n';
return 0;
}


The output of the testing code above:

111 121 131 141 151 161 171 181 191 201
211 221 231 241 251 261 271 281 291 301
311 321 331 341 351 361 371 381 391 401
411 421 431 441 451 461 471 481 491 501
511 521 531 541 551 561 571 581 591 601
611 621 631 641 651 661 671 681 691 701
711 721 731 741 751 761 771 781 791 801
811 821 831 841 851 861 871 881 891 901
911 921 931 941 951 961 971 981 991 1001
1011    1021    1031    1041    1051    1061    1071    1081    1091    1101

*******************
112 122 132 142 152 162 172 182 192 202
212 222 232 242 252 262 272 282 292 302
312 322 332 342 352 362 372 382 392 402
412 422 432 442 452 462 472 482 492 502
512 522 532 542 552 562 572 582 592 602
612 622 632 642 652 662 672 682 692 702
712 722 732 742 752 762 772 782 792 802
812 822 832 842 852 862 872 882 892 902
912 922 932 942 952 962 972 982 992 1002
1012    1022    1032    1042    1052    1062    1072    1082    1092    1102

*******************
113 123 133 143 153 163 173 183 193 203
213 223 233 243 253 263 273 283 293 303
313 323 333 343 353 363 373 383 393 403
413 423 433 443 453 463 473 483 493 503
513 523 533 543 553 563 573 583 593 603
613 623 633 643 653 663 673 683 693 703
713 723 733 743 753 763 773 783 793 803
813 823 833 843 853 863 873 883 893 903
913 923 933 943 953 963 973 983 993 1003
1013    1023    1033    1043    1053    1063    1073    1083    1093    1103

*******************
114 124 134 144 154 164 174 184 194 204
214 224 234 244 254 264 274 284 294 304
314 324 334 344 354 364 374 384 394 404
414 424 434 444 454 464 474 484 494 504
514 524 534 544 554 564 574 584 594 604
614 624 634 644 654 664 674 684 694 704
714 724 734 744 754 764 774 784 794 804
814 824 834 844 854 864 874 884 894 904
914 924 934 944 954 964 974 984 994 1004
1014    1024    1034    1044    1054    1064    1074    1084    1094    1104

*******************
115 125 135 145 155 165 175 185 195 205
215 225 235 245 255 265 275 285 295 305
315 325 335 345 355 365 375 385 395 405
415 425 435 445 455 465 475 485 495 505
515 525 535 545 555 565 575 585 595 605
615 625 635 645 655 665 675 685 695 705
715 725 735 745 755 765 775 785 795 805
815 825 835 845 855 865 875 885 895 905
915 925 935 945 955 965 975 985 995 1005
1015    1025    1035    1045    1055    1065    1075    1085    1095    1105

*******************
116 126 136 146 156 166 176 186 196 206
216 226 236 246 256 266 276 286 296 306
316 326 336 346 356 366 376 386 396 406
416 426 436 446 456 466 476 486 496 506
516 526 536 546 556 566 576 586 596 606
616 626 636 646 656 666 676 686 696 706
716 726 736 746 756 766 776 786 796 806
816 826 836 846 856 866 876 886 896 906
916 926 936 946 956 966 976 986 996 1006
1016    1026    1036    1046    1056    1066    1076    1086    1096    1106

*******************
117 127 137 147 157 167 177 187 197 207
217 227 237 247 257 267 277 287 297 307
317 327 337 347 357 367 377 387 397 407
417 427 437 447 457 467 477 487 497 507
517 527 537 547 557 567 577 587 597 607
617 627 637 647 657 667 677 687 697 707
717 727 737 747 757 767 777 787 797 807
817 827 837 847 857 867 877 887 897 907
917 927 937 947 957 967 977 987 997 1007
1017    1027    1037    1047    1057    1067    1077    1087    1097    1107

*******************
118 128 138 148 158 168 178 188 198 208
218 228 238 248 258 268 278 288 298 308
318 328 338 348 358 368 378 388 398 408
418 428 438 448 458 468 478 488 498 508
518 528 538 548 558 568 578 588 598 608
618 628 638 648 658 668 678 688 698 708
718 728 738 748 758 768 778 788 798 808
818 828 838 848 858 868 878 888 898 908
918 928 938 948 958 968 978 988 998 1008
1018    1028    1038    1048    1058    1068    1078    1088    1098    1108

*******************
119 129 139 149 159 169 179 189 199 209
219 229 239 249 259 269 279 289 299 309
319 329 339 349 359 369 379 389 399 409
419 429 439 449 459 469 479 489 499 509
519 529 539 549 559 569 579 589 599 609
619 629 639 649 659 669 679 689 699 709
719 729 739 749 759 769 779 789 799 809
819 829 839 849 859 869 879 889 899 909
919 929 939 949 959 969 979 989 999 1009
1019    1029    1039    1049    1059    1069    1079    1089    1099    1109

*******************
120 130 140 150 160 170 180 190 200 210
220 230 240 250 260 270 280 290 300 310
320 330 340 350 360 370 380 390 400 410
420 430 440 450 460 470 480 490 500 510
520 530 540 550 560 570 580 590 600 610
620 630 640 650 660 670 680 690 700 710
720 730 740 750 760 770 780 790 800 810
820 830 840 850 860 870 880 890 900 910
920 930 940 950 960 970 980 990 1000    1010
1020    1030    1040    1050    1060    1070    1080    1090    1100    1110

*******************
1726.75 -80.7207    2.5284e-13  -8.64604    -7.77618e-14    -2.82843    1.59162e-13 -1.14371    -1.43473e-13    -0.320717
-807.207    2.57244e-14 -1.24763e-13    -1.00325e-13    4.82332e-14 -7.91025e-14    -9.83958e-14    1.62707e-13 5.91661e-14 6.73658e-14
3.81988e-13 -1.54346e-14    -1.28622e-14    -1.92933e-15    1.41484e-14 -3.85866e-15    -1.28622e-15    3.21555e-15 -5.46643e-15    8.84276e-15
-86.4604    -1.22191e-14    -4.50177e-15    -8.36043e-15    5.14488e-15 1.92933e-15 7.07421e-15 7.71732e-15 2.57244e-15 -2.41166e-15
-5.54792e-14    6.4311e-15  1.92933e-15 -1.28622e-15    -2.57244e-15    3.85866e-15 -4.50177e-15    7.71732e-15 -3.21555e-15    -2.21873e-14
-28.2843    -5.14488e-15    1.28622e-15 -4.50177e-15    3.21555e-15 -5.78799e-15    -7.39576e-15    9.9682e-15  -1.28622e-15    9.64665e-15
1.64164e-13 -7.71732e-15    4.50177e-15 2.57244e-14 -6.4311e-16 3.21555e-16 -1.28622e-15    3.85866e-15 -4.50177e-15    -3.85866e-15
-11.4371    2.18657e-14 -1.57562e-14    -3.21555e-15    1.60777e-15 1.70424e-14 -3.21555e-16    -4.66255e-15    -1.60777e-15    -9.56626e-15
-3.77668e-13    8.68198e-15 1.28622e-15 -5.78799e-15    -1.92933e-15    -1.447e-15  -5.14488e-15    2.73322e-15 2.81361e-15 1.00888e-14
-3.20717    5.30566e-15 8.03887e-16 5.30566e-15 -1.68816e-14    -3.5371e-15 -1.63993e-14    7.39576e-15 1.2822e-14  1.18171e-14

*******************
-8.07207    2.3152e-14  -7.71732e-15    -5.14488e-15    1.92933e-15 -6.4311e-16 -3.85866e-15    1.02898e-14 4.50177e-15 4.82332e-16
-2.4824e-13 5.45697e-15 -5.45697e-15    1.18234e-14 -9.09495e-16    8.18545e-15 0   -4.54747e-15    4.09273e-15 -3.86535e-15
-4.50177e-14    3.63798e-15 -1.81899e-15    6.36646e-15 1.81899e-15 6.36646e-15 3.63798e-15 -4.54747e-15    -4.54747e-15    -2.27374e-16
6.4311e-15  -9.09495e-16    0   1.81899e-15 -9.09495e-16    0   -8.18545e-15    1.81899e-15 -4.54747e-15    4.09273e-15
-1.54346e-14    -8.18545e-15    -1.81899e-15    -1.81899e-15    6.36646e-15 4.54747e-15 7.27596e-15 -1.81899e-15    -2.27374e-15    2.27374e-15
-6.4311e-16 4.54747e-15 -3.63798e-15    -2.72848e-15    8.18545e-15 4.54747e-15 -5.91172e-15    -2.27374e-15    -6.13909e-15    -5.00222e-15
1.02898e-14 6.36646e-15 1.81899e-15 -9.09495e-16    4.54747e-15 5.00222e-15 4.09273e-15 9.09495e-16 -2.27374e-16    8.6402e-15
-1.31838e-14    -9.09495e-16    5.45697e-15 -1.36424e-15    -9.09495e-16    -2.72848e-15    -5.91172e-15    5.68434e-15 6.13909e-15 -8.18545e-15
1.86502e-14 -5.00222e-15    -9.09495e-16    -4.54747e-16    -9.54969e-15    -4.77485e-15    -4.54747e-16    -9.09495e-16    5.00222e-15 3.41061e-15
1.63993e-14 5.68434e-15 -6.82121e-15    -9.09495e-16    -9.77707e-15    -7.04858e-15    2.95586e-15 2.95586e-15 1.02318e-15 1.39266e-15

*******************
3.31056e-13 -5.14488e-15    -1.41484e-14    -9.64665e-15    1.28622e-14 -6.4311e-15 6.4311e-15  3.21555e-15 -4.18021e-15    1.12544e-15
-1.06756e-13    -1.09139e-14    -1.81899e-15    -9.09495e-16    -7.27596e-15    -9.09495e-16    -1.81899e-15    -4.54747e-15    4.54747e-15 -7.50333e-15
-4.1159e-14 5.45697e-15 -3.63798e-15    -1.81899e-15    -5.45697e-15    -4.54747e-15    2.72848e-15 -2.27374e-15    9.09495e-16 -2.72848e-15
-5.14488e-15    -1.45519e-14    -7.27596e-15    -9.09495e-16    -2.72848e-15    9.09495e-15 0   -2.72848e-15    -9.09495e-16    -1.11413e-14
-3.21555e-15    3.63798e-15 1.00044e-14 -3.63798e-15    -4.54747e-15    -8.18545e-15    -4.54747e-15    -3.18323e-15    3.18323e-15 -5.22959e-15
-1.41484e-14    6.36646e-15 0   9.09495e-16 -1.81899e-15    -8.18545e-15    -2.27374e-15    -7.7307e-15 1.59162e-15 -8.6402e-15
1.92933e-14 1.81899e-15 -2.72848e-15    -1.81899e-15    2.27374e-15 -7.27596e-15    -2.72848e-15    2.27374e-15 1.81899e-15 3.18323e-15
1.60777e-15 -9.09495e-16    4.54747e-15 2.27374e-15 0   5.91172e-15 -2.72848e-15    2.04636e-15 3.18323e-15 -5.91172e-15
-1.54346e-14    3.18323e-15 -1.36424e-15    4.54747e-16 1.22782e-14 -4.3201e-15 8.41283e-15 0   -1.13687e-15    6.82121e-16
-6.27032e-15    -4.54747e-16    6.82121e-16 -4.54747e-16    2.27374e-16 -2.38742e-15    -1.13687e-15    2.04636e-15 -6.82121e-16    -6.16751e-15

*******************
-0.864604   -9.64665e-15    6.4311e-16  7.07421e-15 0   -7.07421e-15    1.35053e-14 -3.85866e-15    -6.4311e-16 -2.41166e-15
-2.71392e-13    -8.18545e-15    -3.63798e-15    3.63798e-15 -8.18545e-15    -1.09139e-14    2.72848e-15 -3.63798e-15    -9.09495e-16    7.7307e-15
-5.14488e-15    0   -9.09495e-15    -4.54747e-15    -2.72848e-15    -5.45697e-15    -2.72848e-15    -3.63798e-15    -6.36646e-15    -2.27374e-16
-4.6947e-14 0   5.45697e-15 -4.54747e-15    -9.09495e-16    2.72848e-15 -4.54747e-15    -1.36424e-15    6.36646e-15 -1.18234e-14
1.09329e-14 -8.18545e-15    5.45697e-15 -2.72848e-15    9.09495e-16 8.18545e-15 -1.09139e-14    3.18323e-15 2.27374e-15 5.45697e-15
1.54346e-14 -1.81899e-15    -9.09495e-16    2.72848e-15 4.54747e-15 -2.72848e-15    -6.36646e-15    3.63798e-15 9.09495e-16 3.86535e-15
5.78799e-15 -9.09495e-15    -7.27596e-15    -9.09495e-16    -3.63798e-15    4.54747e-15 1.22782e-14 -1.36424e-15    -2.27374e-16    5.68434e-16
-1.28622e-14    5.45697e-15 5.91172e-15 6.82121e-15 -6.36646e-15    2.27374e-15 -2.27374e-15    1.36424e-15 3.86535e-15 2.6148e-15
-5.46643e-15    2.72848e-15 -4.54747e-15    4.54747e-15 -2.04636e-15    -1.81899e-15    3.86535e-15 2.72848e-15 1.13687e-15 -1.47793e-15
-1.78463e-14    8.6402e-15  -2.04636e-15    -2.27374e-15    1.28466e-14 2.38742e-15 -3.41061e-16    -1.02318e-15    -2.04636e-15    -9.37916e-16

*******************
-4.09273e-14    1.1576e-14  1.09329e-14 3.85866e-15 -3.85866e-15    3.85866e-15 6.4311e-15  1.1576e-14  -5.14488e-15    -9.32509e-15
-5.14488e-15    -6.36646e-15    -1.81899e-15    0   -4.54747e-15    -4.54747e-15    0   -5.45697e-15    -5.00222e-15    2.27374e-15
-2.37951e-14    0   9.09495e-16 -3.63798e-15    -1.00044e-14    -9.09495e-16    3.63798e-15 5.91172e-15 5.00222e-15 -1.25056e-14
6.4311e-16  9.09495e-16 7.27596e-15 -9.09495e-16    4.54747e-15 8.18545e-15 -1.63709e-14    5.00222e-15 -4.54747e-16    1.81899e-15
-1.02898e-14    -4.54747e-15    -5.45697e-15    3.63798e-15 2.72848e-15 -1.81899e-15    1.09139e-14 -8.18545e-15    4.3201e-15  5.22959e-15
-2.25088e-14    -1.81899e-15    -8.18545e-15    0   1.81899e-15 1.81899e-15 1.36424e-15 1.81899e-15 9.09495e-16 1.13687e-15
1.80071e-14 -3.63798e-15    -5.00222e-15    9.09495e-16 5.45697e-15 6.82121e-15 -4.54747e-15    -9.54969e-15    3.86535e-15 -2.27374e-15
-1.38269e-14    9.09495e-16 -3.63798e-15    -2.72848e-15    -1.81899e-15    2.72848e-15 -3.86535e-15    -1.81899e-15    -3.29692e-15    -7.56017e-15
-2.95831e-14    5.45697e-15 -2.27374e-15    2.50111e-15 -1.36424e-15    1.11413e-14 9.54969e-15 5.34328e-15 -1.02318e-15    -1.04023e-14
-1.09329e-14    -1.59162e-15    -9.09495e-16    -1.10276e-14    4.20641e-15 6.36646e-15 3.86535e-15 1.59162e-15 -7.04858e-15    1.29319e-14

*******************
-0.282843   -9.00354e-15    -7.71732e-15    -3.21555e-15    9.64665e-15 -5.78799e-15    2.25088e-15 3.5371e-15  -1.92933e-15    4.50177e-15
-1.08686e-13    2.72848e-15 -5.45697e-15    8.18545e-15 8.18545e-15 9.09495e-16 -4.54747e-16    3.18323e-15 -5.22959e-15    4.54747e-16
4.75901e-14 -4.54747e-15    3.63798e-15 -9.09495e-16    0   -6.36646e-15    8.6402e-15  -1.81899e-15    -2.95586e-15    -8.6402e-15
-3.08693e-14    -1.81899e-15    2.72848e-15 -9.09495e-16    1.00044e-14 0   6.36646e-15 -9.09495e-16    -1.36424e-15    1.29603e-14
3.21555e-15 -1.81899e-15    2.72848e-15 1.81899e-15 3.63798e-15 0   2.27374e-15 -3.63798e-15    0   -6.13909e-15
-3.47279e-14    -1.81899e-15    -5.45697e-15    0   9.09495e-16 5.00222e-15 4.54747e-16 1.13687e-15 -5.22959e-15    2.16005e-15
-2.57244e-15    -2.27374e-15    5.00222e-15 5.45697e-15 2.27374e-15 7.27596e-15 9.09495e-16 -8.6402e-15 -9.09495e-16    -1.20508e-14
-1.41484e-14    3.18323e-15 -5.91172e-15    4.54747e-16 8.18545e-15 -3.18323e-15    -2.27374e-16    -4.3201e-15 3.29692e-15 -1.81899e-15
4.82332e-16 2.72848e-15 -1.13687e-14    5.22959e-15 6.13909e-15 4.09273e-15 5.00222e-15 -1.7053e-15 -2.38742e-15    1.47793e-15
6.4311e-15  7.61702e-15 2.16005e-15 -4.54747e-16    2.27374e-16 -6.36646e-15    2.95586e-15 0   3.9222e-15  9.6918e-15

*******************
1.64164e-13 -1.54346e-14    3.21555e-15 1.41484e-14 -3.21555e-15    5.46643e-15 -5.14488e-15    -1.28622e-15    -7.07421e-15    0
-1.02898e-14    1.00044e-14 5.45697e-15 0   -7.27596e-15    -7.7307e-15 4.09273e-15 -4.54747e-15    2.95586e-15 1.59162e-14
3.85866e-15 1.09139e-14 6.36646e-15 9.09495e-16 4.09273e-15 -3.63798e-15    9.09495e-16 4.54747e-16 9.09495e-16 1.36424e-15
2.12226e-14 -3.63798e-15    1.81899e-15 -3.63798e-15    -1.00044e-14    -2.72848e-15    5.91172e-15 -3.63798e-15    5.68434e-15 1.87583e-14
-4.88764e-14    3.63798e-15 -7.7307e-15 -5.45697e-15    -1.09139e-14    -9.54969e-15    -2.72848e-15    -7.7307e-15 9.32232e-15 -4.09273e-15
-7.71732e-15    8.6402e-15  -9.54969e-15    -3.63798e-15    -8.6402e-15 -9.09495e-16    -9.09495e-16    -2.27374e-16    3.18323e-15 -4.77485e-15
-2.57244e-15    5.91172e-15 -7.7307e-15 2.27374e-15 -5.91172e-15    5.91172e-15 -2.04636e-15    1.81899e-15 4.20641e-15 -5.85487e-15
-2.89399e-15    -1.36424e-15    5.91172e-15 -7.7307e-15 7.50333e-15 -1.13687e-15    4.54747e-16 -1.81899e-15    6.36646e-15 -1.20508e-14
-5.94877e-15    1.13687e-15 -6.82121e-16    9.54969e-15 2.72848e-15 5.91172e-15 -3.86535e-15    7.95808e-16 -2.04636e-15    -4.20641e-15
-2.06599e-14    -6.48015e-15    4.3201e-15  3.97904e-15 2.72848e-15 -3.63798e-15    7.7307e-15  -5.68434e-17    2.70006e-15 6.96332e-15

*******************
-0.114371   2.44382e-14 -3.21555e-16    -8.36043e-15    7.39576e-15 8.68198e-15 -4.50177e-15    -1.60777e-16    -1.92933e-15    -1.21387e-14
2.08046e-13 -1.00044e-14    1.81899e-15 2.27374e-15 1.81899e-15 -4.54747e-15    -4.54747e-15    5.68434e-15 7.04858e-15 1.00044e-14
-4.9841e-14 -2.72848e-15    9.09495e-15 9.09495e-16 -6.36646e-15    4.54747e-16 5.45697e-15 1.81899e-15 4.54747e-16 8.6402e-15
3.85866e-14 5.45697e-15 -4.54747e-16    -4.09273e-15    -3.63798e-15    2.27374e-15 -5.00222e-15    -9.09495e-16    1.13687e-15 -1.37561e-14
6.75265e-15 -2.72848e-15    -3.63798e-15    3.63798e-15 -9.09495e-16    -9.09495e-16    -2.04636e-15    1.13687e-15 2.6148e-15  1.53477e-15
1.1576e-14  -4.09273e-15    -7.7307e-15 5.45697e-15 3.63798e-15 9.09495e-16 -2.95586e-15    2.50111e-15 -2.16005e-15    -1.00044e-14
-1.06113e-14    5.91172e-15 6.82121e-15 -9.54969e-15    3.86535e-15 -1.13687e-15    9.09495e-15 -3.86535e-15    8.6402e-15  5.22959e-15
-8.52121e-15    3.63798e-15 2.04636e-15 4.77485e-15 -3.18323e-15    1.59162e-15 -2.72848e-15    7.04858e-15 -5.22959e-15    1.33014e-14
8.03887e-15 -5.22959e-15    5.00222e-15 -1.13687e-15    -2.72848e-15    1.19371e-14 -4.66116e-15    -7.95808e-16    7.21911e-15 -6.39488e-15
-8.03887e-16    -8.07177e-15    2.27374e-15 5.68434e-16 -3.63798e-15    6.13909e-15 -9.09495e-16    8.69704e-15 -2.84217e-15    8.14282e-15

*******************
-6.29825e-14    8.03887e-15 -8.36043e-15    -3.5371e-15 -6.4311e-15 4.34099e-15 -1.02898e-14    -6.4311e-16 -3.45672e-15    2.23079e-14
6.75265e-14 7.7307e-15  -9.09495e-16    4.54747e-16 5.91172e-15 3.86535e-15 -8.6402e-15 6.82121e-15 -9.09495e-16    6.13909e-15
-4.1159e-14 3.18323e-15 2.27374e-15 -2.27374e-15    9.54969e-15 -1.59162e-15    -4.3201e-15 2.50111e-15 1.36424e-15 5.22959e-15
2.25088e-15 -1.00044e-14    -4.54747e-15    2.72848e-15 5.22959e-15 -7.27596e-15    9.32232e-15 -1.13687e-15    4.54747e-16 -1.47793e-15
-1.67209e-14    -2.72848e-15    1.81899e-15 2.50111e-15 6.82121e-15 5.68434e-15 -1.59162e-15    1.93268e-15 -2.16005e-15    -1.67688e-14
-4.66255e-15    6.36646e-15 4.54747e-15 1.59162e-15 7.95808e-15 -4.54747e-15    4.54747e-16 -2.50111e-15    -3.41061e-16    -9.43601e-15
-7.23499e-15    7.04858e-15 9.32232e-15 -4.54747e-16    6.36646e-15 -1.04592e-14    7.50333e-15 8.07177e-15 -6.59384e-15    1.25056e-15
2.47597e-14 1.38698e-14 -5.45697e-15    2.50111e-15 -6.82121e-16    1.28466e-14 -1.47793e-15    -3.75167e-15    3.12639e-15 -1.84741e-15
1.96952e-14 -5.45697e-15    -3.75167e-15    -1.09139e-14    -3.63798e-15    -3.18323e-15    -1.00613e-14    -4.77485e-15    2.04636e-15 1.6172e-14
9.36529e-15 9.09495e-16 3.52429e-15 4.66116e-15 -1.3415e-14 -1.02318e-14    -3.78009e-15    1.90425e-15 7.24754e-15 -1.09424e-14

*******************
-0.0320717  1.43092e-14 4.34099e-15 6.27032e-15 -9.80743e-15    8.36043e-15 -5.94877e-15    -4.90371e-15    1.25004e-14 1.95345e-14
5.75583e-14 -7.04858e-15    4.54747e-16 4.54747e-15 3.86535e-15 -1.13687e-14    6.59384e-15 2.95586e-15 -9.09495e-16    4.83169e-16
1.94541e-14 -1.81899e-15    -5.22959e-15    0   3.41061e-15 -3.97904e-15    -3.18323e-15    -1.59162e-15    4.54747e-16 -6.62226e-15
-1.12544e-15    -2.27374e-15    -2.95586e-15    5.00222e-15 3.97904e-15 -1.02318e-15    1.7053e-15  -4.20641e-15    -1.36424e-15    3.60956e-15
1.92933e-15 -1.13687e-15    9.09495e-16 -1.93268e-15    9.89075e-15 -1.31877e-14    1.13687e-15 6.82121e-15 -7.50333e-15    -1.02602e-14
-9.00354e-15    4.20641e-15 8.98126e-15 -4.09273e-15    -1.13687e-15    2.04636e-15 1.81899e-15 4.54747e-15 -1.53477e-15    8.78231e-15
-1.80875e-14    -4.20641e-15    5.91172e-15 6.25278e-15 8.18545e-15 6.82121e-16 3.97904e-15 -2.84217e-16    -8.78231e-15    4.23483e-15
-1.17368e-14    4.66116e-15 5.45697e-15 -6.0254e-15 -1.36424e-15    2.50111e-15 6.82121e-16 5.74119e-15 1.93268e-15 4.50484e-15
4.22041e-15 1.93268e-15 -6.82121e-16    1.93268e-15 1.13687e-16 4.54747e-16 1.98952e-16 -5.3717e-15 1.15676e-14 -2.07478e-15
9.305e-15   -1.19371e-15    4.66116e-15 4.74643e-15 3.32534e-15 3.35376e-15 6.79279e-15 3.90799e-15 -1.8332e-15 -3.41061e-16

*******************
1.16529e-12 1.1795e-12  9.66338e-13 1.22213e-12 1.08002e-12 1.19371e-12 1.16529e-12 1.02318e-12 1.25056e-12 1.02318e-12
1.3074e-12  1.27898e-12 1.08002e-12 1.27898e-12 1.10845e-12 1.25056e-12 1.19371e-12 1.02318e-12 1.13687e-12 1.08002e-12
5.11591e-13 5.68434e-13 2.27374e-13 5.11591e-13 3.97904e-13 3.97904e-13 5.11591e-13 4.54747e-13 4.54747e-13 2.84217e-13
1.42109e-12 1.42109e-12 9.66338e-13 1.42109e-12 1.19371e-12 1.3074e-12  1.42109e-12 1.3074e-12  1.3074e-12  1.19371e-12
1.02318e-12 1.13687e-12 4.54747e-13 9.09495e-13 6.82121e-13 6.82121e-13 9.09495e-13 9.09495e-13 9.09495e-13 9.09495e-13
1.13687e-12 1.13687e-12 5.68434e-13 1.02318e-12 7.95808e-13 7.95808e-13 1.02318e-12 1.02318e-12 9.09495e-13 9.09495e-13
1.93268e-12 1.93268e-12 1.25056e-12 1.81899e-12 1.59162e-12 1.59162e-12 1.81899e-12 1.93268e-12 1.7053e-12  1.7053e-12
1.7053e-12  1.59162e-12 7.95808e-13 1.36424e-12 1.13687e-12 9.09495e-13 1.36424e-12 1.47793e-12 1.36424e-12 1.47793e-12
1.81899e-12 1.81899e-12 1.02318e-12 1.59162e-12 1.36424e-12 1.13687e-12 1.59162e-12 1.7053e-12  1.36424e-12 1.59162e-12
1.81899e-12 2.16005e-12 1.13687e-12 1.81899e-12 1.59162e-12 1.36424e-12 1.81899e-12 2.04636e-12 1.81899e-12 2.04636e-12

*******************
1.16529e-12 1.15108e-12 8.81073e-13 1.08002e-12 8.81073e-13 1.10845e-12 1.08002e-12 8.2423e-13  1.08002e-12 9.37916e-13
1.10845e-12 1.08002e-12 9.37916e-13 1.08002e-12 9.37916e-13 1.19371e-12 1.13687e-12 9.66338e-13 1.08002e-12 1.02318e-12
3.41061e-13 3.97904e-13 5.68434e-14 3.41061e-13 2.27374e-13 2.27374e-13 3.41061e-13 2.84217e-13 2.84217e-13 1.13687e-13
1.19371e-12 1.19371e-12 7.38964e-13 1.19371e-12 9.66338e-13 1.08002e-12 1.19371e-12 1.08002e-12 1.08002e-12 9.66338e-13
6.82121e-13 1.02318e-12 3.41061e-13 7.95808e-13 5.68434e-13 5.68434e-13 7.95808e-13 7.95808e-13 7.95808e-13 7.95808e-13
1.02318e-12 1.02318e-12 4.54747e-13 9.09495e-13 6.82121e-13 6.82121e-13 9.09495e-13 9.09495e-13 7.95808e-13 7.95808e-13
1.81899e-12 1.81899e-12 1.13687e-12 1.7053e-12  1.47793e-12 1.47793e-12 1.7053e-12  1.81899e-12 1.59162e-12 1.59162e-12
1.59162e-12 1.47793e-12 6.82121e-13 1.25056e-12 1.02318e-12 7.95808e-13 1.25056e-12 1.36424e-12 1.25056e-12 1.36424e-12
1.7053e-12  1.7053e-12  9.09495e-13 1.47793e-12 1.25056e-12 1.02318e-12 1.47793e-12 1.59162e-12 1.25056e-12 1.47793e-12
1.59162e-12 1.93268e-12 9.09495e-13 1.59162e-12 1.36424e-12 1.13687e-12 1.59162e-12 1.81899e-12 1.59162e-12 1.81899e-12

*******************
1.16529e-12 1.03739e-12 8.2423e-13  1.13687e-12 9.66338e-13 1.02318e-12 1.08002e-12 8.81073e-13 1.02318e-12 9.66338e-13
8.52651e-13 8.81073e-13 6.25278e-13 8.81073e-13 7.38964e-13 9.66338e-13 9.09495e-13 7.38964e-13 8.52651e-13 7.95808e-13
5.68434e-14 1.13687e-13 -2.27374e-13    5.68434e-14 -5.68434e-14    -5.68434e-14    5.68434e-14 0   0   -1.7053e-13
6.25278e-13 6.25278e-13 1.7053e-13  6.25278e-13 3.97904e-13 5.11591e-13 6.25278e-13 5.11591e-13 5.11591e-13 3.97904e-13
2.27374e-13 3.41061e-13 -3.41061e-13    1.13687e-13 -1.13687e-13    -1.13687e-13    1.13687e-13 1.13687e-13 1.13687e-13 1.13687e-13
3.41061e-13 3.41061e-13 -2.27374e-13    2.27374e-13 0   0   2.27374e-13 2.27374e-13 1.13687e-13 1.13687e-13
7.95808e-13 7.95808e-13 1.13687e-13 6.82121e-13 4.54747e-13 4.54747e-13 6.82121e-13 7.95808e-13 5.68434e-13 5.68434e-13
5.68434e-13 4.54747e-13 -3.41061e-13    2.27374e-13 0   -2.27374e-13    2.27374e-13 3.41061e-13 2.27374e-13 3.41061e-13
6.82121e-13 6.82121e-13 -1.13687e-13    4.54747e-13 2.27374e-13 0   4.54747e-13 5.68434e-13 2.27374e-13 4.54747e-13
5.68434e-13 9.09495e-13 -4.54747e-13    2.27374e-13 0   -2.27374e-13    2.27374e-13 4.54747e-13 2.27374e-13 4.54747e-13

*******************
1.20792e-12 1.15108e-12 9.37916e-13 1.27898e-12 1.13687e-12 1.13687e-12 1.16529e-12 1.0516e-12  1.16529e-12 1.10845e-12
1.0516e-12  1.08002e-12 8.2423e-13  1.10845e-12 9.09495e-13 9.66338e-13 9.09495e-13 7.38964e-13 8.52651e-13 7.95808e-13
2.27374e-13 2.84217e-13 -5.68434e-14    2.27374e-13 1.13687e-13 1.13687e-13 2.27374e-13 1.7053e-13  1.7053e-13  0
7.95808e-13 7.95808e-13 3.41061e-13 7.95808e-13 5.68434e-13 6.82121e-13 7.95808e-13 6.82121e-13 6.82121e-13 5.68434e-13
4.54747e-13 5.68434e-13 -1.13687e-13    3.41061e-13 1.13687e-13 1.13687e-13 3.41061e-13 3.41061e-13 3.41061e-13 3.41061e-13
5.68434e-13 5.68434e-13 0   4.54747e-13 2.27374e-13 2.27374e-13 4.54747e-13 4.54747e-13 3.41061e-13 3.41061e-13
1.02318e-12 1.02318e-12 3.41061e-13 9.09495e-13 6.82121e-13 6.82121e-13 9.09495e-13 1.02318e-12 7.95808e-13 7.95808e-13
7.95808e-13 6.82121e-13 -1.13687e-13    4.54747e-13 2.27374e-13 0   4.54747e-13 5.68434e-13 4.54747e-13 5.68434e-13
1.02318e-12 1.02318e-12 2.27374e-13 7.95808e-13 5.68434e-13 3.41061e-13 7.95808e-13 9.09495e-13 5.68434e-13 7.95808e-13
7.95808e-13 1.13687e-12 -4.54747e-13    2.27374e-13 0   -2.27374e-13    2.27374e-13 4.54747e-13 2.27374e-13 4.54747e-13

*******************
7.81597e-13 6.96332e-13 5.40012e-13 7.10543e-13 6.25278e-13 7.38964e-13 7.67386e-13 5.40012e-13 7.38964e-13 6.25278e-13
6.82121e-13 7.10543e-13 4.83169e-13 6.82121e-13 5.40012e-13 7.95808e-13 7.38964e-13 5.68434e-13 6.82121e-13 6.25278e-13
-5.68434e-14    0   -3.41061e-13    -5.68434e-14    -1.7053e-13 -1.7053e-13 -5.68434e-14    -1.13687e-13    -1.13687e-13    -2.84217e-13
6.82121e-13 6.82121e-13 2.27374e-13 6.82121e-13 4.54747e-13 5.68434e-13 6.82121e-13 5.68434e-13 5.68434e-13 4.54747e-13
2.27374e-13 3.41061e-13 -3.41061e-13    1.13687e-13 -1.13687e-13    -1.13687e-13    1.13687e-13 1.13687e-13 1.13687e-13 1.13687e-13
3.41061e-13 3.41061e-13 -2.27374e-13    2.27374e-13 0   0   2.27374e-13 2.27374e-13 1.13687e-13 1.13687e-13
1.02318e-12 1.02318e-12 3.41061e-13 9.09495e-13 6.82121e-13 6.82121e-13 9.09495e-13 1.02318e-12 7.95808e-13 7.95808e-13
6.82121e-13 5.68434e-13 -2.27374e-13    3.41061e-13 1.13687e-13 -1.13687e-13    3.41061e-13 4.54747e-13 3.41061e-13 4.54747e-13
6.82121e-13 6.82121e-13 -1.13687e-13    4.54747e-13 2.27374e-13 0   4.54747e-13 5.68434e-13 2.27374e-13 4.54747e-13
5.68434e-13 9.09495e-13 -2.27374e-13    4.54747e-13 2.27374e-13 0   4.54747e-13 6.82121e-13 4.54747e-13 6.82121e-13

*******************
6.53699e-13 4.54747e-13 3.69482e-13 6.53699e-13 4.54747e-13 5.68434e-13 4.83169e-13 3.69482e-13 5.68434e-13 4.54747e-13
4.54747e-13 4.83169e-13 2.55795e-13 5.11591e-13 2.27374e-13 4.54747e-13 3.97904e-13 2.27374e-13 3.41061e-13 2.84217e-13
-3.97904e-13    -3.41061e-13    -6.82121e-13    -3.97904e-13    -5.11591e-13    -5.11591e-13    -3.97904e-13    -4.54747e-13    -4.54747e-13    -6.25278e-13
2.27374e-13 2.27374e-13 -2.27374e-13    2.27374e-13 0   1.13687e-13 2.27374e-13 1.13687e-13 1.13687e-13 0
-2.27374e-13    -1.13687e-13    -7.95808e-13    -3.41061e-13    -5.68434e-13    -5.68434e-13    -3.41061e-13    -3.41061e-13    -3.41061e-13    -3.41061e-13
-1.13687e-13    -1.13687e-13    -6.82121e-13    -2.27374e-13    -4.54747e-13    -4.54747e-13    -2.27374e-13    -2.27374e-13    -3.41061e-13    -3.41061e-13
3.41061e-13 3.41061e-13 -3.41061e-13    2.27374e-13 0   0   2.27374e-13 3.41061e-13 1.13687e-13 1.13687e-13
2.27374e-13 1.13687e-13 -6.82121e-13    -1.13687e-13    -3.41061e-13    -5.68434e-13    -1.13687e-13    0   -1.13687e-13    0
2.27374e-13 2.27374e-13 -5.68434e-13    0   -2.27374e-13    -4.54747e-13    0   1.13687e-13 -2.27374e-13    0
1.13687e-13 4.54747e-13 -6.82121e-13    0   -2.27374e-13    -4.54747e-13    0   2.27374e-13 0   2.27374e-13

*******************
6.11067e-13 6.39488e-13 2.55795e-13 5.40012e-13 4.54747e-13 4.54747e-13 5.40012e-13 3.69482e-13 4.83169e-13 3.69482e-13
5.40012e-13 5.68434e-13 2.55795e-13 4.83169e-13 3.41061e-13 5.11591e-13 4.54747e-13 2.84217e-13 3.97904e-13 3.41061e-13
-2.27374e-13    -1.7053e-13 -5.11591e-13    -2.27374e-13    -3.41061e-13    -3.41061e-13    -2.27374e-13    -2.84217e-13    -2.84217e-13    -4.54747e-13
3.41061e-13 3.41061e-13 -1.13687e-13    3.41061e-13 1.13687e-13 2.27374e-13 3.41061e-13 2.27374e-13 2.27374e-13 1.13687e-13
-2.27374e-13    -1.13687e-13    -7.95808e-13    -3.41061e-13    -5.68434e-13    -5.68434e-13    -3.41061e-13    -3.41061e-13    -3.41061e-13    -3.41061e-13
-1.13687e-13    -1.13687e-13    -6.82121e-13    -2.27374e-13    -4.54747e-13    -4.54747e-13    -2.27374e-13    -2.27374e-13    -3.41061e-13    -3.41061e-13
5.68434e-13 5.68434e-13 -1.13687e-13    4.54747e-13 2.27374e-13 2.27374e-13 4.54747e-13 5.68434e-13 3.41061e-13 3.41061e-13
3.41061e-13 2.27374e-13 -5.68434e-13    0   -2.27374e-13    -4.54747e-13    0   1.13687e-13 0   1.13687e-13
3.41061e-13 3.41061e-13 -4.54747e-13    1.13687e-13 -1.13687e-13    -3.41061e-13    1.13687e-13 2.27374e-13 -1.13687e-13    1.13687e-13
3.41061e-13 1.13687e-12 0   6.82121e-13 4.54747e-13 2.27374e-13 6.82121e-13 9.09495e-13 6.82121e-13 9.09495e-13

*******************
5.68434e-13 5.96856e-13 3.69482e-13 6.82121e-13 5.11591e-13 5.68434e-13 6.25278e-13 4.26326e-13 5.68434e-13 5.11591e-13
5.11591e-13 5.40012e-13 2.84217e-13 5.40012e-13 3.97904e-13 6.25278e-13 5.68434e-13 3.97904e-13 5.11591e-13 4.54747e-13
-1.7053e-13 -1.13687e-13    -4.54747e-13    -1.7053e-13 -2.84217e-13    -2.84217e-13    -1.7053e-13 -2.27374e-13    -2.27374e-13    -3.97904e-13
5.11591e-13 5.11591e-13 5.68434e-14 5.11591e-13 2.84217e-13 3.97904e-13 5.11591e-13 3.97904e-13 3.97904e-13 2.84217e-13
0   1.13687e-13 -5.68434e-13    -1.13687e-13    -3.41061e-13    -3.41061e-13    -1.13687e-13    -1.13687e-13    -1.13687e-13    -1.13687e-13
1.13687e-13 1.13687e-13 -4.54747e-13    0   -2.27374e-13    -2.27374e-13    0   0   -1.13687e-13    -1.13687e-13
7.95808e-13 7.95808e-13 1.13687e-13 6.82121e-13 4.54747e-13 4.54747e-13 6.82121e-13 7.95808e-13 5.68434e-13 5.68434e-13
5.68434e-13 4.54747e-13 -3.41061e-13    2.27374e-13 0   -2.27374e-13    2.27374e-13 3.41061e-13 2.27374e-13 3.41061e-13
6.82121e-13 6.82121e-13 -1.13687e-13    4.54747e-13 2.27374e-13 0   4.54747e-13 5.68434e-13 2.27374e-13 4.54747e-13
5.68434e-13 1.13687e-12 0   6.82121e-13 4.54747e-13 2.27374e-13 6.82121e-13 9.09495e-13 6.82121e-13 9.09495e-13

*******************
4.83169e-13 4.54747e-13 1.98952e-13 5.11591e-13 2.55795e-13 4.83169e-13 3.97904e-13 2.55795e-13 4.54747e-13 3.12639e-13
3.69482e-13 3.41061e-13 1.42109e-13 3.97904e-13 5.68434e-14 2.84217e-13 2.27374e-13 5.68434e-14 1.7053e-13  1.13687e-13
-5.68434e-13    -5.11591e-13    -8.52651e-13    -5.68434e-13    -6.82121e-13    -6.82121e-13    -5.68434e-13    -6.25278e-13    -6.25278e-13    -7.95808e-13
1.7053e-13  1.7053e-13  -2.84217e-13    1.7053e-13  -5.68434e-14    5.68434e-14 1.7053e-13  5.68434e-14 5.68434e-14 -5.68434e-14
-6.82121e-13    -5.68434e-13    -1.25056e-12    -7.95808e-13    -1.02318e-12    -1.02318e-12    -7.95808e-13    -7.95808e-13    -7.95808e-13    -7.95808e-13
-5.68434e-13    -5.68434e-13    -1.13687e-12    -6.82121e-13    -9.09495e-13    -9.09495e-13    -6.82121e-13    -6.82121e-13    -7.95808e-13    -7.95808e-13
-2.27374e-13    -2.27374e-13    -9.09495e-13    -3.41061e-13    -5.68434e-13    -5.68434e-13    -3.41061e-13    -2.27374e-13    -4.54747e-13    -4.54747e-13
-4.54747e-13    -5.68434e-13    -1.36424e-12    -7.95808e-13    -1.02318e-12    -1.25056e-12    -7.95808e-13    -6.82121e-13    -7.95808e-13    -6.82121e-13
-3.41061e-13    -3.41061e-13    -1.13687e-12    -5.68434e-13    -7.95808e-13    -1.02318e-12    -5.68434e-13    -4.54747e-13    -7.95808e-13    -5.68434e-13
-4.54747e-13    -2.27374e-13    -1.36424e-12    -6.82121e-13    -9.09495e-13    -1.13687e-12    -6.82121e-13    -4.54747e-13    -6.82121e-13    -4.54747e-13

*******************
5.68434e-13 6.25278e-13 2.84217e-13 5.96856e-13 4.54747e-13 5.68434e-13 5.40012e-13 4.54747e-13 6.25278e-13 5.11591e-13
3.97904e-13 4.83169e-13 2.27374e-13 4.83169e-13 4.54747e-13 6.82121e-13 6.25278e-13 4.54747e-13 5.68434e-13 5.11591e-13
-2.84217e-13    -2.27374e-13    -5.68434e-13    -2.84217e-13    -3.97904e-13    -3.97904e-13    -2.84217e-13    -3.41061e-13    -3.41061e-13    -5.11591e-13
1.7053e-13  1.7053e-13  -2.84217e-13    1.7053e-13  -5.68434e-14    5.68434e-14 1.7053e-13  5.68434e-14 5.68434e-14 -5.68434e-14
-5.68434e-13    -4.54747e-13    -1.13687e-12    -6.82121e-13    -9.09495e-13    -9.09495e-13    -6.82121e-13    -6.82121e-13    -6.82121e-13    -6.82121e-13
-4.54747e-13    -4.54747e-13    -1.02318e-12    -5.68434e-13    -7.95808e-13    -7.95808e-13    -5.68434e-13    -5.68434e-13    -6.82121e-13    -6.82121e-13
-1.13687e-13    -1.13687e-13    -7.95808e-13    -2.27374e-13    -4.54747e-13    -4.54747e-13    -2.27374e-13    -1.13687e-13    -3.41061e-13    -3.41061e-13
-3.41061e-13    -4.54747e-13    -1.25056e-12    -6.82121e-13    -9.09495e-13    -1.13687e-12    -6.82121e-13    -5.68434e-13    -6.82121e-13    -5.68434e-13
0   0   -7.95808e-13    -2.27374e-13    -4.54747e-13    -6.82121e-13    -2.27374e-13    -1.13687e-13    -4.54747e-13    -2.27374e-13
-2.27374e-13    0   -1.13687e-12    -4.54747e-13    -6.82121e-13    -9.09495e-13    -4.54747e-13    -2.27374e-13    -4.54747e-13    -2.27374e-13

*******************
Computation finished at Sun Jan  2 08:55:02 2022
elapsed time: 0.0323257


TinyDIP on GitHub

All suggestions are welcome.

The summary information:

• Which question it is a follow-up to?

3D Discrete Cosine Transformation Implementation in C++

• What changes has been made in the code since last question?

• Update dct3_detail and dct3 template functions

• Propose idct3_detail and idct3 template functions implementation.

• Why a new review is being asked for?

If there is any possible improvement, please let me know.

Use std::uint8_t instead of BYTE

I see you are declaring a type alias BYTE. However, I think that is a bad idea. If you really wanted to talk about an opaque byte, then since C++17 there is std::byte. However, you cannot do any arithmetic with std::bytes. Instead, you want to treat bytes as 8 bit integers. There is already a perfect type for this: std::uint8_t. I recommend you use that instead of declaring your own type alias, as it is less confusing for someone who doesn't know what you mean by BYTE without looking up how it is declared.

Avoid macros

There is no reason to make is_size_same() a macro, you can just make it a function:

constexpr void assert_size_same(const auto& x, const auto& y) {
assert(x.getWidth() == y.getWidth());
assert(x.getHeight() == y.getHeight());
}


The drawback of that function is that if it triggers, it will print the line number of the assert() statements inside assert_size_same(), not the line where assert_size_same() itself is called. If you enable coredumps, the backtrace will contain the desired information. Possible workarounds are to make is_size_same() return a boolean:

constexpr bool is_size_same(const auto& x, const auto& y) {
return x.getWidth() == y.getWidth() && x.getHeight() == y.getHeight();
}


Then you can call assert(is_size_same(x, y)). The drawback of that is that it won't tell you if it is the width or height that was not the same, but I don't think that is that important in this case.

throw vs assert()

You are mixing both exceptions and assertions in your code. This can be OK, but you have to be sure when to use exceptions and when to use assertions. assert() helps you as the developer catch bugs early. It will always terminate the program when triggered (except if NDEBUG is defined), so there is no way to recover from assertions. In constrast, exceptions can be handled by the caller by using a try-catch block. Exceptions should be used for things that are exceptional but can happen, assertions should be used for things that are programming bugs and should never happen. It is sometimes hard to decide between the two, but we can look at how things are done in the standard library for some guidance.

For example, STL containers have a member function at(), and this one throws when you are out of bounds. So your checkBoundary() function should also throw a std::out_of_range exception.

Also be consistent. Is it a programmer bug if someone tries to construct an image with a given std::vector<ElementT> that is not exactly newWidth * newHeight long? Is this similar to trying to add two Images with different sizes? I think so, so I would recommend to choose either to throw in both cases or to assert() in both cases.

Avoid passing vectors by value

The constructor of Image that takes a vector as parameter takes it by value. This means a copy will be made by the caller. Consider creating a version instead that takes image as a const reference. This does not avoid a copy, but it does avoid the move (which is very fast, but still not free).

You can also add a constructor that takes image by rvalue reference; this allows the caller to move an existing vector into an Image, avoid all copies.

Printing Images

Your operator<<() duplicates the code of print(). The simple solution is to have operator<<() call print(). As a bonus, it will then also be able to print RGB images using <<.

However, consider whether adding printing functions to an Image class is a good idea. It adds responsibilities to the class. Also, it's not very flexible; what if I want to have the values printed in hex? Or I want an ASCII-art representation of the image? Or I want it in JSON? What if I add a struct RGBA and make an image out of that? Instead, callers can access the image data themselves and print however they like it.

Create a namespace detail

Instead of having a function dct3_detail() in the same namespace as dct3(), consider creating a namespace detail and putting it in there. This is a common technique to indicate that this is a "private" function.

If it is meant to be called by the user, then I wouldn't call it dct3_detail(), but rather dct3_one_plane() or something similar.

Consider using FFTW

You have implemented a naive DCT algorithm that runs in $$\O(N^2)\$$ time, where N = width * height * planes. Similar to the discrete Fourier transform (DFT), you can do a fast cosine transform in only $$\O(N \log N)\$$ time. However, instead of implementing it yourself, I strongly recommend you use the FFTW library to do most of the work for you.

• Thank you for answering. About the question Why not put everything into dct3()?, the answer is mentioned in the question: if the dimension of input DCT cube is large, dct3 may take long execution time. Instead of continuous computation, each plane can be calculated (saved) separately with dct3_detail function (which is called externally). Commented Jan 2, 2022 at 13:20
• About fast cosine transform (FCT) part, is N supposed to be power of 2 for using FCT algorithm? or it can be applied to arbitrary N (positive integer) case? I am checking for this because I assume N in positive integer set can be used for experiments. Commented Jan 2, 2022 at 13:34
• Fast fourier and cosine transforms do not have to be powers of 2, although POT dimensions will be processed the fastest. Commented Jan 2, 2022 at 13:38
• Decide whether to throw or assert() first, as it's kind of part of your API. Whether you want to put this into a function or not is secondary. For throw I think there is no drawback to putting it into a function. Commented Jan 3, 2022 at 19:40
• “Avoid passing vectors by value” — yet the pass by value and move idiom is well established. See for example stackoverflow.com/questions/51705967/… Commented Jan 16, 2022 at 0:35