Skip to content
Snippets Groups Projects
Commit fb47afdd authored by Wuttke, Joachim's avatar Wuttke, Joachim
Browse files

ConvolutionDetectorResolution ditto

parent 5bc83481
No related branches found
No related tags found
1 merge request!2019in convolution context, throw -> assert
Pipeline #115136 passed
......@@ -71,7 +71,7 @@ void ConvolutionDetectorResolution::setResolutionFunction(const IResolutionFunct
void ConvolutionDetectorResolution::apply1dConvolution(Datafield* intensity_map) const
{
ASSERT(m_res_function_1d == nullptr);
ASSERT(m_res_function_1d);
ASSERT(intensity_map->rank() == 1);
const Scale& axis = intensity_map->axis(0);
......@@ -101,7 +101,7 @@ void ConvolutionDetectorResolution::apply1dConvolution(Datafield* intensity_map)
void ConvolutionDetectorResolution::apply2dConvolution(Datafield* intensity_map) const
{
ASSERT(m_res_function_2d);
ASSERT(intensity_map->rank() != 2);
ASSERT(intensity_map->rank() == 2);
const Scale& axis_1 = intensity_map->axis(0);
const Scale& axis_2 = intensity_map->axis(1);
size_t axis_size_1 = axis_1.size();
......
......@@ -13,9 +13,7 @@
// ************************************************************************************************
#include "Device/Resolution/Convolve.h"
#include <iostream>
#include <sstream>
#include <stdexcept> // need overlooked by g++ 5.4
#include "Base/Util/Assert.h"
Convolve::Convolve()
: m_mode(FFTW_UNDEFINED)
......@@ -99,12 +97,11 @@ void Convolve::fftconvolve(const double2d_t& source, const double2d_t& kernel, d
result[i].resize(ws.w_dst, 0);
for (int j = 0; j < ws.w_dst; j++) {
// result[i][j]=ws.dst[i*ws.w_dst+j];
if (m_mode == FFTW_CIRCULAR_SAME_SHIFTED) {
if (m_mode == FFTW_CIRCULAR_SAME_SHIFTED)
result[i][j] = ws.dst_fft[((i + int(ws.h_kernel / 2.0)) % ws.h_fftw) * ws.w_fftw
+ (j + int(ws.w_kernel / 2.0)) % ws.w_fftw];
} else {
else
result[i][j] = ws.dst_fft[(i + ws.h_offset) * ws.w_fftw + j + ws.w_offset];
}
}
}
}
......@@ -121,8 +118,7 @@ void Convolve::fftconvolve(const double1d_t& source, const double1d_t& kernel, d
double2d_t result2d;
fftconvolve(source2d, kernel2d, result2d);
if (result2d.size() != 1)
throw std::runtime_error("Convolve::fftconvolve -> Panic in 1d");
ASSERT(result2d.size() == 1);
result = result2d[0];
}
......@@ -131,12 +127,10 @@ void Convolve::fftconvolve(const double1d_t& source, const double1d_t& kernel, d
/* ************************************************************************* */
void Convolve::init(int h_src, int w_src, int h_kernel, int w_kernel)
{
if (!h_src || !w_src || !h_kernel || !w_kernel) {
std::ostringstream os;
os << "Convolve::init -> Panic! Wrong dimensions " << h_src << " " << w_src << " "
<< h_kernel << " " << w_kernel << std::endl;
throw std::runtime_error(os.str());
}
ASSERT(h_src);
ASSERT(w_src);
ASSERT(h_kernel);
ASSERT(w_kernel);
ws.clear();
ws.h_src = h_src;
......@@ -178,19 +172,12 @@ void Convolve::init(int h_src, int w_src, int h_kernel, int w_kernel)
break;
case FFTW_LINEAR_VALID:
// Valid Linear convolution
if (ws.h_kernel > ws.h_src || ws.w_kernel > ws.w_src) {
ws.h_fftw = 0;
ws.w_fftw = 0;
ws.h_dst = 0;
ws.w_dst = 0;
std::cout << "The 'valid' convolution results in an empty matrix" << std::endl;
throw std::runtime_error("The 'valid' convolution results in an empty matrix");
} else {
ws.h_fftw = find_closest_factor(h_src);
ws.w_fftw = find_closest_factor(w_src);
ws.h_dst = h_src - h_kernel + 1;
ws.w_dst = w_src - w_kernel + 1;
}
ASSERT(ws.h_kernel <= ws.h_src);
ASSERT(ws.w_kernel <= ws.w_src);
ws.h_fftw = find_closest_factor(h_src);
ws.w_fftw = find_closest_factor(w_src);
ws.h_dst = h_src - h_kernel + 1;
ws.w_dst = w_src - w_kernel + 1;
// Here we just take [h_dst x w_dst] elements starting at [h_kernel-1;w_kernel-1]
ws.h_offset = ws.h_kernel - 1;
ws.w_offset = ws.w_kernel - 1;
......@@ -210,11 +197,7 @@ void Convolve::init(int h_src, int w_src, int h_kernel, int w_kernel)
ws.w_offset = 0;
break;
default:
std::cout
<< "Unrecognized convolution mode, possible modes are "
<< "FFTW_LINEAR_FULL, FFTW_LINEAR_SAME, FFTW_LINEAR_SAME_UNPADDED, FFTW_LINEAR_VALID, "
<< "FFTW_CIRCULAR_SAME, FFTW_CIRCULAR_SHIFTED " << std::endl;
break;
ASSERT_NEVER;
}
ws.in_src = new double[ws.h_fftw * ws.w_fftw];
......@@ -228,19 +211,15 @@ void Convolve::init(int h_src, int w_src, int h_kernel, int w_kernel)
// Initialization of the plans
ws.p_forw_src = fftw_plan_dft_r2c_2d(ws.h_fftw, ws.w_fftw, ws.in_src, (fftw_complex*)ws.out_src,
FFTW_ESTIMATE);
if (ws.p_forw_src == nullptr)
throw std::runtime_error("Convolve::init -> Error! Cannot initialise p_forw_src plan.");
ASSERT(ws.p_forw_src);
ws.p_forw_kernel = fftw_plan_dft_r2c_2d(ws.h_fftw, ws.w_fftw, ws.in_kernel,
(fftw_complex*)ws.out_kernel, FFTW_ESTIMATE);
if (ws.p_forw_kernel == nullptr)
throw std::runtime_error("Convolve::init -> Error! Cannot initialise p_forw_kernel plan.");
ASSERT(ws.p_forw_kernel);
// The backward FFT takes ws.out_kernel as input
ws.p_back = fftw_plan_dft_c2r_2d(ws.h_fftw, ws.w_fftw, (fftw_complex*)ws.out_kernel, ws.dst_fft,
FFTW_ESTIMATE);
if (ws.p_back == nullptr)
throw std::runtime_error("Convolve::init -> Error! Cannot initialise p_back plan.");
ASSERT(ws.p_back);
}
/* ************************************************************************* */
......@@ -249,8 +228,8 @@ void Convolve::init(int h_src, int w_src, int h_kernel, int w_kernel)
void Convolve::fftw_circular_convolution(const double2d_t& src, const double2d_t& kernel)
{
if (ws.h_fftw <= 0 || ws.w_fftw <= 0)
throw std::runtime_error("Convolve::fftw_convolve -> Panic! Initialization is missed.");
ASSERT(ws.h_fftw > 0);
ASSERT(ws.w_fftw > 0);
double *ptr, *ptr_end, *ptr2;
......@@ -319,9 +298,8 @@ bool Convolve::is_optimal(int n)
if (n == 1)
return false;
size_t ntest = n;
for (size_t i = 0; i < m_implemented_factors.size(); i++) {
for (size_t i = 0; i < m_implemented_factors.size(); i++)
while ((ntest % m_implemented_factors[i]) == 0)
ntest = ntest / m_implemented_factors[i];
}
return ntest == 1;
}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment